分子模拟教程讲解.ppt

Monte Carlo方法的应用举例： 计算积分： 常用的积分方法求解： 将积分区域[a,b]均匀地划分成N各分区间，则积分结果可近似地表示成： Δx = (b-a)/N 简单的Monte Carlo积分方法求解： 利用均匀分布的随机数发生器，从[a,b]区间产生一系列随机数xi，i=1, 2, ..., N 其中 X为均匀分布，并且 X?[a,b] 近似求解E[g(X)]: 近似求解积分: 随机抽样 当我们用简单Monte Carlo计算积分时，若该函数为常数函数，g(x)=constant，则取样数不管多少，准确度为100％。 如果在积分区间内，g(x)为一平滑函数，则简单Monte Carlo方法较为准确，反之，如果g(x)的变动很剧烈，则简单Monte Carlo方法的误差会变大。 说明： 重要性Monte Carlo抽样方法 在 g(x) 变化剧烈时，如果以Monte Carlo方法取样，最好依据g(x)的大小来决定取样率。 当|g(x)|的值较大时，对∫g(x)dx的贡献也较大，如果没被选中，则结果的误差极大。 解决方式：改变 x被选中的机率，让|g(x)| 值较大的点被选中的机率增加。 采用权重分布函数(Weight distribution function) w(x) :决定每个x被选中的机率。 重要性抽样的定义：根据一定的分布形式进行的随机抽样。 w(x)必须归一化，即在积分区间内∫w(x)dx=1。 由于 x 的选取已被 w(x) 扭曲，所以计算积分时要把这部分[还]回去：若一共取样了N个x，则积分值为： 重要性Monte Carlo抽样方法 Metropolis Monte Carlo方法 我们所模拟的系统最终要达到的平衡分布是Boltzman分布： Boltzmann 概率分布函数： 我们如果能够产生这种分布，我们就能够计算系统的大多数性质，但这是不可能的，因为我们不知道Z的值，但是对于任意两个状态，我们有： 可以在相空间中构造一个马尔科夫链，使相空间中的样本点随着链的增长逐步趋近于Boltzman分布。 一个序列x0, x1, x2, …,xn,如果对任何n都有： 则此序列是一个Markov链。 要求： 任何一次实验的结果依赖于前一次的试验，并且近依赖于前一次的试验。 马尔科夫(Markov)链 可以证明：通过构造Markov链，体系中最终的平衡分布就是Boltzman分布 Metropolis平衡条件 (Detailed balance condition)： 平衡条件： 系统处于状态X的概率正比于其Boltzman因子： 如果?是对称的: Metropolis Monte Carlo方法的算法： 给出一个初始状态，并计算系统的能量Eold 随机产生一个新状态，并计算新系统的能量Enew 如果ΔE＝(Enew?Eold)0, 则接受新状态并回到 b) 如果ΔE＝(Enew?Eold)0, 则计算Boltzman因子： 在(0, 1)区间上产生一个均匀分布的随机数? ； 如果 则接受新状态并回到b) 否则保留原值并回到b) 1. 正则系综蒙特卡罗模拟方法(Canonical MC Simulation) 具有确定的粒子数N、温度T和体积V 对于含N个粒子的系统，位型（构型）的配分函数： 某个特定构型的发生概率为 PNVT(rN) 1. 正则系综蒙特卡罗模拟方法(Canonical MC Simulation) Monte Carlo模拟中任一物理量的计算： 位型积分 概率密度 系统处于位型{rN}的概率密度 给出一个初始状态，并计算系统的能量Uold 随机产生一个新状态，并计算新系统的能量Unew 如果ΔU＝(Unew?Uold)0, 则接受新状态并回到 b) 如果ΔU＝(Unew?Uold)0, 则计算Boltzman因子： 在(0, 1)区间上产生一个均匀分布的随机数? ； 如果 则接受新状态并回到b) 否则保留原值并回到b) 正则系综MC模拟算法的组织： 正则系综MC模拟算法的流程图： 给定每个分子的初始位置，ri(0) 随机选取一个分子，并随机移动到新的位置 计算移动前后的系统能量变化ΔU 拒绝移动 ΔU?0 ? Exp(ΔU)?? (0,1)? 统计系统的热力学性质及其它物理量 统计性质不变？ 打印结果，结束 Yes No No No 接受移动 Yes Yes 大约循环107到108次 Monte Carlo模拟中几个热力学量的计算： N个粒子系统中的总势能： 假设采用截断势能函数： Uc：截断范围内的总势能； Ulrc：截断半径外对势能的长程校正(Long-range correction) 对于LJ流体： 含N个粒子系统中的压力： Wc：截断范围内的总维里项(Viria