[gmx-users] solvation free energy

[Justin Lemkul]

On 7/17/15 12:08 PM, Justin Lemkul wrote:
>
>
> On 7/17/15 12:06 PM, David van der Spoel wrote:
>> On 17/07/15 17:55, Justin Lemkul wrote:
>>>
>>>
>>> On 7/17/15 11:29 AM, Daniele Veclani wrote:
>>>> I'm using gromacs version 5.0.4, so should I upgrade to version 5.0.5 ?
>>>>
>>>
>>> If you need [ pairs_nb ], yes.  But I actually don't think you do.  See
>>> below.
>>>
>>>> I read in the manual (section 5.3.4) that can turn off the 1-4
>>>> interactions
>>>> without ever touching the topology, by using the couple-moltype,
>>>> couple-lambda0, couple-lambda 1, and couple-intramol keywords.
>>>>
>>>
>>> I think you're reading this section backwards.  Admittedly, it's a
>>> confusing section, so let's simplify.
>>>
>>> If you turn off the nonbonded interactions of some molecule in water,
>>> the total contribution to the free energy change is composed of the
>>> solute-solute terms and solute-water terms.  To get the actual solvation
>>> free energy, one needs to correct for the solute-solute terms being
>>> turned off, which is the purpose of the gas-phase simulations.  Turning
>>> off interactions in vacuo tells you how much of the free energy change
>>> is attributed to intramolecular terms.  The difference between the two
>>> transformations is then the solvation free energy.
>>>
>>> But, several years ago, this procedure became much easier in GROMACS
>>> with the introduction of the "couple-intramol" keyword.  It *eliminates*
>>> the need for any gas-phase simulations or topology hacking.  This is
>>> what the manual is trying to tell you.  If you run your transformations
>>> in water with "couple-intramol = no," then the intramolecular
>>> interactions are always full strength, so the free energy difference you
>>> calculate is solely attributed to solute-water interactions, viz. the
>>> solvation free energy.
>>>
>>> If you *want* to make your life more difficult, set "couple-intramol =
>>> yes" and conduct gas-phase calculations, as well, which is what requires
>>> topology hacking.  I wouldn't do this.  Not only is the manual somewhat
>>> unclear, but any time you have to manually modify the topology, there is
>>> a potential for error. The free energy code attempts to make your life
>>> easy and efficient, unless you have some very complex protocol to do
>>> that requires lots of work.  A solvation energy is not one of these
>>> cases; it is straightforward with appropriate use of the .mdp keywords.
>>>
>> If I understand that correctly then the protocol ignores the difference in
>> intramolecular energy between gas-phase and solvation. This is probably OK for
>> small and rigid molecules but not necessarily for larger or flexible molecules.
>> When using polarizable force fields it will probably induce a significant error.
>> Have you tried this yet?
>>
>
> We have not done a comparison.  In CHARMM, we always explicitly do gas phase,
> but that is by necessity.  We don't have an analogous "couple-intramol"-style
> behavior.  It could probably be done, but we have an established protocol that
> has been used for years that is rock solid.
>

I would also agree that there would be massive errors for a polarizable force 
field in this case due to environmental response.  I assume, for now, most 
people are using additive models.  In fact, until I have actually tested the 
free energy code with our Drude model, it is explicitly disabled by grompp.

-Justin

-- 
==================================================

Justin A. Lemkul, Ph.D.
Ruth L. Kirschstein NRSA Postdoctoral Fellow

Department of Pharmaceutical Sciences
School of Pharmacy
Health Sciences Facility II, Room 629
University of Maryland, Baltimore
20 Penn St.
Baltimore, MD 21201

jalemkul at outerbanks.umaryland.edu | (410) 706-7441
http://mackerell.umaryland.edu/~jalemkul

==================================================