[gmx-users] Simulations in NVE ensemble with CHARMM36 - update

[ABEL Stephane 175950]

Hello Justin

Thanks for your reply. I have finally found why I obtained a large drift and decrease of the temp of my system during the NVE : I stupidly mixed different cpt and tpr files. If I use the right files and the following parameters below suggested by M. Shirts in  [1] for my tests:

rlist           = 1.3
;vdW
vdwtype         = Switch
rvdw            = 1.1
rvdw_switch     = 1.0

; Electrostatics
coulombtype     = PME
rcoulomb        = 1.3

The energy and the temp of my system remain stable during the 100 ps of my simulation after  initial small decrease of the energy and the temp to 295.5 K (initial temp was 297.904 K)

[1] http://comments.gmane.org/gmane.science.biology.gromacs.user/42686

Please also note that with the parameters for Verlet options you gave I obtain the following error with 4.6.5 for simulations with NVE, 

 Temperature coupling is required for calculating rlist using the energy
  drift with verlet-buffer-drift > 0. Either use temperature coupling or
  set rlist yourself together with verlet-buffer-drift = -1

If I add the verlet-buffer-drift = -1 option in my mdp (not sure it is correct), grompp gives only two notes 

NOTE 1 [file run_CH_1_NVE_pot_mod.mdp]:
  With Verlet lists the optimal nstlist is >= 10, with GPUs >= 20. Note
  that with the Verlet scheme, nstlist has no effect on the accuracy of
  your simulation.


NOTE 2 [file run_CH_1_NVE_pot_mod.mdp]:
  rlist is equal to rvdw and/or rcoulomb: there is no explicit Verlet
  buffer. The cluster pair list does have a buffering effect, but choosing
  a larger rlist might be necessary for good energy conservation.


Stéphane


On 5/7/14, 9:05 AM, ABEL Stephane 175950 wrote:
> Dear GMX-users
>
> I have done a 90 ns long simulation of a heterogeneous  system (SPCE water, phospholipids and a hydrophobic solvent) in NPT ensemble (298 K at 1 bar) with GROMACS4.6.5. to compute related dynamical properties, I would like to continue this simulation in NVE. But I am puzzled to choose the mdp parameters to have a good energy conservation. Indeed, I have checked the gromacs 2004 JCTC or the following and follow the advices given in the following threads
>
> http://comments.gmane.org/gmane.science.biology.gromacs.user/42686
> https://www.mail-archive.com/gmx-users@gromacs.org/msg60393.html
>
> But I always obtain for my system a large drift and the temp decreases significantly (20 K) after few dozens of ps (the initial temp was around 297.3 K). For this simulation I used the double precision of gromacs with a timestep of 1 fs.
>

Can you provide a full .mdp file so we can look at all your settings?  Are these
drifts during NPT or NVE?

> More, since I have never used the potential-modifier features, so I not sure how to setup  the correct parameters in the context of CHARMM simulations.
>

With the caveat that I have NOT thoroughly tested the following, my
understanding of the new Verlet options indicates that these settings should
perform pretty well using CHARMM36:

cutoff-scheme = Verlet
rlist = 1.2
vdwtype = cutoff
vdw-modifier = Potential-switch
rvdw = 1.2
rvdw-switch = 1.0
coulombtype = PME
rcoulomb = 1.2

Again, barely tested on my part, and never tested for NVE.

-Justin

--
==================================================

Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 601
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

==================================================