SPONGE常用教程：蛋白+配体模拟2

前序课程

目录

1. 应用场景简述；- [ Done ]
2. DSDP:蛋白-配体对接；- [ Done ]
3. XPONGE：蛋白-配体建模，加溶剂；- [ Done ]
4. SPONGE：能量极小化-NVT-NPT-正式模拟；
5. XPONGE：数据简单后处理。

4. SPONGE：能量极小化-NVT-NPT-正式模拟
文件下载
进入建模文件目录，如下图：

SPONGE输入文件采用独特的参数模块化设计，每一个文档对应一种力场项类型及参数（例如1ae5_bond.txt对应键能项，第一行为系统总键数，后面每一行分别是：原子a的id 原子b的id 力常数 平衡位置），方便用户对于体系进行灵活调整及检查，具体文档参见：https://spongemm.cn/zh/%E6%96%87%E6%A1%A3

进行正式的分子动力学模拟之前，体系通常需要3步预平衡：1.能量极小化，调整体系H原子的位置，避免碰撞；2.NVT模拟将体系升温到反应温度；3.NPT模拟将体系压强控制到正确的压强。后两步预平衡主要调整的是体系的溶剂，因此需要对蛋白-配体作用对进行坐标限制，免得在正式模拟前就已经各奔东西。

查看1ae5_pro_lig_sol.pdb文件，溶质共有3388个原子，全盒12651原子

cp -r 1ae5 SPONGE/
#将1ae5文件夹放到SPONGE主程序所在的文件夹中
cd SPONGE
awk -v n=3388 'BEGIN{for(i=1;i<=n;i++)print i-1}' > restrain_id.txt
#标记溶质原子id用于后续NVT，NPT
rm mdin.txt
#删除测试脚本，SPONGE文件会默认读取mdin.txt作为输入
./SPONGE -mode Minimization -dt 1e-9 -minimization_dynamic_dt 1 -step_limit 10000 -default_in_file_prefix 1ae5/1ae5 -coordinate_in_file 1ae5/1ae5_coordinate.txt -rst Hmin -mdout Hmin.out
#能量极小化过程中可以观察到总能量项（potential）降低

接下来的模拟参数输入方式采用文件输入，采用生成的构象Hmin_coordinate.txt进一步控温到常温300K，打开mdin_NVT.txt了解NVT参数配置：

    The first line is the name of the md tast
#本文本文件通过命令行中加入"-mdin 本文本文件名"读入，如果该文件名为"mdin.txt"且在当前命令行的工作路径中，可不在命令行中显式输入
#本文本文件内的所有参数指令均可通过加"-"前缀在命令行里输入，例如 -cutoff 8.0 -skin 2.0

#分子模拟最基本相关设置
mode = NVT #可选Rerun, Minimization, NVE, NVT, NPT

dt = 1e-3 #模拟步长，单位ps
step_limit = 10000 #模拟总步数，通常需要长一点的模拟，本教程用了偏小的值
write_information_interval = 1000 #保存中间信息的间隔
write_restart_file_interval = 10000  #保存重启坐标和速度的间隔

target_temperature = 300.0 #模拟的目标平衡温度，单位为开
thermostat = middle_langevin #热浴方法

restrain 
{
    atom_id = restrain_id.txt
    weight = 10.0
}
#用于限制蛋白-配体初始结构不动，只优化溶剂水

    #输入
        #文件目录
        default_in_file_prefix = 1ae5/1ae5

        #坐标和速度
        coordinate_in_file = Hmin_coordinate.txt #使用能量极小化得到的坐标作为输入坐标
        #velocity_in_file = Hmin_velocity.txt

    #输出
        mdout = mdout.txt   #记录能量轨迹的文件
        mdinfo = mdinfo.txt #记录参数相关的文件
        box = mdbox.txt       #记录盒子轨迹的文件
        crd = mdcrd.dat        #记录坐标轨迹的文件
        rst = Vmin               #记录重启坐标和速度的文件

#其他参数
end_pause = 1 #程序运行结束后是否需要输入任意键退出
device = 0 #使用的GPU设备
dont_check_input = 1 #忽略多余参数

然后开始NVT预平衡

./SPONGE -mdin mdin_NVT.txt
#NVT过程中可以观察到体系温度（temperature）逐步收敛到300K，并在附近扰动。

采用生成的构象Vmin_coordinate.txt进一步控压到标准大气压，打开mdin_NPT.txt了解NPT参数配置：

    The first line is the name of the md tast
#本文本文件通过命令行中加入"-mdin 本文本文件名"读入，如果该文件名为"mdin.txt"且在当前命令行的工作路径中，可不在命令行中显式输入
#本文本文件内的所有参数指令均可通过加"-"前缀在命令行里输入，例如 -cutoff 8.0 -skin 2.0

#分子模拟最基本相关设置
mode = NPT #可选Rerun, Minimization, NVE, NVT, NPT

dt = 1e-3 #模拟步长，单位ps
step_limit = 100000 #模拟总步数，通常需要长一点的模拟，本教程用了偏小的值
write_information_interval = 1000 #保存中间信息的间隔
write_restart_file_interval = 10000  #保存重启坐标和速度的间隔

target_temperature = 300.0 #模拟的目标平衡温度，单位为开
thermostat = middle_langevin #热浴方法
barostat =  andersen_barostat #压浴设置
target_pressure = 1.0 #模拟的目标平衡压强，单位为bar

restrain
{
    atom_id = restrain_id.txt 
    weight = 1.0
}
#用于限制蛋白-配体初始结构不动，只优化溶剂水

    #输入
        #文件目录
        default_in_file_prefix = 1ae5/1ae5

        #坐标和速度
        coordinate_in_file = Vmin_coordinate.txt #使用能量极小化得到的坐标作为输入坐标
        #velocity_in_file = Hmin_velocity.txt

    #输出
        mdout = Pmdout.txt   #记录能量轨迹的文件
        mdinfo = mdinfo.txt #记录参数相关的文件
        box = mdbox.txt       #记录盒子轨迹的文件
        crd = mdcrd.dat        #记录坐标轨迹的文件
        rst = Pmin               #记录重启坐标和速度的文件

#其他参数
end_pause = 1 #程序运行结束后是否需要输入任意键退出
device = 0 #使用的GPU设备
dont_check_input = 1 #忽略多余参数

然后开始NPT预平衡

./SPONGE -mdin mdin_NPT.txt
#NPT过程中可以观察到体系密度（density）逐步收敛到1.0g/cm3，并在附近扰动。

最终成品模拟采用生成的构象Pmin_coordinate.txt，解开restrain原子坐标限制，打开mdin_MD.txt了解MD参数配置：

    The first line is the name of the md tast
#本文本文件通过命令行中加入"-mdin 本文本文件名"读入，如果该文件名为"mdin.txt"且在当前命令行的工作路径中，可不在命令行中显式输入
#本文本文件内的所有参数指令均可通过加"-"前缀在命令行里输入，例如 -cutoff 8.0 -skin 2.0

#分子模拟最基本相关设置
mode = NPT #可选Rerun, Minimization, NVE, NVT, NPT

dt = 1e-3 #模拟步长，单位ps
step_limit = 1000000 #模拟总步数，通常需要长一点的模拟，本教程用了偏小的值
write_information_interval = 1000 #保存中间信息的间隔
write_restart_file_interval = 100000  #保存重启坐标和速度的间隔

target_temperature = 300.0 #模拟的目标平衡温度，单位为开
thermostat = middle_langevin #热浴方法
barostat =  andersen_barostat #压浴设置
target_pressure = 1.0 #模拟的目标平衡压强，单位为bar

##restrain#双#号注释掉整个括号里的内容，放开初始构象限制，进行模拟
{
    atom_id = restrain_id.txt 
    weight = 1.0
}


    #输入
        #文件目录
        default_in_file_prefix = 1ae5/1ae5

        #坐标和速度
        coordinate_in_file = Pmin_coordinate.txt #使用NPT得到的坐标作为输入坐标
        velocity_in_file = Pmin_velocity.txt

    #输出
        mdout = mdout.txt   #记录能量轨迹的文件
        mdinfo = mdinfo.txt #记录参数相关的文件
        box = mdbox.txt       #记录盒子轨迹的文件
        crd = mdcrd.dat        #记录坐标轨迹的文件
        rst = restart              #记录重启坐标和速度的文件

#其他参数
end_pause = 1 #程序运行结束后是否需要输入任意键退出
device = 0 #使用的GPU设备
dont_check_input = 1 #忽略多余参数

运行成品模拟：

./SPONGE -mdin mdin_MD.txt

最终将得到二进制轨迹文件mdcrd.dat和盒子信息mdbox.txt，用于轨迹分析。
【to be continued...】