物理问题的计算机模拟方法(2)—蒙特卡罗方法第三章 随机性模拟方法—蒙特卡罗方法(MC)
§ 3.1 预备知识
例:一个粒子在一个二维正方格点上跳跃运动
随机行走:每一时间步上,粒子可选择跳到四个最近邻格点上的任何一个,而记不得自己来自何方;
自回避行走:粒子记得自己来自什么地方,而回避同它自己的路径交叉。
随机行走的每一步的结果就是系统的一个状态,从一个状态到另一个状态的跃迁只依赖于出发的状态,这些状态形成一个序列,这就是一个马尔可夫链。
状态序列:x0, x1, …, xn, …
已给出状态x0, x1, …, xn+1 的确定值,xn出现...
第三章 随机性模拟方法—蒙特卡罗方法(MC)
§ 3.1 预备知识
例:一个粒子在一个二维正方格点上跳跃运动
随机行走:每一时间步上,粒子可选择跳到四个最近邻格点上的任何一个,而记不得自己来自何方;
自回避行走:粒子记得自己来自什么地方,而回避同它自己的路径交叉。
随机行走的每一步的结果就是系统的一个状态,从一个状态到另一个状态的跃迁只依赖于出发的状态,这些状态形成一个序列,这就是一个马尔可夫链。
状态序列:x0, x1, …, xn, …
已给出状态x0, x1, …, xn+1 的确定值,xn出现的概率叫做条件概率
马尔可夫链的定义:
如果序列x0, x1, …, xn, …对任何n都有
则此序列为一个马尔可夫链(或过程)。
§ 3.2 布朗动力学(BD)
1.郎之万方程
方程右边第一项为随机力,对粒子起加热作用;第二项为摩擦力,避免粒子过热。
将方程变形为:
于是,解可写为:
当随机力R(t)服从高斯分布时,上述方程的解描述的即为布朗运动,于是,布朗运动问题就化为在一些补充条件下求解郎之万方程,即
注:
表示随机力R在t和t’时刻没有关联, q为噪声强度。
为方差(R平均值为零), h为随机力关联时间,一般即为时间步长,在此时间内,粒子受到一个恒定地随机力作用。
由于R(t)服从高斯分布,上述方程的解v(t) 也服从高斯分布。
算法A7 (布朗动力学)
(1) 指定粒子的初始位置r0和初始速度v0;
(2) 从均值为零,方差为的高斯分布中选取一个随机数;
(3) 对速度积分得出
;(
)
(4) 在速度上加上随机分量。(??)
可将上述算法推广到含有规则力的系统。下面考虑粒子系与热浴耦合保持系统温度恒定的方法之一——系统粒子和虚拟粒子发生碰撞(正则细综控制温度恒定的方法之一)。
在系统的哈密顿运动方程中的动量方程里加上一个随机力就可实现上述目的:
算法A8 (随机碰撞)
(1) 按泊松分布
选取粒子发生碰撞的时间间隔(;
(2) 积分运动方程至发生碰撞为止(
);
(3) 若遭受碰撞的粒子为i, 则从温度为T的玻尔兹曼分布中随机选取一个动量pi (
(作为该粒子的新动量);
(4) 回到第2步。
§ 3.3 蒙特卡罗方法
1.什么是蒙特卡罗方法?
建立问题的概率模型
问题的解为模型的数学期望
抽样计算算术平均近似为数学期望
2.蒙特卡罗方法的数学基础
(1)大数定理
算术平均
数学期望I
(2)中心极限定理
收敛程度,估计误差
对给定的置信区间,误差界与方差(成正比,与抽样次数n的平方根成反比。
3.蒙特卡罗方法的基本步骤
算法A9
(1) 规定
;
(2) 产生
;
(3) 计算跃迁概率
;
(4) 产生随机数
;
(5) 若
,则
(不动),回到第二步;
(6) 若
,则接受
,回到第二步。
4.蒙特卡罗方法的局限性
(1) 有限系统大小的限制;
(2) 有限马尔科夫链长短的限制。
5.有关跃迁几率
的问题
系统处于状态x的几率:
( 配分函数
( 物理量的平均值
细致平衡条件:
定义:
对
的限制:
对
的限制:
(1)
, (2)
,(3)
§ 3.4 微正则系综蒙特卡罗方法
1.NVE恒定
配分函数:
系综平均:
中不含动能项,因此是一个由马尔可夫链产生的动力学。
允许的状态:
引入一个妖精(Spirit),具有能量
(允许的能量变化范围)
2.算法A10 NVE—MC方法
(1) 建立一个状态x, 使得
;
(2) 设定初始的ED (例如ED=1);
(3) 选定系统的一部分(例如某个粒子,随机选定);
(4) 改变系统的局部状态:
(例如粒子随机行走一步);
(5) 计算
(6) 若
,则
,
,返回第3步;
(7) 若
,则当
时,
,
,返回第3步;
当
时,返回第3步(不发生跃迁);
当系统达到平衡后,妖精的能量满足玻耳兹曼分布,即
因此,可以通过计算达到平衡后妖精的能量的平均值来确定系统的能量,从这个角度看,妖精像一支温度计,而系统为热库。
3. Ising 模型(d维格子)
广义体积:
, 自旋:
哈密顿量:
J 为交换耦合能,
J > 0, 此哈密顿量为铁磁体模型,自旋倾向于同向排列;
J < 0, 此哈密顿量为反铁磁体模型,自旋倾向于交错反向排列。
(为单个电子的自旋磁矩,B为外磁场。
磁化强度:
—向上自旋数,
—向下自旋数,
—自旋总数,
自发磁化(H=0), 二级相变,临界点Tc
T > Tc, m = 0 (磁化强度)
T < Tc, m ( 0
理论计算结果:
(1) d = 1 (一维)
(2) d = 2 (二维)
(3) d = 2 (三维)
模拟方法:
(1) 设定一个自旋位形
(全部向下);
(2) 随机选择一个格点, 使其自旋反转, 直到总能量达到设定值为止,此时的自旋位形
为模拟开始的初始位形;
(3) 设定一个妖精能量
,选定一个格点,其自旋反转的条件是:
,反转后,新的妖精能量为
;否则位形及妖精能量不变。直到访问完全部N个格点为止,完成一个MC步,产生一个新的位形;
(4) 设经过n0个MC步后系统达到热平衡,然后开始求平均。若进行了n个MC步,则平均磁化强度为:
Si为第步生成的自旋位形(总位形)
妖精的平均能量代表温度
(玻尔兹曼分布下求平均)
求和步长:
(最小能量改变)
因为每个格点有6个最近邻格点(三维立方晶格),当指定格点的自旋与所有最近邻格点的自旋相反时,该格点自旋反向后系统的能量改变为
当有五个最近邻与其相反时,能量改变为
当有四个最近邻与其相反时,能量改变为
当有三个近邻与其相反时,能量改变为
因此,
的最小能量改变4J。
即由
可定出温度T。
边界条件的处理:
周期性边界条件
(
为(方向的格点数)
当
时,
§ 3.5 正则系综蒙特卡罗方法
1. 跃迁几率
Z =
Metropolis 跃迁几率(通常的选择,不是唯一的选择)
2. 算法A11 正则系综MC方法
(1) 规定一个初始位形x;
(2) 产生一个新的位形
;
(3) 计算能量变化
;
(4) 若
,接受
(
)并回到第二步;
(5) 计算
;
(6) 产生一个随机数
;
(7) 若
,则接受
(
)并回到第二步;
(8) 否则,保持x并回到第二步。
若定义
则上述方法中的(3)~(8)步,就可以按算法A9的(3)~6)步进行,但注意
两者在编程上的区别,算法A11的效率更高。
3.Ising模型的模拟(计算程序PL4, 三维Ising模型)
讨论自发磁化,
单个自旋反转的跃迁几率:
Metropolis函数:
也可选择Glauber函数:
两种选择都满足细致平衡条件。
(时间标度),
,
进行一个MC步的Metropolis算法:
(1) 规定一个初始位形
;
(2) 选择一个格点i;(可以随机,也可不随机依次选)
(3) 计算Wi;
(4) 产生一个随机数
;
(5) 若
,则
;
(6) 否则转到第二步,直到作了N次尝试为止,得到一个新位形
,即完成一个MC步。
重复上述过程,直到模拟结束。然后按上一节类似的方法计算磁化强度m。
注意:(1)
的有效尺寸效应;
模拟结果依赖于格点大小L, 如图4.5,
附近,依赖关系明显,
依赖关系不明显;高温下,
,但在热力学极限下不存在自发磁化。
(2)+m(T) 与-m(T) 都有可能存在,因此要用绝对值求平均;
(3)与跃迁几率有关的困难。
当
或
时,
太小,系统运动缓慢,需要采取别的方法克服这一困难。
§ 3.6 等温等压系综MC方法
1. 跃迁几率
令
令
定义:
跃迁几率:
2. 算法A12 NPT蒙特卡罗方法
(1) 规定(、(;(即规定V, x)
(2) 随机产生
; (约化体积)
(3) 产生于
相容的
;(位置不能和体积相矛盾)
(4) 计算
;
(5) 若
,则
,返回(2);
(6) 否则,计算
;
(7) 产生随机数
;
(8) 若
,则
,返回(2);
(9) 否则,返回(2)。
§ 3.7 巨正则系综MC方法
1.跃迁几率
巨正则配份函数:
经典形式:
式中,r为自由度,对三维情况r = 3。
(绝对活度)
系统处在xN状态的几率:
产生一个粒子的几率:
消灭一个粒子的几率:
三个基本步骤:
(1) 粒子坐标方面的位形变化:Move — M
(2) 粒子的产生: Create — C
(3) 粒子的消灭: Destroy — D
2. 算法A13 巨正则系综MC方法I
(1) 指定V内N个粒子的初始位形
;
(2) 以相同的概率选择三个过程(M, C, D);
(3) Move
(3.1) 在V内选择一个粒子并随机移动它到
;
(3.2) 计算
;
(3.3) 若
,
,返回(2);
(3.4) 计算
;
(3.5) 产生随机数
;
(3.6) 若
,
,返回(2);
(3.7) 否则,x不变,返回(2)。
(4) Create
(4.1) 在V内为一个新粒子随机选择坐标;
(4.2) 计算
;
(4.3) 计算
;
(4.4) 若
,则接受
,返回(2);
(4.5) 产生随机数
;
(4.6) 若
,则接受
,返回(2);
(4.7) 否则,不接受,返回(2)。
(5) Destroy
(5.1) 在V内的N个粒子中随机选择一个;
(5.2) 计算
;
(5.3) 计算
;
(5.4) 若
,则接受
,返回(2);
(5.5) 产生随机数
;
(5.6) 若
,则接受
,返回(2);
(5.7) 否则,不接受,返回(2)。
3. 算法A14 巨正则系综MC方法II
基本思路:体积V内保持有衡定的粒子数M(V内能挤得下的数目), 实际系统的粒子数为N, 其余M-N个粒子为虚粒子或鬼粒子,为了保持系统的化学势(恒定, 引进粒子数的变化, 即增加一个粒子就把一个虚粒子变实,减少一个粒子就把一个实粒子变虚。
改写前面的积分:
(???)
1/VM 因子不影响平均值的计算, 可略去, 即上两式可写为
这样处理以后就不需要考虑粒子的产生和湮没。
算法A14
(1) 给体积V中的M个粒子指定坐标;
(2) 选择N个粒子为实粒子;
(3) 从N个实粒子中挑选一个粒子,坐标为x;
(4) 为所挑选的实粒子产生新坐标
;(Move)
(5) 计算能量变化
;
(6) 若
,接受这个位形并转到第11步;
(7) 计算
;
(8) 产生一个随机数
;
(9) 若
,接受这个位形并转到第11步;
(10) 否则,拒绝这个新位形,保持老位形;
(11) 从实粒子的集合中挑选一个准备删去的粒子;(Destroy)
(12) 计算能量变化
和
;
(13) 产生一个随机数
;
(14) 若
,把这个粒子放到虚粒子集合中去,并跳转到第16步;
(15) 否则,此粒子仍保留在实粒子中;
(16) 从虚粒子的集合中挑选一个准备改为实粒子的粒子;(Create)
(17) 计算能量变化
和
;
(18) 产生一个随机数
;
(19) 若
,把这个粒子加到实粒子集合中,并转到第4步;
(20) 否则拒绝这个粒子的产生并转到第4步。
优点:减少了不成功增加粒子试图的次数;
缺点:计算量增大,因为考虑的是M个而不是N个粒子的相互作用。
_1186780533.unknown
_1187512981.unknown
_1188555268.unknown
_1188725371.unknown
_1188728111.unknown
_1191138466.unknown
_1191140913.unknown
_1191143475.unknown
_1207497473.unknown
_1191140958.unknown
_1191140980.unknown
_1191140930.unknown
_1191138533.unknown
_1191139100.unknown
_1191139710.unknown
_1191139895.unknown
_1191138987.unknown
_1191137489.unknown
_1191138412.unknown
_1191138315.unknown
_1191137151.unknown
_1191137481.unknown
_1191136968.unknown
_1188728028.unknown
_1188728068.unknown
_1188728088.unknown
_1188728047.unknown
_1188725432.unknown
_1188728018.unknown
_1188725410.unknown
_1188724791.unknown
_1188724983.unknown
_1188725263.unknown
_1188725331.unknown
_1188725205.unknown
_1188724859.unknown
_1188724931.unknown
_1188724611.unknown
_1188724768.unknown
_1188724671.unknown
_1188555788.unknown
_1188724438.unknown
_1188555755.unknown
_1188462871.unknown
_1188549421.unknown
_1188554244.unknown
_1188554701.unknown
_1188555106.unknown
_1188554428.unknown
_1188550005.unknown
_1188554113.unknown
_1188549550.unknown
_1188463063.unknown
_1188463131.unknown
_1188463161.unknown
_1188462931.unknown
_1188462952.unknown
_1188462990.unknown
_1188462885.unknown
_1187944107.unknown
_1188461896.unknown
_1188462328.unknown
_1188462801.unknown
_1188462079.unknown
_1188461412.unknown
_1188461530.unknown
_1188461308.unknown
_1187942785.unknown
_1187942862.unknown
_1187943768.unknown
_1187942818.unknown
_1187513686.unknown
_1187942484.unknown
_1187513661.unknown
_1187513445.unknown
_1186842017.unknown
_1187509268.unknown
_1187510557.unknown
_1187512393.unknown
_1187512611.unknown
_1187511145.unknown
_1187510466.unknown
_1187510542.unknown
_1187509654.unknown
_1187509837.unknown
_1187510394.unknown
_1187509689.unknown
_1187509436.unknown
_1187507450.unknown
_1187508744.unknown
_1187509058.unknown
_1187509100.unknown
_1187508908.unknown
_1187508433.unknown
_1187508501.unknown
_1187507741.unknown
_1187507179.unknown
_1187507310.unknown
_1187507390.unknown
_1187507105.unknown
_1187507003.unknown
_1187507020.unknown
_1186842307.unknown
_1186781645.unknown
_1186782265.unknown
_1186841799.unknown
_1186841873.unknown
_1186782397.unknown
_1186782056.unknown
_1186782159.unknown
_1186781812.unknown
_1186780931.unknown
_1186781278.unknown
_1186781534.unknown
_1186781066.unknown
_1186780681.unknown
_1186780778.unknown
_1186780602.unknown
_1186733527.unknown
_1186769487.unknown
_1186772804.unknown
_1186774323.unknown
_1186779340.unknown
_1186779726.unknown
_1186774797.unknown
_1186773466.unknown
_1186773636.unknown
_1186774263.unknown
_1186773294.unknown
_1186772082.unknown
_1186772254.unknown
_1186772315.unknown
_1186769843.unknown
_1186769883.unknown
_1186769940.unknown
_1186769556.unknown
_1186735587.unknown
_1186735951.unknown
_1186736237.unknown
_1186736614.unknown
_1186769478.unknown
_1186736249.unknown
_1186736054.unknown
_1186736166.unknown
_1186735907.unknown
_1186735935.unknown
_1186735849.unknown
_1186734354.unknown
_1186734673.unknown
_1186735307.unknown
_1186734595.unknown
_1186734074.unknown
_1186734219.unknown
_1186733983.unknown
_1185371442.unknown
_1185784249.unknown
_1186732664.unknown
_1186733494.unknown
_1186733516.unknown
_1186732826.unknown
_1186733484.unknown
_1186732701.unknown
_1186732496.unknown
_1186732631.unknown
_1186732459.unknown
_1185783378.unknown
_1185783868.unknown
_1185784118.unknown
_1185783465.unknown
_1185783244.unknown
_1185783325.unknown
_1185372095.unknown
_1185783230.unknown
_1185345523.unknown
_1185351951.unknown
_1185352567.unknown
_1185371355.unknown
_1185352255.unknown
_1185348806.unknown
_1185350870.unknown
_1185346851.unknown
_1179843707.unknown
_1179844504.unknown
_1185094354.unknown
_1179844503.unknown
_1179841799.unknown
_1179843447.unknown
_1179841494.unknown
本文档为【物理问题的计算机模拟方法(2)—蒙特卡罗方法】,请使用软件OFFICE或WPS软件打开。作品中的文字与图均可以修改和编辑,
图片更改请在作品中右键图片并更换,文字修改请直接点击文字进行修改,也可以新增和删除文档中的内容。
该文档来自用户分享,如有侵权行为请发邮件ishare@vip.sina.com联系网站客服,我们会及时删除。
[版权声明] 本站所有资料为用户分享产生,若发现您的权利被侵害,请联系客服邮件isharekefu@iask.cn,我们尽快处理。
本作品所展示的图片、画像、字体、音乐的版权可能需版权方额外授权,请谨慎使用。
网站提供的党政主题相关内容(国旗、国徽、党徽..)目的在于配合国家政策宣传,仅限个人学习分享使用,禁止用于任何广告和商用目的。