蒙特卡洛算法

蒙特卡洛算法

一.蒙特卡洛算法的介紹二.隨機點的產(chǎn)生1.偽隨機數(shù)的產(chǎn)生2.準(zhǔn)隨機數(shù)的產(chǎn)生3.隨機點產(chǎn)生程序分析1蒙特卡洛算法的介紹

算法簡介

蒙特卡洛算法，也稱統(tǒng)計模擬方法，是二十世紀四十年代中期由于科學(xué)技術(shù)的發(fā)展和電子計算機的發(fā)明，而被提出的一種以概率統(tǒng)計理論為指導(dǎo)的一類非常重要數(shù)值計算方法，蒙特·卡羅方法在金融工程學(xué)，宏觀經(jīng)濟學(xué)，計算物理學(xué)(如粒子輸運計算、量子熱力學(xué)計算、空氣動力學(xué)計算)等領(lǐng)域應(yīng)用廣泛。2蒙特卡洛算法的介紹

算法描述

以概率和統(tǒng)計理論方法為基礎(chǔ)的一種計算方法。將所求解的問題同一定的概率模型相聯(lián)系，用計算機實現(xiàn)統(tǒng)計模擬或抽樣，以獲得問題的近似解。例如pi的計算：在一個面積為一的正方形里面畫一個圓，半徑0.5，與正方形四邊相切，在正方形中生成一系列隨機點，統(tǒng)計單位圓內(nèi)的點數(shù)與總點數(shù)，（圓面積和正方形面積之比為pi/4=m/n），當(dāng)隨機點取得越多（但即使取10的9次方個隨機點時，其結(jié)果也僅在前4位與圓周率吻合）時，其結(jié)果越準(zhǔn)確。3隨機點的產(chǎn)生偽隨機算法隨機數(shù)生成算法是一類重要的算法，廣泛應(yīng)用于仿真技術(shù)等場合。偽隨機數(shù)是由確定的算法生成的，其分布函數(shù)與相關(guān)性均能通過統(tǒng)計測試。與真實隨機數(shù)的差別在于，它們是由算法產(chǎn)生的，而不是一個真實的隨機過程。一般地，偽隨機數(shù)的生成方法很多，有線性同余法，直接法，逆轉(zhuǎn)法等

4隨機點的產(chǎn)生C/C++語言中偽隨機數(shù)生成算法實際上是采用了“線性同余法“。具體的計算如下：Xi=(Xi-1*A+C)modM

其中A,C,M都是常數(shù)（一般會取質(zhì)數(shù)）。當(dāng)C=0時，叫做乘同余法。引出一個概念叫seed，它會被作為X0被代入上式中，然后每次調(diào)用rand()函數(shù)都會用上一次產(chǎn)生的隨機值來生成新的隨機值?？梢钥闯鰧嶋H上用rand()函數(shù)生成的是一個遞推的序列，一切值都來源于最初的seed。所以當(dāng)初始的seed取一樣的時候，得到的序列都相同。C語言里面有RAND_MAX這樣一個宏，定義了rand()所能得到的隨機值的范圍。在C里可以看到RAND_MAX被定義成0x7fff，也就是32767。rand()函數(shù)里遞推式中M的值就是32767。線性同余法生成的是偽隨機數(shù)，粗略符合均勻分布。5隨機點的產(chǎn)生準(zhǔn)隨機算法

偽隨機算法都存在差異性，不均勻性。因此，不要求新的發(fā)生器模擬真實的均勻分布，而力求任意大小的樣本（尤其是小樣本）都能滿足低差異性。換言之，以犧牲隨機性為代價，換來均勻性的提高，稱其為準(zhǔn)隨機模擬器。目前有3種準(zhǔn)隨機序列可用來輔助生成均勻分布隨機數(shù)，分別是Halton序列、Sobol序列、Latin超立方體序列。

6隨機點的產(chǎn)生Halton序列

以質(zhì)數(shù)為基底產(chǎn)生的序列假設(shè)是以質(zhì)數(shù)2為基底的Halton序列，范圍在0～1之間。則從產(chǎn)生的列為1/21/43/41/85/83/87/81/169/16……..怎么產(chǎn)生的呢？3=1+25=1+47=3+49=1+8………..

假設(shè)以質(zhì)數(shù)3為基底，Halton序列如下1/32/31/94/97/92/95/98/91/2710/2719/27……..7隨機點的產(chǎn)生在二維中，0～1之間產(chǎn)生的點的序列就是（1/2，1/3）（1/4，2/3）（3/4，1/9）（1/8，4/9）（5/8，7/9）（3/8，2/9）（7/8，5/9）（1/16，8/9）（9/16，1/27）….8核心代碼分析privatestaticclassHaltonSequence{//Halton序列的產(chǎn)生staticfinalint[]P={2,3};//以2、3為基底產(chǎn)生序列staticfinalint[]K={63,40};//使2和3產(chǎn)生的序列中元素的個數(shù)相同privatelongindex;//索引項privatedouble[]x;//x數(shù)組定義privatedouble[][]q;//二維數(shù)組q定義privateint[][]d;//二維數(shù)組d定義HaltonSequence(longstartindex){//輸入?yún)?shù)值，得到序列中的元素index=startindex;x=newdouble[K.length];//K.length=2q=newdouble[K.length][];d=newint[K.length][];for(inti=0;i<K.length;i++){//給q，d賦值q[i]=newdouble[K[i]];//q[0]=63,q[1]=40d[i]=newint[K[i]];//d[0]=63,d[1]=40}9for(inti=0;i<K.length;i++){longk=index;//把參數(shù)傳給kx[i]=0;

for(intj=0;j<K[i];j++){//這部分通過運算來演示q[i][j]=(j==0?1.0:q[i][j-1])/P[i];d[i][j]=(int)(k%P[i]);k=(k-d[i][j])/P[i];x[i]+=d[i][j]*q[i][j];}}}10double[]nextPoint(){//根據(jù)上面求得的點來產(chǎn)生第二個點index++;for(inti=0;i<K.length;i++){for(intj=0;j<K[i];j++){d[i][j]++;x[i]+=q[i][j];if(d[i][j]<P[i]){break;}d[i][j]=0;x[i]-=(j==0?1.0:q[i][j-1]);}}returnx;}}11改寫后的程序（C++版）#include<string>#include<iostream>usingnamespacestd;intP[2]={2,3};intK[2]={63,40};longindex;doublex[2];doubleq[2][100];intd[2][100];12voidHaltonSequence(){index=1;for(inti=0;i<2;i++) {longk=index;x[i]=0;

for(intj=0;j<K[i];j++) {q[i][j]=(j==0?1.0:q[i][j-1])/P[i];d[i][j]=(int)(k%P[i]);k=(k-d[i][j])/P[i];x[i]+=d[i][j]*q[i][j];}} }13voidnextPoint() {index++;for(inti=0;i<2;i++) {for(intj=0;j<K[i];j++) {d[i][j]++;x[i]+=q[i][j];if(d[i][j]<P[i]) {break;}d[i][j]=0;x[i]-=(