搜索
分子动力学模拟的基本过程黄敏生基本步骤A.原子位置的初始化•建立分子动力学模拟过程的首要问题和第一步是确定分子体系的初始条件。•两种方式，一是采用实验数据，二是借助各种理论模型得到分子结构的几何参数，如面心立方(FCC)模型等。•1.无论采取哪种方法，给定分子结构的空间坐标都不一定处在分子力场最稳定的位置，即各原子并非处在平衡态，造成体系的能量比较高。•2.要进行一个不施加载荷的弛豫过程，使得系统达到稳定的平衡状态（共轭梯度法）。•3.在这个过程中，系统从人为的初始构形转变成真实初始构形，势能减小并达到稳定。•4.初始条件最好与真实构形类似，FccBCC,固体结果影响较大，气体影响较小。A.原子位置的初始化•1.为使模拟尽快达到平衡，分子初始速度的分布应该尽量接近真实情形。•2.采用近似的Maxwell-Boltzmann统计分布来赋予原子的初始速度是比较合理的。能够使得系统尽快弛豫。A.原子速度的初始化速度的高斯分布x,[0,1]范围内的随机数0，2Pi随机数A初始条件二•1.在满足以上温度的条件下，必须保证系统净总动量为零。•2.另一种获得初始条件的方法是选取模拟过程某一时刻的原子坐标和速度。•3.分子动力学模拟经常分成不同的物理阶段进行，上一个模拟过程结束时的原子位置和速度就可以作为下一次模拟的初始条件。边界条件•分子动力学模拟中，只有足够的粒子数量，才能准确的描述材料的宏观性能。•为了减小计算规模，人们引入了周期性、固定、全反射等边界条件。目前常用的边界条件包括周期性边界条件、对称边界条件和固定边界条件。1.周期性边界(Periodicboundarycondition)像胞元•分子体系的模拟并不是都使用周期性边界条件，在很多情况下，如溶液中沉淀的分子团簇、蛋白质分子、病毒分子、材料的表面等并不需要周期性边界条件。•可以根据分子体系所处的外界环境对非周期边界上的粒子施加一定的限制。•例如，边界上的原子设计为位置固定的，就可以形成刚性边界。（原子始终不动）•对边界上的原子施加一定荷载或考虑边界上原子与外界环境之间的作用力，就形成阻尼边界非周期边界条件力场的截断力场的截断力场的截断力场的截断在分子动力学中在分子动力学中在分子动力学中在分子动力学中，，，，出于计算上的考虑出于计算上的考虑出于计算上的考虑出于计算上的考虑，，，，力场的截力场的截力场的截力场的截断是必须的断是必须的断是必须的断是必须的，，，，即在某一范围内力场是有效的即在某一范围内力场是有效的即在某一范围内力场是有效的即在某一范围内力场是有效的，，，，因因因因此会导致一些计算上的困难此会导致一些计算上的困难此会导致一些计算上的困难此会导致一些计算上的困难。。。。势函数直接截断势函数直接截断势函数直接截断势函数直接截断：：：：≤=cijcijcijijSrrrrrr0-)()(VVV力场的截断力场的截断力场的截断力场的截断力场连续的势函数截断力场连续的势函数截断力场连续的势函数截断力场连续的势函数截断：：：：≤--==cijcijcijrrijijcijijSrrrrrrrrrrcij0)()d)(d(-)()(VVVV近邻表近邻表近邻表近邻表•虽然引入了阶段半径的概念，但然而计算原子间的距离需要耗费大量的CPU时间。•设想研究对象为N个原子构成的粒子系统，由于要计算每个原子与其余原子间的距离，因此需要计算N(N-1)次原子间距，计算量随系统规模的增大成几何次增大。Verlet近邻表近邻表近邻表近邻表网格近邻表网格近邻表网格近邻表网格近邻表•最早由Verlet提出近邻表近邻表近邻表近邻表Verlet近邻表近邻表近邻表近邻表首先通过建立链表记录各原子周围首先通过建立链表记录各原子周围首先通过建立链表记录各原子周围首先通过建立链表记录各原子周围((((截断截断截断截断半径半径半径半径:r:r:r:rcccc))))区域内的所有原子区域内的所有原子区域内的所有原子区域内的所有原子，，，，然后通过链表然后通过链表然后通过链表然后通过链表计算各分子所受到的作用力计算各分子所受到的作用力计算各分子所受到的作用力计算各分子所受到的作用力。。。。为了计算原子为了计算原子为了计算原子为了计算原子1111所受的作用力所受的作用力所受的作用力所受的作用力。。。。假定一个假定一个假定一个假定一个球形的区域球形的区域球形的区域球形的区域，，，，半径为半径为半径为半径为rrrr1111，，，，且且且且rrrr1111rrrrcccc。。。。(r(r(r(rcccc为截为截为截为截断半径断半径断半径断半径))))，，，，并给截断半径并给截断半径并给截断半径并给截断半径rrrrcccc区域内的分子建区域内的分子建区域内的分子建区域内的分子建立一个链表立一个链表立一个链表立一个链表。。。。当截断半径当截断半径当截断半径当截断半径rrrrcccc范围内的分子没有离开范围内的分子没有离开范围内的分子没有离开范围内的分子没有离开rrrr1111球球球球形区域时形区域时形区域时形区域时，，，，只需要根据链表中的分子只需要根据链表中的分子只需要根据链表中的分子只需要根据链表中的分子，，，，即可即可即可即可计算出分子计算出分子计算出分子计算出分子1111所受到的作用力所受到的作用力所受到的作用力所受到的作用力。。。。如果截断半径如果截断半径如果截断半径如果截断半径rrrrcccc范围内的分子离开了范围内的分子离开了范围内的分子离开了范围内的分子离开了rrrr1111，，，，球形区域球形区域球形区域球形区域，，，，则需要建立新的链表则需要建立新的链表则需要建立新的链表则需要建立新的链表，，，，即在计算即在计算即在计算即在计算过程中过程中过程中过程中，，，，每隔一定时间每隔一定时间每隔一定时间每隔一定时间，，，，更新列表更新列表更新列表更新列表。。。。RListUpdateintervalTimeN=256N=500No2.602.702.903.103.433.50-5.7812.5026.3243.4883.33100.003.332.242.172.282.472.89-10.004.934.554.514.79-5.86占用了一定量的计算机内存。分子数较时，此方法具有一定的优势近邻表近邻表近邻表近邻表Verlet近邻表算法近邻表算法近邻表算法近邻表算法近邻表近邻表近邻表近邻表Verlet近邻表近邻表近邻表近邻表—周期性周期性周期性周期性方法（a）方法（b）这种方法的思想是将研究对象看成一个方盒，将这一方盒划分为M×M×M个单元（cell），每个单元的边长必须大于势函数的截断半径。与单元13内的原子间距小于截断半径的其它原子必然在单元13与其邻近单元7，8，9，12，13，14，17，18，19共9个单元内。由于每个单元内原子数为Nc=N/M2，因此对每个原子只要计算9Nc个原子间距，对整个原子系统就要计算9NNc个原子间距。对三维结构则要计算27NNc个原子间距（Nc=N/M3）。关于原子间距的计算量就与微结构的尺度即原子系统的原子数N成正比。近邻表近邻表近邻表近邻表网格近邻表网格近邻表网格近邻表网格近邻表•在计算过程中，每次更新原子位置时，将跨越本网格边界的原子从网格中删除，将其插入到相应的邻近网格中。•与邻域列表法相比，此方法不占用多余内存，在进行大规模分子系统模拟时，此方法可以明显的减少计算量。近邻表近邻表近邻表近邻表网格近邻表网格近邻表网格近邻表网格近邻表近邻表近邻表近邻表近邻表网格近邻表网格近邻表网格近邻表网格近邻表单位的无量纲化•在模拟中涉及很多浮点和指数运算，为提高计算效率，往往将温度、密度、压力等类似量表示成无量纲的形式。以Lennard-Jones势为例势阱常数平衡常数一种合适但不唯一的基本单位取为：系综平衡分子动力学模拟，总是在一定的系综下进行的。微正则系综正则系综平衡系综等温等压系综等压等焓系综系统原子数N，体积V，能量E保持不变。又称为NVE系综。孤立、保守的系统。一般说，给定能量的精确初始条件是无法得到的。能量的调整通过对速度的标度进行，这种标度可能使系统失去平衡，迭代弛豫达到平衡。系统原子数N，体积V，温度T保持不变，且总动量保持不变。又称为NVT系综。保持温度不变虚拟热浴系统动能固定原子速度标度系统原子数N，压强P，温度T保持不变，又称为NPT系综。压强P与体积共轭，控压可以通过标度系统的体积来实现。系统原子数N，压强P，焓值H=E+PV保持不变。在模拟中较少见。系综的控温•温度调控机制可以使系统的温度维持在给定值，也可以根据外界环境的温度使系统温度发生涨落。•一个合理的温控机制能够产生正确的统计系综，即调温后各粒子位形发生的概率可以满足统计力学法则。直接速度标定法Berendsen温控机制Nose-Hoover温控机制Gaussian温控机制•系统温度和粒子的速度直接相关，可以通过调整粒子的速度使系统温度维持在目标值。直接速度标定法∑==NifBiiNktvmtT12)()(•N系统粒子的个数。Nf=3*N-3为体系的自由度数（总动量固定）。波尔兹曼常数kB=1.38×10-23J/K。绝对温度T与体系的总动能密切相关。•温度是一个统计平均的量。对于单个的原子，其温度？？？？•引入速度标定因子直接速度标定法λ•Ndf为系统的总自由度数。Tc为目标温度；K为标度前体系的总动能。标度时实际分子动力学模拟中，并不需要对每一步的速度都进行标定，而是每隔一定的积分步，对速度进行周期性的标定，从而使系统温度在目标值附近小幅波动。或•直接速度标定法的优点是原理简单，易于程序编制.•缺点是模拟系统无法和任何一个统计力学的系综对应起来。•突然的速度标定引起体系能量的突然改变，致使模拟系统和真实结构的平衡态相差较远。直接速度标定法•又称Berendsen外部热浴法。•其其其其基本思想基本思想基本思想基本思想是假设系统和一个恒温的外部是假设系统和一个恒温的外部是假设系统和一个恒温的外部是假设系统和一个恒温的外部热浴耦合在一起热浴耦合在一起热浴耦合在一起热浴耦合在一起，，，，通过热浴吸收和释放能通过热浴吸收和释放能通过热浴吸收和释放能通过热浴吸收和释放能量来调节系统的温度量来调节系统的温度量来调节系统的温度量来调节系统的温度，，，，使之与恒温热浴保使之与恒温热浴保使之与恒温热浴保使之与恒温热浴保持一致持一致持一致持一致。。。。•对速度每一步进行标定对速度每一步进行标定对速度每一步进行标定对速度每一步进行标定，，，，以保持温度的变以保持温度的变以保持温度的变以保持温度的变化率与热浴和系统的温差化率与热浴和系统的温差化率与热浴和系统的温差化率与热浴和系统的温差(Tbath-T(t))成比例成比例成比例成比例。。。。每一步温度变化每一步温度变化每一步温度变化每一步温度变化：：：：Berendsen温控机制Berendsen温控机制•在Gaussian法中，给每个原子施加一个摩擦力，摩擦力的系数依赖于当前系统温度和热浴温度之间的差值。•若摩擦力系数为正，则表示系统温度高于热浴温度，需对系统进行冷却；若系数为负，则表示加热系统。Gaussian温控机制摩擦力项摩擦力系数•保持温度不变就是保持体系的总动能不变。也就是体系内部力所做功（功率为零(F-f).v=0）•根据上式，可以得到摩擦力系数为：Gaussian温控机制•所得期望温度为体系初始温度。程序编码简单的优点，但消除局域的相关运动。•该方法可以把任何数量的原子与一个热浴耦合起来，可以消除局域的相关运动，而且可以模拟宏观系统的温度涨落现象。•是等温系综统计力学研究领域的一个里程碑。Nose-Hoover温控机制基本思想引入一个反映真实系统与热浴相互作用的广义变量s(无量纲量)，将真实系统与热浴作为统一的扩展系统来考虑。扩展系统的哈密顿量为：•根据以上Hamiltian推导真实系统的运动方程：Nose-Hoover温控机制摩擦项摩擦系数率•系统常压条件压力诱导下的相