gromacs模拟退火问题

yangjunfang

各位老师同学大家好，

我想请教一个问题：我想通过退火模拟取不同温度下的结构进行优化，请问可行吗？下面是我的npt系综的mdp文件


title                = NVT with position restraint on heavy atoms of ttvp
;define               = -DPOSRES
integrator           = md               ;integrator
dt                   = 0.002            ;timestep
nsteps               = 25000000         ;50 ns

pbc                  = xyz              ;Periodic Boundaries in XYZ

nstlist              = 10               ;Neighb actualization steps
rlist                = 1.2              ;Cut for near neigh
ns_type              = grid             ;Type of neighb mapping

coulombtype          = PME              ;Coul type
pme_order            = 4                ;Order of cubic interpolation
rcoulomb             = 1.2              ;Cut for coulomb int
fourierspacing       = 0.16             ;Gird spaccing for FFTW
optimize_fft         = yes              ;Do optimize FFTW at the beggining of the MD
ewald_rtol           = 1e-5             ;Strength of Displaced edwald potential
DispCorr             = EnerPres

vdwtype              = switch           ;VDW type
rvdw_switch          = 1.0              ;where to start switching the LJ potential
rvdw                 = 1.2              ;Cut for LJ or Bukingham

; Temperature coupling is on
tcoupl               = V-rescale        ;Thermostat on (T)
tc-grps              = system
tau_t                = 0.1              ;Temperature coupling dT in ps
ref_t                = 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 x-y box vectors, independent z
tau_p           = 2.0                   ; time constant, in ps
ref_p           = 1.0                   ; reference pressure, x-y, z (in bar)
compressibility = 4.5e-5                ; isothermal compressibility, bar^-1               ;

continuation           = yes              ;continuing the simulation from the previous NVT equilibration phase
;gen_vel               = yes               ;Do gen initial velocities
;gen_temp              = 50               ;Do generate the initial temperatures at 50K
;gen_seed              = -1               ;Generate every time

;annealing simulation
annealing           =  single    ;Type of annealing
annealing_npoints   = 9        ;Number of annealing points
annealing_time      =  0 5000 10000 15000 20000 25000 30000 35000 40000   ;Time to reach the target
annealing_temp      =  300 310 320 330 340 350 360 370 380  ;Initial and final temperature

constraints          = hbonds           ;Bonds restrictions
constraint_algorithm = lincs            ;Algorithm to restrict
lincs_iter           = 1                ;Number of iterations of LINCS
lincs_order          = 4                ;Max exp order

nstcomm              = 10               ;Num of steps to remove the COM
comm-mode            = Linear           ;Type of COM removal
comm-grps            = system           ;Remove COM to what?

nstxout              = 20000          ;Frecuency to write coordinates to the TRR
nstvout              = 20000          ;Frecuency to write velocities to the TRR
nstxtcout            = 5000           ;Frequency to write to the XTC

nstlog               = 2000             ;Frecuency to write to the LOG
nstenergy            = 2000             ;Frecuency to write energies to the LOG

sobereva

退火这种非平衡过程甭用Parrinello-Rahman压浴

程序上可行不可行是一码事，有没有实际意义是另一码事，后者没细节没法说

yangjunfang

好的 老师  那我如果用已经平衡好的构型进行退火的话是不是就可以用Parrinello-Rahman压浴了？

还有我的思路是我想模拟聚集体对温度依赖性，所以想通过退火得到不同温度下的构型，然后将这些构型进行优化，想知道的是这种方法可不可行？

yangjunfang

sobereva 发表于 2022-2-22 21:02
退火这种非平衡过程甭用Parrinello-Rahman压浴

程序上可行不可行是一码事，有没有实际意义是另一码事， ...

老师，您好，我按照您说的把压浴换成了Berendsen，然后按照上述的mdp文件在NPT系综尝试对聚集体水溶液进行退火，结果出现了问题，如下：

Source file: src/gromacs/ewald/pme-redistribute.cpp (line 282)
MPI rank:    35 (out of 48)
147 particles communicated to PME rank 2 are more than 2/3 times the cut-off
out of the domain decomposition cell of their charge group in dimension x.
This usually means that your system is not well equilibrated.

For more information and tips for troubleshooting, please check the GROMACS
website at http://www.gromacs.org/Documentation/Errors
-------------------------------------------------------

Fatal error:
126 particles communicated to PME rank 5 are more than 2/3 times the cut-off
out of the domain decomposition cell of their charge group in dimension x.
This usually means that your system is not well equilibrated.

For more information and tips for troubleshooting, please check the GROMACS
website at http://www.gromacs.org/Documentation/Errors
-------------------------------------------------------
--------------------------------------------------------------------------
MPI_ABORT was invoked on rank 35 in communicator MPI_COMM_WORLD
with errorcode 1.

NOTE: invoking MPI_ABORT causes Open MPI to kill all MPI processes.
You may or may not see output from other processes, depending on
exactly when Open MPI kills them.

请问这是什么原因，是mdp文件有问题吗？应该怎么做呢？

sobereva

yangjunfang 发表于 2022-2-23 08:52
好的 老师  那我如果用已经平衡好的构型进行退火的话是不是就可以用Parrinello-Rahman压浴了？

还有我的 ...

退火是变温的过程，期间体系必然处于非平衡状态，当然不能用PR压浴

即便同一种温度下模拟，模拟过程中也会大产生大量不同构型，能量极小化完了也会得到不同极小点结构（除非体系真的极度简单）。做能量极小化没什么意义。直接取MD的一些帧说明当前温度下有代表性的结构是什么样就完了，必要时可以做簇分析取最大的簇中的代表性结构。

sobereva

yangjunfang 发表于 2022-2-24 10:41
老师，您好，我按照您说的把压浴换成了Berendsen，然后按照上述的mdp文件在NPT系综尝试对聚集体水溶液进 ...

mdrun接上-ntmpi 1再试，如果还这样，输入文件本身就有问题，导致体系严重不稳定

yangjunfang

sobereva 发表于 2022-2-24 11:28
退火是变温的过程，期间体系必然处于非平衡状态，当然不能用PR压浴

即便同一种温度下模拟，模拟过程中 ...

嗯嗯，好的，谢谢您的解答，我按照您说的把压浴换成了Berendsen，然后在NPT系综尝试对聚集体水溶液进行退火，结果出现了问题，如下：
Fatal error:
126 particles communicated to PME rank 5 are more than 2/3 times the cut-off
out of the domain decomposition cell of their charge group in dimension x.
This usually means that your system is not well equilibrated.
能麻烦您看一下该怎么解决呢？

yangjunfang

sobereva 发表于 2022-2-24 11:31
mdrun接上-ntmpi 1再试，如果还这样，输入文件本身就有问题，导致体系严重不稳定

好的，谢谢您

hanmuhua

yangjunfang 发表于 2022-2-24 10:41
老师，您好，我按照您说的把压浴换成了Berendsen，然后按照上述的mdp文件在NPT系综尝试对聚集体水溶液进 ...

您好，我也遇到了这样的问题，想请教一下。您解决了吗？