定積分的Monte-Carlo計(jì)算方法的實(shí)現(xiàn)省公開課金獎(jiǎng)全國賽課一等獎(jiǎng)微課獲獎(jiǎng)?wù)n件

《計(jì)算機(jī)高級(jí)語言》

認(rèn)知實(shí)習(xí)匯報(bào)

課題名稱：定積分MonteCarlo計(jì)算方法實(shí)現(xiàn)

指導(dǎo)老師：王玉蘭小組組員：07020302曾穎超07020301李海全0701向慧

實(shí)習(xí)時(shí)間：.12.27——.01.08摘要：此次實(shí)習(xí)設(shè)計(jì)并實(shí)現(xiàn)定積分MonteCarlo計(jì)算方法，采取VC++6.0開發(fā)完成。1/26

一、設(shè)計(jì)任務(wù)與要求1、設(shè)計(jì)內(nèi)容2、設(shè)計(jì)要求二、MonteCarlo算法引導(dǎo)1、MonteCarlo算法概念2、MonteCarlo算法例子三、MonteCarlo思想計(jì)算定積分實(shí)現(xiàn)1、定積分定義2、MonteCarlo算法積分原理四、詳細(xì)設(shè)計(jì)及實(shí)現(xiàn)五、總結(jié)六、參考資料附錄：源程序2/26一、設(shè)計(jì)任務(wù)與要求1、設(shè)計(jì)內(nèi)容針對(duì)某詳細(xì)函數(shù)f(x)形式，利用MonteCarlo算法實(shí)現(xiàn)積分。3/262、設(shè)計(jì)要求1)學(xué)習(xí)MonteCarlo計(jì)算數(shù)值積分方法；2)就某個(gè)詳細(xì)被積函數(shù)f(x)定積分，設(shè)計(jì)算法步驟;3)編寫程序?qū)崿F(xiàn)算法；4）對(duì)程序進(jìn)行調(diào)試和測(cè)試，使用有解析解定積分方便計(jì)算對(duì)比；4/26二、MonteCarlo算法引導(dǎo)

1、MonteCarlo算法概念

非形式化地說，MonteCarlo算法泛指一類算法。在這些算法中，要求解問題是某隨機(jī)事件概率或某隨機(jī)變量期望。這時(shí)，經(jīng)過“試驗(yàn)”方法，用頻率代替概率或得到隨機(jī)變量一些數(shù)字特征，以此作為問題解。5/262、MonteCarlo算法例子如圖(1)所表示，在一個(gè)1平方米正方形木板上，隨意畫一個(gè)圈，求返個(gè)圈面積。圖(1)假設(shè)我手里有一支飛鏢，我將飛鏢擲向木板。木板，我們假定每一次都能擲在木板上，不會(huì)偏出木板，但每一次擲在木板什么地方，是完全隨機(jī)。即，每一次擲飛鏢，飛鏢扎進(jìn)木板任何一點(diǎn)概率相等。返樣，我們投擲屢次，比如100次，然后我們統(tǒng)計(jì)返100次中，扎入不規(guī)則圖形內(nèi)部次數(shù)，假設(shè)為k，那么，我們就能夠用k/100*1近似預(yù)計(jì)不規(guī)則圖形面積，比如100次有32次擲入圖形內(nèi)，我們就能夠預(yù)計(jì)圖形面積為0.32平方米從上述能夠看出，Monte-Carlo算法區(qū)分不確定性算法，它解不一定是準(zhǔn)確或正確，其準(zhǔn)確或正確性依賴于概率和統(tǒng)計(jì)，但在一些問題上，當(dāng)重復(fù)試驗(yàn)次數(shù)越夠大時(shí)，可從很大約率上（返個(gè)概率是能夠在數(shù)學(xué)上證實(shí)，但依賴于詳細(xì)問題）確保解準(zhǔn)確或正確性，所以，我們能夠依據(jù)詳細(xì)概率分析，設(shè)定試驗(yàn)次數(shù)，從而將誤差或錯(cuò)誤率降到一個(gè)可容忍程度。6/26三、Monte-Carlo思想計(jì)算定積分實(shí)現(xiàn)1、定積分定義2、MonteCarlo積分法原理

7/261.定積分定義如圖2，定積分就是求函數(shù)在區(qū)間[a,b]中圖線下包圍面積，即y=0,x=a,x=b,y=f(x)所包圍面積。普通地，設(shè)函數(shù)在區(qū)間上連續(xù)，用分點(diǎn)將區(qū)間等分成個(gè)小區(qū)間，每個(gè)小區(qū)間長度為（），在每個(gè)小區(qū)間上取一點(diǎn)，作和式：8/26假如無限靠近于（亦即）時(shí)，上述和式無限趨近于常數(shù)，那么稱該常數(shù)為函數(shù)在區(qū)間上定積分。記為：其中成為被積函數(shù)，叫做積分變量，為積分區(qū)間，積分上限，積分下限。9/262、MonteCarlo積分法原理

求定積分直觀意義就是求面積，所以，用MonteCarlo求積分原理就是經(jīng)過模擬統(tǒng)計(jì)方法求解面積。即經(jīng)過向特定區(qū)域隨機(jī)產(chǎn)生大量點(diǎn)，然后統(tǒng)計(jì)點(diǎn)落在函數(shù)區(qū)域內(nèi)頻率，以此頻率預(yù)計(jì)面積，從而得到積分值。10/26四、詳細(xì)設(shè)計(jì)及實(shí)現(xiàn)依據(jù)蒙特卡洛算法概念，將落在數(shù)學(xué)函數(shù)f(x)里面隨機(jī)點(diǎn)搜集起來并計(jì)算出數(shù)量，這么能夠計(jì)算出蒙特卡洛因子。然后就可計(jì)算出函數(shù)f(x)積分。我們把積分區(qū)間較長函數(shù)進(jìn)行分段求解，然后再求和。希望以此來提升準(zhǔn)確性和縮短時(shí)間。11/26五測(cè)試被積函數(shù)：f(x)=sin(x),x∈[0,20]積分準(zhǔn)確解為：0.591917把區(qū)間分為N段，隨機(jī)點(diǎn)個(gè)數(shù)為n個(gè)以下是測(cè)試數(shù)據(jù)12/26nN積分誤差（取正）時(shí)間/毫秒10610.566960.02495723810620.58220.00971746510640.592270.00035393110680.5909450.0009721,813106160.595310.0033933,566106320.5951810.0032647,17913/2610710.5986440.0067272,30810720.5960520.0041354,60710740.5977180.0058019,01110780.59310.00118317,847107160.5903590.00155835,35214/26六總結(jié)用蒙特卡洛方法計(jì)算定積分含有普遍意義。依據(jù)以上試驗(yàn)我們了解到，對(duì)于求解一些非常規(guī)問題我們有了新思緒，即一個(gè)不具隨機(jī)性事件能夠經(jīng)過一定方法用隨機(jī)事件來模擬或迫近。15/26參考文件WilliamH.Press等著，傅祖蕓等譯，C數(shù)值算法（第二版）電子工業(yè)出版社，.01同濟(jì)大學(xué)數(shù)學(xué)系編，高等數(shù)學(xué)（第六版）高等教育出版社，.06孫淑霞李思明劉煥君編著，C/C++程序設(shè)計(jì)（第三版）電子工業(yè)出版社，.0916/26附錄：源程序#defineIM12147483563#defineIM22147483399#defineAM(1.0/IM1)#defineIMM1(IM1-1)#defineIA140014#defineIA240692#defineIQ153668#defineIQ252774#defineIR112211#defineIR23791#defineNTAB32#defineNDIV(1+IMM1/NTAB)17/26#defineEPS1.2e-7#defineRNMX(1.0-EPS)/*以上是為ran2函數(shù)做準(zhǔn)備*/#include<stdlib.h>#include<stdio.h>#include<math.h>#include<time.h>18/26doubleran2(long*idum)/*產(chǎn)生隨機(jī)數(shù)*/{intj;longk;staticlongidum2=123456789;staticlongiy=0;staticlongiv[NTAB];doubletemp;

if(*idum<=0)19/26

{if(-(*idum)<1)*idum=1;else*idum=-(*idum);for(j=NTAB+7;j>=0;j--){k=(*idum)/IQ1;*idum=IA1*(*idum-k*IQ1)-k*IR1;if(*idum<0)*idum+=IM1;if(j<NTAB)iv[j]=*idum;}iy=iv[0];}20/26k=(*idum)/IQ1;*idum=IA1*(*idum-k*IQ1)-k*IR1;if(*idum<0)*idum+=IM1;k=idum2/IQ2;idum2=IA2*(idum2-k*IQ2)-k*IR2;if(idum2<0)idum2+=IM2;j=iy/NDIV;iy=iv[j]-idum2;iv[j]=*idum;if(iy<1)iy+=IMM1;if((temp=AM*iy)>RNMX)returnRNMX;elsereturntemp;}21/26doublefun1(doublex)/*被積函數(shù)*/{return(sin(x));}doubleMonte(doublen,doublea,doubleb,double(*p)(),doublefMAX)/*計(jì)算蒙特卡洛因子*/{doublex,y;doublesx=0,sy=0,vx=0,vy=0;doubles,t;doublei;longko1,ko2;22/26ko1=ko2=1;for(i=0;i<n;i++){t=ran2(&ko1);s=ran2(&ko2);x=a+(b-a)*t;y=(-1)*fMAX+2*fMAX*s;if(y>=0){if(y==0)sy++;if(y<=(*p)(x))sx++;}23/26else{if(y>=(*p)(x))vx++;}}return((sx-vx)/n);}main(){doublen,fMAX,a,b,Mont,A,B,h,i=0,N;longtime0,time1;printf("Entera=:");scanf("%lf",&a);24/26printf("\nEnterb=:");scanf("%lf",&b); printf("\nEnterN=:");scanf("%lf",&N);/*N表示區(qū)間個(gè)數(shù)*/printf("\nEntern=:");scanf("%lf",&n);printf("\nEnterfMAX(此函數(shù)界)=:");scanf("%lf",&fMAX);

h=(b-a)/N;time0=clock();25/26fo