第六章 分子动力学模拟

分子动力学模拟
三、模拟细节
3.1 初始化位型
分子动力学模拟
模拟开始时首先要初始化位型分布，既首先要给定微 元中分子的初始位臵和初始速度。
①、简单立方晶格
②、体心立方晶格
c
③、面心立方晶格
a b
c
分子动力学模拟
a
b c c c c
c c a a a ba c
b c b
b
b a a a ba b
分子动力学模拟
6
1
5
4
D
2
6
1
5
4
C
2
6
1
5
4
B
2
3 6
1
3 6
1
3
6
1
5
E
4
2
基本元胞 5
2
5
A
4
2
3 6
1
4
3 6
1
3
6
1
5
F
4
2
5
4
G
2
5
H
4
2
3
3
L
3
3.3 Calculation of interactions
i , j 两个离子的相互作用势 体系总势能可表示为：
分子动力学模拟
U u (rij )
3/ 2
mv 2 2 exp 2kT v dv
3.2 Periodic boundary conditions
分子动力学模拟
通常选取小数体系（几十个到数千个分子）作为研究 对象，但由于位于表面和内部的分子受力差别较大，将会 产生表面效应。
为了消除此效应，同时建造 出一个准无穷大体积，使其可以 更精确地代表宏观系统，必须引 入周期性边界条件。 周期性边界条件： 把小数体系看作一个中心元胞，此中心元胞被与中心 元胞相同的其它元胞所包围，当某分子离开中心元胞时， 该分子将在整个点格上移动以致它从中心元胞对面重新进 入这个元胞。
两式相加得
d 2 ri (t ) dt 2 1 1 [ri (t h) 2ri (t ) ri (t h)] f i (t ) 2 m h
分子动力学模拟 上式提供了一个方法，从粒子在前两步（t和t-h）时刻 的位臵以及t时刻的作用力来得到粒子在t+h时刻的位臵。
1 2 ri (t h) 2ri (t ) ri (t h) h f i (t ) m
b
初始速度
分子动力学模拟
模拟时，各粒子的初始速度按麦克斯韦速度分布取样。
Maxwell’s distribution law of velocity 指平衡状态下理想气体分子 速度分布的统计规律。在平衡状 态下，分布在任一速率区间内的 分子数与总分子数的比率为
m f (v)dv 4 2kT
分子动力学模拟
V
HS
(r ) 0
r r
3. 软球势（Soft-Sphere）
V ( r ) ( ) r r
SS

通常，v 是为整数的参数。
4. 方阱势（Square-Well）
SW V ( r ) 0
r 1
1 r 2 r 2
h i (i ) r (t h) r (t ) r (t ) Rn i 1 i!
n 1
1、有限差分方法-预测校正法
r p (t t ) r (t ) tv (t ) t 2a(t )/2 t 3b(t )/6 v p (t t ) v (t ) ta(t ) t 2b(t )/2 a p (t t ) a(t ) tb(t ) b p (t t ) b(t )
1
ÆÜ ÊÄ
0.5
rij 6 2
0

-0.5
-1 0.8
1
1.2
1.4
1.6
1.8 ¼¾ äà
2
2.2
2.4

2.6
对势能的最大贡献来自于粒子的近邻区域，位势截断
常用的方法是球形截断，截断半径一般取2.5σ或3.6 σ，对
截断距离之外分子间相互作用能按平均密度近似的方法进 行校正。
分子动力学模拟
方法首先是由Alder和Wainwright提出的，现已逐渐成为预测系统特性、
验证理论和改进模型的计算工具。 MD方法的基本思想是把物质看成由原子和分子组成的粒子系统， 从该体系的某一假定的位能模型出发，并假定体系粒子的运动遵循经典 力学或量子力学描述的规律，若已知粒子的所有受力作用，则可以求解
出运动方程而得到系统中全体粒子在相空间中的轨道，然后统计得到系
2.2 Potential energy functions
分子动力学模拟
对一个由N个原子构成的简单系统，其势能项由下式给出
V V1 ( ri ) V2 (ri , r j ) V3 (ri , r j , rk )
i i j i i j i k j i
' ' ' ri ri vi t ri vi h
The MD simulation can describe systems that evolve in time. The new positions are derived from the Newtonian law of motions and therefore deterministic.
i j
在计算势能时必须考虑相互作用的力程，对粒子数为
N的模拟系统，原则上任何两个粒子间都有相互作用，在
计算体系势能时须进行N（N-1）/2次运算。 为了提高计算效率，实际模拟 过程中通常在一个方便的力程上 将位势截断，用以减少计算势能 所消耗的时间。
分子动力学模拟
Lennard-Jones Ê Ä Ä Ð ÆÜ£Í 2 1.5
The force causes an acceleration
分子动力学模拟
fi mai
Which in turn modifies the initial velocity vi as
' vi vi ai t vi ai h
And modifies the initial position ri as
vi (t 1 1 1 h) vi (t h) hfi (t ) 2 2 m 1 ri (t h) ri (t ) hvi (t h) 2
1 1 vi (t h) [ri (t ) ri (t h)] 2 h 1 1 vi (t ) [vi (t h) vi (t h)] / 2 2 2
1 3 2 * Ekin m vi NkT 2 i 2
3NkT T /T * m vi2
i
v'i vi
分子动力学模拟
2、有限差分方法-Verlet算法
①、Verlet算法的一般形式
为了用数值方法求解微分方程，
fi mai
采用有限差分方法离散化得到差分格式
1 2 ri (t h) ri (t ) hvi (t ) h f i (t ) 2m 1 2 ri (t h) ri (t ) hvi (t ) h f i (t ) 2m
uri uij ri , rj
j
In that potential the particle feels a force
fi uri
According to Newton’s second law
2 vi ri mai m m 2 fi t t
ri
n 1
2ri ri
n
n 1
1 2 n h fi m
ri n 1 ri n 1 v in 2h
速度对轨道计算没有关系，但对 估算动能及速度相关函数（用来研 究粒子的输运特性）非常有用。速 度形式为
②、Verlet 蛙跳格式（leap-frog）
分子动力学模拟
The velocity varies so one has to choose a reasonable average value to be used. The velocity at the middle of the step ought to be a good compromise,
分子动力学模拟 下面仅对简单系统的相互作用模型给予简介
1、Lennard-Joans势
Lennard-Jones Ê Ä Ä Ð ÆÜ£Í 2
1.5
1
ÆÜ ÊÄ
0.5
rij 6 2
0百度文库

-0.5
-1 0.8
1
1.2
1.4
1.6
1.8 ¼¾ äà
2
2.2
2.4

2.6
2. 硬球势（Hard-Sphere）
分子动力学模拟
2.3 Calculations of force, velocity, position
The initial distribution of the Molecular dynamics simulation is generated in a random distribution. Each particle is also assigned an initial velocity vi
统的热力学参数、结构和输运特性等。也就是由体系的微观性质来求算 其宏观性质。属于微观尺度的模拟技术。
二、分子动力学方法
2.1 Newtonian mechanics
分子动力学模拟
In the MD method, every new distribution is derived from the previous one by using the interactions between the particles. The interactions depends on the position of the particles.
式中右端第一项是外场（如电场、 磁场、声场等）对系统的作用；第二项 是两体势即系统中每两个粒子间的相互 作用；第三项是三体势，表示系统中每 三个粒子间的相互作用……
有效两体势
V V1 ( ri ) V
i i j i
eff 2
( rij )
V ( rij )
它包含多体效应，可很好地反映系统粒子间的相互作用。
③、Verlet速度相加形式
分子动力学模拟
1 2 ri (t h) ri (t ) hvi (t ) h f i (t ) 2m
1 vi (t h) vi (t ) h[ f i (t ) f i (t h)] 2m
1 1 vi (t h) vi (t ) hfi (t ) 2 2m v (t h) v (t 1 h) 1 hf (t h) i i i 2 2m
2.4 Equations of motion
分子动力学模拟
为了在计算机上解运动方程，必须为微分方程建立一个 有限差分格式，从差分方程中再导出位臵和速度的递推关系 式。这些算法是一步一步执行的，先算t 时刻的位臵和速度， 然后在此基础上计算t+1时刻的位臵和速度。 微分方程最为直接的离散化格式来自泰勒展开:
f x 48 f y 48 f z 48
1 ( x i x j )[( )14 ( ) 8 ] rij 2 rij 2
1 ( y i y j )[( )14 ( ) 8 ] rij 2 rij 2
In simulation:
1 ( z i z j )[( )14 ( ) 8 ] rij 2 rij 2
★Cell method ★Neighbor list method
1、近邻表- Verlet近邻表
2、近邻表- 网格近邻表
分子动力学模拟
3.4 标度与趋衡
标度：正则系综的MD模拟，模拟过程中必须保证体系温度 恒定，所以在每一时间步要对分子的速度进行标度使其满 2 足正则系综的条件，标度因子为 (T / T *)1 /，其中T是 指定的体系温度，T*是每一时间步的瞬时动能温度。
分子动力学模拟
中国石油大学
分子动力学模拟
本章主要内容
一、系综理论 二、分子动力学方法
三、模拟细节
四、参量的计算 五、液态水的MD模拟 六、误差分析 七、分子动力学模拟方法的应用
分子动力学模拟
一、系综理论
分子动力学模拟（molecular dynamics simulation，简称MD）