GROMACS 从轨迹逐帧分析距离约束

下载输入、真实结果与绘图资料包，解压进入资料包目录。inputs 可修改，reference 是集群已计算的参考结果；脚本每次在 work 下建立独立目录。计算在分配好的资源内执行，可用 sbatch run.sbatch 跑整组，结果路径写入 LAST-RUN.txt。

什么时候用这篇#

EDR 与 NMR 输出教程说明如何在 MD 时记录距离约束。本篇从已经保存的坐标轨迹逐帧分析拓扑中的距离约束，适合看什么时候违反、违反多少。

本例使用同一份柔性水实际 MD：40 fs、201 帧。参考包另有宽阈值与多距离对分析控制，它们是对同一轨迹的后处理，不是重新施加这些约束后运行的新 MD。

bash run.sh disre

1. 选择分析入口#

module purge
module load gromacs/2026.4-cpu-double
gmx-disre-compat -h

集群提供独立的 gmx-disre-compat，参数沿用 disre。它由固定 v2026.4 工具源码、GNU13.3 和 CPU 双精度库构建；程序显示 2026.4-dev。原 gmx、gmx_d 和库保持原入口；此兼容命令用于本篇的距离约束后处理。源码、差异与 LGPL 许可在安装目录 share/source/disre-compat，也随本包提供源码。

2. 准备匹配的拓扑与轨迹#

mkdir my-disre
cp inputs/restraint-ordinary/* my-disre/
cd my-disre
gmx-disre-compat -s analysis.tpr -f md.trr -ds sum.xvg -da average.xvg -dn count.xvg -dm maximum.xvg -l details.log

analysis.tpr 包含实际原子顺序与距离约束参数，TRR 的原子顺序必须相同。只用 GRO 不足以恢复约束的上下限、分组和类型；也不要把另一个体系的 TPR 配过来。

原单距离对参数为：

; ai aj funct index type low up1 up2 fac
2 3 1 0 1 0.0 0.149 0.18 1.0

3. 读取四类逐帧输出#

head -25 sum.xvg
head -25 average.xvg
head -25 count.xvg
head -25 maximum.xvg
cat details.log

分别对应总违反量、平均违反量、违反的约束数与最大违反量。单距离对本例的瞬时违反是 max(rHH−0.149,0)，单位 nm。不要把瞬时平均和 <r⁻³>⁻¹/³ 得到的有效平均距离混为一谈。

4. 四组可复现控制#

• ordinary：普通周期体系中的单个 H—H 距离约束。
• periodic：相同坐标，以 periodic-molecules 对应拓扑处理周期边界。
• wide：把 H—H 上限改成 0.2 nm，同一轨迹所有帧都不违反。
• multi：两个 O—H 距离组成一个含等价对的分组，再加单独 H—H 约束。

分组例的拓扑条目是：

1 2 1 0 1 0.084 0.086 0.18 1.0
1 3 1 0 1 0.084 0.086 0.18 1.0
2 3 1 1 1 0.0   0.149 0.18 1.0

这一分组中两个距离的有效量为 (rOH1⁻⁶+rOH2⁻⁶)⁻¹/⁶，不是两条键长的算术平均。四组每组 201 帧，均与直接从 TRR 坐标计算的违反量一致，最大差小于 1e-7 nm。

修改分析上下限可以研究轨迹对候选约束的满足程度，但不能推断“如果 MD 当时用了这些新约束，轨迹也会保持一样”。要研究改变约束后的动力学，须用相应拓扑重新运行。

选项含义见 GROMACS disre；集群命令使用本页注明的独立兼容入口。