|
用周期扰动法对甘油三酯的粘度进行模拟,我目前使用的流程是:
建模:建立分子模型-结构优化-packmol建立盒子-用LigParGen得到itp文件(电荷用Multiwfn计算的1.2*CM5原子电荷取代),然后include进top文件。
模拟:em-nvt退火-npt_eq平衡相-nvt_eq平衡相-nvt_vis生产相
细节:npt平衡先用的使用V-rescale方法和Berendsen方法在1000bar下进行模拟(1bar收敛太慢),然后换成1bar,等温度和压强达到设定值后换用Nose-hoover 方法和Parrinello-Rahman进行平衡。nvt平衡相使用的是Nose-hoover 方法控温。平衡相是否平衡参考sob老师的博文《谈谈怎么判断分子动力学模拟是否达到了平衡 - 思想家公社的门口:量子化学·分子模拟·二次元》,检查的rmsd曲线均已收敛。
我的猜测:在npt平衡相后使用gmx energy 从平衡开始对密度采样,平均值和实验值差不多,系统能量、体积和box的z方向长度也已经收敛了。所以我怀疑是后面nvt模拟的问题,这里贴上nvt_vis.mdp
另:附上itp文件和npt_eq.mdp(模拟是用 gmx convert-tpr 一段一段跑的,最后加到rmsd收敛为止,所以mdp文件里的模拟时长不够)。我浏览了所有含有粘度的帖子,还是不知道自己哪里有问题,解决不了。请各位老师帮忙看看,非常感谢。
integrator = md
nsteps = 5000000
dt = 0.001
comm-grps = system
; Output control
nstxout = 1000
nstvout = 1000
nstenergy = 10
nstlog = 1000
nstxout-compressed = 10
compressed-x-grps = system
nstcheckpoint = 1000
; Neighbor searching and cutoffs
pbc = xyz
cutoff-scheme = Verlet
coulombtype = PME
rcoulomb = 1.5
vdwtype = cut-off
rvdw = 1.5
Dispcorr = Enerpres
cos-acceleration = 0.025
; Temperature coupling with annealing
Tcoupl = nose-hoover
tc-grps = System
tau_t = 0.2
ref_t = 273.15
; Pressure coupling (disabled for NVT)
Pcoupl = no
Pcoupltype = isotropic
tau_p = 2.0
ref_p = 1.0
compressibility = 4.5e-5
gen_vel = no
gen_temp =273.15
gen_seed = -1
|
|