[gmx-users] Constraint and Temperature Drop in NVE

[Mark Abraham]

Hi,

On Wed, Sep 10, 2014 at 5:40 PM, Johnny Lu <johnny.lu128 at gmail.com> wrote:

> hmm. If that is the case, I have two questions.
> 1. Why gromacs doesn't adopt the method of measuring temperature as
> described in the JCTC paper?
>     (
>
> https://www.deshawresearch.com/publications/Equipartition%20and%20the%20Calculation%20of%20Temperature%20in%20Biomolecular%20Simulations.pdf
> )
>

As they note in the discussion, this is not a practical approach for use in
a *thermostat* because you need to estimate the time derivative of the
position. During the simulation you can only do a one-sided interpolation
from previous points, so using such an estimator makes the integration
scheme irreversible, which has its own drawbacks. So, an error exists in
the momentum-squared-based estimator of the temperature (from
finite-time-step size effects), and feeds back into the NVT integration
scheme. Whether this matters for the quality of the observables from a
given simulation is an open question, but their NVT results suggest that
time steps of the usual sizes are not good in combination with multiple
thermostats on small-ish systems. They used a Berendsen thermostat, which
is known to sample the wrong ensemble, so their NVT results may not be
quantitatively accurate, but offhand I'd expect to observe the trend they
see.

Naturally, one could write code to do two-sided interpolation as they do
there and so get a reliable temperature estimate in post-processing, but
nobody has bothered to do that. BTW, GROMACS does report the conserved
quantity for all integration schemes that have one.

Mark

2. As the paper mentioned that using two thermostats can cause hot solvent
> cold solute problem for the timesteps that we commonly use, why do your
> lysozyme tutorial still use two thermostats?
>     (
>
> http://www.bevanlab.biochem.vt.edu/Pages/Personal/justin/gmx-tutorials/lysozyme/Files/nvt.mdp
> )
>
> Thank you.
>
>
> On Wed, Sep 3, 2014 at 4:07 PM, Mark Abraham <mark.j.abraham at gmail.com>
> wrote:
>
> > Hi,
> >
> > It's impossible to say anything in the absence of a full description of
> the
> > system and .mdp file - expecting good conservation can be a fool's
> errand.
> > You could
> > a) investigate the effect of time step on the difference in the observed
> > temperature
> > b) investigate the effect of the restraint strength constant on the
> > temperature difference (it could be too large for the time step)
> > c) try double precision (to see if the drift is intrinsic or due to
> > accumulated rounding imprecision)
> > d) assess which of the approximations you are making (drift from finite
> > time step/model/implementation, static point charges, broken ergodicity,
> > etc.) dominates the quality of your model physics.
> >
> > And consider that the measurement of temperature from momentum (as in
> > GROMACS) is intrinsically flawed - see
> >
> >
> https://www.deshawresearch.com/publications/Equipartition%20and%20the%20Calculation%20of%20Temperature%20in%20Biomolecular%20Simulations.pdf
> >
> > Mark
> >
> >
> > On Wed, Sep 3, 2014 at 6:21 PM, Johnny Lu <johnny.lu128 at gmail.com>
> wrote:
> >
> > > I try the following:
> > >
> > > 0. Minimization
> > > 1. NVT (100,000 steps of 2fs, 300K)
> > > 2. NPT (100,000 steps of 2fs, 300K, V-rescale thermostat)
> > > 3. NPT (5000 steps of 2fs, 308K, V-rescale thermostat)
> > > 4. 9 x NVE (each 32768 steps of 0.75 fs)
> > > 5. go to step 3.
> > >
> > > I plan to get rid of the first of nine NVE steps, as an NVE
> equilibration
> > > step.
> > >
> > > At 0.75 fs, the total energy drop was about 0.03% per NVE step.
> > >
> > > Will that be fine? I am a bit worry that step 3 is too short.
> > > I pick 308K because last time when I choose 300K for step 3, step 4
> drops
> > > to 292 K.
> > > I hope that the NVE will be about 300K.
> > >
> > > May be I should look at the .edr file of step 3.
> > >
> > > Thanks.
> > >
> > >
> > > On Wed, Sep 3, 2014 at 11:54 AM, Johnny Lu <johnny.lu128 at gmail.com>
> > wrote:
> > >
> > > > reading this post in maillist again:
> > > >
> > >
> >
> https://mailman-1.sys.kth.se/pipermail/gromacs.org_gmx-users/2012-May/071225.html
> > > >
> > > > I think that is exactly the same problem.
> > > >
> > > >
> > > > On Wed, Sep 3, 2014 at 10:24 AM, Johnny Lu <johnny.lu128 at gmail.com>
> > > wrote:
> > > >
> > > >> Hi.
> > > >>
> > > >> After equilibrating with NVT (for 10,000 steps of 2fs) and then NPT
> > > >> (10,000 steps of 2fs with V-rescale thermostat at 300K), when I run
> > NVE
> > > >> (0.5fs timesteps) with a harmonic restraint with force constant
> > > 4184kJ/mol
> > > >> rad^2, the temperature then fluctuate around 292 K (
> > > >> http://oi60.tinypic.com/sni5at.jpg) in gromacs 4.6.6.
> > > >>
> > > >> When I instead run NVE with same setting except without the
> constraint
> > > >> after the two equilibration steps, the temperature stays around
> 300K.
> > > >>
> > > >> What is the cause of this drop of temperature to 292K?
> > > >>
> > > >> Thanks again.
> > > >>
> > > >>
> > > >>
> > > >
> > > --
> > > Gromacs Users mailing list
> > >
> > > * Please search the archive at
> > > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > > posting!
> > >
> > > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> > >
> > > * For (un)subscribe requests visit
> > > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > > send a mail to gmx-users-request at gromacs.org.
> > >
> > --
> > Gromacs Users mailing list
> >
> > * Please search the archive at
> > http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> > posting!
> >
> > * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
> >
> > * For (un)subscribe requests visit
> > https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> > send a mail to gmx-users-request at gromacs.org.
> >
> --
> Gromacs Users mailing list
>
> * Please search the archive at
> http://www.gromacs.org/Support/Mailing_Lists/GMX-Users_List before
> posting!
>
> * Can't post? Read http://www.gromacs.org/Support/Mailing_Lists
>
> * For (un)subscribe requests visit
> https://maillist.sys.kth.se/mailman/listinfo/gromacs.org_gmx-users or
> send a mail to gmx-users-request at gromacs.org.
>