3.1 mdp 参数逐项精读

3.1 mdp 参数逐项精读

本节摘要：mdp 文件是模拟的"导演剧本"：每一个键值对都是一次物理决策——步长多长、哪个系综、静电怎么算、多久存一次盘。本节给一份带逐行注释的 NVT 平衡 mdp 与一份成品 mdp，按"积分、系综、邻居与截断、静电、输出"五组拆解；再配一个差异对比脚本，让两份文件之间的每个不同都无处藏身。读完你应当能把任何一份来路不明的 mdp 讲清楚，而不是复制粘贴了事。

上一章把体系装进了盒子，本章从剧本开始走台。mdp 之所以值得用一整节精读，是因为它是 GROMACS 里唯一"不写就取默认、写了就必须兑现"的文件：grompp 把它原样编译进 tpr，运行期不再解释。一个 dt 写错数量级的 mdp，程序不会提醒你——它只会忠实地把体系炸掉。本节的读法是先分组后逐项，最后用脚本做"差异审计"。

一份带注释的平衡 mdp

以下是 NVT 平衡阶段的全参数文件，注释即讲解（生产 mdp 只列差异）：

; ===== 平衡段 nvt.mdp =====
define                   = -DPOSRES     ; 打开拓扑里 #ifdef POSRES 的位置限制（3.6 详述）

; ---- 运行长度与积分 ----
integrator               = md           ; 蛙跳/Velocity Verlet 家族，3.2 讲透
dt                       = 0.002        ; ps；2 fs 的前提：全体系含氢键已约束
nsteps                   = 250000       ; 250000 × 0.002 ps = 500 ps 平衡
comm-mode                = Linear       ; 每步扣除质心整体平动，防"自由飘移"

; ---- 邻居与截断 ----
cutoff-scheme            = Verlet       ; 现代默认，配 nstlist 自适应
nstlist                  = 20           ; 每 20 步（40 fs）重建邻居表
rlist                    = 1.0          ; 邻居表半径 nm；Verlet 方案下程序自动加缓冲
vdwtype                  = Cut-off      ; LJ 用简单截断（LJ 衰减快，误差小）
rcoulomb                 = 1.0          ; 静电实空间截断 nm，与 PME 配套
rvdw                     = 1.0          ; LJ 截断 nm
DispCorr                 = EnerPres     ; 长程色散尾校正：修能量与压强（略改体积涨落）

; ---- 静电 ----
coulombtype              = PME          ; 3.4 的主角
pme_order                = 4            ; 四阶插值（三次样条级精度）
fourierspacing           = 0.12         ; 倒易格点间距 nm，程序据此定网格

; ---- 温度耦合 ----
tcoupl                   = V-rescale    ; 随机速度重标定，平衡与生产都可用
tc-grps                  = PEO SOL NA CL ; 分组耦合：聚合物与溶剂/离子分开控温
tau_t                    = 0.5 0.5      ; 弛豫时间 ps；太小=过阻尼，太大=控不住
ref_t                    = 350 350      ; 目标温度 K；PEO 熔体研究常取 350 以上

; ---- 压力耦合（NVT 阶段关掉）----
pcoupl                   = no

; ---- 约束 ----
constraints              = h-bonds      ; 所有含氢键全约束，2 fs 步长的门票
constraint-algorithm     = LINCS        ; 3.2 末尾与 4.3 都会再见到它

; ---- 输出频率 ----
nstxout                  = 5000         ; trr 坐标每 10 ps 一帧（平衡段够用）
nstvout                  = 5000
nstenergy                = 500          ; 能量每 1 ps 一条，涨落分析靠它
nstlog                   = 1000
nstxout-compressed       = 1000         ; xtc 每 2 ps 一帧
compressed-x-precision   = 1000         ; xtc 精度 1/1000 nm
​

成品 md.mdp 与它的差异通常只有五处：define 行删掉（撤掉位置限制）、pcoupl 换成 C-rescale 或 Parrinello-Rahman、nsteps 拉长到亿步级、nstxout-compressed 加密到每 5 ps 或更细、nstcalcenergy 与 nstlist 视性能微调。差异越少越好审——这正是下面脚本的用武之地。

差异审计脚本

两份 mdp 摆在一起人眼比对容易漏行，写个脚本让差异自己站出来：

def parse_mdp(path):
    params = {}
    for raw in open(path, encoding="utf-8"):
        line = raw.split(";")[0].strip()      # 分号后是注释
        if not line or "=" not in line:
            continue
        key, val = line.split("=", 1)
        params[key.strip().lower()] = val.strip()
    return params

def diff_mdp(path_a, path_b):
    a, b = parse_mdp(path_a), parse_mdp(path_b)
    all_keys = sorted(set(a) | set(b))
    print(f"{'参数':<22}{'A 文件':<18}{'B 文件':<18}备注")
    risky = {"dt", "integrator", "tcoupl", "pcoupl", "rcoulomb", "constraints"}
    for k in all_keys:
        va, vb = a.get(k, "(无)"), b.get(k, "(无)")
        if va == vb:
            continue
        tag = "  <-- 高风险项，务必人工确认" if k in risky else ""
        print(f"{k:<22}{va:<18}{vb:<18}{tag}")
    if a == b:
        print("两份文件完全一致")

diff_mdp("nvt.mdp", "md.mdp")
​

对上面那对文件，脚本会揪出 define、pcoupl、nsteps、nstxout-compressed 四行，其中 pcoupl 属于高风险项——它一变，系综就从 NVT 换成了 NPT，平衡判据与涨落语义全部随之改变。把 diff 输出贴进实验记录，是"可复现科研"最便宜的一步。

参数背后的三个"为什么"

为什么 dt 是 0.002？ 不是习惯，而是约束与振动的联合决定：不约束含氢键时最快振动周期约 10 fs（3.2 推导），稳定步长上限约 1 fs；全约束后最快振动变成 O–C–O 角振动（周期约 20 fs），2 fs 才有安全余量。步长翻倍对产能是翻倍，所以"全约束 + 2 fs"是全行业的默认交易。

为什么 tc-grps 要分组？ 单组控温时溶剂浴温由恒温器直接保证，溶质（尤其远离溶剂的大链段）可能通过弱耦合"偷凉"或"偷热"，产生所谓溶剂-溶质温差。分组各配 tau_t 是便宜且有效的防波堤；代价是组间能量流的核算要靠 gmx energy 的 ΔH 项，分析时别漏。

为什么 nstlist 能到 20 步？ Verlet 截断方案给邻居表加了动态缓冲（rlist 比 rvdw 大出一截），原子在两次重建之间跑不出缓冲带就不会漏对。缓冲策略由 verlet-buffer-tolerance 控制误差预算（默认每原子每步 0.005 kJ/mol），程序自己把 rlist 加到够用——你只管把 nstlist 定在 20 到 40，别再手调 rlist 数值（那是 Verlet 之前的旧习惯）。

易错点与本节速记

三处高频翻车：tc-grps 写了组名却没写对应的 tau_t 与 ref_t（组数必须三处一致，grompp 会拦）；NVT 的 mdp 忘了 pcoupl = no，程序会当你要 NPT 平衡； dispersion 校正 DispCorr 在压强敏感的体系（脂膜！）里一开一关直接改变平衡面积——7.2 的案例里它是第一个被审问的参数。

速记：mdp 是物理承诺书，每个默认值都有来历；五组记忆法——积分（integrator/dt/nsteps）、系综（tcoupl/pcoupl/tau/ref）、邻居截断（nstlist/rlist/rcoulomb/rvdw）、静电（PME 三件套）、输出（nstxout 系列）；改任何参数前先问它改的是什么承诺。下一节把 integrator 与 dt 这组承诺讲透。