写点学术-获取小分子力场参数

04 Sep, 2025

参考

GROMACS中文手册：https://jerkwin.github.io/9999/12/31/GROMACS%E4%B8%AD%E6%96%87%E6%89%8B%E5%86%8C/
GROMACS英文手册： https://manual.gromacs.org/documentation/current/reference-manual/index.html

适用于GROMACS 不定时更新

这个东西应该只会对我自己有用吧，防止自己后面忘了所以记录一下
所以这篇略去了不少东西，如果没接触过可能会看不懂吧

↓先从做MD（molecular dynamics）模拟的必备常识说起

分子力学方法利用牛顿力学原理,通过解力学方程,得到分子体系的运动轨迹,从而研究分子体系的能量和各种性质.由于使用的势函数和参数都是由经验和实验值拟合获得,故亦称经验力场方法

1 分子力场

对含 $N$ 个相互作用原子的系统, MD模拟求解其牛顿运动方程:
$m_{i} \ddot{r_{i}} (t) = F_{i}, i = 1, \dots, N$ (1.1)
原子所受的力是一个势能函数导数的负值:
$𝐅_{i} = - \frac{\partial V}{\partial 𝐫_{i}}$ (1.2)
模拟时, 会以很小的时间步长同时求解这些方程, 使系统演化一段时间, 同时小心地使温度和压力维持在设定值, 并以一定的时间间隔将原子坐标写入输出文件中。随着时间变化的坐标代表了体系的一个轨迹. 经过初始的一些变化之后, 体系通常会达到平衡态. 通过对平衡后的轨迹进行平均, 就可以从输出文件中提取出许多宏观性质.
一个力场由两个明显不同的两部分组成:用于产生势能及其导数（也即力）的方程（或者叫势函数）组+用于方程组的参数.
例如：AMBER力场总势能可以写成 $V = \sum_{bonds} k_{r} (r - r_{0})^{2} + \sum_{angles} k_{θ} (θ - θ_{0})^{2} + \sum_{torsions} k_{ϕ} [1 + \cos (n ϕ - γ)] & + \sum_{i < j} [4 ϵ_{i j} ((\frac{σ_{i j}}{r_{i j}})^{12} - (\frac{σ_{i j}}{r_{i j}})^{6}) + \frac{q_{i} q_{j}}{4 π ε_{0} r_{i j}}]$ （1.3）
E(V) = 成键相互作用能 + 非键相互作用能 + 附加能
E(V）= 伸缩能 + 弯曲能 + 扭曲能 + VDW + ELE + 交叉项能和氢键能
力场参数就是方程里的这些键参数、角度参数等
局限性:
1 模拟是经典力学的：
当共振频率 $v$ 接近或超过 $k_{B} T / h$ 时, 经典谐振子与真正的量子谐振子的统计力学明显不同,在经典模拟中的行为都可能不正确. 这意味着所有键和键角的振动实际上都是可疑的, 甚至氢键键合基团的运动, 比如平动或氢键的摆动都超出了经典的极限.除了进行真正的量子动力学模拟（显然目前对于大体系做不到）, 我们可以采取下面两种方法中的一种:

如果进行MD模拟时对键采用了谐振子近似, 我们应该对总内能和比热容 (以及熵、自由能如果要计算它们) 进行校正.
我们可以将键(和键角)视为运动方程的约束. 其理由是, 与经典振子相比, 处于基态的量子振子与被约束的键更相似. 作这种选择有一个很实用的理由: 当最高频率被移除后, 可以在算法中使用更大的时间步长. 实际中, 相比把键当作谐振子, 当键被约束时的时间步长可以取为原先的4倍. GROMACS对键和键角提供了这种方案. 灵活的键角约束是相当重要的, 它使得运动更加真实, 并能覆盖构型空间，例如.mdp文件里的：constraints 决定了哪些键被约束，这里没有特殊需要就取constraints = h-bonds。

2 电子处于基态：
B-O近似，认为原子核是不运动的

3 力场是近似的：\

4 力场是对势累加的：\

5 长程相互作用被截断：\

6 边界条件是不自然的：\