MonteCarlo方法.ppt

《MonteCarlo方法.ppt》由会员分享，可在线阅读，更多相关《MonteCarlo方法.ppt（48页珍藏版）》请在三一办公上搜索。

1、Monte Carlo 方法,黄世萍,起源,这一方法源于美国在第二次世界大战进研制原子弹的曼哈顿计划。Monte Carlo方法创始人主要是这四位：Stanislaw Marcin Ulam,Enrico Fermi,John von Neumann和 Nicholas Metropolis。,Monte Carlo 方法基础,Monte Carlo（MC）方法是在简单的理论准则基础上（如简单的物质与物质以及物质与环境相互作用），采用反复随机抽样的手段，解决复杂系统的问题。该方法采用随机抽样的手法，可以模拟对象的概率与统计的问题。通过设计适当的概率模型，该方法还可以解决确定性问题，如定积分等。

2、随着计算机的迅速发展，MC 方法已在应用物理、原子能、固体物理、化学、材料、生物、生态学、社会学以及经济学等领域得到了广泛的应用。,MC meothds,classical Monte Carlo:samples are drawn from a probability distribution,often the classical Boltzmann distribution,to obtain thermodynamic properties or minimum-energy structures;quantum Monte Carlo:random walks are used to

3、 compute quantum-mechanical energies and wavefunctions,often to solve electronic structure problems,using Schrdingers equation as a formal starting point;volumetric Monte Carlo:random number generators are used to generate volumes per atom or to perform other types of geometrical analysis;kinetic Mo

4、nte Carlo:simulate processes using scaling arguments to establish timescales or by introducing stochastic effects into molecular dynamics.,Monte Carlo方法的基本思想,Monte Carlo 方法的基本思想是:为了求解某个问题，建立一个恰当的概率模型或随机过程，使得其参量（如事件的概率、随机变量的数学期望等）等于所求问题的解，然后对模型或过程进行反复多次的随机抽样试验，并对结果进行统计分析，最后计算所求参量，得到问题的近似解。,Monte Carl

5、o 方法是随机模拟方法；它不仅限于模拟随机性问题，还可以解决确定性的数学问题。对随机性问题，可以根据实际问题的概率法则，直接进行随机抽样试验，即直接模拟方法。对于确定性问题采用间接模拟方法，即通过统计分析随机抽样的结果获得确定性问题的解。用 Monte Carlo 方法解决确定性的问题主要是在数学领域，如计算重积分、求逆矩阵、解线性代数方程组、解积分方程、解偏微分方程边界问题和计算微分算子的特征值等。用 Monte Carlo 方法解决随机性问题则在众多的科学及应用技术领域得到广泛的应用，如中子在介质中的扩散问题、库存问题、随机服务系统中的排队问题、动物的生态竞争、传染病的蔓延等。,简单的例子

6、,对积分进行变换，构造新的被积函数 g（x），使得该函数满足下列条件：,g（x）是连续随机变量 的概率密度函数。,定积分是概率积分，其积分值等于概率 Pr（a b），即,这个步骤就是将一个积分转化为一个概率模型的过程；然后，反复多次的随机抽样试验，以抽样结果的统计平均作为索求概率的近似值，从而求得该积分。,具体试验步骤如下：,（1）产生服从给定分布函数 g（x）的随机变量值 xi,（2）检查 xi是否落入积分区域（a x b），如果满足条件，则记录一次。,反复进行上述试验。假设在 N 次试验后，xi落入积分区域的总次数为 m，那么，积分值近似表示为,对于随机性问题，可直接将实际的随机问题抽象为

7、概率数学模型，然后与求解确定性问题一样进行抽样试验和统计计算。,Monte Carlo 方法解决实际问题的过程中，主要有以下几个内容,建立简单而又便于实现的概率统计模型，使所求的解正是该模型的某 一事件的概率或数学期望，或该模型能够直接描述实际的随机过程。根据概率统计模型的特点和计算的需求，改进模型，以便减小方差和减低费用，提高计算效率。建立随机变量的抽样方法，包括伪随机数和服从特定分布的随机变量的产生方法。给出统计估计值及其方差或标准误差。,One-Dimensional Integrals,Methodical approachesrectangle rule,trapezoid rule

8、,Simpsons rule Quadrature formula,n uniformly separated points,Sum areas of shapes approximating shape of curve,Evaluating the general integral,假设对于下列积分,首先，按照均匀分布，随即采样区域a,b，记采样点为Xi，每个采样点的值p(x)=1/(b-a)。利用Monte Carlo算法，上述积分可以模拟为：,当N趋近于无穷大的时候，我们认为两者的值是一致的，所以可以利用上述离散的方式模拟这个积分,Monte Carlo Integration,Sto

9、chastic approachSame quadrature formula,different selection of points,n points selected from uniform distribution p(x),p(x),我们在对于a,b进行采样的时候，完全没有必要进行均匀采样，这样做只是简单而已，我们可以根据一种概率分布来进行采样，从而使得某些区域的采样密度更大。加入采样点的概率分布函数为P(x)，那么我们可以利用如下公式计算定积分的值：,证明：,3 Monte Carlo方法的收敛性和基本特点,设所求的量 x是随机变量的数学期望 E（x），那么，Monte Car

10、lo 方法通常使用随机变量的简单子样，N 的算术平均值，即,作为所求量 x 的近似值。由柯尔莫哥罗夫（Kolmogorov）大数定理可知，,即当 N 充分大时，有,成立的概率等于，亦即可以用 N 作为所求量 x 的估计值。,根据中心极限定理，如果随机变量 的标准差 不为零，那么 Monte Carlo 方法的误差为,式中，为正态差，是与置信水平有关的常量。Monte Carlo 方法的收敛速度的阶为 o（N-1/2），误差是由随机变量的标准差 s和抽样次数 N 决定的。,提高精度一位数，抽样次数要增加 倍；减小随机变量的标准差，可以减小误差。,Monte Carlo 方法具有以下四个重要特征：

11、,由于 Monte Carlo 方法是通过大量简单的重复抽样来实现的，因此，方法和程序的结构十分简单。收敛速度比较慢，因此，较适用于求解精度要求不高的问题。收敛速度与问题的维数无关，因此，较适用于求解多维问题。问题的求解过程取决于所构造的概率模型，而受问题条件限制的影响较小，因此，对各种问题的适应性很强。,随机数的产生,1 随机数与伪随机数,Monte Carlo 方法的核心是随机抽样。在该过程中往往需要各种各样分布的随机变量其中最简单、最基本的是在，区间上均匀分布的随机变量。在该随机变量总体中抽取的子样，N 称为随机数序列，其中每个个体称为随机数。,用数学的方法产生随机数是目前广泛使用的方法

12、。该方法的基本思想是利用一种递推公式：,对于给定的初始值，逐个地产生，,数学方法产生的随机数存在两个问题：,整个随机数序列是完全由递推函数形式和初始值唯一确定的，严格地说不满足随机数的相互独立的要求。存在周期现象。,基于这两个原因，将用数学方法所产生的随机数称为伪随机数。伪随机数的优点是适用于计算机，产生速度快，费用低廉。目前，多数计算机均附带有“随机数发生器”。,选择递推函数必须注意以下几点：,随机性好；在计算机上容易实现；省时；伪随机数的周期长。,2 伪随机数的产生方法,最基本的伪随机数是均匀分布的伪随机数。该方法是首先给一个 r位的数，取其中间的 r位数码作为第一个伪随机数，然后将这个数

13、平方，构成一个新的 r位的数，再取中间的 r位数作为第二个伪随机数。如此循环可得到一个伪随机数序列。该方法的递推公式为,x表示对 x 取整，运算 B（Mod M）表示 B被 M 整除后的余数。数列i 是分布在，上的。该方法由于效率较低，有时周期较短，甚至会出现零。,同余法,该方法也是由选定的初始值出发，通过递推产生伪随机数序列。由于该递推公式可写成数论中的同余式，故称同余法。该方法的递推公式为,a，c，m分别称作倍数（multiplier）、增值（Increment）和模（modulus），均为正整数。x 称为种子或初值，也为正整数。该方法所产生伪随机数的质量，如周期的长度、独立性和均匀性都与

14、式中三个参数有关。,3 伪随机数的统计检验,伪随机数的特性好坏将直接影响到 Monte Carlo 的计算结果，因此要对所产生的伪随机数序列进行随机性检验。随机性检验主要包括均匀性检验、独立性检验、组合规律性检验和无连贯性检验。检验是伪随机数检验最常用的方法。,1.均匀性就是伪随机数列的 N 个数是否均匀分布在，区间上。若将，区间分成 k个相等的子区间（一般 k，），若所得伪随机数在，区间上是均匀分布的，则虚假设H 应为“每个伪随机数属于第 i组的概率为,频率检验 检验每组观测频数 ni与理论频数mi N/k之间相差的显著性,3.独立性,按先后顺序排列的 N 个伪随机数中，每个数的出现是否与其

15、前后各个数独立无关。对于两组伪随机数来说，独立性就是指它们不相关。,4.组合规律性,将 N 个伪随机数按一定的规律组合起来，则各种组合的出现具有一定的概率。,5.无连贯性,将一次出现的 N 个伪随机数，按其大小分为两类或 k类，则各类数的出现是否没有连贯现象。,确定性问题的 Monte Carlo 方法求解,Monte Carlo 方法所能解决的问题可以分为两大类,确定性问题随机性问题,确定性问题主要包括：求解线性和非线性方程组、逆矩阵、椭圆型差分方程的边值、积分方程以及多重积分等。,用 Monte Carlo 方法求解确定性问题的基本思想是：对于给定的确定性问题，设计一个概率统计模型（或随机

16、过程模型），然后采用一定的抽样方法，按照所设计的概率统计模型进行抽样，最后把这个模型产生的一个数字特征作为该确定性问题的近似解。,蒲丰试验,年法国著名学者蒲丰（Buffon）就发表了采用随机抽样法计算圆周率的论文。这就是蒲丰随机投针试验，即著名的蒲丰问题。,蒲丰概率模型是:在平面上画相距均为a的平行线束，向该平面上随机投置一枚长为 l 的针。为了避免针同时与两条平行线相交的复杂情况，设定 l a。如图所示，M 为针的中点，y为针的中点M 到与之最近的平行线的距离，为针与平行线的夹角。显然，y a，。该随机试验所有可能的集为,针与平行线相交的充分必要条件是 y lsin，针与平行线相交事件的集为

17、,则针与平行线相交的概率为,由于 l和 a均为已知常数，只要通过大量抽样试验求得该概率 p，由上式即可算出圆周率。设投针总次数为 N，其中针与平行线相交次数为 v，由贝努里（Bernoulli）定理可知，当 N充分大时，该事件出现的频率接近于其概率，即,这就是蒲丰试验求圆周率的过程。虽然，该方法要想获得较高的精确度，需要数以百万次的抽样试验，效率很低，但蒲丰试验具有 Monte Carlo 方法解决确定性问题的基本思想。,A Classic Example The Calculation of,随机性问题的 Monte Carlo 模拟,该过程是先建立一个随机过程模型，使得该过程的随机变量的数

18、学期望等于所要求解确定问题的解。这种计算方法叫做间接模拟方法。另一方面，采用随机试验的方法直接模拟随机过程，解决随机问题，这就是所谓直接 Monte Carlo 模拟。如模拟布朗运动、扩散过程、有机高分子形态、晶粒生长等随机问题的模拟。,1 随机行走（random walk）模拟,随机行走是一种典型的简单抽样方法，可用以模拟扩散、溶液中长而柔性的大分子的性质等.,随机行走主要有三种类型:,无限制的随机行走（Unrestricted random walk，RW）不退行走（Nonreversal random walk，NRRW）自回避行走（Self-avoiding walk，SAW）,无限制

19、随机行走就是指，某一个质点的每一次行走没有任何限制，既与前一次行走无关，也与以前任何一步所到的位置无关。这种模型可以用于模拟质点的扩散等过程，但是，不能用于模拟高分子的位形。因为，用随机行走方法模拟高分子位形是用随机行走的轨迹代表高分子的位形，行走过的位置代表的是构成分子的原子或官能团，因此，无限制随机行走忽略了体斥效应。不退行走就是禁止在每一步行走后立即倒退，可以解决刚走的一步与上一步重叠的问题。但不退行走没有完全解决高分子的体斥效应问题。自回避行走就是所有已走过的位置不能再走，这样就完全解决了体斥效应问题。,二维方格子上的三种随机行走的示意图。可以用四个矢量记述从某个节点向邻近节点的行走方

20、向（设格子间距为）：,N 步无限制随机行走的算法如下：,取 r（，）（坐标原点），并令 k；取一个在 和 之间的随机整数 vk；k k，rk rk（vk）；若 k N，行走结束，否则回到第 步。,对于不退行走，可选择的行走方向不再是，而是。禁止在每一步行走后立即倒退，即第k步的方向矢量不能与第 k 步的方向矢量相逆。由式（）可以看出方向矢量（），（）互为逆方向，（），（）互为逆方向。,2 Markov链,对于简单抽样，每一次的抽样都是独立的。如上述的随机行走过程，每行走一步都与前一步无关，更与初始位置无关。,Sampling,Integrate Over a Simple Shape?,Sta

21、tistical-mechanics integrals typically have significant contributions from miniscule regions of the integration space contributions come only when no spheres overlap e.g.,100 spheres at freezing the fraction is 10-260 Evaluation of integral is possible only if restricted to region important to integ

22、ralmust contend with complex shape of regionMC methods highly suited to“importance sampling”,Importance Sampling,Put more quadrature points in regions where integral receives its greatest contributionsReturn to 1-dimensional exampleMost contribution from region near x=1Choose quadrature pointsnot un

23、iformly,but accordingto distribution(x)linear form is one possibilityHow to revise the integral to remove the bias?,The Importance-Sampled Integral,Consider a rectangle-rule quadrature with unevenly spaced abscissasSpacing between pointsreciprocal of local number of points per unit lengthImportance-

24、sampled rectangle ruleSame formula for MC sampling,Greater p more points smaller spacing,choose x points according to p,在正则系综中，任意观察量 A（x）的热平均为,重要性抽样：在随机过程中，选取一个随机抽样的分布，使生成的随机数满足选取分布形式。根据一定的分布形式进行的随机抽样称为重要性抽样。,重要抽样 Monte Carlo 方法的实质是每次抽样试验不是完全独立的，而是与前一次或者与以前的所有抽样结果具有一定的概率关系，如不退随机行走和自回避随机行走。,设一个系统的状态序

25、列（随机变量序列）为 x，x，xn，如果对于任何一个状态 xn只与前一个状态 xn 有关，而与初始状态无关，即状态 n的概率为,则称此序列为 Markov 链。,Markov 链是一种随机行走状态，从状态 i单步行走到状态 j 的概率叫做转移概率，或跃迁概率，即,设所有可能的状态数为 N，由 pij 构成的 N N 矩阵叫做转移矩阵 p，该矩阵的每一行元素的和等于。Markov 链的重要性质是，无论初始状态如何，最终状态（足够多的时间步长次数）会遵从某一个唯一的分布，该分布叫做极限分布 xlim，即,也就是说，极限状态乘转移概率后状态不再发生变化，即系统达到一个平衡状态。因此，Markov 链

26、在平衡态 Monte Carlo 模拟中具有重要的意义。,3 Metropolis Monte Carlo法,Monte Carlo 方法主要分为简单随机抽样方法和重要随机抽样方法。简单抽样就是以平均分布进行抽样，每次抽样是完全独立的。正如前面关于积分问题中所述，很多问题难以用简单抽样方法解决，而重要随机抽样能够获得很好结果。,Metropolis Monte Carlo 方法是一种重要随机抽样方法。Metropolis 等人提出了一种基于建立一个 Markov 过程的方法。该方法的实质是，系统的各状态不是彼此独立无关地选取，而是建造一个 Markov 过程，过程中每一个状态 xi 是由前一个

27、状态 xi 通过一个适当的跃迁概率W（xi xi）得到的，并且，该概率能够使得在 M 的极限下，Markov 过程产生的状态的分布函数 P（xi）趋于所要的平衡分布，即,满足上述要求的充分条件为,也就是说，两个状态正向与反向的跃迁概率之比只依赖于两者的能量差 dH H（xj）H（xi）。但是满足该条件的跃迁概率 W 的形式并不是唯一的。通常采用以下两种形式：,ts 可取为,Metropolis Monte Carlo 方法的具体步骤如下:,Metropolis Monte Carlo（）建立体系状态与能量的关系模型。（）由初始状态出发，通过简单抽样设立新状态。（）根据新旧状态的哈密顿量 dH，

28、判断新状态的舍选，判断舍选有以下 种情况：dH，接受新状态，并在该状态基础上，进一步进 行步骤（）；dH，不是直接否决，而是进一步判断，抽取一个随机数 如果新状态被拒绝，则把原来的状态作为新状态，重复进行步骤（），并记录一次。,如果系统的粒子数为 M，每次新状态的抽样均随机抽选一个粒子，并不是每个粒子逐一地进行。只要伪随机数的质量足够高，各粒子被抽样的概率是均等的。当抽样次数达到系统粒子总数 M 时 该过程叫做一个 Monte Carlo step MCS,Metropolis 方法在状态抽样时虽然采用的是简单抽样，但是通过新旧状态的能量判断，实现新状态的舍选，建立 Markov 过程。对于恒定组成的正则和微正则系综，系统的能量用哈密顿量表示。对于变化组成的巨正则系综，随机选取一个粒子并通过改变粒子的种类得到一个新的组态，系统的能量用混合能及混合物化学位之和表示。,