1. 把上一步生成的Enzyme.gro和Enzyme.top文件改个名字

cp Enzyme.top topol.top
cp Enzyme.gro complex.gro

2. 能量最小化

首先准备能量最小化mdp文件em_real.mdp,示例如下

; minim.mdp - used as input into grompp to generate real_em.tpr
; Parameters describing what to do, when to stop and what to save
integrator  = steep ; Algorithm (steep = steepest descent minimization)
emtol       = 1000.0        ; Stop minimization when the maximum force < 1000.0 kJ/mol/nm
emstep      = 0.01          ; Minimization step size
nsteps      = 50000         ; Maximum number of (minimization) steps to perform; Parameters describing how to find the neighbors of each atom and how to calculate the interactions
nstlist         = 1         ; Frequency to update the neighbor list and long range forces
cutoff-scheme   = Verlet ; Buffered neighbor searching
ns_type         = grid ; Method to determine neighbor list (simple, grid)
coulombtype     = PME ; Treatment of long range electrostatic interactions
rcoulomb        = 1.0       ; Short-range electrostatic cut-off
rvdw            = 1.0       ; Short-range Van der Waals cut-off
pbc             = xyz ; Periodic Boundary Conditions in all 3 dimensionsgmx grompp -f em_real.mdp -c complex.gro -r complex.gro -p topol.top -o em.tpr
gmx mdrun -v -deffnm em

em_real.mdp文件中的各项参数说明:

(1)integrator:

=steep:用最陡下降法做能量极小化。

=cg:用共现桃度法做能量极小化(不支持约束)。

=l-bigs:用L-BFGS法做能量极小化。

=nm:做正则振动分析( 必须双精度版)。

(2)emtol: 能量极小化时最大受力小于多少就认为收敛(kJ/mol/nm).默认为10.0。

当最大作用力小于此值,认为最小化过程收敛。一般蛋白-配体中是1000、500。

(3)emstep: 最陡下降法最大步长(nm)。默认0.01。

(4)nsteps: 动力学或能量极小化的步数上限,指定最大迭代次数。

之后运行能量最小化

gmx grompp -f em_real.mdp -c complex.gro -r complex.gro -p topol.top -o em.tpr
gmx mdrun -v -deffnm em

3. NVT平衡

首先准备NVT所需mdp文件nvt.mdp,示例如下:

title                   = NVT
define                  = -DPOSRES ; position restrain the protein; Run parameters
integrator              = md ; leap-frog integrator
nsteps                  = 50000     ; 2 * 50000 = 100 ps
dt                      = 0.002     ; 2 fs; Output control
nstxout                 = 500       ; save coordinates every 1.0 ps
nstvout                 = 500       ; save velocities every 1.0 ps
nstenergy               = 500       ; save energies every 1.0 ps
nstlog                  = 500       ; update log file every 1.0 ps; Bond parameters
continuation            = no        ; first dynamics run
constraint_algorithm    = shake ; holonomic constraints
constraints             = h-bonds ; bonds involving H are constrained
shake-SOR               = no
shake-tol               = 0.0001 ; Nonbonded settings
cutoff-scheme           = Verlet ; Buffered neighbor searching
ns_type                 = grid ; search neighboring grid cells
nstlist                 = 10        ; 20 fs, largely irrelevant with Verlet
rcoulomb                = 1.0       ; short-range electrostatic cutoff (in nm)
rvdw                    = 1.0       ; short-range van der Waals cutoff (in nm)
DispCorr                = EnerPres ; account for cut-off vdW scheme; Electrostatics
coulombtype             = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order               = 4         ; cubic interpolation
fourierspacing          = 0.16      ; grid spacing for FFT; Temperature coupling is on
tcoupl                  = V-rescale ; modified Berendsen thermostat
tc-grps                 = Protein Non-Protein ; two coupling groups - more accurate
tau_t                   = 0.1     0.1           ; time constant, in ps
ref_t                   = 300     300           ; reference temperature, one for each group, in K; Pressure coupling is off
pcoupl                  = no        ; no pressure coupling in NVT; Periodic boundary conditions
pbc                     = xyz ; 3-D PBC; Velocity generation
gen_vel                 = yes       ; assign velocities from Maxwell distribution
gen_temp                = 300       ; temperature for Maxwell distribution
gen_seed                = -1        ; generate a random seed

关于各参数的含义,参见后面md.mdp文件中的参数解释。

需要注意的是,前面我们利用的AmberTools来产生top文件,因此如果我们这里用POSRES来做位置约束的话,需要将topol.top文件中的ifdef flexible修改为ifdef POSRES,否则此步中会出现warning。

之后运行NVT平衡

gmx grompp -f nvt.mdp -c em.gro -r em.gro -p topol.top -o nvt.tpr
gmx mdrun -deffnm nvt  

4. NPT平衡

同样,首先准备NPT所需的npt.mdp文件

title                   = NPT
define                  = -DPOSRES ; position restrain the protein; Run parameters
integrator              = md ; leap-frog integrator
nsteps                  = 50000     ; 2 * 50000 = 100 ps
dt                      = 0.002     ; 2 fs; Output control
nstxout                 = 500       ; save coordinates every 1.0 ps
nstvout                 = 500       ; save velocities every 1.0 ps
nstenergy               = 500       ; save energies every 1.0 ps
nstlog                  = 500       ; update log file every 1.0 ps; Bond parameters
continuation            = yes       ; Restarting after NVT
constraint_algorithm    = shake ; holonomic constraints
constraints             = h-bonds ; bonds involving H are constrained
shake-SOR               = no
shake-tol               = 0.0001 ; Nonbonded settings
cutoff-scheme           = Verlet ; Buffered neighbor searching
ns_type                 = grid ; search neighboring grid cells
nstlist                 = 10        ; 20 fs, largely irrelevant with Verlet scheme
rcoulomb                = 1.0       ; short-range electrostatic cutoff (in nm)
rvdw                    = 1.0       ; short-range van der Waals cutoff (in nm)
DispCorr                = EnerPres ; account for cut-off vdW scheme; Electrostatics
coulombtype             = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order               = 4         ; cubic interpolation
fourierspacing          = 0.16      ; grid spacing for FFT; Temperature coupling is on
tcoupl                  = V-rescale ; modified Berendsen thermostat
tc-grps                 = Protein Non-Protein ; two coupling groups - more accurate
tau_t                   = 0.1     0.1           ; time constant, in ps
ref_t                   = 300     300           ; reference temperature, one for each group, in K; Pressure coupling is on
pcoupl                  = Parrinello-Rahman ; Pressure coupling on in NPT
pcoupltype              = isotropic ; uniform scaling of box vectors
tau_p                   = 2.0                   ; time constant, in ps
ref_p                   = 1.0                   ; reference pressure, in bar
compressibility         = 4.5e-5                ; isothermal compressibility of water, bar^-1
refcoord_scaling        = com; Periodic boundary conditions
pbc                     = xyz ; 3-D PBC; Velocity generation
gen_vel                 = no        ; Velocity generation is off

之后运行NPT

gmx grompp -f npt.mdp -c nvt.gro -r nvt.gro -t nvt.cpt -p topol.top -o npt.tpr
gmx mdrun -deffnm npt 

5. 成品模拟

md.mdp文件示例如下

title                   = MD
; Run parameters
integrator              = md ; leap-frog integrator
nsteps                  = 50000000    ; 2 * 50000000 = 100000 ps (100 ns)
dt                      = 0.002     ; 2 fs
; Output control
nstxout                 = 0         ; suppress bulky .trr file by specifying
nstvout                 = 0         ; 0 for output frequency of nstxout,
nstfout                 = 0         ; nstvout, and nstfout
nstenergy               = 5000      ; save energies every 10.0 ps
nstlog                  = 5000      ; update log file every 10.0 ps
nstxout-compressed      = 5000      ; save compressed coordinates every 10.0 ps
compressed-x-grps       = System ; save the whole system
; Bond parameters
continuation            = yes       ; Restarting after NPT
constraint_algorithm    = shake ; holonomic constraints
constraints             = h-bonds ; bonds involving H are constrained
shake-SOR               = no
shake-tol               = 0.0001
; Neighborsearching
cutoff-scheme           = Verlet ; Buffered neighbor searching
ns_type                 = grid ; search neighboring grid cells
nstlist                 = 10        ; 20 fs, largely irrelevant with Verlet scheme
rcoulomb                = 1.0       ; short-range electrostatic cutoff (in nm)
rvdw                    = 1.0       ; short-range van der Waals cutoff (in nm)
; Electrostatics
coulombtype             = PME ; Particle Mesh Ewald for long-range electrostatics
pme_order               = 4         ; cubic interpolation
fourierspacing          = 0.16      ; grid spacing for FFT
; Temperature coupling is on
tcoupl                  = V-rescale ; modified Berendsen thermostat
tc-grps                 = Protein Non-Protein ; two coupling groups - more accurate
tau_t                   = 0.1     0.1           ; time constant, in ps
ref_t                   = 300     300           ; reference temperature, one for each group, in K
; Pressure coupling is on
pcoupl                  = Parrinello-Rahman ; Pressure coupling on in NPT
pcoupltype              = isotropic ; uniform scaling of box vectors
tau_p                   = 2.0                   ; time constant, in ps
ref_p                   = 1.0                   ; reference pressure, in bar
compressibility         = 4.5e-5                ; isothermal compressibility of water, bar^-1
; Periodic boundary conditions
pbc                     = xyz ; 3-D PBC
; Dispersion correction
DispCorr                = EnerPres ; account for cut-off vdW scheme
; Velocity generation
gen_vel                 = no        ; Velocity generation is off

md.mdp文件中各项参数含义:

(1)title=MD:标题,可任意定义(最长64个字)

(2)define:常涉及到的是-DFLEXIBLE,会使水分子不被SETTLES算法约束成刚性。当利用posre.itp对位置进行限制时可用-DPOSRES使之生效.

(3)integrator: 决定当前任务干什么。

=md (默认) :leap-frog( md蛙跳算法)方式做动力学。

=md-vv:用Velocity-Verlet方式做动力学,在功能上有一些限制。

=sd:做随机动力学(Langevin动力学)。

=bd:做布朗动力学。

(4) tinit:动力学开始的时间(ps)。    默认0

(5) dt:动力学步长(ps)。默认0.001 (1fs),一般我们都使用0.002。

(6)nsteps:动力学或能量极小化的步数上限,指定最大迭代次数。

(7)init-step: 起始步数, 对非平衡模拟, 精确重启或重做某部分模拟时, 设定为重启步编号。

(8)nstxout = 0:输出trr坐标10000

nstvout = 0:输出trr速度 10000

nstfout = 0:输出trr受力 10000

nstxtcout = 5000:如果这个数比前三者小的多,那么其就会很大,比trr都大。

compressed-x-precision =5000:写入xtc轨迹的精确度。默认为 (1000)

compressed-x-grps:选择输出到xtc文件里的group, 默认为system。该参数只影响xtc轨迹,不影响trr轨迹。如果想在做动力学的时候就只把蛋白写入轨迹文件中,可以mdp中设compressed-x-grps = non-water。一般我们用体系即可。之后再根据需要用trjconv来处理获得包含水的轨迹文件。

nstlog=5000:每多少步输出一次各种能量信息到log文件,默认为1000。

nstenergy:每多少步输出一次能量信息到.edr文件,默认为1000。

(9)constraini-algorithm            ;设定约束算法

=LINCS(默认):使用LINCS方法约束,但不能约束键角。

=SHAKE:使用SHAKE方法约束,没LINCS可靠不能用于能量极小化。

之后运行md.mdp

gmx grompp -f md.mdp -c npt.gro -r npt.gro -t npt.cpt -p topol.top -o md.tpr
gmx mdrun -deffnm md 

分子动力学模拟学习2-Gromacs运行分子动力学模拟相关推荐

  1. linux的gromacs模拟分子运动,GROMACS分子动力学模拟技术与应用

    该楼层疑似违规已被系统折叠 隐藏此楼查看此楼 "GROMACS分子动力学模拟技术与应用" 课 程 内 容 分子动力学模拟基础 1 分子模拟入门理论--对分子动力学模拟原理介绍1.1 ...

  2. AMBER:运行分子动力学中vmd的使用教程

    标题转载来源于 san's note 来自于大佬的个人博客 写得很详细.对初学者的我帮助很大 本教程旨在向您介绍如何将VMD与AMBER一起使用,介绍如何加载AMBER轨迹文件和inpcrd文件以及如 ...

  3. Windows进程与线程学习笔记(五)—— 模拟线程切换

    Windows进程与线程学习笔记(五)-- 模拟线程切换 ThreadSwitch代码分析 ThreadSwitch.cpp ThreadCore.h ThreadCore.cpp 总结 Thread ...

  4. 操作系统学习之用C语言模拟伙伴(Buddy)算法

    前言 学到了操作系统的的虚拟内存部分,硬件不太好的我学起来有些吃力,概念性知识点太多,所以我决定用软件的方式,实现一下虚拟内存常用的算法,因为用到了指针,暂时用C语言写一下Buddy算法.FIFO算法 ...

  5. 在我方某前沿防守地域 matlab,[matlab]Monte Carlo模拟学习笔记

    理论基础:大数定理,当频数足够多时,频率可以逼近概率,从而依靠概率与$\pi$的关系,求出$\pi$ 所以,rand在Monte Carlo中是必不可少的,必须保证测试数据的随机性. 用蒙特卡洛方法进 ...

  6. 2020制冷与空调设备运行操作模拟考试题库及制冷与空调设备运行操作复审模拟考试

    题库来源:安全生产模拟考试一点通公众号小程序 2020制冷与空调设备运行操作模拟考试题库及制冷与空调设备运行操作复审模拟考试,包含制冷与空调设备运行操作模拟考试题库答案解析及制冷与空调设备运行操作复审 ...

  7. 2021年制冷与空调设备运行操作考试内容及制冷与空调设备运行操作模拟试题

    题库来源:安全生产模拟考试一点通公众号小程序 安全生产模拟考试一点通:2021年制冷与空调设备运行操作考试内容为正在备考制冷与空调设备运行操作操作证的学员准备的理论考试专题,每个月更新的制冷与空调设备 ...

  8. 非静压模型NHWAVE学习(6)——波浪模拟算例学习(Periodic wave over a submerged bar)

    波浪模拟算例学习(Periodic wave over a submerged bar) 算例简介 网格与地形 运行参数配置 模型的编译运行及模拟结果 本blog介绍了NHWAVE模型自带算例Peri ...

  9. Web前端学习仿拼多多软件——模拟练习

    Web前端学习仿拼多多软件模拟练习,作者声明:该项目只是学习模拟软件前端的界面,不允许用于商业行为,否则后果自负! 代码界面截图(因为只是学习分享,所以已经多处避免雷同): 下载资源后运行的时候要记得 ...

  10. 操作系统学习之用C语言模拟FIFO算法

    前言 FIFO算法比较简单,没什么好说的,就是先进先出.不过我添加了3状态,不过也只有堵塞,没有将阻塞进程唤醒的过程. Buddy算法:操作系统学习之用C语言模拟伙伴(Buddy)算法 FIFO算法: ...

最新文章

  1. boost::geometry::densify用法的测试程序
  2. Linux如何在栈上分配内存,Java堆与栈内存分配及String小记
  3. java cygwin 乱码_windows10乱码怎么解决
  4. android box2d运行小例子,配置eclipse【转】
  5. pmo与敏捷团队可以共存吗_匿名与问责制可以共存吗?
  6. Linux设置node的process.env.NODE_ENV
  7. 国内银行卡BIN号速查简表(2016)
  8. 跌疯了,理性投资人会这样做 附国内外免费股票数据源
  9. java速算24,java枚举速算24正解集合
  10. 算法:什么是宠物收养所问题?
  11. chatgpt智能提效职场办公-ppt怎么设置背景图片
  12. cgroup driver: cgroupfs还是systemd
  13. ftl和html的区别?
  14. BCD工艺、CMOS工艺、BiCMOS工艺
  15. mysql 表 视图 索引吗_Mysql多表查询, 视图,事务,索引,函数,go连接数据库
  16. Python是什么及小龙虾的面向对象
  17. 前端加密解密 crypto-js
  18. 不用学Python的自动化办公 确定不来看看?
  19. 【演示工具】Focusky教程 | 删除录音
  20. SAMSON阀门定位器3730-1的操作指南

热门文章

  1. java多分支流程图_Java 流程控制 之 分支结构(条件判断)
  2. 雾里探花之计算机硬件及组成原理(前言)
  3. u8系统计算机上启动不了,用友erp u8装好后为何启动不了
  4. Python数据分析之理论知识
  5. 发几个iphone助手相关代码,供大家参考
  6. 达梦共享存储集群DMDSC-2节点部署手册
  7. 视频教程_干货视频教程分享
  8. c# 中崎_C# 打开钱箱支持北洋、佳博、爱普生
  9. 系统分析师视频教程-张友生
  10. 数据库系统概论(第5版)