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

随机数产生的办法
关于随机数的几点注意
注1 由于均匀分布的随机数的产生总是采用某个确定 的模型进行的，从理论上讲，总会有周期现象出现的。 初值确定后，所有随机数也随之确定，并不满足真正 随机数的要求。因此通常把由数学方法产生的随机数 成为伪随机数。 但其周期又相当长，在实际应用中几乎不可能出 现。因此，这种由计算机产生的伪随机数可以当作真 正的随机数来处理。 注2 应对所产生的伪随机数作各种统计检验，如独 立性检验，分布检验，功率谱检验等等。
考虑平面上的一个边长为1的正方形及其内部的一个形状不规 则的“图形”，如何求出这个“图形”的面积呢？Monte Carlo方法是这样一种“随机化”的方法：向该正方形“随机 地”投掷N个点，若有M个点落于“图形”内，则该“图形” 的面积近似为M/N。
用该方法计算π的基本思路是： 1 、根据圆面积的公式： s=πR^2 ,当R=1时，
11
面积的计算
辛普逊方法
蒙特-卡洛方法
在长方形中均匀投N0组(x,y) 如 y<f(x), 则 N=N+1
I = ΣSn
f (x)
Hale Waihona Puke I =(N/N0)×S0f (x)
S0
S
x x
MC 的优点 MC与传统数学方法相比，具有直观性强，简便易行的优点，该方
法能处理一些其他方法无法解决的负责问题，并且容易在计算机 上实现，在很大程度上可以代替许多大型、难以实现的复杂实验 和社会行为。无污染、无危险、能摆脱实验误差。
Monte Carlo方法之随机数的产生
许多计算机系统都有随机数生成函数 F90: call random_seed
call random_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。
Monte Carlo方法之典型算法与应用
实例一、计算π值
计算过程： 1、构造或描述问题的概率过程 2、从概率密度函数出发进行随机抽样，实 现从已知概率分布的抽样，得到特征量的一 些模拟结果——计算均值
实验者 年代 投掷次数 相交次数 圆周率估计值 沃尔夫 1850 5000 2531 3.1596 史密斯 1855 3204 1219 3.1554 德摩根 1880 600 383 3.137 福克斯 1884 1030 489 3.1595 拉泽里尼 1901 3408 1808 3.1415929 赖纳 1925 2520 859 3.1795 布丰投针实验是第一个用几何形式表达概率问题的 例子，他首次使用随机实验处理确定性数学问题， 为概率论和蒙特卡罗方法的发展起到一定的推动作 用。
S=π。 2、由于圆的方程是：x^2+y^2=1(x^2为x的平方的 意思),因此1/4圆面积为x轴、y轴和上述方程所包围 的部分。 3、如果在1*1的正方形中均匀地落入随机点，则落 入1/4圆中的点的概率就是1/4圆的面积。其4倍， 就是圆面积。由于半径为1，该面积的值为π的值。
REAL R,R1,R2,PI ISEED=RTC() N0=0 N=300000 DO I=1,N R1=RAN(ISEED) R2=RAN(ISEED) R=SQRT(R1*R1+R2*R2) IF(R<1.0)N0=N0+1 END DO PI=4.0*N0/N WRITE(*,*)PI END
y xi 2
i 1 n
REAL Y Y=0 N=300000 ISEED=RTC() DO I=1,N X=RAN(ISEED) Y=Y+X**2/N END DO WRITE(*,*)Y END
lim∑x2dx (dx 0)
Monte Carlo方法另一个重要问题：随机数
随机数：由单位矩阵分布中所产生的简单子样称为随 机数序列，其中的每一个个体称为随机数。 但真正的随机数的不适合电子计算机上使用，因 为它需要很大的存储量。利用某些物理现象可以在 电子计算机上产生随机数，且其产生的序列无法重 复实现，使程序无法复算，结果无法验证，同时需 要增添随机数发生器和电路联系等附加设备。
注意以下两点： Monte Carlo方法利用随机抽样的方法来求解物理问题; 数值解法:从一个物理系统的数学模型出发,通过求解一 系列的微分方程来的导出系统的未知状态; Monte Carlo方法并非只能用来解决包含随机过程的问题: 例如:用Monte Carlo方法计算定积分. 对这样的问题可将其转换成相关的随机过程, 然后用Monte Carlo方法进行求解
MC的应用
自然现象的模拟： 宇宙射线在地球大气中的传输过程； 高能物理实验中的核相互作用过程； 数值分析： 数学问题，求积分，求逆矩阵，解线性代数等 经济学模拟：库存问题，随机服务系统中排队问题 人口问题：人口的出生，传染病的蔓延；乃至动物的生态竞争
金属学：扩散、组织长大、相 变过程
蒙特-卡洛模拟的意义
1) 取一张白纸，在上面画上许多条间距为d的平行线。 2) 取一根长度为l（l<d） 的针，随机地向画有平行直线的纸上 掷n次，观察针与直线相交的次数，记为m 3）计算针与直线相交的概率． 布丰本人证明了，这个概率是: p=2l/(πd) ，π为圆周率 ： 利用这个公式可以用概率的方法得到圆周率的近似值。下面是 一些资料
Monte-Carlo模拟，即随机模拟（重复“试验”）
计算机模拟
Monte Carlo方法： 亦称统计模拟方法，statistical simulation method 利用随机数进行数值模拟的方法 Monte Carlo名字的由来： • 是由Metropolis在二次世界大战期间提出的：Manhattan 计划，研究与原子弹有关的中子输运过程； • Monte Carlo是摩纳哥（monaco)的首都，该城以赌博闻名
计算材料学概述
第三章
蒙特卡罗方法 （Monte Carlo）
主要内容
Monte Carlo模拟发展简介
Monte Carlo模拟基本原理
Monte Carlo模拟典型算法 Monte Carlo模拟典型应用
蒙特卡洛法是什么？
蒙特卡洛(Monte Carlo)方法，是在简单的理论 准则基础上，采用反复随即抽样的方法，解决复杂 系统的问题。其实质是一种概率和统计的问题。 蒙特· 卡罗方法（Monte Carlo method），也 称统计模拟方法，是二十世纪四十年代中期由于科 学技术的发展和电子计算机的发明，而被提出的一 种以概率统计理论为指导的一类非常重要的数值计 算方法。是指使用随机数（或更常见的伪随机数） 来解决很多计算问题的方法。
Nicholas Metropolis (1915-1999)
Monte-Carlo, Monaco
Monte Carlo方法简史
简单地介绍一下Monte Carlo方法的发展历史 1、Buffon投针实验： 18世纪，法国数学家Comte de Buffon利用投针实验估计的值
L
d
p
d
2L
s= rand(‘state’,0) 重设该生成函数到初始状态 注意：上述随机数序列均具周期性，如上页 random子程序的周期 约230。
Finite difference approximation of differential equations
A differential equation can be approximated by a finite difference scheme. For example t t t
能研究不同边界、不同材料的影响 理论不可能、实验耗费太大 用于实验设计 无污染 反应堆防护 核弹爆炸 能摆脱实验误差 作理论和实验的桥梁
Monte Carlo模拟的步骤： 1. 根据欲研究的物理系统的性质，建立能够描述该 系统特性的理论模型，导出该模型的某些特征量 的概率密度函数； （即构造或描述问题的概率过程） 2. 从概率密度函数出发进行随机抽样，实现从已知 概率分布的抽样，得到特征量的一些模拟结果； 有了明确的概率过程后，为了实现过程的数 值模拟，必须实现从已知概率分布的随机数的抽 样，进行大量的随机模拟实验，从中获得随机变 量的大量试验值。产生已知概率分布的随机变量， 是实现MC方法的关键步骤，其中最基本的是（0， 1）均匀分布。 4. 模拟结果的检验
3. 对模拟结果进行分析总结，预言物理系统的某些特性。
Monte Carlo算法的主要组成部分
概率密度函数(pdf)— 必须给出描述一个物理系统的一组概 率密度函数; 随机数产生器—能够产生在区间[0,1]上（均匀）分布的随机数 抽样规则—如何从在区间[0,1]上均匀分布的随机数出发,随 机抽取服从给定的pdf的随机变量; 模拟结果记录—记录一些感兴趣的量的模拟结果
du u u dx t
Forward time
du u tx u txx dx x d 2u u txx 2u tx u txx 2 dx x 2
Backward space
Central space
t 1 t t t t ux ux 2 u u u 2u u x x x x 2 x t x t x 2
Monte Carlo方法之随机数的产生
许多计算机系统都有随机数生成函数 F90: call random_seed
call random_number(a) 2、ISEED=RTC()
X=RAN(ISEED) Y=RAN(ISEED）
Matlab: x=rand(N) 产生元素在(0, 1)间随机分布的N*N矩阵
f x dx
b a
! 计算x**2在(0,1)上积分 计算过程： 1、构造或描述问题的概率过程:产生服从分 布f（x）的随机变量Xi( n ) （i=1,2, ···,N） 2 xi 2、从概率密度函数出发进行随机抽样，实现 i 1 从已知概率分布的抽样，得到特征量的一些 模拟结果——计算均值( )
伪随机数：
是有数学递推公式所产生的随机数。（近似 的具备随机数的性质。） An+1=T(A);An+1= An+k+1 伪随机的优点和缺点： 判断伪随机数好坏的方法： 1、它能够有较好的均匀性和独立性； 2、它的费用大小，即指所消耗计算机的时间； 3、容量要求尽可能大。
产生均匀分布随机数的几种方法； （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 random
MC的基本思想
MC基本思想很早以前就被人们所发现和利用。 17世纪，人们就知道用事件发生的“频率”来决定 事件的“概率”。但要真正实现随机抽样是很困难 的，甚至几乎是不可能的。 高速计算机的出现，使得用数学方法在计算机 上大量、快速地模拟这样的试验成为可能。
确定性系统
模拟 随机性系统 自然界
重复试验