GROMACS 蛋白质建模与溶剂化

本文以鸢算科技集群 GROMACS 2026.4 实际结果为例。下载输入、参考轨迹、分析与绘图资料包，资料包约 529 MB，包含参考轨迹。解压后进入 gromacs-examples。参考输出在 cases；新计算由 new-case.py 建立独立目录。

本篇任务与操作要点#

例子使用公开结构 1UBQ 泛素，保留蛋白 ATOM 记录，去掉晶体水，重新加溶剂。力场为 AMBER99SB-ILDN，水为 TIP3P。此处没有配体、金属中心、缺失环或非标准残基；遇到这些情况，应先补全结构并准备与蛋白力场兼容的参数，不能把 pdb2gmx 当作任意配体自动参数化程序。

-ignh 让程序按力场重建氢；真实课题需要先确定 pH、质子化态、端基、二硫键与组装体。本例采用该输入的默认端基处理，完整拓扑随附件提供。-d 1.0 -bt cubic 给蛋白边缘留出至少 1 nm 的初始溶剂层。

pdb2gmx 还会生成位置约束文件；平衡时通过 define = -DPOSRES 启用，生产阶段移除。图片来自本次计算的实际坐标，只画 Cα 骨架，便于看清折叠。

1 准备与提交#

python new-case.py ubiquitin
cd work/ubiquitin
sbatch run.sh

run.sh 申请资源后执行 workflow.sh。不要在登录节点直接执行完整模拟。下列每个算例可单独新建副本并提交。

2 实际终端命令#

以下是在对应算例目录中实际执行并生成结果的命令。分析命令在该算例的 analysis 子目录执行；交互回答已用 printf 写入，组名需随自己的体系核对。

ubiquitin

gmx pdb2gmx -f protein.pdb -o processed.gro -p topol.top -ff amber99sb-ildn -water tip3p -ignh
gmx editconf -f processed.gro -o box.gro -c -d 1.0 -bt cubic
gmx solvate -cp box.gro -cs spc216.gro -o solv.gro -p topol.top
gmx grompp -f em.mdp -c solv.gro -r solv.gro -p topol.top -o ions.tpr
printf %s 'SOL
' | gmx genion -s ions.tpr -o ions.gro -p topol.top -pname NA -nname CL -neutral -conc 0.15 -seed 76266
gmx grompp -f em.mdp -c ions.gro -r ions.gro -p topol.top -o em.tpr -po em-processed.mdp
gmx grompp -f nvt.mdp -c em.gro -r em.gro -p topol.top -o nvt.tpr -po nvt-processed.mdp
gmx grompp -f npt.mdp -c nvt.gro -r nvt.gro -p topol.top -o npt.tpr -po npt-processed.mdp -t nvt.cpt
gmx grompp -f md.mdp -c npt.gro -r npt.gro -p topol.top -o md.tpr -po md-processed.mdp -t npt.cpt

3 输入文件与关键设置#

ubiquitin/protein.pdb 完整文件见资料包，保留实际坐标、原子顺序与索引。用 head -n 20 protein.pdb 查看开头，修改前先复制备份。

cutoff-scheme = Verlet
nstlist = 20
coulombtype = PME
rcoulomb = 1.0
rvdw = 1.0
pbc = xyz
constraints = h-bonds
integrator = md
dt = 0.002
nsteps = 10000
nstenergy = 50
nstcalcenergy = 50
nstlog = 1000
nstxout-compressed = 50
tcoupl = v-rescale
tc-grps = System
tau-t = 0.5
ref-t = 300
pcoupl = no
gen-vel = yes
gen-temp = 300
gen-seed = 76266
define = -DPOSRES

4 实际结果与画图#

下表为本次 100 ps 生产轨迹后半段的简单均值，不是收敛物性或误差分析。

用附件中的 Python 脚本重画所有参考图：

# 回到 gromacs-examples 根目录
module load ase/3.29.0
python plot.py --root .

plot.py 使用 NumPy 和 Matplotlib，直接读取 XVG/COLVAR/坐标数值，图片保存在 figures。对单项新结果，可照着脚本对应段替换文件路径；不要把不同模型或单位的数据直接拼在一起。

参数查阅#

以本集群 gmx 命令名 -h 为准；GROMACS 2026.4 命令手册和MDP 参数手册可查完整选项。