计算物理学_蒙特卡罗方法
- 格式：doc
- 大小：2.31 MB
- 文档页数：43
h
h
第八讲 蒙特卡罗方法
蒙特卡罗（Monte Carlo简称MC）方法又称随机抽样法（Random Sampling）、随机模拟（Random Simulation）或统计试验法(Statistic Testing）。这个方法的起源可以追溯到十七世纪或更早的年代。Monte Carlo 是摩纳哥（Monaco）的一个著名城市，位于地中海之滨，以旅游赌博闻名。Von Neumann 等人把计算机随机模拟方法定名为Monte Carlo方法，显然反映了这种方法带有随机的性质。
简单地说，MC方法是一种利用随机统计规律，进行计算和模拟的方法。它可用于数值计算，也可用于数字仿真。在数值计算方面，可用于多重积分、线性代数求解、矩阵求逆以及用于方程求解，包括常微分方程、偏微分方程、本征方程、非齐次线性积分方程和非线性方程等。在数字仿真方面，常用于核系统临界条件模拟、反应堆模拟以及实验核物理、高能物理、统计物理、真空、地震、生物物理和信息物理等领域。
§8.l蒙特卡罗方法的基础知识
8.1.1 基本概念
为了对MC方法有一点初步的认识，请先看使用MC方法的几个例子。
蒲丰投针问题：
蒲丰（Buffon-法国著名数学家）在1777年发现随机投针的概率与无理数之间的关系．这个问题是说，若在平面上画有距离为a的平行线束，向平面上投掷长为()lla的针，试求针与一平行线相交的概率。
这个问题的解法如下：以M表示落下后针的中点，x表M与最近一平行线的距离，表针与此线的交角，见上图。
可见，h
h 02 0/,xa
这两式决定x平面上一矩形R；为了使针与一平行线（这线必定是与针中点M最近的平行线）相交，充分而且必要条件是
2sinlx
这个不等式决定R中一个子集G。因此，我们的问题等价于向R中均匀分布地掷点而求点落于G中的概率P.根据概率的几何意义，得
0 222sin()ldlPaa
此式提供了求值的一个方法：可以通过投针事件求得针与平行线相交概率P，求得值：
2/()lPa （8.1）
若投掷次数为m，针与平行线相交的次数为n，那么
/pnm
即 2/()lman
于是，可用投针试验来求无理数的近似值．下表列举了历史上若干学者投针试验计算值的结果：
试验者 年份 投针次数 的计算值
Wolf 1850 5000 3.1596
Smith 1855 3204 3.1553
Fox 1894 1120 3.1419
Lazzarini 1901 3408 3.1415929
射击问题（打靶游戏）：
设r表示射击运动员的弹着点到靶心的距离，()gr表示击中r处相应的得分数（环数），分布密度函数()fr表示该运动员的弹着点分布，它反映运动员射击水平。积分
0()()ggrfrdr （8.2）
表示这个运动员的射击成绩。用概率语言说，g就是随机变量()gr的数学期望值，记为{}gEg。现在，假设这个射击运动员射击N次，弹着点依次是环数分别为h
h 12,,,Nggg，则自然地认为N次射击得分的平均值
11NiiggN （8.3）
这个平均值相当好地代表了这个射击运动员的成绩。换句话说，g是积分（8.2）式的一个估计值（或近似值）。这个例子通常称为打靶游戏，它直观地说明了蒙特卡罗方法计算定积分的基本思想。为进一步阐明这个思想，我们再举个例子：计算积分
10 1(),max()Ifxdxfx （8.4）
直观上，就是在边长为1的正方形里随机投点，当点落在()yfx曲线下面，对积分值有“贡献”，否则对积分值无“贡献”。为此，假设向这个边长为1的正方形里随机投点N次，点落在()yfx曲线下面n次，则（8.4）式积分值近似为/InN来近似。
从上述几个例子可以看到，当所要求解的问题是某种事件出现的概率，或者是某个随机变量的期望值时，它们可以通过某种“试验”的方法，得到这种事件出现的频率，或者这个随机变量的平均值，并用它们作为问题的解。这就是蒙特卡罗方法的基本思想。
所以用蒙特卡罗方法求解问题时，首先要建立一个随机模型，然后要构造一系列的随机数用以模拟这个过程，最后要作统计性的处理。关于建立随机模型，因问题而异。
8.1.2 随机变量及其分布函数
在一定条件下发生的事件分为必然事件（必然发生）、不可能事件（恒不发生）和随机事件（可能发生也可能不发生），事件发生的可能性大小用概率p表示。必然事件的发生概率为1，不可能事件的概率为零；随机事件发生的概率01p。由于测量的随机误差和物理现象本身的随机性。一次测量得到的某个值是随机的，因此，实验观测的物理量是随机变量，被研究的物理问题是一个随机事件．通常，描写随机事件A发生的概率用()pA来表示，显然，0()1pA．
经常碰到的随机变量有两类：一类是离散型随机变量，这种随机变量只能取有限个数值，能够—一列举出来；另一类是连续型随机变量，这种随机变量的可能值连续的分布某个区间．h
h 离散型随机变量可用分布列：
1212,,,,,,,,nnxxxppp （8.5）
来描写，它表示取值ix的概率为ip（1,2,,,in）．即
{}iippx （8.6）
分布列描写了离散型随机变量的概率分布．
连续型随机变量可能取的值是不可列的，因此不能象离散型随机变量那样，用分布列来描写它的概率分布，而要用概率分布密度函数来描写．
考虑连续型随机变量落在区间[,]xxx内的概率()pxxx，如果极限
{}lim()xpxxxfxx （8.7）
存在，则函数()fx描写了在点x的概率密度，把()fx叫做随机变量的概率分布密度，简称分布密度或密度函数．于是，随机变量落在[,]ab内的概率
()pab可写为
{}()bapabfxdx （8.8）
显然，上式只有当()fx可积时才有意义。
此外，还可定义随机变量的分布函数来描写随机变量的概率分布规律．
分布函数()Fx一定义为
(){}Fxpx （8.9）
即分布函数()Fx在x处的值等于随机变量取值小于或等于x这样一个随机事件的概率．显然，()0,()1FxFx．
对于离散型随机变量，分布函数()Fx为阶梯函数
(){}iixxFxpxp （8.10）
对于连续型随机变量，分布函数()Fx与分布密度()fx满足如下关系：
()()()xFxpxfxdx
所以分布密度函数()fx和分布函数()Fx的性质：
① ()0fx
② ()1fxdxh
h ③ 对于任意实数1212,()xxxx：211221{}()()()xxpxXxFxFxfxdx
④ 若()fx在x点连续，则有'()()Fxfx
为了描述随机变量的数学特征还需引进数学期望和方差的概念．对于离散型随机变量其数学期望定义为
1niiiExp （8.11）
方差定义为
21()niiiDpxE （8.12）
对连续型随机变量，若已知分布密度()fx，则数学期望是
()Exfxdx （8.13）
方差
2()()DxEfxdx （8.14）
在实际问题中，对于一组N个实验观测数据，相应于某一个随机变量的一个样本，可以用直方图来形象地表示样本中数据分布的规律性．先将随机变量的取值范围划分为若干个区间，将落入每一区间的数据的个数m（称为频数）与随机变量的取值区间之间的关系画成阶跃曲线，即构成了直方图．数 m／N是观测数据落入这个小区间的频率．直方图的根坐标为随机变量的取值范围，纵坐标为频数．
图表示的是某一个随机事件的直方图。
8.1.3 大数定理和中心极限定理
概率论中的大数定理和中心极限定是蒙特卡罗方法的数学基础。
大数定理： 设12,,,,n为一随机变量序列，独立同分布，数学期望值Ea存在，则对任意0，有h
h lim1npa （8.15）
11niin是算术平均。大数定理指出，当n时，算水平均收敛到数学期望（或统计平均）a．也就是说，可以用算术平均代替数学期望。对于收敛的程度，和误差估计，要用到下面的中心极限定理．
中心极限定理：设12,,,,n为一随机变量序列，独立同分布，数学期望为Ea，方差2D，则当n时，
2/212XxapXedxn （8.16）
利用中心极限定理，当n很大时，可得到误差
Xan （8.17）
成立的概率为1。通常将称为置信度，1为置信水平。的定义为
2/20212Xxedx （8.18）
和X的关系可以计算得到，下面给出常用的几组和X的值：
0.5 0.05 0.02 0.01
X
0.6745 1.9600 203863 2.5758
从（8.17）式可以看到，算术平均值收敛到a的阶为(1/)On，可见，蒙特卡罗方法收敛的阶很低，收敛速度很慢。当0.05，误差1.96/n称为概率误差。误差由和1/n决定．在固定的情况下，要提高 1位精确度，就要增加 100倍试验次数．相反，若减小10倍，就可以减少100倍工作量．因此，控制方差是应用蒙特卡罗方法中很重要的一点．如果要求置信水平为0.95（95%），计算得1.96X，则说明误差估计
111.96niiann
的概率是0.95。
§8.2随机数和随机抽样
在计算机上用蒙特卡罗方法模拟某过程时，需要产生具有各种概率分布的随机变量．最简单和最基本的随机变量就是[0,1]区间上均匀分布的随机变量．这h
h 些随机变量的抽样值就称为随机数．以后谈到随机数，如果不加特别说明，专指
