[gmx-users] Unphysical conformations in decoupled free energy simulation

[Justin Lemkul]

On 8/5/13 9:52 AM, Joerg Sauter wrote:
> On 08/05/2013 02:15 PM, Justin Lemkul wrote:
>>
>>
>> On 8/5/13 8:12 AM, Joerg Sauter wrote:
>>> On 08/05/2013 01:54 PM, Justin Lemkul wrote:
>>>>
>>>>
>>>> On 8/5/13 7:44 AM, Joerg Sauter wrote:
>>>>> Dear all,
>>>>>
>>>>> I am trying to compute the free energy of hydration for cellobiose (a beta
>>>>> (1-4)
>>>>> glucose dimer) using BAR in Gromacs 4.6.1. However, I encounter a problem.
>>>>>
>>>>> I find that the vacuum conformations of the molecule in a regular vacuum
>>>>> simulation differ from the conformations in the decoupled simulation in the
>>>>> free
>>>>> energy case i.e., with the additional mdp entries:
>>>>>
>>>>> free-energy              = yes
>>>>> init-lambda                 = 0
>>>>> couple-lambda0           = none
>>>>> couple-moltype           = solute
>>>>> couple-intramol          = no
>>>>>
>>>>> Here is a histogram of the dihedral angles of the glycosidic linkage in
>>>>> the vacuum case
>>>>> https://dl.dropboxusercontent.com/u/70358077/reg.pdf (it stays in the global
>>>>> free energy minimum)
>>>>> and this is the decoupled free energy case (same starting conformation in the
>>>>> global free energy minimum)
>>>>> https://dl.dropboxusercontent.com/u/70358077/fe.pdf
>>>>>
>>>>> I understand that Gromacs replaces the intramolecular interactions with
>>>>> explicit
>>>>> pair
>>>>> interactions. Therefore, I had to increase table-extension but that did not
>>>>> change much. I was thinking that maybe this could be a problem specific to
>>>>> this
>>>>> topology (a GLYCAM06h conversion from Amber), however, I do not see how
>>>>> this can
>>>>> occur.
>>>>> I hope someone has an idea what is going wrong.
>>>>>
>>>>> mdp file: https://dl.dropboxusercontent.com/u/70358077/md.mdp (same behaviour
>>>>> for longer runs)
>>>>> top file: https://dl.dropboxusercontent.com/u/70358077/cellob.top
>>>>> gro file: https://dl.dropboxusercontent.com/u/70358077/cellob.gro
>>>>>
>>>>> This is sort of a minimal example. The same problem occurs
>>>>> in a regular simulation when decoupling from water using multiple lambdas.
>>>>>
>>>>
>>>> Have you tried running with rlist=rcoulomb=rvdw = 0 rather than 99? If nothing
>>>> else, I would try again with version 4.6.3.  There have been tons of fixes to
>>>> the free energy code since 4.6.1, though I don't know offhand which, if any,
>>>> would be affecting your results. Regardless, the current version should always
>>>> be the "best" version ;)
>>>>
>>>> -Justin
>>>>
>>> Hi Justin,
>>>
>>> I just tried it with version 4.6.3. and rlist=rcoulomb=rvdw = 0. Unfortunately,
>>> I get the same results.
>>>
>>
>> OK, good to know.  We can rule out 4.6.1-era bugs, then.  Does visualization
>> or an analysis of energy terms stored in the .edr file provide any clues as to
>> the source of the problem?
>>
>> -Justin
>>
> Visualization does not show any abnormalities beyond odd confomations in the
> free energy coupling case.
>
> Regarding the energies: I am not sure whether this is normal behaviour, but, in
> the free energy coupling case I do not see energy values for LJ 14 or Coulomb 14
> interactions. I thought that by selecting couple-intramol=0 all intramolecular
> interactions would be replaced by pairtype (14) interactions.
>
> In addition, for Coulomb (SR) I get in the free energy case
> LJ (SR)                    -22.8084      0.076    5.90879 -0.379238 (kJ/mol)
> in comparisson to
> LJ (SR)                     16.8295       0.81    13.2303 4.51371  (kJ/mol)
> in the regular vacuum simulation. So this looks like the 14 interactions get
> shifted to the LJ/Coulomb terms in the coupling simulation?! (However, the
> numbers do not add up)
>

This is smelling buggy to me, and there have been issues with 1-4 interactions 
in the free energy code before (http://redmine.gromacs.org/issues/1225), so 
there could be a lingering problem.  Hopefully Michael will chime in.  It might 
be worthwhile to open up a bug report on redmine.gromacs.org with all of the 
information you have provided (input files, output data, and the very nice 
description of the issue you have posted).

-Justin

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

Justin A. Lemkul, Ph.D.
Postdoctoral Fellow

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

jalemkul at outerbanks.umaryland.edu | (410) 706-7441

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