3. 键合相互作用：键伸缩能、键角弯曲能、二面角扭转能、improper 二面角
好，咱们今天聊聊力场里最基础、也最绕不开的一块——键合相互作用。说白了，就是分子里那些原子之间通过化学键连起来，然后它们怎么动、怎么晃、怎么扭，我们得用数学公式描述出来。
我个人习惯把这块内容叫做「分子骨架的弹簧模型」。你想想看，分子不是铁板一块，原子之间像用弹簧连着一样，可以拉长、压短、弯折、扭转。力场要做的，就是给这些运动定个「规矩」——能量怎么算，力怎么来。
3.1 键伸缩能（Bond Stretching）
先讲最简单的：两个原子之间的键，像一根弹簧。你拉它，它想缩回去；你压它，它想弹开。这个能量怎么描述？
最经典的公式是谐振子模型：
E_bond = ½ * k_b * (r - r₀)²这里：
- k_b 是力常数，弹簧的「硬」度
- r₀ 是平衡键长，也就是原子最舒服的距离
- r 是当前键长
嗯，这里要注意：这个公式只在平衡位置附近好用。如果拉得太远，比如快到断键了，谐振子模型就不准了。我在项目中遇到过模拟聚合物拉伸，结果键长偏离平衡值超过 0.3 Å，能量就开始离谱地飙升。后来我换成了 Morse 势，才把断键行为描述得更合理。
3.2 键角弯曲能（Angle Bending）
三个原子组成一个角，比如水分子里的 H-O-H。这个角也不是固定死的，它会弯来弯去。能量公式和键伸缩很像：
E_angle = ½ * k_θ * (θ - θ₀)²参数含义类似：
- k_θ 是角度力常数
- θ₀ 是平衡键角
- θ 是当前键角
为什么也是二次项？因为在小变形范围内，大多数化学键的角弯曲行为确实接近简谐振动。但如果你把角度掰到 180° 附近，比如线性分子，这个模型就会出问题——因为角度从 179° 变到 181°，能量应该连续，但余弦函数会跳变。
我曾经踩过这个坑：模拟一个线性分子（比如 CO₂），用的力场里角度项是标准的谐振子，结果 O-C-O 角度在 180° 附近震荡时，能量曲线出现不连续，导致体系温度失控。后来我换成了余弦型角度势（比如 Urey-Bradley 项），才把问题解决。
3.3 二面角扭转能（Dihedral Torsion）
这个就更有意思了。四个原子连成一条链，比如乙烷里的 H-C-C-H，它们可以绕着中间的 C-C 键旋转。这个旋转不是自由的，因为有空间位阻和电子效应。
二面角能量的标准形式是傅里叶级数：
E_dihedral = Σ (V_n / 2) * [1 + cos(n * φ - γ)]其中：
- V_n 是势垒高度
- n 是周期性（比如乙烷是 3 重周期）
- φ 是二面角
- γ 是相位偏移
为什么用傅里叶级数？因为二面角旋转一圈，能量变化是周期性的。比如乙烷，每转 120° 能量重复一次，所以 n=3。但像丁烷这样的分子，有 gauche 和 anti 两种构象，就需要多个傅里叶项叠加才能准确描述。
我记得有一次做药物分子构象搜索，力场里二面角参数是从类似分子 transfer 过来的，结果预测的扭转势垒和量子化学计算差了 2 kcal/mol。后来我手动拟合了 V₁、V₂、V₃ 三个项，才把能量曲线对齐。所以，二面角参数千万别偷懒，一定要验证。
3.4 Improper 二面角（Out-of-Plane Bending）
这个名词有点绕，但概念很简单。Improper 二面角不是描述绕键旋转，而是描述一个原子偏离平面的运动。比如苯环上的碳原子，如果被拉出环平面，能量会急剧上升。
常见的 improper 势有两种：
- 谐振子型： E_improper = ½ * k_ψ * (ψ - ψ₀)²
- 余弦型： E_improper = k_ψ * (1 + cos(n * ψ - ψ₀))
这里 ψ 是 improper 二面角，通常定义为四个原子中，中心原子与另外三个原子构成的平面之间的偏离角。
为什么要单独搞一个 improper 项？因为普通的二面角（proper dihedral）只能描述绕键旋转，不能描述平面外弯曲。比如 sp² 杂化的碳（像乙烯、苯），它们倾向于保持平面构型，一旦被拉出平面，能量惩罚必须很大。Improper 项就是干这个的。
我在做石墨烯力学性质模拟时，深刻体会到了 improper 项的重要性。如果力场里没有 improper 项，石墨烯的弯曲刚度会严重偏低，模拟出来的杨氏模量只有实验值的一半。加上 improper 项后，结果才和实验吻合。
3.5 知识体系总览
为了让你更直观地理解这四种键合相互作用的关系，我画了一张图：
3.6 参数验证的实用建议
讲了这么多理论，最后给点实操建议。怎么验证你选的键合参数对不对？
- 第一步： 拿一个小分子（比如你体系里的一个片段），用量子化学（DFT 或 MP2）优化结构，得到平衡键长、键角、二面角。和力场参数对比，偏差超过 0.05 Å 或 5° 就要警惕。
- 第二步： 做一次约束扫描。固定一个键长或二面角，逐步改变它，计算能量变化。把力场结果和量子化学结果画在一起，看曲线是否吻合。
- 第三步： 跑一段短 MD（比如 100 ps），检查键长、键角的分布是否合理。如果出现异常大的波动，说明力常数可能偏小。
好了，键合相互作用这块就聊到这儿。下一节我们会讲非键相互作用——范德华力和静电作用，那又是另一片天地了。