《化学统计热力学》课件Molecular Dynamics

1、分子动力学方法Molecular Dynamics Method7/10/20221Molecular DynamicsMolecular Dynamics与Monte Carlo模拟的区别：Monte Carlo模拟中，分子坐标的改变是完全随机的，跟分子所受的力毫无关系分子力学中，根据力的大小和方向来决定移动的方向7/10/20222Molecular DynamicsMolecular Dynamics动力学与热力学的关系：热力学研究的是宏观平衡态的性质动力学研究的是态与态之间的转化的过程7/10/20223Molecular DynamicsMolecular Dynamics中心问题：

2、解牛顿运动方程，得到分子坐标和速度是如何随时间而改变的轨迹：trajectoryN个粒子，这样的二阶微分方程就有N个7/10/20224Molecular Dynamics1. MD的简单模型硬球模型(Alder and Wainwright in 1957)：粒子作匀速直线运动，直到它碰撞到其它粒子粒子间的相互作用只有在粒子的相互距离小于其半径之和时才存在：Udistance7/10/20225Molecular DynamicsMD的简单模型，continue硬球模型粒子间的相互作用的小改进：Udistance7/10/20226Molecular DynamicsMD的简单模型，cont

3、inue硬球模型的模拟步骤：按速率分布随机给出初始的位置和速度计算碰撞时的位置根据动量守恒定律计算碰撞后的新速度Analyze Output是否继续计算下一次碰撞时间7/10/20227Molecular Dynamics2. 使用连续势能函数的MD作用在分子或粒子上的力其实是时刻在改变的无法求解N个粒子的运动方程，因此采用其它法求解运动方程F(t)的分析形式很复杂物理学到目前为止还为求出三体问题的精确分析解i=1,.,N，共N个二阶微分方程!7/10/20228Molecular Dynamics解决办法：有限差分方法时间t的状态（速度，位置）是通过把0-t的积分分为很小的时间间隔t的运动的

4、加和，假定粒子在这段很短的时间间隔之内受到的力是常数ttr(t)r(t)r(t+t)7/10/20229Molecular Dynamics有限差分方法(Finite Difference Methods)假定粒子在某时刻附近的位置和速度可以用Taylor阶数展开：7/10/202210Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法：利用t-t的位置r(t-t)和t时刻的位置r(t)和加速度a(t)计算t+t的位置r(t+t)7/10/202211Molecular Dynamics有限差分方法(Finite Differ

5、ence Methods)Verlet算法速度r(t)r(t+t)r(t-t)t-tt+ttr(t+t)-r(t-t)7/10/202212Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法速度也可以计算half-step (t+1/2t)的速度7/10/202213Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法的优点存储量小：r(t-t)，r(t)，a(t)Verlet算法的缺点精度初始位置r(0)和r(-t)：7/10/202214Molecular

6、Dynamics有限差分方法(Finite Difference Methods)Verlet算法的改进1leap-frog (Hockney, 1970)，蛙跳蛙跳的由来7/10/202215Molecular Dynamics有限差分方法(Finite Difference Methods)蛙跳算法的优点直接包括速度，不再对较大的量(r(t)做差分7/10/202216Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法的改进2velocity Verlet (Swope, 1982)速度也可以这样算7/10/202217M

7、olecular Dynamics三种算法的比较a. Verlet; b. Leap flog; c. Velocity Verlet7/10/202218Molecular Dynamics有限差分方法(Finite Difference Methods)Verlet算法的改进3Beemans 算法(Beeman, 1976)7/10/202219Molecular DynamicsPredictor-Corrector integration methods (Gear, 1971)利用Taylor展开式计算下一时刻的位置，速度，加速度及其它高级项把实际计算ac(t+t)的与预测的a(t+

8、t)作比较7/10/202220Molecular DynamicsPredictor-Corrector integration methods (Gear, 1971)Gear 给出了一套他认为最好的系数这些系数Taylor阶数的截断阶数有关，这里（三阶后截断）7/10/202221Molecular DynamicsMD模拟的一些实际考虑该选择什么方法上述算法的计算时间并不是最需要考虑的，最耗时间的步骤其实是计算能量，力和加速度需要保持能量和动量守恒不要期望MD得到的运动轨迹和真实的轨迹非常一致甚至很接近7/10/202222Molecular Dynamics实际的一些考虑，conti

9、nue时间的步长(t)太小：只探索了相空间的一小部分太大：能量的涨落很大，能量和动量不守恒时间模拟-真实10fs50fs100fs两个Ar原子间的碰撞7/10/202223Molecular Dynamics实际的一些考虑，continue时间的步长(t)的设定SystemTypes of MotionSuggested time step (fs)Atomstranslation10Rigid moleculestranslation, rotation5Flexible molecules, rigid bondstranslation, rotation, torsion2Flexibl

10、e molecules, flexible bondstranslation, rotation, torsion, vibration1 or 0.57/10/202224Molecular Dynamics7.2-3 如何开始和进行分子动力学模拟初始位置随机确定初始动量Maxwell Boltzmann分布这是一个Gaussian分布，可以用7/10/202225Molecular Dynamics7.2-3 如何开始和进行分子动力学模拟初始动量x使用生成随机数的方法产生：这样产生的x，它的平均值为零，且2=1，其它平均值不为零， 21的x这样得到：注意动量守恒和能量守恒7/10/2022

11、26Molecular Dynamics7.2-3 如何开始和进行分子动力学模拟计算相互作用，力和加速度例：Lennard-Jones potential运行一段时间，达到平衡能量涨落变小温度基本恒定速度和动量分布符合Maxwell-Boltzmann分布7/10/202227Molecular Dynamics7.2-3 如何开始和进行分子动力学模拟计算热力学函数温度Nc：对体系的限制数目对孤立体系，1. 体系总的平动为常数（守恒）；2. 总的角动量也守恒，所以Nc = 6.对于周期性的体系，角动量不守恒， Nc = 3.3N-Nc：总的平动自由度7/10/202228Molecular D

12、ynamics7.2-3 如何开始和进行分子动力学模拟计算热力学函数其它热力学量E(pi,qi)：每一步的体系能量，包括动能和势能G(pi,qi)：要求平均值的物理量N：MD模拟所进行的步数7/10/202229Molecular Dynamics7.2-4 受约束条件限制的动力学模拟模拟实际体系的困难：刚性分子相互作用于分子的取向有关柔性分子，质心随分子的构象改变；相互作用同构象有关；要考虑振动，模拟的时间步长取决于数值最大的振动频率。SHAKE Procedure动力学模拟中固定键长一个完全限制(holonomic)的方法：7/10/202230Molecular DynamicsSHAK

13、E方法Constraint的实质就是对分子中的原子加上了额外的力对键长固定，即两原子间的距离在任意时刻满足7/10/202231Molecular DynamicsSHAKE方法代入Verlet算法中下一时刻位置矢量的表达式如何求出k呢？7/10/202232Molecular DynamicsSHAKE方法根据限制条件，共k个SHAKE Procedure:迭代的方法求解，先猜测k-1个的值，根据两个原子间距离(|r12|)的限制求出一个的值(12)，再假定一个未知(13)，代入12和其它k-2个猜测的值，求出13，如此循环往复，直至收敛。7/10/202233Molecular Dynam

14、icsSHAKE方法SHAKE Procedure:猜测13, 14,.,1k根据|r13(t)|=d13，代入12 , 14,.,1k，求13求出所有其它步骤收敛与否？根据|r12(t)|=d12求127/10/202234Molecular DynamicsSHAKE方法对键角的固定：转化为对距离的固定7/10/202235Molecular Dynamics7.2-5 恒温和恒压下的动力学模拟上面的MD没有提温度和压力，所以是在NVE(NVEP)系综下进行的。固定温度7/10/202236Molecular Dynamics控制温度因此只要在每一时刻（模拟的每一步）将这点的动量（或速度）

15、乘以一个因子缺点：hot solvent, cold solute7/10/202237Molecular Dynamics控制温度的改进方法：将体系与外部恒定温度的热源接触，因此允许体系的温度有波动。参考：Aderson, H. C. Molecular Dynamics Simulations at Constant Pressure and/or Temperature: J. Chem. Phys. 1980, 72, 2384-93.Nos, S. Molecular Dynamics Method for Simulations in the Canonical Ensemble:

16、 Mol. Phys. 1984, 53, 255-268.Hoover, W. G. Canonical Dynamics: Equilibrium Phase-Space Distributions: Phys. Rev. A 1985, 31, 1695-1697.7/10/202238Molecular Dynamics固定压力保持压力恒定的方式通过改变体积来实现的。把体系的体积乘以一个恰当的因子把体系同一个外部的体系相作用，这个外部体系被作为“Pressure Bath”p：耦合因子把体积乘以一个因子（相应的坐标乘以1/3）：物质的压缩率7/10/202239Molecular Dy

17、namics7.2-5 MD和MC的应用构象分析系统搜索（手动，工作量大，陷入局部最低点）普通分子动力学（自动，陷入局部最低点）模拟退伙算法（自动，全局最优）Kirkpatrick, S.; Gelatt, C. D.; Vecchi, M. P. Optimization by Simulated Annealing: Science 1983, 220, 671-680.遗传算法（自动，全局最优）7/10/202240Molecular Dynamics模拟退火算法在高温进行Monte-Carlo模拟一定的步数改变（降低）温度进行Monte-Carlo模拟7/10/202241Molecular Dynamics遗传算法这已经不是动力学模拟了产生一些代表分子不同构象的population010111010011123000110111011123基因7/10/202242Molecular Dynamics遗传算法将基因解码为构型参数，并计算其适应度(fitnes