计算材料学概述-之-蒙特卡洛方法.详解

计算材料学概述第三章蒙特卡罗方法（MonteCarlo）主要内容MonteCarlo模拟发展简介MonteCarlo模拟基本原理MonteCarlo模拟典型算法MonteCarlo模拟典型应用蒙特卡洛法是什么？蒙特卡洛(MonteCarlo)方法，是在简洁的理论准则基础上，接受反复随即抽样的方法，解决困难系统的问题。其实质是一种概率和统计的问题。蒙特·卡罗方法（MonteCarlomethod），也称统计模拟方法，是二十世纪四十年头中期由于科学技术的发展和电子计算机的独创，而被提出的一种以概率统计理论为指导的一类特别重要的数值计算方法。是指运用随机数（或更常见的伪随机数）来解决很多计算问题的方法。MC的基本思想MC基本思想很早以前就被人们所发觉和利用。17世纪，人们就知道用事务发生的“频率”来确定事务的“概率”。但要真正实现随机抽样是很困难的，甚至几乎是不行能的。高速计算机的出现，使得用数学方法在计算机上大量、快速地模拟这样的试验成为可能。确定性系统随机性系统模拟自然界Monte-Carlo模拟，即随机模拟（重复“试验”）重复试验计算机模拟MonteCarlo方法：亦称统计模拟方法，statisticalsimulationmethod

利用随机数进行数值模拟的方法MonteCarlo名字的由来：是由Metropolis在二次世界大战期间提出的：Manhattan支配，探讨与原子弹有关的中子输运过程；MonteCarlo是摩纳哥（monaco)的首都，该城以赌博著名NicholasMetropolis(1915-1999)Monte-Carlo,MonacoMonteCarlo方法简史简洁地介绍一下MonteCarlo方法的发展历史1、Buffon投针试验：18世纪，法国数学家ComtedeBuffon利用投针试验估计的值dL1777年法国科学家布丰提出的一种计算圆周率的方法——随机投针法，即著名的布丰投针问题。这一方法的步骤是：

1)取一张白纸，在上面画上很多条间距为d的平行线。

2)取一根长度为l（l<d）的针，随机地向画有平行直线的纸上掷n次，视察针与直线相交的次数，记为m

3）计算针与直线相交的概率．

布丰本人证明白，这个概率是:

p=2l/(πd)，π为圆周率：

利用这个公式可以用概率的方法得到圆周率的近似值。下面是一些资料试验者年头投掷次数相交次数圆周率估计值沃尔夫1850500025313.1596史密斯1855320412193.1554德摩137福克斯188410304893.1595拉泽里尼1901340818083.1415929赖纳192525208593.1795布丰投针试验是第一个用几何形式表达概率问题的例子，他首次运用随机试验处理确定性数学问题，为概率论和蒙特卡罗方法的发展起到确定的推动作用。

MonteCarlo方法之随机数的产生很多计算机系统都有随机数生成函数F90:callrandom_seedcallrandom_number(a)2、ISEED=RTC()X=RAN(ISEED)Y=RAN(ISEED）Matlab:x=rand(N)产生元素在(0,1)间随机分布的N*N矩阵

s=rand(‘state’,0)重设该生成函数到初始状态留意：上述随机数序列均具周期性，如上页random子程序的周期约230。实例一、计算π值计算过程：1、构造或描述问题的概率过程2、从概率密度函数动身进行随机抽样，实现从已知概率分布的抽样，得到特征量的一些模拟结果——计算均值MonteCarlo方法之典型算法与应用考虑平面上的一个边长为1的正方形及其内部的一个形态不规则的“图形”，如何求出这个“图形”的面积呢？MonteCarlo方法是这样一种“随机化”的方法：向该正方形“随机地”投掷N个点，若有M个点落于“图形”内，则该“图形”的面积近似为M/N。用该方法计算π的基本思路是：

1

、依据圆面积的公式：

s=πR^2

,当R=1时，S=π。

2、由于圆的方程是：x^2+y^2=1(x^2为x的平方的意思),因此1/4圆面积为x轴、y轴和上述方程所包围的部分。

3、假如在1*1的正方形中匀整地落入随机点，则落入1/4圆中的点的概率就是1/4圆的面积。其4倍，就是圆面积。由于半径为1，该面积的值为π的值。

REALR,R1,R2,PI ISEED=RTC() N0=0 N=300000 DOI=1,N R1=RAN(ISEED) R2=RAN(ISEED) R=SQRT(R1*R1+R2*R2) IF(R<1.0)N0=N0+1 ENDDO PI=4.0*N0/N WRITE(*,*)PI END面积的计算f(x)x辛普逊方法I=ΣSn蒙特-卡洛方法f(x)x在长方形中均匀投N0组(x,y)如y<f(x),则N=N+1I=(N/N0)×S0SS011MC的优点MC与传统数学方法相比，具有直观性强，简便易行的优点，该方法能处理一些其他方法无法解决的负责问题，并且简洁在计算机上实现，在很大程度上可以代替很多大型、难以实现的困难试验和社会行为。无污染、无危急、能摆脱试验误差。MonteCarlo方法利用随机抽样的方法来求解物理问题;数值解法:从一个物理系统的数学模型动身,通过求解一系列的微分方程来的导出系统的未知状态;MonteCarlo方法并非只能用来解决包含随机过程的问题:例如:用MonteCarlo方法计算定积分.

对这样的问题可将其转换成相关的随机过程,然后用MonteCarlo方法进行求解留意以下两点：MC的应用自然现象的模拟：宇宙射线在地球大气中的传输过程；高能物理实验中的核相互作用过程；数值分析：数学问题，求积分，求逆矩阵，解线性代数等经济学模拟：库存问题，随机服务系统中排队问题人口问题：人口的诞生，传染病的扩散；乃至动物的生态竞争金属学：扩散、组织长大、相变过程蒙特-卡洛模拟的意义能探讨不同边界、不同材料的影响理论不行能、试验耗费太大用于试验设计无污染反应堆防护核弹爆炸能摆脱试验误差作理论和试验的桥梁MonteCarlo模拟的步骤：依据欲探讨的物理系统的性质，建立能够描述该系统特性的理论模型，导出该模型的某些特征量的概率密度函数；（即构造或描述问题的概率过程）从概率密度函数动身进行随机抽样，实现从已知概率分布的抽样，得到特征量的一些模拟结果；有了明确的概率过程后，为了实现过程的数值模拟，必需实现从已知概率分布的随机数的抽样，进行大量的随机模拟试验，从中获得随机变量的大量试验值。产生已知概率分布的随机变量，是实现MC方法的关键步骤，其中最基本的是（0，1）匀整分布。对模拟结果进行分析总结，预言物理系统的某些特性。4.模拟结果的检验MonteCarlo算法的主要组成部分概率密度函数(pdf)—必需给出描述一个物理系统的一组概率密度函数;随机数产生器—能够产生在区间[0,1]上（匀整）分布的随机数抽样规则—如何从在区间[0,1]上匀整分布的随机数动身,随机抽取听从给定的pdf的随机变量;模拟结果记录—记录一些感爱好的量的模拟结果误差估计—必需确定统计误差（或方差）随模拟次数以及其它一些量的变更；削减方差的技术—利用该技术可削减模拟过程中计算的次数；并行和矢量化—可以在先进的并行计算机上运行的有效算法实例二定积分计算事实上，不少的统计问题，如计算概率、各阶距等，最终都归结为定积分的近似计算问题。下面考虑一个简洁的定积分!计算x**2在(0,1)上积分计算过程：1、构造或描述问题的概率过程:产生听从分布f（x）的随机变量Xi()（i=1,2,···,N）2、从概率密度函数动身进行随机抽样，实现从已知概率分布的抽样，得到特征量的一些模拟结果——计算均值()

REALY Y=0 N=300000 ISEED=RTC() DOI=1,N X=RAN(ISEED) Y=Y+X**2/NENDDO WRITE(*,*)Y END

lim∑x2dx(dx0)MonteCarlo方法另一个重要问题：随机数随机数：由单位矩阵分布中所产生的简洁子样称为随机数序列，其中的每一个个体称为随机数。但真正的随机数的不适合电子计算机上运用，因为它须要很大的存储量。利用某些物理现象可以在电子计算机上产生随机数，且其产生的序列无法重复实现，使程序无法复算，结果无法验证，同时须要增加随机数发生器和电路联系等附加设备。伪随机数：是有数学递推公式所产生的随机数。（近似的具备随机数的性质。）An+1=T(A);An+1=An+k+1伪随机的优点和缺点：推断伪随机数好坏的方法：1、它能够有较好的匀整性和独立性；2、它的费用大小，即指所消耗计算机的时间；3、容量要求尽可能大。随机数产生的方法产生匀整分布随机数的几种方法；（1）物理方法；（2）数学方法。伪随机数产生方法：加同余法乘同余法乘加同余法取中方法逆变换法合成法筛选法。。。。。。关于随机数的几点留意注1由于匀整分布的随机数的产生总是接受某个确定的模型进行的，从理论上讲，总会有周期现象出现的。初值确定后，全部随机数也随之确定，并不满足真正随机数的要求。因此通常把由数学方法产生的随机数成为伪随机数。但其周期又相当长，在实际应用中几乎不行能出现。因此，这种由计算机产生的伪随机数可以当作真正的随机数来处理。注2应对所产生的伪随机数作各种统计检验，如独立性检验，分布检验，功率谱检验等等。FORTRAN语言产生随机数的实例random_number(x)

产生一个0到1之间的随机数（x可以是向量），但是每次总是那几个数。

用了random_seed

()后，系统依据日期和时间随机地供应种子，使得随机数更随机了。

program

random

real

::

x

call

random_seed

()

!

系统依据日期和时间随机地供应种子

call

random_number

(x)

!

每次的随机数就都不一样了

write(*,*)

x

stop

end

program

randomMonteCarlo方法之随机数的产生很多计算机系统都有随机数生成函数F90:callrandom_seedcallrandom_number(a)2、ISEED=RTC()X=RAN(ISEED)Y=RAN(ISEED）Matlab:x=rand(N)产生元素在(0,1)间随机分布的N*N矩阵

s=rand(‘state’,0)重设该生成函数到初始状态留意：上述随机数序列均具周期性，如上页random子程序的周期约230。FinitedifferenceapproximationofdifferentialequationsAdifferentialequationcanbeapproximatedbyafinitedifferencescheme.ForexampleForwardtimeBackwardspaceCentralspaceForwardtaime-CentralspaceFICK其次定律Fick其次定律——稳态扩散解

REALC0(0:1000+1),C(0:1000+1)!C(DISTANCE)

C0=0.1 C0(0)=0.8!BOUNDARYCONDITION C0(1001)=0.1!BOUNDARYCONDITION C=C0OPEN(1,FILE="F:\DIF.dat") DOJT=1,50000！Time DOIX=1,1000!distance C(IX)=C0(IX)+0.45*(C0(IX+1)+C0(IX-1)-2.0*C0(IX))! C(0)=0.8! C(1001)=0.1 IF(JT==50000)WRITE(1,*)IX,C(IX) ENDDO C0=CENDDO END应用之二生日问题MC模拟假设有n个人在一起，各自的生日为365天之一，依据概率理论，与很多人的直觉相反，只需23个人便有大于50％的几率人群中至少有2个人生日相同。n理论几率模拟几率0.1170.1100.4110.4120.5270.5200.7060.6920.9410.936500.9860.987INTEGERM(1:10000),

NUMBER1(0:364),

NUMBER2REALX,Y

ISEED=RTC() DOJ=1,10000 NUMBER1=0 X=RAN(ISEED) NUMBER1(0)=INT(365*X+1)

JJJ=1 DOI=1,365 Y=RAN(ISEED) NUMBER2=INT(365*Y+1) ETR=COUNT(NUMBER1.EQ.NUMBER2) IF(ETR==1)THEN EXIT ELSE JJJ=JJJ+1 M(J)=JJJ NUMBER1(I)=NUMBER2ENDIFENDDO ENDDO

DOI=1,10000IF(M(I).LE.23)SUM=SUM+1 ENDDOPRINT*,SUM/10000ENDMC在材料学领域的应用

——随机行走背景如，布朗运动－－－最简洁、无限制随机行走(Unrestrictedrandonwalk,RW)startend平均平方端-端位移:,自然科学和社会生活中很多现象都与随机运动有关可以模拟的内容？扩散；分子运动；。。。。。。如图所示，第i个分子在经过N步随机行走后距原点距离为R，对n个分子每步的位移平方求和后取平均值就得到了全部分子距原点的方均距离<R2>：

！MonteCarloSimulationofOneDimensionalDiffusion

INTEGERX,XX(1:1000,1:1000) REALXXM(1:1000)! X:INSTANTANEOUSPOSITIONOFATOM! XX(J,I)：X*X,J:第几天试验,I:第几步跳动! XXM(I):THEMEANOFXX WRITE(*,*)"试验天数JMAX，试验次数IMAX" READ(*,*)JMAX,IMAX ISEED=RTC() DOJ=1,JMAX!第几天试验 X=0!!! DOI=1,IMAX!第几步跳动 RN=RAN(ISEED) IF(RN<0.5)THEN X=X+1 ELSE X=X-1 ENDIF XX(J,I)=X*X ENDDO ENDDO OPEN(1,FILE=“f:\DIF1.DAT") DOI=1,IMAX XXM=0.0 XXM(I)=1.0*SUM(XX(1:JMAX,I))/JMAX!! WRITE(1,*)I,XXM(I) ENDDO CLOSE(1) END! MonteCarloSimulationofTwoDimensionalDiffusion

INTEGERX,Y,XY(1:1000,1:1000) REALXYM(1:1000)! X:INSTANTANEOUSPOSITIONOFATOM! XY(J,I)：X*Y,J:第几天试验,I:第几步跳动! XYM(I):THEMEANOFXY WRITE(*,*)"试验天数JMAX，试验次数IMAX" READ(*,*)JMAX,IMAX ISEED=RTC() DOJ=1,JMAX!第几天试验 X=0!!! Y=0!!! DOI=1,IMAX!第几步跳动 RN=RAN(ISEED) IF(RN.LT.0.25)THEN x=x y=y-1 ENDIF IF(RN.LT.0.5.AND.RN.GE.0.25)THEN x=x y=y+1 ENDIF IF(RN.LT.0.75.AND.RN.GE.0.5)THEN x=x-1 y=y ENDIF IF(RN.GE.0.75)THEN x=x+1 y=y ENDIF XY(J,I)=X*X+Y*Y ENDDO ENDDO OPEN(1,FILE=“f:\DIF2.DAT") DOI=1,IMAX XYM=0.0 XYM(I)=1.0*SUM(XY(1:JMAX,I))/JMAX!! WRITE(1,*)I,XYM(I) ENDDO CLOSE(1) END! MonteCarloSimulationofTwoDimensionalDiffusion

INTEGERX,XY(1:1000,1:1000),y,XN(1:4),YN(1:4),RN REALXYM(1:1000)! X:INSTANTANEOUSPOSITIONOFATOM! XY(J,I)：X*Y,J:第几天试验,I:第几步跳动! XYM(I):THEMEANOFXY WRITE(*,*)"试验天数JMAX，试验次数IMAX" READ(*,*)JMAX,IMAX XN=(/0,0,-1,1/) YN=(/-1,1,0,0/) ISEED=RTC() DOJ=1,JMAX!第几天试验 X=0!!! Y=0!!! DOI=1,IMAX!第几步跳动 RN=4*RAN(ISEED)+1 X=X+XN(RN) Y=Y+YN(RN) XY(J,I)=X*X+Y*Y ENDDO ENDDO OPEN(1,FILE="C:\DIF2.DAT") DOI=1,IMAX XYM=0.0 XYM(I)=1.0*SUM(XY(1:JMAX,I))/JMAX!! WRITE(1,*)I,XYM(I) ENDDO CLOSE(1)！做三维空间随机行走？ END实际环境是困难的:可变，缺陷，风，季节，电场等。

空间中有随机性缺陷不退行走自回避行走例如模拟高分子的位形用随机行走方法模拟高分子位形是用随机行走的轨迹代表高分子的位形，行走过的位置代表的是构成分子的原子或官能团，因此，无限制随机行走忽视了体斥效应。不退行走就是禁止在每一步行走后马上倒退，可以解决刚走的一步与上一步重叠的问题。但不退行走没有完全解决高分子的体斥效应。自回避行走就是全部已走过的位置不能再走，这样就完全解决了体斥效应问题。

气体分子的随机行走

假设在一个空旷封闭的房间中心滴上一滴香水，挥发的香水分子随机的与空气中的粒子发生碰撞，最终会扩散到整个房间。这里我们将通过计算机模拟400个分子在二维平面内的随机运动探讨熵变现在来探讨一下熵在我们这个气体扩散模型中是什么意思。在刚起先的时候，我们的系统处在一个特别有序的状态：全部的分子都处于原点。这时系统的熵为零，系统不存在任何的无序度。随着时间的推移，分子起先不断的向外扩散，这时系统出现了无序状态，熵起先渐渐增加。系统的熵和我们预料的一样在随时间增大，但当全部分子布满整个边界以内区域时，系统的熵起先趋于稳定。从这一结果中我们可以了解：当全部的分子随机的布满整个区域时，虽然当我们跟踪某一个确定分子时，它还是在区域内到处乱窜，但每个小区域内分子的密度却不会再变更了。所以，一旦气体分子扩散到整个区域以后，不管我们再等上多少时间，系统的熵都不会再有太大的起伏。换句话说，让系统自动回到起先的状态，即全部分子都在原点的状态，已经不行能了。画一个圆晶粒

USEMSFLIB! INTEGERXR,YR!在的区域中画一个圆

PARAMETERXR=400,YR=400 INTEGERR,S(1:XR,1:YR) X0=XR/2!圆心位置X0，YO Y0=YR/2 R=MIN(X0-10,Y0-10)!圆半径

S=0!像素的初始状态（颜色）

DOI=1,XR DOJ=1,YR IF((I-X0)**2+(J-Y0)**2<=R**2)S(I,J)=10 IER=SETCOLOR(S(I,J)+3) IER=SETPIXEL(I,J) ENDDO ENDDO END画晶界! 画一个圆 USEMSFLIB INTEGERXR，YR!在的区域中画一个圆 PARAMETERXR=400,YR=400 INTEGERR,S(0:XR+1,0:YR+1),XN(1:4),YN(1:4),SNS XN=(/0,0,-1,1/) YN=(/-1,1,0,0/) X0=XR/2!圆心位置X0，YO Y0=YR/2 R=MIN(X0-10,Y0-10)!圆半径 S=0!像素的初始状态（颜色） DOI=1,XR DOJ=1,YR IF((I-X0)**2+(J-Y0)**2<=R**2)S(I,J)=10 IER=SETCOLOR(S(I,J)) IER=SETPIXEL(I,J) ENDDO ENDDO DOI=1,XR!画晶界 DOJ=1,YR NDS=0 DOK=1,4 IF(S(I,J).NE.S(I+XN(K),J+YN(K)))NDS=NDS+1 ENDDO IF(NDS>0)THEN IER=SETCOLOR(9) ELSE IER=SETCOLOR(8) ENDIF IER=SETPIXEL(I,J) ENDDO ENDDO END如何画有确定宽度的晶界？例题例求前N个自然数倒数和。即求程序须要几个变量：累加的项数I(整型)，存放通项值的T(实型)，存放总和的S(实型)，待输入的总项数N(整型)。PROGRAMMAINIMPLICITNONEINTEGERN,IREALT,SREAD*,NS=0.0;I=1DOT=1.0/IS=S+TI=I+1IF(I>N)EXITENDDOPRINT*,"SUM=",SENDPROGRAMMAIN程序的运行结果如下15

SUM=3.182290用DO结构处理循环问题，最关键的是要处理好初值、循环终止条件和循环体语句的关系。FinitedifferenceapproximationofdifferentialequationsAdifferentialequationcanbeapproximatedbyafinitedifferencescheme.ForexampleForwardtimeBackwardspaceCentralspaceForwardtime-CentralspaceFICK其次定律Fick其次定律——稳态扩散解

REALC0(0:1000+1),C(0:1000+1)!C(DISTANCE)

C0=0.1 C0(0)=0.8!BOUNDARYCONDITION C0(1001)=0.1!BOUNDARYCONDITION C=C0OPEN(1,FILE="F:\DIF.dat")

DOJT=1,50000！Time DOIX=1,1000!distance C(IX)=C0(IX)+0.45*(C0(IX+1)+C0(IX-1)-2.0*C0(IX)) C(0)=0.8 C(1001)=0.1 IF(JT==50000)WRITE(1,*)IX,C(IX) ENDDO C0=CENDDO PRINT*,"PROGRAMMEISFINISHED"CLOSE(1)END蒙特卡罗(MonteCarlo)方法1、简介2、相关几个例子MonteCarlo方法：亦称统计模拟方法，statisticalsimulationmethod

利用随机数进行数值模拟的方法MonteCarlo名字的由来：是由Metropolis在二次世界大战期间提出的：Manhattan支配，探讨与原子弹有关的中子输运过程；MonteCarlo是摩纳哥（monaco)的首都，该城以赌博著名NicholasMetropolis(1915-1999)Monte-Carlo,MonacoMonteCarlo方法的基本思想很早以前就被人们所发觉和利用。早在17世纪，人们就知道用事务产生的“频率”来近似事务的“概率”。18世纪人们用投针试验的方法来确定圆周率π。本世纪40年头电子计算机的出现，特殊是近年来高速电子计算机的出现，使得用数学方法在计算机上大量、快速地模拟这样的试验成为可能。1930年，EnricoFermi利用MonteCarlo方法探讨中子的扩散，并设计了一个MonteCarlo机械装置，Fermiac,用于计算核反应堆的临界状态VonNeumann是MonteCarlo方法的正式奠基者,他与StanislawUlam合作建立了概率密度函数、反累积分布函数的数学基础，以及伪随机数产生器。在这些工作中，StanislawUlam意识到了数字计算机的重要性合作起源于Manhattan工程：利用ENIAC(ElectronicNumericalIntegratorandComputer)计算产额MonteCarlo模拟的应用：自然现象的模拟：宇宙射线在地球大气中的传输过程；高能物理实验中的核相互作用过程；试验探测器的模拟数值分析：利用MonteCarlo方法求积分材料学中的应用：扩散、晶粒的长大、再结晶等MonteCarlo模拟在物理探讨中的作用MonteCarlo模拟的步骤：依据欲探讨的物理系统的性质，建立能够描述该系统特性的理论模型，导出该模型的某些特征量的概率密度函数；从概率密度函数动身进行随机抽样，得到特征量的一些模拟结果；对模拟结果进行分析总结，预言物理系统的某些特性。留意以下两点：MonteCarlo方法与数值解法的不同:MonteCarlo方法利用随机抽样的方法来求解物理问题;数值解法:从一个物理系统的数学模型动身,通过求解一系列的微分方程来的导出系统的未知状态;MonteCarlo方法并非只能用来解决包含随机的过程的问题:很多利用MonteCarlo方法进行求解的问题中并不包含随机过程例如:用MonteCarlo方法计算定积分.对这样的问题可将其转换成相关的随机过程,然后用MonteCarlo方法进行求解MonteCarlo算法的主要组成部分概率密度函数—必需给出描述一个物理系统的一组概率密度函数;随机数产生器—能够产生在区间[0,1]上匀整分布的随机数抽样规则—如何从在区间[0,1]上匀整分布的随机数动身,随机抽取听从给定的pdf的随机变量;模拟结果记录—记录一些感爱好的量的模拟结果误差估计—必需确定统计误差（或方差）随模拟次数以及其它一些量的变更；削减方差的技术—利用该技术可削减模拟过程中计算的次数；并行和矢量化—可以在先进的并行计算机上运行的有效算法确定性系统随机性系统模拟自然界Monte-Carlo模拟，即随机模拟（重复“试验”）重复试验计算机模拟蒙特卡罗方法的基本原理及思想当所要求解的问题是某种事务出现的概率，或者是某个随机变量的期望值时，它们可以通过某种“试验”的方法，得到这种事务出现的频率，或者这个随机变数的平均值，并用它们作为问题的解。把蒙特卡罗解题归结为三个主要步骤： 构造或描述概率过程； 实现从已知概率分布抽样； 建立各种估计量。对于本身就具有随机性质的问题，主要是正确描述和模拟这个概率过程；对于不是随机性质的确定性问题，比如计算定积分，就必需事先构造一个人为的概率过程。最简洁、最基本、最重要的一个概率分布是(0,1)上的匀整分布（或称矩形分布）。实现从已知概率分布抽样(模拟试验)：构造了概率模型以后，产生已知概率分布的随机变量（或随机向量），就成为实现蒙特卡罗方法模拟试验的基本手段，这也是蒙特卡罗方法被称为随机抽样的缘由。建立各种估计量：构造了概率模型并能从中抽样后，即实现模拟试验后，就要确定一个随机变量，作为所要求的问题的解，我们称它为无偏估计。收敛性:

由大数定律,Monte-Carlo模拟的收敛是以概率而言的．误差:

用频率估计概率时误差的估计，可由中心极限定理，给定置信水平的条件下，有：

模拟次数:由误差公式得随机数的产生原理-计算机抽样试验随机数的产生是实现MC计算的先决条件。而大多数概率分布的随机数的产生都是基于匀整分布U(0,1]的随机数。首先，介绍听从匀整分布U(0,1]的随机数的产生方法。其次，介绍听从其他各种分布的随机数的产生方法。以及听从正态分布的随机数的产生方法。最终，关于随机数的几点注。(0,1]匀整分布随机数的产生最常用的是在(0,1]区间内匀整分布的随机数，其他分布的随机数可利用匀整分布随机数来产生。匀整分布随机数的产生：乘同余法产生区间(0,1]内匀整分布的随机数的递推公式是：

其中，是乘因子，M是模数。第一式称做以M为模数（modulus）的同余式，即以M除后得到的余数记为。当给定了一个初值（称为种子），计算出的即在(0,1]上匀整分布的随机数。

例1：R=RAND()ISEED=RTC()R=RAN(ISEED)其他各种分布的随机数的产生基本方法有如下三种：逆变换法合成法筛选法逆变换法设随机变量的分布函数为，定义定理设随机变量听从上的匀整分布，则的分布函数为。因此，要产生来自的随机数，只要先产生来自的随机数，然后计算即可。其步骤为

合成法

合成法的应用最早见于Butlter的书中。构思如下：假如的密度函数难于抽样，而关于的条件密度函数以及的密度函数均易于抽样，则的随机数可如下产生：可以证明由此得到的听从。筛选抽样

假设我们要从抽样，假如可以将表示成，其中是一个密度函数且易于抽样，而，是常数，则的抽样可如下进行：定理设的密度函数，且，其中，，是一个密度函数。令和分别听从和，则在的条件下，的条件密度为

标准正态分布的随机数的随机数产生方法很多。简要介绍三种。法1、变换法（Box和Muller1958）设，是独立同分布的变量，令则与独立，均听从标准正态分布。法2、结合合成法与筛选法。（略）法3、近似方法（利用中心极限定理）即用个变量产生一个变量。其中是抽自的随机数，可近似为一个变量。由于正态随机数的独立性，当很小时，按(a)、(b)所构成的过程的相关时间很短，从而具有很高的截止频率。当然，过小，样本的计算量过大。因此，选取适当小即可。关于随机数的几点注注1由于匀整分布的随机数的产生总是接受某个确定的模型进行的，从理论上讲，总会有周期现象出现的。初值确定后，全部随机数也随之确定，并不满足真正随机数的要求。因此通常把由数学方法产生的随机数成为伪随机数。但其周期又相当长，在实际应用中几乎不行能出现。因此，这种由计算机产生的伪随机数可以当作真正的随机数来处理。注2应对所产生的伪随机数作各种统计检验，如独立性检验，分布检验，功率谱检验等等。MonteCarlo方法相关引例:1、Buffon投针试验：18世纪，法国数学家ComtedeBuffon利用投针试验估计的值dLProblemofBuffon’sneedle:Ifaneedleoflengthlisdroppedatrandomonthemiddleofahorizontalsurfaceruledwithparallellinesadistanced>lapart,whatistheprobabilitythattheneedlewillcrossoneofthelines?Solution:Thepositioningoftheneedlerelativetonearbylinescanbedescribedwitharandomvectorwhichhascomponents:Therandomvectorisuniformlydistributedontheregion[0,d)×[0,).Accordingly,ithasprobabilitydensityfunction1/d.Theprobabilitythattheneedlewillcrossoneofthelinesisgivenbytheintegral圆周率的值π=3.

14159265358979323846264338327950288419716939937510

58209749445923078164062862089986280348253421170679

82148086513282306647093844609550582231725359408128

48111745028410270193852110555964462294895493038196

44288109756659334461284756482337867831652712019091

45648566923460348610454326648213393607260249141273

72458700660631558817488152092096282925409171536436

78925903600113305305488204665213841469519415116094

33057270365759591953092186117381932611793105118548

07446237996274956735188575272489122793818301194912

98336733624406566430860213949463952247371907021798

60943702770539217176293176752384674818467669405132

00056812714526356082778577134275778960917363717872

14684409012249534301465495853710507922796892589235

420199561121290219608640344181598136297747713.....什么是蒙特-卡洛模拟？计算机试验物理试验——用探测器取数据5RHIC-STAR实验探测器蒙特-卡洛模拟

——

用计算机取数据一个例子——道尔顿板理论——分布曲线试验——丢铁球模拟——计算机丢球6模拟的动身点——模型依据球和钉子碰撞的规律追踪每个球的运动依据落入每个格子的几率给出这一几率的每次实现89令n为铁球数，m为格子数。p(i)是铁球落入第i个格子中的概率再产生n个匀整分布的随机数： ra(j),j=1,2,...,n将p(i)累加s=0.r(1)=0.do2i=1,ms=s+p(i)r(i+1)=senddodo3k=1,mnp(k)=0do3j=1,nif(ra(j).ge.r(k).and.ra(j).lt.r(k+1))[r(k)≤ra(j)<r(k+1)]np(k)=np(k)+1Enddonp(k)是落入第k个盒子中的粒子数。r(m)Σp(i)=1r(3)p(1)+p(2)r(2)p(1)r(1)0···i=1,m···通过比较第j个铁球的ra(j)和r(i)来确定这个球落入那一格子中：例一：道尔顿板的蒙特-卡洛模拟p3p4p5p2p1p6p7例二面积的计算f(x)x辛普逊方法I=ΣSn蒙特-卡洛方法f(x)x在长方形中均匀投N0组(x,y)如y<f(x),则N=N+1I=(N/N0)×S0SS011考虑平面上的一个边长为1的正方形及其内部的一个形态不规则的“图形”，如何求出这个“图形”的面积呢？MonteCarlo方法是这样一种“随机化”的方法：向该正方形“随机地”投掷N个点，若有M个点落于