lj/cut力场公式：
E
=
4
ϵ
[
(
σ
r
)
12
−
(
σ
r
)
6
]
r
<
r
c
E=4\epsilon\left[\left(\frac{\sigma}{r}\right)^{12}-\left(\frac{\sigma}{r}\right)^{6}\right]\quad\quad r<r_{c}
E=4ϵ[(rσ)12−(rσ)6]r<rc 在lammps使用lj/cut只需要设置3个参数：
ϵ
\epsilon
ϵ、
σ
\sigma
σ、
r
c
r_{c}
rc lj/cut有两种写法，第一种是写到in文件中，另一种是写到data文件中，两种写法稍有不同。
lj/cut/coul/cut力场公式：
E
=
4
ϵ
[
(
σ
r
)
12
−
(
σ
r
)
6
]
r
<
r
c
E=4\epsilon\left[\left(\frac{\sigma}{r}\right)^{12}-\left(\frac{\sigma}{r}\right)^{6}\right]\quad\quad r<r_{c}
E=4ϵ[(rσ)12−(rσ)6]r<rc
E
=
C
q
i
q
j
ϵ
r
r
<
r
c
E=\frac{Cq_{i}q_{j}}{\epsilon r}\quad\quad r<r_{c}
E=ϵrCqiqjr<rc lj/cut/coul/cut的截断半径有两个：lj势的截断半径与库仑力的截断半径，若只写一个截断半径，默认cutoff2 = cutoff。
长程库仑力的计算需要特定算法，kspace_style命令。但是这个算法经常提示错误，最常见的错误，“Out of range atoms - cannot compute PPPM”。出现原因分为两种情况： 1.原子移动速度相对较低，在一个时间步内移动距离超过1/2的neighbor bin，同时中心原子的邻居列表更新不及时。解决方法：提高邻居列表的更新频率，减小时间步长。 2.原子移动速度相对较高，当原子重叠度较高或者力场参数设置不对，原子之间的受力过大，导致原子移动速度过快，在一个时间步长内跑出盒子。解决方法：重新建模，若用MS建模，可以在MS跑平衡，获得一个近似平衡稳定的体系。
3 morse
morse力场与lj力场相似，用来描述没有键相互连接的原子之间的作用力。 morse力场公式：
E
=
D
0
[
e
−
2
α
(
r
−
r
0
)
−
2
e
−
α
(
r
−
r
0
)
]
r
<
r
c
E=D_{0}\left[e^{-2\alpha\left(r-r_{0}\right)}-2e^{-\alpha\left(r-r_{0}\right)}\right]\quad\quad\quad\quad r<r_{c}
E=D0[e−2α(r−r0)−2e−α(r−r0)]r<rc 在使用morse力场时需要设置4个参数，分别为：
D
0
D_{0}
D0、
α
\alpha
α、
r
0
r_{0}
r0、
r
c
r_{c}
rc。 morse力场使用方法：
pair_style morse cutoff
pair_coeff ** D0 a r0 [cutoff]
pair_style tersoff
pair_coeff ** SiC.tersoff Si Si Si Si C
5 EAM
金属原子之间没有键连接，因此，在lammps模拟中，金属体系的势函数类型为pair_style，而不是bond_style。 模拟金属体系时，可以用lj势描述金属原子之间的受力，不过更精确的是嵌入原子势(EAM)，eam势函数公式为：
E
i
=
F
α
(
∑
j
≠
i
ρ
β
(
r
i
j
)
)
+
1
2
∑
j
≠
i
ϕ
α
β
(
r
i
j
)
E_{i}=F_{\alpha}\left(\sum_{j\ne i}\rho_{\beta}\left(r_{ij}\right)\right)+\frac{1}{2}\sum_{j\ne i}\phi_{\alpha\beta}\left(r_{ij}\right)
Ei=Fα⎝⎛j=i∑ρβ(rij)⎠⎞+21j=i∑ϕαβ(rij) eam势由两部分组成，在原子对势(pair)的基础上添加了电子云密度相关项，比单纯的对势精确度更高。 eam势函数写在一个以“.eam”为后缀的文件中，lammps自带的势函数包含一部分eam势文件，也可以到网站下载eam文件，下载eam文件后，保存到in文件所在的文件夹。 常用的势函数下载网站有： http://www.ctcms.nist.gov/potentials https://openkim.org
#elt lat z iel atwt# alpha b0 b1 b2 b3 alat esub asub# t0 t1 t2 t3 rhozero ibar'Zr''hcp'12191.22004.45019083282.4501.0003.0002.0003.20000000006.3600.6801.006.300-3.300-10.001.0003'Cu''fcc'12163.54605.15483008303.8302.2006.0002.2003.61331565193.5400.9401.002.7203.0401.9501.0003
TIP4P水在正常水分子的基础上增加了一个虚原子，具体可查阅lammps官方手册(https://docs.lammps.org/Howto_tip4p.html)，虚原子不是实体原子，因此在水分子建模中不需要建出虚原子，仅在pair势中设置即可。 这两种势对应的公式一样的，由两部分组成：LJ势和库仑势。 LJ势对应公式：
E
=
4
ϵ
[
(
σ
r
)
12
−
(
σ
r
)
6
]
r
<
r
c
E=4\epsilon\left[\left(\frac{\sigma}{r}\right)^{12}-\left(\frac{\sigma}{r}\right)^{6}\right]\quad\quad r<r_{c}
E=4ϵ[(rσ)12−(rσ)6]r<rc 库仑势对应公式：
E
=
C
q
i
q
j
ϵ
r
r
<
r
c
E=\frac{Cq_{i}q_{j}}{\epsilon r}\quad\quad r<r_{c}
E=ϵrCqiqjr<rc
10.2 TIP4P设置案例
units real #real单位体系
atom_style full #full原子类型
boundary p p p #周期性边界条件
pair_style lj/cut/tip4p/long12110.154612.0#pair_style 参数设置：氧原子类型，氢原子类型，键类型，角类型，虚原子距离，截断半径
bond_style harmonic #键类型
angle_style harmonic #角类型
kspace_style pppm/tip4p 1.0e-4#长程库仑力算法
pair_modify shift yes mix arithmetic
read_data H20.data
#设置力场参数
pair_coeff 110.18523.1589#O-O原子
pair_coeff 2200#H-H原子
bond_coeff 110000.9572#H-O键
angle_coeff 11000104.52#H-O-H角
fix 1all shake 0.0001500 b 1 a 1#通过SHAKE算法对键角进行约束，迭代精度差0.0001，迭代步数50，不输出运行结果，对类型1的键和类型2的键约束
比较常见的一个错误是：All pair coeffs are not set，错误提示为“不是所有的原子都设置了势参数”，表示有的原子可能没有设置势参数。 出现这种错误主要有两种原因： 1.原子数量多，的确少写了某种原子的势参数。 2.设置了全部的原子势函数，但是设置方法出错。 下面一个简单的Cu-Al模型为例，给出这种错误的解决方法。 Cu原子类型为1，Al原子类型为2，Cu和Al均使用meam势，Cu-Al之间使用LJ势。 势函数设置如下：
pair_style hybrid meam/c lj/cut 10
pair_coeff ** meam/c library.Cu.meam Cu Cu.meam Cu NULL #Cu
pair_coeff ** meam/c library-Al.meam Al Al.meam NULL Al #Al
pair_coeff 12 lj/cut 0.42.47#Cu-Al
如果按照以上代码设置，就会提示“All pair coeffs are not set”。 出现这种错位，主要是因为同一in文件中使用了两个meam势，第二个Al的meam势会覆盖掉Cu的meam势，导致Cu的meam势参数丢失。 解决方法： 如果在同一in文件中使用多个同种类型的势，**为防止覆盖，需要对相同类型的势进行编号区分，**所以正确写法：
pair_style hybrid meam/c meam/c lj/cut 10
pair_coeff ** meam/c 1 library.Cu.meam Cu Cu.meam Cu NULL #Cu
pair_coeff ** meam/c 2 library-Al.meam Al Al.meam NULL Al #Al
pair_coeff 12 lj/cut 0.42.47#Cu-Al
OPLS力场不需要力场文件，但需要设置力场参数。OPLS包含了键(bond)、角(angle)、非正常二面角(dihedral或torison)以及非键接势(non-bonded)。在lammps设置中，这些势需要单独设置。 OPLS对应公式：
B
o
n
d
s
t
r
e
t
c
h
i
n
g
:
E
b
o
n
d
=
∑
b
o
n
d
s
K
r
(
r
−
r
e
q
)
2
Bond stretching:\quad \quad E_{bond}=\sum_{bonds}K_{r}\left(r-r_{eq}\right)^{2}
Bondstretching:Ebond=bonds∑Kr(r−req)2
A
n
g
l
e
b
e
n
d
i
n
g
:
E
a
n
g
l
e
=
∑
a
n
g
l
e
s
K
θ
(
θ
−
θ
e
q
)
2
Angle bending:\quad \quad E_{angle}=\sum_{angles}K_{\theta}\left(\theta-\theta_{eq}\right)^{2}
Anglebending:Eangle=angles∑Kθ(θ−θeq)2
T
o
r
s
i
o
n
:
E
(
ϕ
)
=
V
1
2
[
1
+
c
o
s
(
ϕ
+
f
1
)
]
+
V
2
2
[
1
−
c
o
s
(
2
ϕ
+
f
2
)
]
+
V
3
2
[
1
+
c
o
s
(
3
ϕ
+
f
3
)
]
Torsion:\quad E\left(\phi \right) =\frac{V_{1}}{2}\left[1+cos\left(\phi+f1\right)\right]+\frac{V_{2}}{2}\left[1-cos\left(2\phi+f2\right)\right]+\frac{V_{3}}{2}\left[1+cos\left(3\phi+f3\right)\right]
Torsion:E(ϕ)=2V1[1+cos(ϕ+f1)]+2V2[1−cos(2ϕ+f2)]+2V3[1+cos(3ϕ+f3)]
N
o
n
−
b
o
n
d
e
d
:
E
a
b
=
∑
i
o
n
a
∑
j
o
n
b
[
q
i
q
j
e
2
r
i
j
+
4
ε
i
j
(
σ
i
j
12
r
i
j
12
−
σ
i
j
6
r
i
j
6
)
]
f
i
j
Non-bonded: \quad E_{ab}=\sum^{ona}_{i}\sum^{onb}_{j}\left[\frac{q_{i}q_{j}e^{2}}{r_{ij}}+4\varepsilon_{ij}\left(\frac{\sigma^{12}_{ij}}{r^{12}_{ij}}-\frac{\sigma^{6}_{ij}}{r^{6}_{ij}}\right)\right]f_{ij}
Non−bonded:Eab=i∑onaj∑onb[rijqiqje2+4εij(rij12σij12−rij6σij6)]fij
f
i
j
=
0.5
i
f
i
,
j
a
r
e
1
,
4
;
o
t
h
e
r
w
i
s
e
,
f
i
j
=
1.0
\quad \quad \quad \quad f_{ij}=0.5\quad if\quad i,j\quad are \quad1,4;\quad otherwise,\quad f_{ij}=1.0
fij=0.5ifi,jare1,4;otherwise,fij=1.0 (1)bond
born势对应公式：
E
=
A
e
x
p
(
σ
−
r
ρ
)
−
C
r
6
+
D
r
8
r
<
r
c
E=Aexp\left(\frac{\sigma -r}{\rho}\right)-\frac{C}{r^{6}}+\frac{D}{r^{8}}\quad \quad r<r_{c}
E=Aexp(ρσ−r)−r6C+r8Dr<rc &enap;&enap;&enap;&enap;born势不需要势文件，仅需要设置势类型和势参数。
# 不计算库仑力
pair_style born [cutoff]# 计算库仑力
pair_style born/coul/long[cutoff]
15 Buckingham势
Buckingham势基本公式：
E
=
A
e
−
r
ρ
−
C
r
6
r
<
r
c
E=Ae^{-\frac{r}{\rho}}-\frac{C}{r^{6}}\quad \quad r<r_{c}
E=Ae−ρr−r6Cr<rc pair_style buck命令包括多个类型，pair_style buck/coul/cut(包括长程库仑力)，pair_style buck/coul/long(包括长程库仑力)，pair_style buck/coul/msm(包括长程MSM库仑力)。 举例：
