OpenMM：分子动力学模拟开源引擎

OpenMM 是 GPU 加速的开源分子动力学引擎，最大特点是用 Python 直接描述模拟体系，不需要写 GROMACS/AMBER 那样的输入文件格式。这让它极易与化学信息学和机器学习代码集成，也是今天做蛋白-配体 MD 最推荐的入门选择。

安装

mamba install -c conda-forge openmm openmmforcefields openff-toolkit pdbfixer
python -m openmm.testInstallation      # 检查可用平台与速度

testInstallation 会报告 CUDA / OpenCL / CPU 各平台的可用性与相对速度，装完先跑一遍确认 GPU 被识别。

一个完整的蛋白模拟脚本

from openmm.app import *
from openmm import *
from openmm.unit import *
from sys import stdout

# 1) 读结构并补齐（先用 PDBFixer 修复缺失原子）
pdb = PDBFile("protein_fixed.pdb")
ff  = ForceField("amber14-all.xml", "amber14/tip3pfb.xml")

# 2) 加溶剂盒子与抗衡离子
modeller = Modeller(pdb.topology, pdb.positions)
modeller.addSolvent(ff, model="tip3p", padding=1.0*nanometer,
                    ionicStrength=0.15*molar)

# 3) 建体系
system = ff.createSystem(modeller.topology,
                         nonbondedMethod=PME,
                         nonbondedCutoff=1.0*nanometer,
                         constraints=HBonds)
system.addForce(MonteCarloBarostat(1*bar, 300*kelvin))

# 4) 积分器与模拟对象
integrator = LangevinMiddleIntegrator(300*kelvin, 1/picosecond, 2*femtoseconds)
sim = Simulation(modeller.topology, system, integrator)
sim.context.setPositions(modeller.positions)

# 5) 能量最小化 → 平衡 → 产出
sim.minimizeEnergy()
sim.reporters.append(DCDReporter("traj.dcd", 5000))
sim.reporters.append(StateDataReporter(stdout, 5000, step=True,
    potentialEnergy=True, temperature=True, volume=True, speed=True))
sim.step(50000)          # 100 ps 平衡
sim.step(5000000)        # 10 ns 产出

参数为什么这么选

设置	取值	理由
时间步长	2 fs	配合 constraints=HBonds 固定含氢键长；不加约束只能用 1 fs
非键截断	1.0 nm + PME	PME 处理长程静电，截断太小会有伪影
盒子 padding	1.0 nm	保证蛋白与其周期镜像不直接接触
离子强度	0.15 M	近似生理条件
积分器	LangevinMiddle	控温稳定，比老式 Langevin 精度更好
系综	NPT（加 Barostat）	平衡阶段让密度自然收敛

带小分子配体时

蛋白有现成力场，小分子没有，必须先做参数化。标准做法是用 OpenFF（见 202《OpenFF Toolkit》）生成配体参数再与蛋白力场合并：

from openff.toolkit import Molecule
from openmmforcefields.generators import SMIRNOFFTemplateGenerator

lig = Molecule.from_file("ligand.sdf")        # 需含正确的键级与电荷
gen = SMIRNOFFTemplateGenerator(molecules=lig)
ff  = ForceField("amber14-all.xml", "amber14/tip3pfb.xml")
ff.registerTemplateGenerator(gen.generator)

配体的质子化态、互变异构体和电荷必须先确定好——参数化不会替你纠正化学上的错误输入，这是配体 MD 最常见的失败源。

常见问题

• 体系炸了（NaN / 粒子飞出）：几乎总是初始结构有原子重叠或缺失。先用 PDBFixer（见 215《PDBFixer》）修复，再充分最小化，必要时用位置约束分阶段升温。
• 模拟时长不够：10 ns 只能看局部稳定性；构象变化、配体解离需要百纳秒到微秒量级。用太短的轨迹下构象结论是常见错误。
• 只跑一条轨迹：单条轨迹的采样极其有限，重要结论应跑 3～5 条独立重复（不同初速度种子）。
• 忘了平衡：直接从最小化后开始收集数据，前段属于非平衡态，分析时要丢弃。

延伸资源

• 教程：015《OpenMM 教程》；力场参数化：202《OpenFF Toolkit》、014《OpenFF 教程》；
• 结构修复：215《PDBFixer》；轨迹分析：205《MDAnalysis》、206《MDTraj》；
• 自由能计算：201《OpenFE》、203《Perses》、204《YANK》。