SPONGE模拟丙氨酸二肽：分子运动轨迹与分析

目录

1. 引言
2. 超算互联网软件与环境
3. 理论目标与判据
4. 准备项目文件
5. 编译环境配置
6. 编译与安装 SPONGE
7. 能量最小化作业
8. NPT 短程模拟作业
9. 能量与密度结果分析
10. 使用 Xponge 和 MDAnalysis 计算 RMSD/RMSF
11. 结果文件同步与复现检查
12. 小结

1.引言

分子动力学模拟是研究分子体系微观运动行为的重要方法。通过数值积分原子间相互作用产生的运动方程，可以获得分子在一定时间内的结构变化轨迹，从而观察分子的构象变化、热运动特征以及局部结构波动。对于生物分子体系而言，分子动力学模拟不仅能够补充实验观测中难以直接获得的动态信息，也为理解分子结构稳定性和构象变化提供了直观工具。
丙氨酸二肽是分子模拟中常用的模型体系。由于其结构简单、原子数较少，同时又包含典型的肽键和主链二面角，因此常被用于演示和测试分子动力学模拟流程。通过对丙氨酸二肽进行短时间模拟，可以观察其在溶液或真空环境中的分子运动，并进一步分析结构随时间的偏移和不同原子的柔性差异。
SPONGE是一款由北京大学高毅勤教授课题组开发的分子动力学模拟引擎，结合GPU加速与人工智能技术，支持传统力场和神经网络力场建模和模拟，广泛应用于生物体系、材料体系的动力学模拟。
本次实践将基于超算互联网环境，从 GitHub 源码编译 SPONGE2，并对丙氨酸二肽体系进行一次短时间分子动力学模拟。实践重点不是获得充分采样的生产模拟结果，而是完整走通“源码获取、环境配置、编译安装、能量最小化、NPT 短程模拟、RMSD/RMSF 后处理”的基础流程。

2.超算互联网软件与环境配置

A目前，超算互联网已上线SPONGE，用户可以通过软件服务或智能体方式快速进入环境，这种方式适合直接运行已有示例、熟悉软件功能和完成常规教学演示。
由于分子模拟软件和相关前沿算法迭代很快，GitHub 上的开发分支往往会先加入新的后端适配、编译修复或功能更新。当需要使用最新版本、指定分支，或希望完整掌握软件从源码到可执行文件的部署过程时，仅依赖平台预装版本就不够灵活。
因此，本实践我们不直接使用平台已有的 SPONGE 软件，而是从 GitHub 获取 SPONGE 源码开始，在超算互联网环境中用 pixi、DTK 模块和 Slurm 作业脚本完成一次可复现的编译、运行和分析流程。
在实际科研中，常常需要自己管理环境，科研领域中常见的环境管理工具是conda，但是Anaconda本身因为开发年代久远效率较低，我这里推荐使用pixi，使用rust重写的conda包管理器，存储、包管理效率远超Anaconda。在命令行中使用curl -fsSL https://pixi.sh/install.sh | bash 安装pixi。

3.理论目标与判据

SPONGE 可以在给定的丙氨酸二肽体系上进行分子动力学模拟，得到体系在短时间内的运动轨迹。本实践的目标不是进行长时间生产模拟，而是检查一次短程模拟是否能够稳定运行，并从轨迹中提取结构偏移和局部柔性指标。 判断标准包括：能量最小化阶段势能降低并生成下一步运行的坐标文件；NPT 阶段温度围绕 300 K 波动、密度进入相对稳定区间；溶质 RMSD 在短时间内没有持续发散；溶质 RMSF 能够反映不同原子的柔性差异。
因此，本教程定位为“超算互联网 SPONGE 短程模拟与轨迹分析流程”，不是长时间采样或构象统计收敛分析。若要进一步研究丙氨酸二肽构象空间，需要延长模拟并分析 phi/psi 二面角的 Ramachandran 图分布。

4.准备项目文件

本实践使用两个主要目录：SPONGE 源码目录和 ala2 体系目录。SPONGE 源码从 GitHub 仓库获取并压缩为 SPONGE.zip，上传到超算互联网后解压得到 ~/SPONGE；丙氨酸二肽体系则通过 pixi 管理的 Xponge 建模环境在远端生成。
ala2 项目目录中的 pixi.toml 指定 python=3.10.*，并通过 PyPI 依赖安装 xponge>=1.4.4,<2。进入 ala2 目录后执行 pixi install，即可得到独立的建模环境。

建模脚本 build.py 载入 Xponge 以及 AMBER ff14SB、TIP3P 力场，使用 ACE + ALA + NME 构造丙氨酸二肽，并调用 add_solvent_box 加入水盒。随后输出 ala2.pdb，并用 save_sponge_input 生成 SPONGE 所需的 top/ala2_* 输入文件。
建模环境和体系生成命令如下。

mkdir -p ~/ala2
cd ~/ala2
pixi init
pixi add python=3.10
pixi add --pypi xponge

cat > build.py <<'PY'
import Xponge
import Xponge.forcefield.amber.ff14sb
import Xponge.forcefield.amber.tip3p
mol = ACE + ALA + NME
sys = Xponge.add_solvent_box(mol, WAT, 30)
Xponge.save_pdb(sys, "ala2.pdb")
Xponge.save_sponge_input(sys, "ala2", dirname="top")
PY

pixi run python build.py

源码压缩包可以从本地上传；丙氨酸二肽体系则在远端通过 pixi + Xponge 生成。实际操作中建议把可复现脚本集中放在 ~/scripts 下，避免在多次调试后远端目录过于混乱。

scp SPONGE.zip hpccube-wuzh02:~/SPONGE.zip
ssh hpccube-wuzh02
unzip -q ~/SPONGE.zip -d ~
mkdir -p ~/scripts

5.编译环境配置

本次实践使用的是华东三区，系统glibc为2.17。采用与集群系统兼容的 pixi 环境，显式指定 sysroot_linux-64=2.17，并加载集群提供的 compiler/dtk/24.04.3 模块。华东三区的 DCU 架构对应 gfx906，CMake 配置时只编译 gfx906 单架构。

module load compiler/dtk/24.04.3
cd ~/SPONGE
~/.pixi/bin/pixi add --feature hip --platform linux-64   'gcc_linux-64=11.*' 'gxx_linux-64=11.*' 'sysroot_linux-64=2.17'
~/.pixi/bin/pixi install -e hip

DTK=/public/software/compiler/rocm/dtk-24.04.3
ENV=$HOME/SPONGE/.pixi/envs/hip
export PATH="$ENV/bin:$DTK/bin:$DTK/hip/bin:$DTK/llvm/bin:$PATH"
export LD_LIBRARY_PATH="$ENV/lib:$DTK/lib:$DTK/lib64:$DTK/hip/lib:$DTK/hipfft/lib:$DTK/rocfft/lib:$DTK/llvm/lib:${LD_LIBRARY_PATH:-}"
export ROCM_PATH="$DTK" HIP_PATH="$DTK/hip"
export CC="$ENV/bin/x86_64-conda-linux-gnu-cc"
export CXX="$ENV/bin/x86_64-conda-linux-gnu-c++"

cmake -S . -B build-hip -G Ninja -DPARALLEL=hip   -DCMAKE_INSTALL_PREFIX="$ENV"   -DCMAKE_HIP_ARCHITECTURES=gfx906   -DAMDGPU_TARGETS=gfx906 -DGPU_TARGETS=gfx906

6.编译与安装 SPONGE

配置完成后，使用 Ninja 编译并安装到 pixi 环境中。实际复现实验中，可执行文件安装到 ~/SPONGE/.pixi/envs/hip/bin/SPONGE。

bash ~/scripts/configure_SPONGE.sh 2>&1 | tee ~/scripts/repro_configure.log
bash ~/scripts/compile_SPONGE.sh 2>&1 | tee ~/scripts/repro_compile.log
cmake --build build-hip --config Release --target install --parallel 4

7.能量最小化作业

能量最小化用于消除初始结构中明显的局部不合理接触，并为后续 NPT 模拟提供 restart_coordinate.txt。作业脚本保存为 ~/scripts/min.sh，通过 Slurm 提交到 国产异构加速卡 队列。

#!/usr/bin/env bash
#SBATCH -J ala2-min
#SBATCH -p wzhdtest
#SBATCH --gres=dcu:1
#SBATCH -N 1
#SBATCH -n 1
#SBATCH -c 4

module load compiler/dtk/24.04.3
export PATH="$HOME/SPONGE/.pixi/envs/hip/bin:/public/software/compiler/rocm/dtk-24.04.3/bin:$PATH"
mkdir -p ~/ala2/min
cd ~/ala2/min
SPONGE -mode minimization -step_limit 5000 -default_in_file_prefix ../top/ala2

本次运行完成后，min/mdout.txt 中势能从 -111950.16 下降到 -116731.59，并生成 restart_coordinate.txt，说明最小化阶段正常结束。

8.NPT 短程模拟作业

完成最小化后，使用最小化输出坐标作为初始构型，进行 50000 步 NPT 模拟。时间步长为 2 fs，因此总模拟时间为 100 ps。

#!/usr/bin/env bash
#SBATCH -J ala2-npt
#SBATCH -p wzhdtest
#SBATCH --gres=dcu:1
#SBATCH -N 1
#SBATCH -n 1
#SBATCH -c 4

module load compiler/dtk/24.04.3
export PATH="$HOME/SPONGE/.pixi/envs/hip/bin:/public/software/compiler/rocm/dtk-24.04.3/bin:$PATH"
mkdir -p ~/ala2/npt
cd ~/ala2/npt
ln -sfn ../top top
ln -sfn ../min min
SPONGE -mode npt -step_limit 50000 -dt 2e-3 -constrain_mode SHAKE   -default_in_file_prefix top/ala2   -coordinate_in_file min/restart_coordinate.txt   -thermostat middle_langevin -langevin_gamma 10.0   -target_temperature 300.0   -barostat andersen_barostat -barostat_tau 0.1   -barostat_update_interval 10 -target_pressure 1.0

本次 NPT 作业正常结束，核心计算墙钟时间为 281.63 s，输出速度约 30.68 ns/day。最后一步为 50000 步，对应 100 ps，温度为 301.53 K，密度为 0.9843。

9.能量与密度结果分析

NPT 阶段的 mdout.txt 记录了温度、势能、密度、压力和各能量项。短程模拟结果如下。
输出帧数：50
时间范围：2-100 ps
末步温度：301.53 K
末步密度：0.9843
后半段势能均值：-89359.74
后半段密度均值：0.984356
后半段温度均值：301.04 K
后半段压力均值：-56.89 bar
温度在初始升温后围绕 300 K 波动，密度在约 0.98-0.99 区间内变化。压力在 100 ps 的短模拟中仍有较大瞬时涨落，这是小体系短时间 NPT 模拟中常见的现象。
NPT 阶段势能、密度、温度和压力随时间变化

10.使用 Xponge 和 MDAnalysis 计算 RMSD/RMSF

SPONGE 的 mdcrd.dat 可以通过 Xponge.analysis.md_analysis 接入 MDAnalysis。拓扑前缀使用 top/ala2，轨迹文件使用 npt/mdcrd.dat，盒子文件使用 npt/mdbox.txt。

import MDAnalysis as mda
import Xponge.analysis.md_analysis as xmda
from MDAnalysis.analysis import rms, align

u = mda.Universe(
    "top/ala2", "npt/mdcrd.dat",
    topology_format=xmda.SpongeInputReader,
    format=xmda.SpongeTrajectoryReader.with_arguments(box="npt/mdbox.txt"),
)
solute = u.select_atoms("resname ACE ALA NME")

ref = mda.Universe(
    "top/ala2", "npt/mdcrd.dat",
    topology_format=xmda.SpongeInputReader,
    format=xmda.SpongeTrajectoryReader.with_arguments(box="npt/mdbox.txt"),
)
R = rms.RMSD(u, ref, select="resname ACE ALA NME")
R.run()
rmsd = R.results.rmsd[:, 2]

ref.trajectory[0]
ref_atoms = ref.select_atoms("resname ACE ALA NME")
mobile_atoms = u.select_atoms("resname ACE ALA NME")
aligned_positions = []
for ts in u.trajectory:
    align.alignto(mobile_atoms, ref_atoms, weights=None)
    aligned_positions.append(mobile_atoms.positions.copy())

RMSD 以第一帧为参考，选择 ACE、ALA、NME 三个残基作为溶质；RMSF 在逐帧对齐溶质后计算每个原子的平均位置涨落。由于 SPONGE reader 不自动携带时间步长信息，作图时按 mdout 输出间隔手动设置为每帧 2 ps。
体系总原子数：28054
体系总残基数：9347
轨迹帧数：50
溶质选择：resname ACE ALA NME
溶质原子数：22
RMSD 末值：1.128835 Å
RMSD 均值：0.918388 Å
RMSD 最大值：1.399309 Å
RMSF 均值：0.566397 Å
RMSF 最大值：1.094953 Å，3-H2
RMSD 在 100 ps 内主要分布在 1 Å 左右，最大值约 1.40 Å，没有出现持续单调增大的发散趋势。RMSF 显示端基氢原子的波动较大，而肽键和主链重原子波动较小，符合短肽溶质端基更灵活的直观预期。
基于 Xponge 与 MDAnalysis 计算的 RMSD 和 RMSF

11.结果文件同步与复现检查

计算结束后，可以将关键输出文件同步回本地，用于绘图、写报告或进一步分析。建议至少保留 min/mdout.txt、min/restart_coordinate.txt、npt/mdout.txt、npt/mdbox.txt、npt/mdcrd.dat、npt/restart_coordinate.txt 以及 Slurm 输出日志。

scp -r hpccube-wuzh02:~/ala2/min E:/超算互联网/ala2_repro_final/
scp -r hpccube-wuzh02:~/ala2/npt E:/超算互联网/ala2_repro_final/
scp -r hpccube-wuzh02:~/ala2/top E:/超算互联网/ala2_repro_final/

12.小结

本实践完成了 SPONGE2 在超算互联网国产算力环境中的一次端到端运行：通过 pixi 和 DTK 24.04.3 建立兼容 glibc 2.17 的 HIP 编译环境，并针对 gfx906 单架构配置和编译 SPONGE；随后对丙氨酸二肽体系进行 5000 步能量最小化和 50000 步 NPT 短程模拟；最后使用 Xponge 与 MDAnalysis 读取 SPONGE 轨迹，计算 RMSD 和 RMSF。 从结果看，最小化阶段势能降低，NPT 阶段温度和密度进入合理波动范围，溶质 RMSD 未出现发散，RMSF 能够反映端基和局部原子的柔性差异。