[gmx-users] Effect of time step size

Faezeh Pousaneh fpoosaneh at gmail.com
Tue Feb 10 21:41:45 CET 2015


Dear Mark,

Thank you so much for the reply.
I had sent the .mdp file in my previous email (below). In more details I
had equilibriuted 1728 methanol (NPT ensemble at T=290 and P=1atm, using
all-bond constraint) for long enough time (10 ns) and I am sure I have
reached the equilibrium according to RMSD analysis.

1- You are right about time step with constraint, however I have been
running for dt=0.002 as well, and I get density 790 kg/m^3. it means that
still there is a big difference in result of dt=0.001 and dt=0.002:
>  783 kg/m^3 for dt=0.001
>  790 kg/m^3 for dt=0.002
isn't big difference? then how it is said that dt can be at most 0.002?

2- I see that in the experiment the density of ethanol is 791 Kg/m^3 at
T=290. That is not in agreement with my result for dt=0.001. What does it
mean? does it mean that the given topology of ethanol (I take it from
Gromacs homepage) is not well defined?

 integrator                = md
 dt                            = 0.001
 nsteps                    = 10000000
 nstxout                   = 0                 ; save coordinates every 0 ps
 nstvout                   = 0                 ; save velocities every 0 ps
 nstlog                     = 0                 ; update log file every 2 ps
 nstenergy               = 1000           ; save energies every 2 ps
 nstxtcout                = 100000        ; Output frequency for xtc file
 xtc-precision          = 100000        ; precision for xtc file
 ; Neighborsearching
 ns_type                   = grid        ; search neighboring grid cells
 nstlist                      = 5             ;
 pbc                          = xyz         ; 3-D PBC
 rlist                          = 1.0         ; short-range neighborlist
cutoff (in nm)
 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
 vdw-type                   = Cut-off
 ; T Coupling
 Tcoupl                       = v-rescale   ; modified Berendsen thermostat
 tau_t                          = 0.1            ; time constant, in ps
 ref_t                           = 290.          ; reference temperature,
one for each group, in K
 tc-grps                       = system
 ; P Coupling
 Pcoupl                       = Berendsen; Parrinello-Rahman
 Pcoupltype                = Isotropic
 tau_p                         = 1.0
 compressibility           = 5e-5
 ref_p                          = 1.0
 ; Velocity generation
 gen_vel                      = no
 ; Dispersion correction
 DispCorr                     = EnerPres   ; account for cut-off vdW scheme
; OPTIONS FOR BONDS
 constraints                    = all-bonds   ; all bonds constrained
(fixed length)
 continuation                  = yes            ; Restarting after NPT
 constraint-algorithm      = lincs          ; holonomic constraints
 lincs_iter                       = 1                ; accuracy of LINCS
 lincs_order                    = 4                ; also related to
accuracy


Many thanks in advance.
------------------
Best regards
*Faezeh Pousaneh*

On Tue, Feb 10, 2015 at 6:33 PM, Mark Abraham <mark.j.abraham at gmail.com>
wrote:

> On Tue, Feb 10, 2015 at 12:16 PM, Faezeh Pousaneh <fpoosaneh at gmail.com>
> wrote:
>
> > Hi,
> >
> > I am simulating methanol (using its topology given in Gromacs homepage,
> for
> > 1728 molecules). I run the simulation (npt, P=1atm and T=290) for
> different
> > time step size, dt=0.001, 0.003 fs. But I obtain different results for
> each
> > of them:
> >
> >  783 kg/m^3 for dt=0.001
> >  790 kg/m^3 for dt=0.003
> >
> > Why is it so?
> >
>
> One needs to have sampled long enough to have reached equilibrium and then
> made enough observations for a sensible average. You haven't told us what
> you've done there, so we'll have to assume that. The size of the time step
> is closely tied to the model physics - you want to have a healthy number of
> time steps per fastest motion in the system, which for atomistic MD with
> some kind of bond constraints normally is at most 2fs. See manual section
> on removing fast motions for some details. So if dt=3fs is reproducibly
> wrong, then it looks like a case of garbage in = garbage out.
>
> Mark
>
>
> > I have also checked the point for lutidine molecule and I see similar
> > behavior. Below is my production.mdb file:
> >
> >
> >  ; RUN CONTROL PARAMETERS
> >  integrator                = md
> >  dt                        = 0.001       ; 2 fs
> >  nsteps                    = 30000000     ; 0.002 * 50000 = 1000 ps = 1ns
> >
> >  ; OUTPUT CONTROL OPTIONS
> >  nstxout                   = 0           ; save coordinates every 0 ps
> >  nstvout                   = 0           ; save velocities every 0 ps
> >  nstlog                    = 0        ; update log file every 2 ps
> >  nstenergy                 = 1000        ; save energies every 2 ps
> >  nstxtcout                 = 100000        ; Output frequency for xtc
> file
> >  xtc-precision             = 100000        ; precision for xtc file
> >
> >  ; Neighborsearching
> >  ns_type                   = grid        ; search neighboring grid cells
> >  nstlist                   = 5           ; 10 fs
> >  pbc                       = xyz         ; 3-D PBC
> >  rlist                     = 1.0         ; short-range neighborlist
> cutoff
> > (in nm)
> >  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
> >  vdw-type                  = Cut-off
> >
> >  ; T Coupling
> >  Tcoupl                    = v-rescale   ; modified Berendsen thermostat
> >  tau_t                     = 0.1         ; time constant, in ps
> >  ref_t                     = 290.        ; reference temperature, one for
> > each group, in K
> >  tc-grps                   = system
> >
> >  ; P Coupling
> >  Pcoupl                    = Berendsen; Parrinello-Rahman
> >  Pcoupltype                = Isotropic
> >  tau_p                     = 1.0
> >  compressibility           = 5e-5
> >  ref_p                     = 1.0
> >
> >  ; Velocity generation
> >  gen_vel                   = no
> >
> >  ; Dispersion correction
> >  DispCorr                  = EnerPres   ; account for cut-off vdW scheme
> >
> >  ; OPTIONS FOR BONDS
> >  constraints               = all-bonds   ; all bonds constrained (fixed
> > length)
> >  continuation              = yes         ; Restarting after NPT
> >  constraint-algorithm      = lincs       ; holonomic constraints
> >  lincs_iter                = 1           ; accuracy of LINCS
> >  lincs_order               = 4           ; also related to accuracy
> >
> > Best regards
> > *Faezeh Pousaneh*
> > Best regards
> > *Faezeh Pousaneh*
> > --
> > 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.
>


More information about the gromacs.org_gmx-users mailing list