[gmx-users] Free energy calculation problems (bug?)

[Олег Титов]

Thank you for your reply, Dr Shirts.

I'll recalculate my energies with vdw-type = shift and DispCorr = EnerPres
and enabled softcore with both old and new version of GROMACS. However, I'm
not sure if softcore is necessary since aromatic hydrogens have VdW
potentials in GAFF so no atom disappears in this particular system.

If nothing changes, i'll file a readmine asap.


2014-02-08 18:49 GMT+04:00 Michael Shirts <mrshirts at gmail.com>:

> Two thoughts that probably don't explain the problem, since they
> should have similar behavior between versions:
>
> * vdw-type                 = Cut-off
>
> is almost never a good idea.  If you are doing this to be
> reproducible, it's fine, but I'd avoid it in the future.
>
> * You almost certainly want DispCorr = EnerPres, or else the density
> will change with cutoff, and dispersion interactions beyond the cutoff
> will be neglected (also changing the energy with cutoff).
>
> Can you file a redmine (http://redmine.gromacs.org) with .top  and
> .gro files as well?  One never knows where a problem might be coming
> from (maybe only with certain top files?), so having those helps as
> well.  It would be better to have the softcore version, since that is
> likely to be better behaved anyway, and it has the same difference.
>
>
> On Sat, Feb 8, 2014 at 12:00 AM, Олег Титов <titovoi at qsar.chem.msu.ru>
> wrote:
> > Dear GROMACS users!
> >
> > I've encountered a problem in thermodynamic integration with GROMACS. I'm
> > interested in calculating of solvation free energy differences between
> > holbenzenes and benzene. To calculate this values I use thermodynamic
> > integration approach (5 point formula).
> >
> > Half a year ago I did my calculation with GROMACS 4.6.1. The results
> looked
> > good and were consistent with literature data. Recently, I've discovered
> > that versions 4.6.0 and 4.6.1 had a bug in free energy code, so I decided
> > to recalculate everything with 4.6.5 version of GROMACS. New results look
> > not so good and, moreover, are not consistent with literature.
> >
> > I use GAFF forcefield with TIP3P water models. To obtain GROMACS topology
> > and coordinate files I convert AMBER topologies (generated with tleap)
> with
> > amb2gmx.pl script and this topologies look fine.
> >
> >  For phenylbromide molecule calculated solvation free difference should
> be
> > near -0.1 - -0.2 kcal/mol (Experimental data is 0.59 kcal/mol).
> > 4.6.1 version gave me -0.1+-0.13 kcal/mol - Fine!
> > I've performed the same calculation with AMBER11  - -0.27 +-0.16
> kcal/mol -
> > also fine.
> > 4.6.5 version for the same system gives me 1.1 - 1.2 kcal/mol
> >
> > I've played with some options (tried adding soft-core potentials and
> turned
> > on/off couple-intramol option). The results are:
> > no_sc no_intramol 1.27 +- 0.18
> > no_sc intramol 1.20 +- 0.20
> > sc no_intramol 1.35 +- 0.21
> > sc intramol 1.22 +- 0.19
> >
> > So now I don't know what to do. Old buggy version produced good results
> > consistent with literature and other MD packages and new fixed version
> > produces different results. The dH/dl curve for two versions of GROMACS
> > looks similar for lambdas 0-0.7 (mutation is PhBr -> PhH ) but becomes
> > different at lambda = 0.9 So something changed there... Moreover
> different
> > elecrostatic approaches for PhBr molecule are consistent with each other
> > (dddG are fine), but new version of GROMACS gives me overestimated ddG as
> > shown above.
> >
> > Can anyone give me some advice what to do with it? May be there is some
> bug
> > in new version? May be I'm doing something wrong?
> >
> > Attaching sample .mdp (for solvated and vacuum system) and .top files.
> >
> > Thanks for the help,
> > Oleg
> >
> > --
> > 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.
>