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

Fri Feb 14 10:24:46 CET 2014

Dear Prof. Shirts

We also have a few questions on that issue and we will really  
appreciate your help.

If we understand correctly:

  coupl-intramol = yes - remove all interactions of the decoupled atoms
  coupl-intramol = no - remove only solute solvent interactions pf the  
decoupled atoms.

You previuosly wrote:
" if one sets coupl-intramol = yes, LJ and Coul 1-4's
are left on (uncoupled) at all lambda. The actual bug is fixed in 5.0
(actually, in 4.6.4). "

The questions are the following:
1. Is there a way to remove also these 1-4's?
2. Is there a tutorial for relative free energy calculation so we can  
verify what we do?
3. Can we use charge at the end state B (rather than zero) ? is the  
same behavior expected with the coupl-intramol parameter when if we  
change the charge to a value other than zero?

4. Is this the proper way to turn charge off (through the topology file)?

  nr type resnr resid atom cgnr charge mass typeB chargeB  massB
  1   HC   1    PCH3   H1   1   0.123 1.0080 HC     0.0    1.0080

5. What is the proper way to turn VDW off and could the VDW parameters  
also be changed to another value rather than zero?

Thanks in advance,
Best regards,
Asaf

Quoting Олег Титов <titovoi at qsar.chem.msu.ru>:

> What should I do to be able to file a redmine? I've registered at
> redmine.gromacs.org, but the account is waiting for confirmation for 2
> days. Should I do anything else, or it's ok?
>
>
> 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.
>>
> --
> 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.