[gmx-users] Very large fluctuations in dg/dl

Tsjerk Wassenaar tsjerkw at gmail.com
Tue May 8 06:39:50 CEST 2007


Hi Berk,

Just a question relating to your reply (should be a new thread
maybe...). The Gromos force fields were parameterized (for use) with
cutoff/RF. PME is of course better for energy conservation, but the
inclusion of the interactions beyond the cut-off may give rise to
deviant behaviour from the force field. Could you enlighten me on this
issue?

Thanks,

Tsjerk

On 5/7/07, Berk Hess <gmx3 at hotmail.com> wrote:
>
>
>
> >From: Patrick Fuchs <Patrick.Fuchs at ebgm.jussieu.fr>
> >Reply-To: Discussion list for GROMACS users <gmx-users at gromacs.org>
> >To: gmx-users at gromacs.org
> >Subject: [gmx-users] Very large fluctuations in dg/dl
> >Date: Mon, 07 May 2007 16:41:26 +0200
> >
> >Hi Gromacs users,
> >I have a few questions related to solvation free energy calculation via
> >thermodynamic integration.
> >I'm trying to reproduce some literature data (on e.g. methane,
> >methanol...) using the GROMOS G53a6 force field. I followed the tutorial
> >of David Mobley (thanks to him BTW), but I used the standard non bonded
> >options of the G53a6 force field (instead of OPLS). For each lambda
> >value I do a minimization, a 10 ps NVT followed by a 20 ps NPT
> >equilibration, and a 1 ns NVT production using the sd integrator. I used
> >21 lambda values (0.00, 0.05...1.00).
> >Here's my topology file:
> >----------------begining of methane.top------------------------
> >; topology for a methane molecule
> >
> >; include GROMOS53a6 force field
> >#include "ffG53a6.itp"
> >
> >;;;;;;; begin methane definition ;;;;;;;
> >[ moleculetype ]
> >; Name           nrexcl
> >METH             3
> >
> >[ atoms ]
> >;nr type resnr residue atom cgnr charge mass    typeB chargeB massB
> >  1  CH4  1     METH    C1   0    0.0000 16.0430 DUM   0.0000  16.04300
> >;;;;;; end methane definition ;;;;;;;;
> >
> >; include water topology
> >#ifdef FLEX_SPC
> >#include "flexspc.itp"
> >#else
> >#include "spc.itp"
> >#endif
> >
> >[ system ]
> >; name
> >1 methane molecule in water
> >
> >[ molecules ]
> >; name  number
> >METH    1
> >SOL     893
> >-----------------end of methane.top------------------------
> >
> >And here is my mdp file for lambda=0:
> >---------------begining of prod.mdp---------------------
> >title               = production NVT methane/water
> >cpp                 = /lib/cpp
> >; OPTIONS FOR BOND CONSTRAINTS
> >constraints         = all-bonds
> >; RUN CONTROL PARAMETERS
> >integrator          = sd
> >tinit               = 0
> >dt                  = 0.002
> >nsteps              = 500000 ; 1000 ps
> >; NUMBER OF STEPS FOR CENTER OF MASS MOTION REMOVAL
> >nstcomm             = 100
> >; OUTPUT CONTROL OPTIONS
> >nstxout                  = 500
> >nstvout                  = 500
> >nstfout                  = 0
> >nstlog                   = 500
> >nstenergy                = 100
> >nstxtcout                = 5000
> >xtc-precision            = 1000
> >; NON BONDED STUFF
> >ns_type             =  grid
> >nstlist                = 5
> >rlist                  = 0.8
> >coulombtype            = generalized-reaction-field
> >rcoulomb               = 1.4
> >rvdw                   = 1.4
> >epsilon_rf             = 54.0
> >;OPTIONS FOR TEMPERATURE COUPLING
> >tc_grps                  = system
> >tau_t                    = 0.1 ; inverse langevin friction cst
> >ref_t                    = 300
> >;OPTIONS FOR PRESSURE COUPLING
> >Pcoupl                   = no
> >tau_p                    = 0.5
> >compressibility          = 4.5e-5
> >ref_p                    = 1.0
> >; FREE ENERGY CONTROL STUFF
> >free_energy              = yes
> >init_lambda              = 0.00
> >delta_lambda             = 0
> >sc_alpha                 = 0.5
> >sc-power                 = 1.0
> >sc-sigma                 = 0.3
> >; VELOCITY GENERATION
> >gen_vel                  = yes
> >gen_temp                 = 300
> >gen_seed                 = -1
> >-----------------end of prod.mdp------------------------
> >I find a reasonable DeltaGsol value of 8.6 kJ/mol for methane (compared
> >to 8.7 in Geerke & van Gunsteren, ChemPhysChem 2006, 7, 671 – 678) but
> >I get really huge fluctuations in the values of dg/dl:
> >lambda=0.00: 5.0 +/- 10.8 (mean +/- standard deviation)
> >lambda=0.05: 4.3 +/- 11.2
> >...
> >lambda=1.00: -0.3 +/- 4.0
> >Furthermore, each of these mean value is very slow at converging (1 ns
> >seems a minimum for certain lambda values...).
> >I can't get reasonable fluctuations even if I sample more. In addition,
> >there are very frequent warnings in the log file such as:
> >----
> >Large VCM(group rest):      0.01363,      0.00818,      0.01147, ekin-cm:
> >3.09490e+00
> >----
> >Here are my questions:
> >1) Has someone an idea of what could be the cause of these [very] large
> >fluctuations? Does it come from my setup, or is this a normal behavior?
> >2) Are these 'Large VCM(group rest)' warnings related to the use of sd
> >integrator (when I switch to md integrator, I no longer get these warnings)
> >?
> >Thanks for your answer,
> >
> >Patrick
>
> Large fluctuations are very common, depending on your system
> and how your choose your intermediate states.
>
> But sampling would improve a lot if you use a larger tau_t, for instance
> 1 ps, this will also solve your large vcm warnings.
> Note that with reaction-field your will have bad energy conservation.
> With a tau_t of 1 ps your temperature might go up somewhat.
> Using PME is always better.
>
> Berk.
>
> _________________________________________________________________
> Watch free concerts with Pink, Rod Stewart, Oasis and more. Visit MSN In
> Concert today. http://music.msn.com/presents?icid=ncmsnpresentstagline
>
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/search before posting!
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>


-- 
Tsjerk A. Wassenaar, Ph.D.
Junior UD (post-doc)
Biomolecular NMR, Bijvoet Center
Utrecht University
Padualaan 8
3584 CH Utrecht
The Netherlands
P: +31-30-2539931
F: +31-30-2537623



More information about the gromacs.org_gmx-users mailing list