如果让我用一个词概括 GROMACS 模拟里最容易翻车的地方我会选文件而不是力场。.tpr、.xtc、.edr、.cpt这四类文件贯穿了从 grompp 生成输入、mdrun 生产轨迹、再到期后分析的完整链条。很多人跑模拟只盯着 log 里的温度和压力却不知道轨迹碎片化、续跑失败、能量统计跑偏根子都在这几个文件的生成和使用细节上。这篇文章不打算展开讲分子模拟理论只讲文件四种文件各管什么、怎么生成、怎么校验、怎么在分析和续跑时避开最常见的坑。适合刚接触 GROMACS 的学生也适合跑了好几年但偶尔还在文件上栽跟头的老人。1. 先把四种文件的职能盘清楚模拟流程才不会乱1.1 一次模拟会产出哪些文件很多新手第一次跑完 mdrun看到目录里冒出来一大堆扩展名直接懵了。我先给一张全家福把每个文件的角色说清楚。文件全称/内容生成者典型用途.tprrun input拓扑、力场参数、坐标、速度、盒子、全部 mdp 参数gromppmdrun 的输入也是后续所有分析工具对齐系统的标准参照.xtc压缩轨迹只含坐标有损压缩mdrun结构分析、RMSD、距离、氢键等.trr全精度轨迹坐标速度力无损mdrun需要速度/力的分析文件极大.edr能量数据能量项、温度、压力、密度、体积等mdrun热力学量分析、平衡判断.cpt检查点完整模拟状态mdrun 周期写入断点续跑、扩展模拟.gro坐标文件末尾帧坐标mdrun下一阶段 grompp 的输入.log运行日志与性能统计mdrun排查崩溃、查看步数.mdp模拟参数用户编写grompp 输入.top拓扑用户编写或 pdb2gmx 生成grompp 输入这张表里最容易被忽略的是.tpr和.cpt的锚点属性。.tpr是空间上的锚点所有分析都要靠它知道这个体系由哪些原子组成、原子叫什么名字、力场参数是什么.cpt是时间上的锚点它记录了模拟进行到哪一步、当时体系处于什么状态。这两个文件一旦丢了或弄混后面所有操作都会出问题。1.2 文件之间的依赖顺序正确的数据流是这样的用户准备.mdp参数、.top拓扑、.gro坐标。执行gmx grompp把三者打包成.tpr。执行gmx mdrun读入.tpr产出.xtc、.edr、.cpt、.log和最终坐标.gro。分析阶段用.tpr .xtc做结构分析用.edr做热力学量分析必要时用.cpt续跑。我见过最典型的错误是把.tpr当成一次性用品。grompp跑完生成md.tprmdrun跑完有人嫌文件乱直接把md.tpr删了。等到要分析 RMSD 时发现gmx rms -s没有参照物只能重新 grompp 一个新 tpr。如果记忆中的力场版本、加氢方式、盐浓度有偏差新 tpr 和旧 xtc 虽然是同一套原子细节上已经有微妙差别分析结果的说服力就打了折扣。所以我的第一条铁律是模拟一旦开始生成这个模拟的.tpr就要永久保留和.xtc、.edr放在一起归档。2. .tpr模拟参数的封存现场grompp 的两个细节别跳过2.1 grompp 到底把什么装进了 tpr从功能上说grompp 干的事相当于把散落的零件组装成一台只能执行固定动作的机器。你写的.mdp里每一项参数.top里的每一根键、每一个电荷.gro里的坐标和盒子信息全都会被写进二进制的.tpr文件里。此后mdrun运行期间GROMACS 完全按照.tpr里固化下来的参数执行.mdp文件后续改了什么它一概不认。实际命令通常是gmx grompp -f md.mdp -c npt.gro -t npt.cpt -p topol.top -o md.tpr -r npt.gro几个容易忽视的细节-c是输入坐标-r是位置限制的参考结构。如果体系里有位置限制比如把蛋白质重原子限制在初始位置而-r没给或给错了文件grompp 会把你-c提供的坐标当作参考。如果-c是续跑后的npt.gro那问题不大但如果你把-c换成别的结构比如想换一个初始朝向位置限制的参考也跟着换了平衡阶段蛋白就会被慢慢推往你给定的新结构这个偏差很隐蔽跑完看 RMSD 才会发现。-t是从 checkpoint 读取速度。从 NVT 接到 NPT 时用-t nvt.cpt能保留已有的原子速度和耦合状态避免重新从零初始化速度。很多人从 NPT 开始跑生产-t不带结果速度场重新随机化前面 NVT/NPT 的平衡白做了一半。-maxwarn只是忽略警告计数不是解决问题。grompp 遇到警告时默认会停下来。加-maxwarn 1或-maxwarn 10能跳过一些非致命警告但很多人把它当成万能钥匙遇到任何 warning 都-maxwarn 10压过去。这里头最常见的隐患是温度耦合组和索引组不匹配盒子尺寸与坐标文件不一致总电荷不为零这类警告它们往往意味着你的体系设置有问题压掉警告等于带着装错零件的机器出厂。2.2 换参数忘了重新 grompp是最隐蔽的翻车点我辅导过的学生里几乎每个人都犯过同一个错改完.mdp里的nsteps直接执行gmx mdrun -s md.tpr以为延长模拟只需要改步数就行。错得离谱。mdrun读的是.tpr不是.mdp。你改参数的一瞬间旧的md.tpr里固化的仍然是旧参数。你以为跑了 20 ns实际上mdrun还是按原来nsteps跑到时间就自然退出log 里也看不出异常因为 GROMACS 认为你本来就想跑那么长。正确的做法是修改.mdp后重新执行 grompp生成新的.tpr再用新的.tpr运行。但这里有个配套问题——如果你只是单纯想延长一段已经跑完的模拟直接重新 grompp 会生成一个全新的.tpr它和原来的.cpt在 GROMACS 的校验体系里对不上续跑就会报错。这个情况我在第 5 节专门讲正确工具是gmx convert-tpr不是 grompp。2.3 从 tpr 反查参数没留 mdp 也能验尸排错时最常遇到的问题是这个模拟当初到底用的什么参数如果 mdp 没进版本管理别慌.tpr里全都有。用gmx dump -s md.tpr | less能翻出完整参数表包括温度、压力、步长、每一步的时间、力场类型、非键作用参数、约束算法等等。我排查过几次为什么两组模拟结果差异巨大的案例最后都是靠gmx dump -s对比发现其中一组 grompp 时用错了 mdp 文件或者用了老版本的力场。所以排查的第一步永远是 dump 出 tpr 看参数而不是盯着 xtc 里的奇怪构象猜原因。另外gmx check -s md.tpr会做基础一致性检查能发现拓扑和坐标之间的原子数量不一致等问题。我习惯在归档前跑一遍这个命令比把文件放三个月后再找问题省事得多。3. .xtc压缩轨迹用起来爽分析前这几道手续不能省3.1 xtc 和 trr 的取舍.xtc是 GROMACS 默认输出的坐标轨迹它做了有损压缩只保留坐标不保存速度和力。精度通常在千分之一纳米量级对绝大多数结构分析RMSD、距离、氢键、接触面积完全够用。.trr则是什么都存坐标、速度、力全都有而且是无损全精度文件体积轻松比.xtc大两个数量级。我的建议很实际默认跑.xtc就够只有两种情况必须开.trr——一是你要用gmx trajectory做需要速度协方差的分析二是你要做精确的动能/温度相关分析或者需要精确还原坐标。.trr的开法是在.mdp里设置nstxout、nstvout、nstfout的步长但注意频率不要设太高否则一份几微秒的轨迹能把硬盘写爆。如果只是偶尔需要一帧的完整信息trjconv -dump单帧导出可比全程开.trr划算得多。3.2 PBC 处理为什么分析的轨迹经常是碎的这是分析阶段最高频的问题。MD 模拟开了周期边界PBC水盒子的原子流出左边界就会从右边界流回来坐标文件里记录的是折叠回盒子后的位置。于是在.xtc里直接画图或算距离会看到一条完整的蛋白链被拦腰截断或者两个原子明明在物理上靠得很近坐标却差出半个盒子。解决办法是在分析前统一做一次 PBC 处理printf Protein\nProtein\n | gmx trjconv -s md.tpr -f md.xtc -o md_pbc.xtc -pbc mol -center-pbc mol会把每个分子恢复完整后再决定怎么回卷-center让蛋白一直待在盒子中央方便可视化。执行时 trjconv 会要求你输入两次组名一次是居中组一次是输出组上面命令里的两个Protein就对应这两次输入。不同分析场景要选不同的-pbc选项我整理成速查表你要做的分析推荐选项原因RMSD、RMSF、距离、氢键-pbc mol分子恢复完整直接可比扩散系数MSD-pbc nojump消除原子跨盒子的跳变否则 MSD 被极大高估可视化出图-pbc whole只补完整不改写折叠逻辑适合渲染膜体系-pbc res或-pbc cluster防止脂分子和蛋白被盒子切断、甩散这里我踩过一个实打实的坑做扩散系数时没用-pbc nojump结果水分子每跨一次边界MSD 就多出几乎半个盒子平方的贡献扩散系数算出来比文献值高了一个数量级。后来用-pbc nojump处理后再算数值立刻回到正常范围。所以别嫌多一步脏活这一步直接决定定量分析对不对。3.3 轨迹拼接和完整性校验跑长模拟被集群作业超时中断是常态你可能会得到md_part1.xtc、md_part2.xtc好几段轨迹。拼接用gmx trjcat -f md_part1.xtc md_part2.xtc -o md_all.xtc如果两段轨迹在时间上有重叠比如因为续跑时没正确 append直接拼接会出现重复帧后面算时间平均时会莫名其妙地偏向重叠区。所以拼接前先各跑一遍gmx check -f看每段轨迹的起始和结束时间确认没有重叠再拼。gmx check还有另一个妙用核对轨迹和能量文件时间轴是否对齐。gmx check -f md.xtc -f2 md.edr它会打印两个文件各自覆盖的时间范围如果轨迹产出 10000 步而能量文件只记到 9000 步说明模拟在最后阶段被异常中断能量文件缺了尾巴后续算温度平均时就得小心。这条检查我建议每次分析前必跑成本极低收益极高。4. .edr热力学量的体检表别只会看温度一条线4.1 先澄清一个命名误会如果你是搜到这篇文章的可能查过.edr相关的检测卸载占用资源这些词。这里必须说明那些词指的是终端安全防护工具和 GROMACS 的.edr没有任何关系。GROMACS 里.edr是 energy data file 的缩写一个二进制格式的能量数据库每次 mdrun 都会自动产出记录每个输出步的能量项、温度、压力、密度、体积等状态量。它只能被 GROMACS 自己的工具读取不要试图用文本编辑器打开。4.2 gmx energy 的正确打开方式分析能量文件的主工具是gmx energy。基本用法gmx energy -f md.edr -o temperature.xvg -b 1000 -e 5000 -xvg none执行后工具会列出所有能量项每个项前面有编号。输入编号可以同时选多个用空格分隔最后输入0结束选择。结束后它不只导出数据还会在屏幕上打印一张统计表包含每个能量项的平均值、标准误差、RMSD 涨落和总漂移。这张表很多人不看其实它是判断平衡是否到位最直接的证据。选中 Temperature 时要注意如果你的 mdp 里设了多个温度耦合组比如tc-grps Protein Non-Protein能量项列表里会有Temperature-Protein、Temperature-non-Protein和总的Temperature好几个相似项。别随手选第一个要想清楚你关心的是哪个对象的温度。我就见过有人把蛋白质组的温度当体系温度汇报差了好几 K 还没发现。-b和-e参数的单位是皮秒ps用来跳过平衡段。生产模拟一般是前面几百 ps 平衡后面的几十 ns 才能用于统计。不加-b直接全段平均平衡期的温度弛豫会把平均值拉偏等温线看起来达不到目标温度很多人误以为温控坏了其实只是统计范围错了。4.3 最容易骗到自己的三个数值细节第一压力只看瞬时值或短区间平均没有意义。水盒子里的瞬时压力波动经常在 ±200 bar 量级这是正常的机械涨落。要判断 NPT 平衡是否达到目标压力 1 bar至少取几千步以上的平均看gmx energy统计表里的平均压力和总漂移。如果平均值在 1 bar 附近且漂移很小才算真正平衡。第二能量单位是 kJ/mol 不是 kcal/mol。同样一个数数值上差 4.184 倍。写成文章时如果想用 kcal/mol记得除以 4.184并注明换算关系。压强单位是 bar密度单位是 kg/m³换算成常用的 g/cm³ 要除以 1000。第三看能量图时重点不是单个时刻的势能绝对值而是它的漂移趋势。平衡良好的体系势能应该围绕一个稳定值上下小幅波动如果势能一路下行不带回头说明体系还在缓慢结构调整这时做的任何时间平均都不可靠。用gmx energy选 Potential 导出后直接看曲线的后 1/3 是否平比看平均值更直观。5. .cpt续跑和扩产的正确姿势别让几千步白跑5.1 checkpoint 里到底存了什么.cpt是模拟的完整快照所有原子的坐标、速度、力盒子向量每一步的能量累加器还有热浴和压浴的耦合器内部状态连随机数生成器的状态都记在里面。这就是为什么它比.gro大得多也是为什么只有它能做真正的无缝续跑。如果没有.cpt你只能从md.gro最终坐标重新开始但原子速度需要重新随机初始化温度耦合器的历史状态也没了。后续几万步体系会重新经历一段遗忘旧状态的过程这段轨迹的统计价值和连续模拟相比要大打折扣。所以在集群上跑任务.cpt比.xtc还娇贵丢了它前面几百度 CPU 小时基本白烧。mdrun 默认每隔一段时间自动写一次 checkpoint默认约 15 分钟可用 mdrun 的-cpt参数修改或在.mdp里设置nstcheckpoint按步数控制。作业被中断后找到最新的.cpt就能续跑。5.2 续跑命令与 append 的讲究最常见的续跑命令gmx mdrun -deffnm md -cpi md.cpt -append -v-cpi指定 checkpoint 文件-append让新产生的轨迹和能量直接追加到已有的.xtc、.edr后面时间轴保持连续。在 GROMACS 较新版本里-append是默认行为但我习惯显式写出来因为这样命令本身就能说明意图半年后再回看命令历史也一目了然。如果你忘了写-cpimdrun 会闷头从 step 0 重新跑一遍而且日志里不会主动提醒你你本来应该续跑。等你看时间线才发现不对已经又跑掉几万步。另一个常见失误是在有旧输出文件的情况下直接重跑GROMACS 可能会拒绝覆盖现有轨迹文件。这时候不要急着删文件先用gmx check看看旧轨迹的时间范围确认没有价值再清理。5.3 扩展模拟用 convert-tpr而不是重新 grompp想延长一段已经跑完或正在跑的模拟新手会改.mdp里的nsteps再 grompp然后用新.tpr加旧.cpt续跑。GROMACS 会报错告诉你 tpr 和 checkpoint 不匹配——因为两个 tpr 的参数指纹对不上系统拒绝把旧状态硬塞进新参数模板里。正确做法是用gmx convert-tpr只修改停止条件gmx convert-tpr -s md.tpr -extend 100000 -o md_ext.tpr gmx mdrun -s md_ext.tpr -cpi md.cpt -deffnm md -append -v-extend的单位是 ps含义是在原有停止时间基础上延长这么长时间。如果不喜欢相对延长可以用-until直接指定绝对停止时间。这个工具只改停止时间其他所有参数都保持原样所以生成的md_ext.tpr能被md.cpt接受。同理如果你希望跑完的体系在相同条件下多跑几段每次都该走convert-tpr路线而不是重新 grompp。重新 grompp 生成的 tpr 和旧 checkpoint 没有血缘关系续跑必然失败。5.4 续跑后核对时间轴续跑完成后别急着分析。先做两件事一是打开md.log看开头有没有 Restarting from checkpoint 之类的记录确认读对了 checkpoint二是跑gmx check -f md.xtc -f2 md.edr看轨迹和能量的时间轴是否都延伸到预期终点。如果 xtc 只到一半而 edr 到了终点说明续跑过程中轨迹写入出了状况可能需要回去找更早的 checkpoint 重跑。另外提醒一点-append续跑会把新轨迹追加进同一个文件这个文件在续跑前后的完整性依赖 GROMACS 内部的写入逻辑正常情况下安全。但如果你的集群环境在续跑过程中又崩了一次那就以最新的.cpt继续下一次续跑即可不要手动去改 xtc 文件。6. 一条龙实例从 grompp 到能量统计的完整命令链6.1 一套可直接套用的三阶段流程假设你已经从 pdb2gmx 得到了protein.gro和topol.top以下是我常用的完整流程# NVT 平衡 gmx grompp -f nvt.mdp -c protein.gro -p topol.top -o nvt.tpr gmx mdrun -deffnm nvt -v # NPT 平衡继承 NVT 的坐标和速度 gmx grompp -f npt.mdp -c nvt.gro -t nvt.cpt -p topol.top -o npt.tpr -r nvt.gro gmx mdrun -deffnm npt -v # 生产模拟继承 NPT 的坐标和速度 gmx grompp -f md.mdp -c npt.gro -t npt.cpt -p topol.top -o md.tpr -r npt.gro gmx mdrun -deffnm md -v生产跑完后标准分析链# 1. PBC 处理 printf Protein\nProtein\n | gmx trjconv -s md.tpr -f md.xtc -o md_pbc.xtc -pbc mol -center # 2. 骨架 RMSD printf Backbone\nBackbone\n | gmx rms -s md.tpr -f md_pbc.xtc -o rmsd.xvg # 3. 只统计生产段的温度假设前 1000 ps 是平衡段 printf Temperature\n | gmx energy -f md.edr -o temp_prod.xvg -b 1000 -e 10000 -xvg none这套流程我把每个阶段的-t都显式带上了就是为了让速度场和耦合器状态一路继承下去。很多人省掉-t也能跑但模拟前段会有一段重新平衡的尾巴等于每换一次系综就浪费一部分计算资源。6.2 常见报错速查表为了让你排查时不抓瞎我把这些年见过的高频问题整理成了一张表现象/报错根本原因处理办法Mismatch between checkpoint and tpr续跑用的 tpr 不是产生这个 cpt 的那个 tpr找回原 tpr想延长时间用gmx convert-tpr不要重新 grompplog 里显示 From step 0而你以为在续跑忘了加-cpi中断后立即用-cpi 最新的.cpt -append重启Group Backbone not found索引文件里没有这个组用gmx make_ndx -f md.gro生成 index.ndx分析命令加-n index.ndx轨迹里分子碎成几段没做 PBC 处理trjconv 加-pbc mol扩散分析用nojump温度涨落几百 K怎么都压不下平衡段没截掉或选错了温度耦合组gmx energy加-b-e只统计生产段核对温度组名grompp 警告一堆-maxwarn 10压掉后跑完结果离谱maxwarn 只是忽略警告不解决问题逐条读警告尤其是盒子尺寸、温度耦合组、电荷总量这几类6.3 我自己的归档与自检习惯最后分享几个已经形成肌肉记忆的习惯。我在集群上跑生产模拟提交脚本里会写一段自动续跑逻辑如果存在md.cpt就执行mdrun -s md.tpr -cpi md.cpt -append如果不存在才从头开始跑。这样作业排队超时、节点被杀下次重新提交时自动从断点接着跑不用人肉盯。每次生产模拟结束后我会把md.tpr、md.cpt、md.gro、md.xtc、md.edr五件套放在一个独立目录并设为只读然后跑一遍gmx check -f md.xtc -f2 md.edr确认时间轴完整。归档前我会顺手用gmx dump -s md.tpr | grep -i nsteps之类的方式再确认一次步数和时间设定。这些小动作看起来琐碎但能挡住绝大多数模拟白跑的惨剧。写完这些回头想想这四种文件其实对应着模拟的四个关键词.tpr是确定性.xtc是采样.edr是验证.cpt是可持续。能把这四个关键词吃透GROMACS 这条路你会走得比大多数人稳。