分子模拟第七、八章资料下载.pdf
- 文档编号:16120239
- 上传时间:2022-11-20
- 格式:PDF
- 页数:62
- 大小:3.51MB
分子模拟第七、八章资料下载.pdf
《分子模拟第七、八章资料下载.pdf》由会员分享,可在线阅读,更多相关《分子模拟第七、八章资料下载.pdf(62页珍藏版)》请在冰豆网上搜索。
即可以由t+t和t-t的位置得到t时刻的速度。
缺点:
此式含有1/t,因为t很小,因此非常容易导致误差,为了矫正这个缺点,Verlet又发展出了Leapfrog法(跳蛙法).对两式分别求导得到:
22)(21)()()()(21)()()(ttrttatvttvttrttrtrttr+=+=+22)(21)()()()(21)()()(ttrttatvttvttrttrtrttr+=+=ttattvttv)
(2)()(+=+tt替换用21对位置矢量的表达式重新整理得到:
)21()(21)()()(21)()()(22ttvttrtatvttrttrdtdtrdtdttrttr+=+=+=+以在x轴上运动的粒子为例:
其初始位置为x1(0)=0,x2(0)=2,初始速度为1=110-3,2=-0.510-3,设粒子质量均为1,积分间隔t=0.01,其VDW势能参数为:
=1,=1.任一瞬间,两个粒子间的距离为:
所受的力分别为:
运行20步后粒子的位置、速度和能量。
还有另外一种方法,即预测矫正法:
以一维运动为例,位置矢量的二阶导数可以表示为:
预测t+h时刻时的位置的表达式P(x)为:
与泰勒展开式对比:
对于给定的k值,这个式子给出和泰勒展开式同样的效果。
其中的系数满足下式:
预测x+h时刻的速率表达式为:
在计算x(t+h)和x(t+h)的预测值时为了修正误差,又引入修正函数:
其中的系数和满足:
大多数的动力学问题可以被表达为二阶微分方程,但也有一阶微分方程:
其预测函数和修正函数为:
系数为几种数值方法的比较:
1.Velert:
最快2.Leapfrog:
更少的工作时间和内存,能量收敛快3.Predictorcorrector:
精度高,变量多,占内存实例:
观察不同积分方法下能量的变化情况7.3积分步长的选取分子动力学中,最重要的工作为如何选取合适的积分步长,在节省时间的同时也保证计算的精确性。
原则:
积分步长小于系统中最快运动周期的十分之一。
htt以水分子为例,较快的运动为分子内部的振动。
通过红外光谱可以得到其最大的振动频率约为1.81014/s,因此其振动周期为0.9210-14s,因此对水分子计算的积分步长约为0.9210-15s即0.92fs。
以氩原子的分子动力计算为例,氩原子的质量为39.95原子质量单位,其LJ12-6势的参数为=3.504,=0.24kcal/mol,求二阶导数得到:
势能最低点的二次微分相当于力常数k,因为势能为极小值时r=,带入后得到:
此时根据简谐振子的振动频率公式知道:
其中u为简化质量:
可求出频率为:
其周期为:
因此积分步长可取为:
即1.3ns和前面的水分子动力学计算相比,氩原子的积分步长为它的1000倍,因此可以研究更长时间内的氩原子运动。
如果研究10-9s内的运动,对于氩原子需要约800步,而水分子则需要1.05106步,而每步需要的时间大体类似,因此氩的整体模拟时间要远小于水分子体系。
利用LJ12-6势能估计原子系统积分步长的方法也可应用于其它多原子系统,具体可以取势能参数中最大的来估计原子间的最快运动。
例如水分子中的LJ作用有O-O、O-H和H-H,因为O-O之间LJ势能的最大,所以可以用来估计最大积分步长。
实例:
不同步长对结果(能量)的影响采用Verlet算法,分别在0.0001ps、0.0005ps、0.001ps、0.005ps、0.01ps的步长下计算晶体的能量。
7.4简化单位研究分子或原子系统时,如果采用国际单位制,原子质量以g为单位,则通常的原子质量约为10-22g级别;
若位置以cm为单位,则通常的量纲为10-8cm;
同样积分步长用s做单位通常在10-1310-16s。
这些量纲非常小,实验中很容易引起误差,因此实际计算时通常采用简化单位。
其中,和为LJ函数的参数值。
原子的质量单位简化为原子的质量数,长度单位简化为。
例如:
氩(39.5)的熔点为84K,=0.24kcal/mol,=3.504,系统的密度为1.784g/cm3,设计算时的步长10-12s,压力1atm,则各物理量可以简化为:
Note:
如果系统中含有不止一种原子,则选取其中一种原子的LJ参数作为简化单位的标准,例如水分子中有H和O两种原子,可以选取O的LJ参数进行简化。
7.5分子动力计算流程系统容量在进行分子动力学计算前,必须先估计运算的可行性,一般分子动力学计算的容量为数千个原子级别,如果系统过大则超出计算的范围。
对于N个原子的系统,每计算一步需要计算1/2N(N+1)组远程作用力,这是运算过程中最耗时的部分,若原子增加一倍,计算时间则为1/22N(2N+1)为原来的四倍多。
研究对象所处的时间范围通常分子动力研究所选取的积分步长为飞秒fs(1fs=10-15s),若以目前一般的个人工作站或较强的个人电脑从事1000个原子系统的计算,累积100万步即研究10-9s(1ns)的时间范围,需要两星期的时间?
。
因此从实际的角度来讲,分子动力学适合研究反应或运动时间小于1ns的体系,而不适合较慢的反应或运动。
例如蛋白质折叠在10-3s(1ms)级别,则需要非常长的时间。
计算过程执行分子动力学计算时,将一定数目的分子放在一定形状的盒子中,并使它的密度和实验密度相符合,再选定实验的温度,即可以着手计算。
通常是将分子随机的放置在盒子中,或是按照其结晶位置排列,这将作为分子动力模拟的初始位置。
计算开始时,还必须知道原子的初始速度。
系统中所有原子速度产生的动能总和应当满足玻尔兹曼分布:
式中N为原子数目,T为热力学温度,kB为玻尔兹曼常数。
原子速度可以按此式随机产生,而总动能符合3/2NTkB即可。
由初始位置和速度开始,计算每一步新的速度和坐标,再由新的速度计算系统的温度:
如果得到的速度过大或过小,则会使得到的温度大于或小于设定的温度很多,这时需要校正速度,一般所允许的温度范围为:
如果计算的温度超过允许范围,则将所有原子的速度乘以一个校正因子:
使系统温度重新调整为设定温度:
实际的计算中,开始时每隔几步就需要校正一次,随着计算的进行,校正间隔逐渐增加到几千步,一直到原子速度不需要校正,而系统动能在3/2NTkB上下10%浮动,此时系统已经达到了热平衡。
在到达热平衡之前的轨迹信息和速度不需要保存,因为那时的物理意义不严格,真正对我们有意义的是达到平衡后的系统信息。
通常的分子动力学模拟需要累积数百万步以供分析,因而会产生极大的存储。
以1000个原子为例,储存每步的位置和速度约要10kB空间,则10000步需要100MB。
因为分子动力学计算的步长很短,每一步移动的距离也很小,通常每隔1020步存储一次来节省硬盘空间。
NVE系综的一般性分子动力计算7.6分子动力计算的初始设定执行分子动力学计算时,必须选择恰当的初始条件,如恰当的初始温度、初始速度、初始位置和积分步长等,如果初始条件选择不当则会浪费很长的时间才能达到热平衡,或者说根本达不到平衡。
初始结构初始结构的选取以越接近模拟体系的理想结构越好。
如要模拟体心立方的固体系统时选取了面心立方作为初始结构,便是不明智之举。
通常在模拟前作一些实验,降低系统的局部高能区,使其接近于能量最低点的理想结构,避免计算中由于偏离理想结构太多而导致失败。
对于均匀相的液体,通常选用它的固相晶体结构作为初始模型,如果有XRD实验得到的结构数据,则可以直接采用;
如果没有,通常采用均匀分布的面心立方结构排列这些分子或原子。
对于分子系统,通常将分子的质心放在面心立方的格点上,分子的取向视系统而定。
对于小分子可以随意指定其方向,如水和二氧化碳分子等;
对于大分子则考虑取向时以不发生重叠为宜。
对于复杂液体,通常参考XRD或NMR实验得到的结构位置,或者采用分子力学能量最小化的结果作为起点。
初始速度通常系统初始速度的产生,是先在-1+1之间选取一个呈高斯分布的随机数,再将其乘以系统所有粒子的平均速度:
即可以得到符合麦克斯韦-玻尔兹曼分布的粒子速度:
通常在实验前还要检查各方向动量是否为0,避免计算过程中系统移动而导致能量不稳定。
THEEND分子模拟牛继南2011.3第八章分子动力学原理-II回顾下我们上节课所讲的内容:
一般性分子动力学计算的步骤8.1控温方法Nos-Hoover控温法牛顿运动方程被修改为:
其中(t)为摩擦系数,由下式定义:
为实时温度,为指定的时间常数(通常在0.52ps之间),为系统的自由度。
如果采用LeapFrog时间积分Berendsen控温法对每一步缩放速率,使得即时温度被拉向目标温度:
采用LF时的时间积分:
限制控温法引入温度限制:
通过下式的最小平方实现动量的恒定:
对于LF时间积分:
8.2控压方法Hoover控压对于各向同性体系的运动方程:
温控有效常数压控有效常数体系自由度压控摩擦系数体系质心位置温控时间常数压控时间常数实时压力对LF时间积分:
仅盒子尺寸变化:
盒子尺寸和形状同时发生变化:
Berendsen控压系统遵循移动方程:
盒子体积乘以因子,而坐标和晶胞矢量乘以因子:
为等温压缩因子。
压控的同时并应用Berendsen温控。
8.3各种系综的分子动力计算方法均以LF法为例:
NVE系综NVT系综
(1)速度校正(rescalingvelocity)为速度校正因子,Tv为设定的温度,kB为玻尔兹曼因子。
(2)阻尼力方法(dampedforce)运动方程修改为:
速度计算公式为为常数NPH系综速度定义为:
运动方程为:
NPT系综结合NVT中阻尼力计算法和NPH中算法:
其中8.4各种系综计算结果的对比Brown和Clark等以256个氩原子为例,对各种系综进行了测试,其中均方位移可由下式得到:
由此可求出扩散系数:
系统的速度相关函数为:
也可以此求出扩散系数:
8.5系综的选择通常根据实际情况,选择最贴近现实的系综。
研究相变时,要选择NPT系综研究扩散系数,一般选择NVE系综
- 配套讲稿:
如PPT文件的首页显示word图标,表示该PPT已包含配套word讲稿。双击word图标可打开word文档。
- 特殊限制:
部分文档作品中含有的国旗、国徽等图片,仅作为作品整体效果示例展示,禁止商用。设计者仅对作品中独创性部分享有著作权。
- 关 键 词:
- 分子 模拟 第七