分子动力学模拟中的力场参数基础
做分子动力学模拟，你首先得回答一个问题：力到底是怎么算出来的？ 答案很直接——靠力场（Force Field）。力场说白了就是一套数学公式，加上一堆参数。没有它，你的原子就是一堆自由粒子，到处乱飞。我个人习惯在项目开始时，先花点时间把力场构成过一遍，这事真不能省。
2.1 力场的组成：势能面的拆解
一个典型的全原子力场，通常把体系的总势能拆成两大部分：键合相互作用 和 非键合相互作用。你想想看，这个拆分的逻辑其实很自然——原子通过化学键直接相连是一种相互作用；不相连的原子隔着老远也能感受到静电和范德华力，那是另一种相互作用。
我遇到过不少新手，总觉得力场参数特别神秘。其实它就是一个数学拆解问题。以下是常见的几个能量项：
- 键伸缩（Bond Stretching）：两个成键原子像弹簧一样，拉伸或压缩会带来能量变化。
- 键角弯曲（Angle Bending）：三个连续原子形成的角度变形，同样会产生能量。
- 二面角扭转（Torsion/Dihedral）：四个连续原子绕中间键的旋转，这是分子构象变化的核心。
- 范德华相互作用（van der Waals）：非键合原子之间的短程吸引力与斥力。
- 静电相互作用（Electrostatic）：携带部分电荷的原子之间的库仑作用。
2.2 函数形式：参数藏在数学公式里
有了上面这些势能项，接下来就看它们怎么用数学表达。力场参数就是这些公式里的系数和常数。以最经典的 AMBER 和 CHARMM 力场为例，它们的函数形式大同小异。我自己的习惯是，先把这些函数形式写下来，再对着看参数，会清楚很多。
2.2.1 键伸缩项
最简单的模型——谐振子（Harmonic Oscillator）。想象一下，原子之间就像用弹簧连着。
E_bond = Σ (1/2) * k_r * (r - r₀)²k_r：力常数，单位通常为 kcal/(mol·Å²) 或 kJ/(mol·nm²)。它决定了弹簧的“硬度”。r₀：平衡键长，单位 Å 或 nm。r：当前键长。
嗯，这里要注意：你可能会看到有的力场输出的是 k_r/2 还是 k_r，这个细节我曾经在迁移力场参数时吃过亏。不同力场的参数定义方式略有不同，搞混了结果就跑偏了。
2.2.2 键角弯曲项
同样是谐振子形式，只不过变量变成了角度。
E_angle = Σ (1/2) * k_θ * (θ - θ₀)²k_θ：角度力常数，单位 kcal/(mol·rad²)。θ₀：平衡键角，单位度。θ：当前键角。
我遇到过一些初学者，把角度直接代入弧度制，但参数却是按度定义的，跑出来的能量直接爆表。这种事情，犯一次就记住了。
2.2.3 二面角扭转项
这一项稍微复杂点，用傅里叶级数来拟合周期性旋转能垒。
E_torsion = Σ (V_n/2) * [1 + cos(n * φ - γ)]V_n：旋转能垒，单位 kcal/mol，决定了这个旋转的“力”有多大。n：多重度，表示旋转能垒的波峰数量。γ：相位角，控制余弦曲线的偏移。φ：当前二面角。
我在项目中遇到过一种情况：同一个分子，不同力场给的二面角参数有时能差好几倍。这种时候，一定要核对参数来源。你可能会问，这会影响什么？直接决定你分子构象的占有率，影响大了去了。
2.2.4 非键合相互作用
非键合项分成两大块，一块是范德华，一块是库仑。我一直觉得，非键合项才是力场参数里最难搞的部分，因为它的截断策略、组合规则到处是坑。
范德华（Lennard-Jones 12-6 势）
E_LJ = 4 * ε * [ (σ/r)^12 - (σ/r)^6 ]ε：势阱深度，单位 kcal/mol，代表相互作用的强弱。σ：碰撞直径，单位 Å，代表原子的有效尺寸。r：两原子间的距离。
静电（库仑势）
E_coulomb = Σ (q_i * q_j) / (ε_r * r_ij)q_i, q_j：原子所带的部分电荷，单位 e。ε_r：相对介电常数。r_ij：原子间距。
2.3 参数意义：一个参数，一段故事
参数不是数字拼凑，每一个都有明确的物理意义。我们来总结一下。
| 参数 | 所属项 | 物理意义 | 常见单位 |
|---|---|---|---|
k_r |
Bond | 键伸缩的刚性（力常数） | kcal/(mol·Å²) |
r₀ |
Bond | 键的平衡长度 | Å |
k_θ |
Angle | 键角弯曲的刚性 | kcal/(mol·rad²) |
θ₀ |
Angle | 键角的平衡值 | ° |
V_n |
Dihedral | 旋转能垒高度 | kcal/mol |
n |
Dihedral | 旋转周期内的波峰数 | — |
γ |
Dihedral | 相位偏移 | ° |
ε |
vdW (LJ) | 势阱深度，吸引力强度 | kcal/mol |
σ |
vdW (LJ) | 碰撞直径，原子有效尺寸 | Å |
q |
Coulomb | 原子部分电荷 | e |
2.4 参数背后的“版本锁”影子
你可能注意到了，上面这堆参数——无论是 k_r、ε 还是 q——它们都是带有“版本归属”的。比如，同一根 C-C 键，在 CHARMM36 和 AMBER ff14SB 里，r₀ 可能分别为 1.335 Å 和 1.340 Å。这个微小的差异，决定了你模拟出来的结构正确与否。
所以，每次拿到新版本的力场，我用 diff 或者 验证工具 去比对参数文件。你想想看，一个参数变了，整条势能曲线就变了。如果没记录下来，后续的论文结果、合作者的复现，都会变成一团乱麻。
我曾在一次对外的协助项目中，因为疏忽，在跑了一个月之后才发现用于二面角的 V_n 参数写的是旧版本。结果重跑了一整轮。从那以后，我在项目开始前一定会把力场参数版本锁起来，用配置文件固定住，并且在模拟日志里打印每个势能项的参数来源。这一步，真不能省。
嗯，以上就是力场参数的基础回顾。你掌握了这些，后面章节讲版本锁定的具体操作，就非常容易理解了——因为你已经知道，我们到底在锁什么。