蛋白结构准备这件事，说它是整个分子动力学模拟里最不性感的一环，大概没什么人会反对——没有轨迹图，没有收敛曲线，只有一堆文本文件和满屏 warning。但我见过的绝大多数“跑飞”案例，根子都埋在这一步。上个月 A 同学拿着一套跑了两周的轨迹来找我，RMSD 从三百帧往后一路往上爬，配体整个翻出结合口袋，他第一反应是力场参数不对、积分步长太大、恒温器设置有问题。我让他把 prmtop 打开，查电荷那一栏的总和，结果是 -0.99983。tleap 当时只丢了一句 warning，谁也没当回事。这类问题在 Amber 环境下的蛋白结构准备阶段极其典型：工具链本身很成熟，但它默认你会读日志、会核对中间结果，一旦你把它当成“一键转换”，它就会安静地把错误带进后面所有的计算里。

这篇内容围绕Amber 分子动力学模拟 5.1 版本的蛋白结构准备流程展开，核心是讲清楚从一份原始坐标文件到一份能进动力学的拓扑与坐标，中间到底要过几道关，每一关的工具在做什么、为什么这么做、哪里最容易翻车。适合刚接手体系搭建的新手照着走一遍，也适合跑了很久但总在排查“莫名其妙的能量爆炸”的老手对照检查。文中涉及的体系、残基编号、文件命名都是示意性的，实际操作请以自己手上的体系为准。

1. 结构准备在模拟流程里的真实分量

1.1 一个体系从拿到坐标到能跑起来，中间隔了多少道工序

很多人对流程的想象是一条直线：下载坐标、跑 tleap、开始模拟。实际上的链条要长得多，而且每一环都有回退的可能。完整走一遍大概是：确定结构来源、清洗原始记录、处理交替构象与无序区、补全缺失重原子、判定质子化状态、处理末端残基与二硫键、给配体和辅因子做参数化、检查并修正原子命名、选择力场与水模型、构建溶剂盒子、中和电荷并加盐、生成拓扑与坐标、逐项自检、最后才是能量最小化。这条链上有十几个决策点，任何一个点选错，后面的计算都要推倒重来。

我习惯把这个阶段的工时占比说得直白一点：一个中等复杂度的蛋白-配体体系，结构准备加上反复验证，通常要吃掉整个项目 40% 到 60% 的时间。听起来夸张，但如果你经历过一次“体系跑到第 3 纳秒突然出现 vlimit exceeded，回头发现是某个赖氨酸侧链缺了两个重原子被 tleap 随手补在了一个不可能的位置”，你就会认同这个比例。更麻烦的是，这类错误的潜伏期很长，前面几十步全都能正常通过，直到某个局部应力累积到阈值才爆出来，那时候你已经没有精力回溯了。

所以结构准备阶段的第一原则不是“快”，而是“每一步都留下可核查的证据”。我自己的做法是，凡是被修改过的结构都另存一份，凡是被工具自动补全的东西都在日志里标出来单独看一眼，凡是涉及电荷和原子数的数字都记在同一个文本文件里。听起来琐碎，但它把“事后追溯”的成本从几个小时压缩到几分钟。

1.2 Amber 5.1 里这套工具各自负责哪一段

先把分工说清楚，后面就不会乱。这一版套件里，结构准备主要靠几个命令行工具协同：pdb4amber负责格式清洗和命名规范化，reduce或者外部工具负责加氢和判定部分可滴定基团，antechamber与parmchk2负责小分子的原子类型与缺失参数，tleap（也就是 LEaP 的命令行前端）负责把所有零件组装成一个单元、补全缺失原子、溶剂化、加离子并输出拓扑，cpptraj和parmed则用来做生成后的核查。

这里有个容易被忽略的点：这些工具之间的“契约”是文件格式和命名规范，不是语义。tleap 看到HIS这个残基名，它不会去猜你想表达哪种互变异构体，它只会按内部映射表处理或者直接报错；antechamber 看到 mol2 里的原子名重复，它不会替你改名，它会把两个原子当成一个或者直接输出错误结果。也就是说，工具链假设你已经在输入文件里把该表达的语义表达清楚了。这个假设在教程数据上成立，在真实数据上经常不成立。

另外提醒一句，不同小版本之间，leaprc 的文件名、默认水模型、参数文件的检查严格程度都可能有调整。我见过有人直接复制网上的 tleap 脚本，结果因为加载顺序问题导致 gaff2 的参数没生效，配体所有原子都找不到类型，报错信息里却只写了“does not have a type”，排查了半天。所以脚本要逐行核对，不要整段照搬。

2. 整体思路拆解：先想清楚体系要什么，再动手改坐标

2.1 四条常见的结构来源与各自的坑

拿到一份坐标，先问它是从哪来的，这个信息决定了你后面要在哪里多花时间。

晶体衍射结构是最常见的一类。它的优点是重原子位置可信，缺点是几乎一定缺氢，侧链末端和表面 loop 经常有缺失，还会有交替构象（altloc）和多条链共用一套编号的情况。处理这类结构的重点是老老实实读一遍原始记录的 REMARK 段，看看作者标注了哪些残基没解析出来，不要指望工具帮你发现。

冷冻电镜结构这几年越来越多。骨架通常还行，但侧链的摆放质量参差不齐，尤其是表面残基。用这类结构做模拟，我的建议是先看局部分辨率或者模型质量评估数据，把置信度低的区域在后续分析里单独对待，别把那些区域的 RMSD 当作有物理意义的指标来解读。

预测模型类的结构，比如各类折叠预测工具给出的结果，优势是重原子齐全、几乎不缺残基，缺点是柔性区域和低置信度区段的构象基本不可信，而且这类模型往往不带氢、残基编号也可能与你的参考序列对不上。用它做模拟之前，先按置信度分数把结构切一刀，把分数极低的尾巴截掉或者用封端基团处理，比硬着头皮全保留要靠谱。

同源建模的结构则要注意 loop 区的构象是从模板或者从头生成的，建模工具本身不做能量优化，直接拿来跑动力学，前几百皮秒基本都在驰豫那些不合理的局部接触。

不管哪种来源，第一步都是把原始文件原封不动地存一份到raw/目录，所有后续操作都在副本上做。这条习惯帮我躲过好几次“改坏了想回头但已经覆盖”的窘境。

2.2 力场与水模型的搭配逻辑，这一步选错后面全白搭

力场选型是结构准备阶段最有“结构性”的一个决策，因为它一旦定下来，后面所有参数、所有对比实验都要跟着走。蛋白部分常用的两代选择是 ff14SB 和 ff19SB，前者成熟、文献基数大，后者对骨架二面角做了修正，在长时间尺度上的构象分布更合理。

但这里有个必须一起考虑的约束：ff19SB 的官方搭配是 OPC 水模型，而 ff14SB 的经典搭配是 TIP3P。很多人只换了蛋白力场，水还是原来那套三点模型，跑出来的结果当然也能收敛，但蛋白-水界面的描述和力场作者的验证条件已经不一致了，遇到跟实验对不上的时候，你连问题出在哪都说不清。所以这两个选择要绑定考虑，不要拆开。

配体的参数化路线也要同步定：GAFF2 是老牌方案，覆盖面广、参数生成快，配合 AM1-BCC 电荷在多数结合自由能计算里够用；如果做的是精细的构象分析或者电荷转移明显的体系，走 RESP 电荷会更稳，代价是要跑量子化学计算，一个中等大小的配体从构象搜索到电荷拟合，一两天是常态。

水模型的选择还会间接影响算力。四点模型比三点模型多了虚拟位点，积分时的约束处理稍复杂一些，相同原子数下的单步耗时略有增加；而盒子里的水分子数取决于你给的缓冲距离和盒子形状，这个影响比水模型本身大得多，后面讲盒子估算的时候会具体算。

2.3 决策顺序：质子化、末端、配体、溶剂、离子，一步都不能跳

顺序搞反是新手最容易犯的错。比如先溶剂化再加配体，那配体的参数就得在已经生成的水盒子里改，几乎等于重做；比如先加离子再改质子化状态，体系净电荷变了，离子数就得重新算。

我自己的固定顺序是这样：先确定要保留哪些链、哪些残基、哪些辅因子和金属，把无关的水分子、重复的链、结晶添加剂全部清掉；然后处理缺失重原子和末端残基，因为末端是否封端会影响后面算出来的净电荷；接着判定质子化状态和组氨酸异构体，这一步直接决定每个残基带多少电荷；再给配体、辅因子做参数化并单独验证；最后才是回 tleap 里组装、溶剂化、按净电荷中和并加盐、输出拓扑。

这个顺序的内在逻辑是“电荷从上游往下游传递”。质子化状态决定了蛋白的净电荷，配体参数化决定了配体的净电荷，两者加起来才知道要加多少反离子，而反离子数又必须在溶剂化之后加，因为加离子的实现方式是替换水分子。你把这个链条理顺了，tleap 里那些看起来莫名其妙的电荷报错基本都能自己对上号。

3. 从原始 PDB 到 prmtop：完整实操链路

3.1 第一步清洗：pdb4amber 能做什么、不能做什么

清洗这一步我用得最多的是 pdb4amber，典型调用是这样：

pdb4amber -i raw/complex_raw.pdb -o clean/complex_clean.pdb --dry --nohyd

它做的事情主要是几类：把原子名重排成 Amber 认得的写法（比如把甘氨酸的第二个氢改成 HA3 这类约定），处理交替构象只保留占据率最高的那一份，把非标准残基名标出来，去掉指定不要的内容。--dry是去掉结晶水，--nohyd是去掉原有氢原子。

关于结晶水要不要保留，我的立场很明确：一般不保留。晶体里的水氧位置来自衍射数据，但氢原子位置根本没有实验信息，而且力场里的水模型是自带特定几何的刚性或半刚性结构，你把它硬塞到实验位置上，等于给体系引入一堆几何不合理的初始接触。真要保留某个关键的结构水，做法是只留那个氧的位置，然后在 tleap 里用单独的水分子单元去对齐，而不是直接把整段结晶水搬过来。

--reduce这个选项要谨慎。它会在清洗阶段顺手做加氢和部分可滴定基团的判定，看起来很省事，但它的判定逻辑和你后面可能用的 pKa 计算工具（比如各种基于连续介质的经验算法）不是一套，两者给出的组氨酸互变异构体结论经常打架。我踩过这个坑：reduce 认为某个组氨酸是 HIE，我按 pKa 工具的结果手动改成了 HID，结果文件里两套氢原子名混在一起，tleap 加载时直接崩。所以要么全程以 reduce 为准，要么全程自己判定，不要混着来。

这里必须强调一点：pdb4amber 的--add-missing-atoms类选项在补全能力上是有边界的，它能处理一部分简单情况，但对缺失整段侧链、缺失多个连续残基这类问题无能为力，而且它补出来的坐标位置往往只是几何上的“合理猜测”，未必符合局部的氢键网络。我的做法是把 pdb4amber 的输出当成“清洗后的输入”，缺失部分单独用建模工具处理，处理完再回来走后面的流程。

3.2 缺失原子与末端残基：tleap 补出来的坐标为什么不能直接用

tleap 在loadPdb的时候如果发现残基缺原子，它会自动补上，日志里留一句 “Added missing heavy atom” 之类的提示。这个机制方便，但有两个致命问题。

第一，它补的只是坐标，不是构象。一个缺失的赖氨酸侧链，tleap 会沿着现有键长键角方向往外接，接出来的位置可能正对着另一个残基的侧链，形成严重的原子重叠。你在拓扑层面完全看不出问题，跑最小化的时候初始能量能到 10 的七八次方，然后要么报错退出，要么在几步之内把结构炸散。我教新人的第一课就是：最小化开始前先看初始能量，正常体系一般在几千到几万 kcal/mol 的量级，如果一上来就是百万级，别急着调参数，先回去查重叠原子。

第二，它不会告诉你补得对不对。日志里那一行很容易被淹没，尤其是体系大的时候。我习惯在 tleap 跑完之后固定做一次检索：

grep -n -i "added missing\|created a new atom\|does not have a type" leap.log

只要有任何一条 “added missing”，就回到结构上人工确认那个位置合不合理。

末端残基是另一个高频问题。tleap 默认给肽链的 N 端加上带电的氨基（三个氢），C 端加上带电的羧基（含 OXT）。如果你手上的结构是从一个完整蛋白中间截出来的片段，那么这个默认处理可能完全不符合真实情况——真实体系里那个位置可能是被别的结构域挡住、或者干脆是经过翻译后修饰的。这时候要么用乙酰基和甲酰胺基团封端，要么在建模时保留足够的上下文残基。

还有一个细节：当你的结构里出现残基编号不连续（比如 45 后面直接跳到 52），tleap 默认会认为这是一个断点，直接在两段之间不做连接，各加各的末端。如果你的本意是“这里缺了 6 个残基，应该补上”，那就必须在进 tleap 之前用建模工具处理掉；如果本意是“这就是两条独立的链”，那就要确认两条链的相对位置在空间上没有冲突。这两种情况的后续处理完全不同，别让它默认替你决定。

3.3 质子化状态与组氨酸异构体：一套工具说了不算

质子化状态是结构准备里最需要判断力的一环。默认 pH 7 附近，天冬氨酸和谷氨酸去质子化、赖氨酸和精氨酸质子化，这些没什么争议，真正麻烦的是组氨酸和半胱氨酸。

组氨酸有三种常见状态：HID（质子化在 δ 位）、HIE（质子化在 ε 位）、HIP（双质子化、带正电）。这三种状态的差别不只是多一个氢少一个氢，而是整个咪唑环的氢键给体/受体方向变了，直接影响它跟周围残基的相互作用模式。判定方法我一般走两步：先用经验 pKa 计算工具跑一遍全蛋白，拿到整体判断；然后对活性位点附近、金属配位附近、以及处在强氢键网络里的组氨酸逐个手工检查。

为什么要手工检查，因为经验算法对金属离子、强极性环境和埋藏残基的感知能力有限，它们本质上是在算溶剂可及性和静电贡献，对“这个氮到底有没有可能形成氢键”这种几何信息用得不够。我遇到过一个相当典型的案例：某个组氨酸在结构上明显跟相邻的天冬氨酸形成了一条短氢键，从几何距离看只能是一个作给体一个作受体，那么它的质子化位置就是确定的了，而经验算法给出的却是另一种状态。

半胱氨酸的处理同样要留神。游离的半胱氨酸需要加氢变成硫醇形式，而参与二硫键的半胱氨酸必须去氢并在 tleap 里显式成键。这一步的关键是显式：tleap 不会因为两个硫原子距离在 2.05 Å 左右就自动帮你连上，你要在脚本里写bond命令。反过来，如果你没做这个处理，两个硫原子会以非键形式相互作用，局部结构会在动力学里慢慢分开，几十纳秒后你才发现二硫键“断了”。

提醒：质子化状态定完之后，一定要重新核一遍体系的总电荷。后面所有关于加多少反离子的计算都建立在那个数字上，这个数字错了，后面全错。

3.4 配体参数化：antechamber + parmchk2 的手工细节

配体的参数化是另一个高发问题区，因为它涉及的信息传链条最长：三维坐标 → 原子类型 → 电荷 → 缺失参数 → 拓扑。

先从坐标说起。配体的三维坐标必须从复合物结构中提取出来，保留原始结合位姿，绝对不要单独建模。原因很简单：配体在溶液里的低能构象和在蛋白口袋里的构象经常不是同一个，你如果用一个真空优化的构象去对接位置，等于一开始就把配体摆错了。提取之后的配体文件要单独检查原子数、键连关系、有没有缺失的氢，特别是可旋转键上的氢。

参数化的典型流程是这样：

antechamber -i params/lig.pdb -fi pdb -o params/lig.mol2 -fo mol2 -c bcc -nc -1 -rn LIG -at gaff2 -s 2 parmchk2 -i params/lig.mol2 -f mol2 -o params/lig.frcmod -s gaff2

-c bcc指定用 AM1-BCC 方法生成电荷，速度快、质量够用；-nc -1必须填对，它告诉程序这个分子整体带多少电荷；-rn LIG是输出里的残基名，要和你 tleap 脚本里用的一致；-at gaff2指定原子类型体系。后面 parmchk2 的-s gaff2必须跟前面保持一致，否则生成的缺失参数文件里原子类型对不上，等于白做。

这里有三个实操点值得单独说：

第一，-nc填错是最常见的翻车原因，而且错得很有迷惑性。程序不会报错，它会按你给的电荷去拟合，出来的电荷总和确实接近你填的值，但每个原子的电荷分布是错的。所以流程上应该反过来：先自己数清楚分子的形式电荷（哪几个基团在这个 pH 下带电），再填这个数字，而不是让程序猜。

第二，电荷总和的核对。BCC 方法算出来的电荷总和不会精确等于你指定的整数，通常会差小数点后三四位，这个是可以接受的。但如果差了 0.02 以上，说明原子类型识别或者电荷拟合过程中出了问题，比如某个环被误判成了芳香性不同的形式。我一般用一个简单的核对脚本去读 mol2 的电荷列求和，不达标就回去查。

第三，配体里如果包含需要特殊处理的基团（比如羧酸根、磷酸基、磺酸基），GAFF2 的通用参数可能描述得不够好，这时候要专门查一下有没有现成的参数可以借用，或者考虑换用另一套覆盖更广的通用力场。不要在这上面省事，一个磷酸基的参数不准确，能让整条结合曲线的形状都变形。

配体如果带多个可滴定基团，还要考虑在不同 pH 下是否要准备多套参数，以及模拟的 pH 条件跟实验条件是否一致。这件事在设计阶段就要想清楚，不要等到跑完才发现选错了质子化态。

3.5 tleap 脚本逐行拆解与盒子尺寸估算

组装这一步把所有零件拼起来。给一份结构清晰、可以直接对照修改的脚本：

# 加载力场 source leaprc.protein.ff19SB source leaprc.water.opc source leaprc.gaff2 # 载入配体参数与坐标 loadAmberParams params/lig.frcmod LIG = loadMol2 params/lig.mol2 # 载入体系 sys = loadPdb clean/complex_fixed.pdb # 二硫键（按实际残基编号填） bond sys.45.SG sys.88.SG # 隐式溶剂用的半径集（显式溶剂可省） set default PBRadii mbondi3 # 检查净电荷 charge sys # 溶剂化：缓冲距离 12 埃，截角八面体盒子 solvateOct sys OPC 12.0 # 按净电荷中和（下面两行二选一，先看 charge 的输出再决定） addIons sys Na+ 0 # addIons sys Cl- 0 check sys saveAmberParm sys leap/complex.prmtop leap/complex.inpcrd savePdb sys leap/complex_leap.pdb quit

几个关键点解释一下。

加载顺序有讲究：先蛋白力场、再水模型、最后小分子力场，参数文件要在使用它的分子之前加载。顺序错了，某些参数可能被后加载的文件覆盖。

charge sys这一步千万别省。它会打印整个单元的净电荷，你根据这个数字决定加钠还是加氯。如果体系净电荷是负的，就加addIons sys Na+ 0；是正的，就加addIons sys Cl- 0。把两行都写上看起来更保险，实际上会让 tleap 在第二个命令里做一件你没预期的事，所以我习惯二选一，把另一行注释掉。

盒子尺寸的估算值得单独讲。solvateOct的第三个参数是溶质到盒子边界的最小距离，单位埃。这个值的选取逻辑是：周期性边界下，溶质不能通过最小镜像看到自己的周期副本，所以边界厚度至少要覆盖非键截断半径加上一定余量。如果你用的是 10 Å 的截断，边界给 10 到 12 Å 是通行做法，我一般固定用 12 Å，留一点余量。

盒子的形状影响算力，这点常被低估。截角八面体在同样缓冲距离下，里面的水分子数比立方盒子少两三成。举个具体的估算：一个最大跨度约 60 Å 的蛋白-配体复合物，缓冲 12 Å 之后，立方盒子的边长约 84 Å，体积约 59 万立方埃，按每个水分子占约 30 立方埃折算，约两万个水分子，也就是六万个水原子；换成截角八面体，水分子数通常落在 1.4 万左右，四万多个原子。再加上蛋白本身的原子，整个体系的规模差了两万原子，在同样的硬件上，这个差别直接反映在每天的采样量上。所以除非有特殊理由，我默认用截角八面体。

还有一点：溶剂化之后一定要用savePdb把 tleap 眼里的最终结构导出来看一眼。我见过太多次坐标文件在文本层面看着正常，但用可视化工具打开发现配体跑到了盒子外面、或者某个残基穿到了另一个周期副本里。这种问题在拓扑文件里是看不出来的，必须看图。

3.6 出拓扑之后的四道自检

拓扑生成了不等于可以开跑。我固定做四道检查，缺一道心里都不踏实。

第一道，电荷核对。prmtop 里的电荷是以特定单位存储的，直接把那一栏的数值加起来，再除以约定的换算系数，应该得到一个非常接近整数的值。可以写个小脚本，读 prmtop 里的电荷标志段，把数值抽出来求和：

import sys def read_charges(path): flag, vals = False, [] with open(path) as f: for line in f: if line.startswith('%FLAG'): flag = (line.split()[1] == 'CHARGE') continue if line.startswith('%FORMAT'): continue if flag: vals.extend(float(x) for x in line.split()) return vals q = read_charges(sys.argv[1]) # Amber 电荷单位换算系数 total = sum(q) / 18.2223 print(f"原子数: {len(q)} 净电荷: {total:.4f}")

跑出来如果离整数超过 0.01，就回去查质子化状态和配体电荷。

第二道，原子与残基计数。把 prmtop 里的原子数、残基数和你在可视化工具里看到的对上，特别是配体是不是完整、有没有哪个残基被拆成了两段或者多出几个原子。

第三道，初始几何检查。用 cpptraj 加载拓扑和坐标，跑一次原子间距检查，把距离过近的原子对找出来：

cpptraj -p leap/complex.prmtop -y leap/complex.inpcrd <<EOF check EOF

它会把距离异常近的原子对报出来。正常的晶体结构经过清洗之后，除了氢原子之外的重原子间距基本都在合理范围，出现 0.8 Å 以下的接触基本可以确定是补全或者坐标转换出了问题。

第四道，跑一个极短的最小化做验证。不用等完整流程，只要几百步最速下降，看看初始能量和能量下降的形态。如果第一步能量就在百万级，或者能量在下降过程中出现剧烈的非单调跳变，说明结构里有严重冲突，回去查。

4. 常见问题排查速查与避坑实录

4.1 电荷类报错：从 “unperturbed charge” 说起

tleap 最常被无视的一句提示是类似 “unperturbed charge of the unit: -0.500000” 这样的输出。它的意思是，这个单元的净电荷不是一个整数。对于要以周期性边界加 PME 做显式溶剂模拟的体系来说，非整数净电荷是个硬问题，因为你需要加反离子来中和，而反离子只能加整数个，剩下的 0.5 无处安放。

出现这个数字的原因通常有三类：质子化状态判错（某个本该带电的残基被当成中性，或者反过来）、配体电荷填错、以及某些残基的末端处理不当（比如链被截断在中间却没有封端）。排查的顺序我一般是：先按残基逐个核对可滴定基团的带电情况，再看配体的形式电荷，最后检查链的末端。

还有一种情况是新人不常遇到的：体系中存在多个配体或者多个金属离子，它们各自的电荷都是整数，但加起来跟蛋白的电荷凑成了半整数。这种时候不要急着改某个配体的电荷去凑整，那等于用一个错误掩盖另一个错误，应该回去确认每个组分的真实质子化状态。

4.2 原子名与残基名类报错

这类报错的信息通常很直白，但理解起来需要一点背景。比如某个原子 “does not have a type”，意思是在当前加载的所有参数文件里都找不到这个原子的类型定义。原因可能是：残基名不在力场库里（比如一个非标准残基被当成了标准残基）、配体的 frcmod 没加载、或者加载顺序让参数被覆盖了。

原子名冲突是另一种形态。antechamber 生成的 mol2 里如果出现重名原子，tleap 加载时可能只认其中一个，导致配体少原子，而报错信息里完全不会提到“少了一个原子”，只会说某个键的参数找不到。所以配体 mol2 生成之后，我固定会检查一遍原子名的唯一性。

还有一个细节是组氨酸的命名。Amber 体系里它必须写成 HID、HIE 或 HIP，不能是 HIS。如果你的结构文件里还留着 HIS，tleap 在部分版本里会报未识别、在另一部分版本里会按默认映射处理，两种情况的后果都不可控。所以清理阶段就应该把所有的 HIS 改掉，一个不留。

4.3 拓扑通过但动力学炸掉的隐性错误

最让人头疼的就是这一类：tleap 一声不响地通过了，最小化也收敛了，加热阶段也过了，跑到几十皮秒甚至几纳秒才出问题。这类错误的共同特征是“在拓扑层面合法，在几何层面不合理”。

典型的有这么几种。缺失侧链被补在了不合理的位置，导致局部密度异常，这种问题在能量最小化阶段可能被强行压下去，但在恒温恒压条件下会反弹。二硫键没显式成键，两个硫原子相距 2.05 Å 却按非键相互作用处理，初始时还好，跑着跑着会因为排斥而分开，把周围结构带歪。配体参数里的二面角参数缺失，被 parmchk2 用通用估计值补上，能量面在某个角度附近出现不合理的势垒，配体在口袋里转不动或者转过某个角度后突然被弹开。末端残基的带电形式在真实的埋藏环境中不合理，导致末端附近的局部结构持续变形。

排查这类问题的方法我觉得只有一个有效的：把结构准备阶段的每一个“自动”决策都变成“显式”决策。tleap 补了原子，你要看它在哪补的；parmchk2 补了参数，你要知道补了哪些项；pdb4amber 重命名了原子，你要确认语义没变。听起来费事，但它把不可见的问题变成了可见的清单。

4.4 一张速查表

下面这张表是我自己常用的排查对照，按现象查原因比按工具查要快。

注意：表里每一条都对应一次真实的返工。如果你的体系同时命中两条以上，建议直接从结构清洗那一步重来，别在拓扑层面打补丁。

5. 我自己的工程化习惯

5.1 目录、命名与脚本化

体系一多，靠记性管理文件必然出事。我固定用一套目录结构：raw/放原始文件且只读不改，clean/放清洗后的结构，params/放配体的 mol2 和 frcmod，leap/放 tleap 脚本与输出，md/放动力学输入输出，ana/放分析脚本与结果。文件名带上体系和日期，比如complex_v3_clean.pdb，改一版就升一个号，绝不覆盖。

更重要的是把所有操作脚本化。tleap 脚本、清洗脚本、自检脚本都写成文件，用一个顶层脚本顺序调用。这么做的好处不只是省事，而是当你三个月后要复现某个体系时，整个流程是一条命令，而不是“我记得当时好像还手动改了一个残基”。

5.2 版本与复现

我现在的习惯是，每完成一个体系的准备工作，就在leap/目录下留一个README式的记录：结构来源和下载日期、用的力场和水模型、质子化状态的判定依据、配体电荷的计算方法、最终的原子数和净电荷、以及任何一个手工修改的说明。这些记录加起来不到一页纸，但它能在半年后救你一次。

还有一个我坚持的小技巧：每次生成 prmtop 之后，把leap.log完整保存一份，并且用固定的一组关键词扫一遍。日志里的 warning 数量是可以作为质量指标的——同一个体系的准备流程，如果这次 warning 比上次多了五条，那多半是某个输入文件出了变化。

跑动力学之前，最后再看一眼导出的 pdb：蛋白在盒子中央，配体在口袋里，没有哪个残基的一半在盒子外面，没有两个原子叠在一起。这三十秒的检查，比跑起来之后花三个小时排查要划算得多。我在这上面栽过的次数不算少，现在基本上形成肌肉记忆了，你也值得养成。