Hi Berk and Tsjerk,
sounds reassuring since the last version of the GROMOS force field (G53a5
and a6) has been parametrized to reproduce DGsol using reaction field.
But, you agree Berk that when you want to calculate the solvation free
energy in the "GROMOS philosophy" (e.g. for reproducing or computing new
parameters), you have to use RF even if energy conservation might not
be as good as with PME?
Thanks for your answers and your interest in the thread.
Cheers,
Patrick
On Tue, 8 May 2007, Berk Hess wrote:
From: "Tsjerk Wassenaar" <[EMAIL PROTECTED]>
Reply-To: Discussion list for GROMACS users <gmx-users@gromacs.org>
To: "Discussion list for GROMACS users" <gmx-users@gromacs.org>
Subject: Re: [gmx-users] Very large fluctuations in dg/dl
Date: Tue, 8 May 2007 06:39:50 +0200
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
Luckily it does not.
And PME is not only better for energy conservation, but also captures
long-range interaction, which for instance leads to more stable proteins.
I have compared solvation free energies for G53A6 with reaction-field
and PME and found no differences:
B. Hess and N. van der Vegt
Hydration Thermodynamic Properties of Amino Acid Analogues: A Systematic
Comparison of Biomolecular Force Fields and Water Models
J. Phys. Chem. B 110, 17616 (2006)
http://dx.doi.org/10.1021/jp0641029
Berk.
On 5/7/07, Berk Hess <[EMAIL PROTECTED]> wrote:
>From: Patrick Fuchs <[EMAIL PROTECTED]>
>Reply-To: Discussion list for GROMACS users <gmx-users@gromacs.org>
>To: gmx-users@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@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 [EMAIL PROTECTED]
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
_______________________________________________
gmx-users mailing list gmx-users@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 [EMAIL PROTECTED]
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
_________________________________________________________________
Live Search, for accurate results! http://www.live.nl
_______________________________________________
gmx-users mailing list gmx-users@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 [EMAIL PROTECTED]
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
_______________________________________________
gmx-users mailing list gmx-users@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 [EMAIL PROTECTED]
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
