[gmx-users] Using the md integrator for calculating free energy of solvation
Justin Lemkul
jalemkul at vt.edu
Wed May 3 16:46:31 CEST 2017
On 5/3/17 10:43 AM, Dan Gil wrote:
> Hi Dr. Kausar,
>
> If I want the molecule to have all of its nonbonded interactions with the
> rest of the system (solvent) at each lambda, don't I want to have vdwq
> option ON for both couple-lambda0 and couple-lambda1?
>
> What I am attempting to do is mutate molecule A to molecule B through the
> simulation.
>
Relative free energy calculations are a little bit of a strange beast; you do
need to set the B-state (couple-lambda1) to "none" and the B-state types and
charges are read from [atoms].
-Justin
> Best Regards,
>
> Dan
>
> On Wed, May 3, 2017 at 1:04 AM, Tasneem Kausar <tasneemkausar12 at gmail.com>
> wrote:
>
>> In the last line of mdp file couple-lambda0 and couple-lambda1 have the
>> same vdwq variable. It indicates A and B state of system. So grompp is
>> complaining about the identical states.
>>
>> You can also see the alchemistry.org for detail information about the free
>> energy calculations.
>>
>> On Wed, May 3, 2017 at 2:14 AM, Dan Gil <dan.gil9973 at gmail.com> wrote:
>>
>>> Thank you Dr. Lemkul,
>>>
>>> Following your advice, I was able to calculate the solvation free energy
>> of
>>> two molecules that I am interested in. The results, unfortunately, were
>> not
>>> what I expected. I want to try an alternative method, that is, the
>>> alchemical thermodynamic integration method.
>>>
>>> If I understand this correctly, I should be able to calculate the free
>>> energy difference between two systems A and B as I change the topology of
>>> the molecule from A to B in small steps.
>>>
>>> I wasn't able to find a tutorial online, but I attempted to try the
>> method
>>> and I obtained this error from grompp:
>>>
>>> WARNING 1 [file grompp.mdp, line 43]:
>>> The lambda=0 and lambda=1 states for coupling are identical
>>>
>>> I was hoping if you could help me understand what I am missing in my
>>> topology or mdp file.
>>>
>>> Topology:
>>> [ moleculetype ]
>>> HEPT 3
>>>
>>> [ atoms ]
>>> 1 opls_135 1 HEPT CH3 1 -0.180 12.011
>>> opls_961 0.360 12.011
>>> 2 opls_136 1 HEPT CH2 2 -0.120 12.011
>>> opls_962 0.240 12.011
>>> 3 opls_136 1 HEPT CH2 3 -0.120 12.011
>>> opls_962 0.240 12.011
>>> 4 opls_136 1 HEPT CH2 4 -0.120 12.011
>>> opls_962 0.240 12.011
>>> 5 opls_136 1 HEPT CH2 5 -0.120 12.011
>>> opls_962 0.240 12.011
>>> 6 opls_136 1 HEPT CH2 6 -0.120 12.011
>>> opls_962 0.240 12.011
>>> 7 opls_135 1 HEPT CH3 7 -0.180 12.011
>>> opls_961 0.360 12.011
>>> 8 opls_140 1 HEPT H 1 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 9 opls_140 1 HEPT H 1 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 10 opls_140 1 HEPT H 1 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 11 opls_140 1 HEPT H 2 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 12 opls_140 1 HEPT H 2 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 13 opls_140 1 HEPT H 3 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 14 opls_140 1 HEPT H 3 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 15 opls_140 1 HEPT H 4 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 16 opls_140 1 HEPT H 4 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 17 opls_140 1 HEPT H 5 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 18 opls_140 1 HEPT H 5 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 19 opls_140 1 HEPT H 6 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 20 opls_140 1 HEPT H 6 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 21 opls_140 1 HEPT H 7 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 22 opls_140 1 HEPT H 7 0.060 1.008
>>> opls_965 -0.120 18.998
>>> 23 opls_140 1 HEPT H 7 0.060 1.008
>>> opls_965 -0.120 18.998
>>> ... And so on with the bonds, pairs, angles, and dihedrals directive.
>>>
>>> MDP File:
>>>
>>> ;Integration Method and Parameters
>>> integrator = md
>>> nsteps = 10000
>>> dt = 0.002
>>> nstenergy = 1000
>>> nstlog = 5000
>>>
>>> ;Output Control
>>> nstxout = 0
>>> nstvout = 0
>>>
>>> ;Cutoff Schemes
>>> cutoff-scheme = group
>>> rlist = 1.0
>>> vdw-type = cut-off
>>> rvdw = 2.0
>>>
>>> ;Coulomb interactions
>>> coulombtype = pme
>>> rcoulomb = 1.0
>>> fourierspacing = 0.4
>>>
>>> ;Constraints
>>> constraints = all-bonds
>>>
>>> ;Temperature coupling
>>> tcoupl = v-rescale
>>> tc-grps = system
>>> tau-t = 0.1
>>> ref-t = 300
>>>
>>> ;Pressure coupling
>>> pcoupl = parrinello-rahman
>>> ref-p = 1.01325
>>> compressibility = 4.5e-5
>>> tau-p = 5
>>>
>>> ;Free energy calculation
>>> free-energy = yes
>>> init-lambda = 0
>>> delta-lambda = 0.00001
>>> couple-moltype = HEPT
>>> couple-lambda0 = vdwq
>>> couple-lambda1 = vdwq
>>>
>>> On Mon, Apr 3, 2017 at 2:05 PM, Justin Lemkul <jalemkul at vt.edu> wrote:
>>>
>>>>
>>>>
>>>> On 4/3/17 2:01 PM, Dan Gil wrote:
>>>>
>>>>> Thank you Dr. Lemkul,
>>>>>
>>>>> I am trying to run grompp with the md integrator, but I am getting
>> this
>>>>> error:
>>>>> "For proper sampling of the (nearly) decoupled state, stochastic
>>> dynamics
>>>>> should be used"
>>>>>
>>>>> Should I ignore this warning with the -maxwarn option and try running
>>> it?
>>>>> I
>>>>> will see if I obtain comparable values for ethanol.
>>>>>
>>>>>
>>>> People have already done such a comparison and that's why grompp is
>>>> telling you this - it's better to use the Langevin integrator. You'll
>>> get
>>>> better sampling, particularly towards the end states. Using -maxwarn
>>> tells
>>>> grompp "you're trying to prevent me from making a mistake, but I want
>> to
>>> do
>>>> it anyway" :) You'd better have a really, really good reason to try to
>>>> override it.
>>>>
>>>> -Justin
>>>>
>>>>
>>>> Best Regards,
>>>>>
>>>>> Dan
>>>>>
>>>>>
>>>>> On Mon, Mar 27, 2017 at 12:51 PM, Justin Lemkul <jalemkul at vt.edu>
>>> wrote:
>>>>>
>>>>>
>>>>>>
>>>>>> On 3/26/17 9:40 PM, Dan Gil wrote:
>>>>>>
>>>>>> Hi,
>>>>>>>
>>>>>>> I am following Dr. Sander Pronk's and Dr. Justin Lemkul's tutorial
>> on
>>>>>>> calculating free energy of solvation. Is it possible and
>> theoretically
>>>>>>> sound to use the md integrator instead of the sd integrator for
>> these
>>>>>>> calculations?
>>>>>>>
>>>>>>>
>>>>>>> Langevin dynamics gives better sampling so it is frequently used for
>>>>>> free
>>>>>> energy calculations. You may get comparable results with the
>> leap-frog
>>>>>> integrator, but I haven never done a side-by-side comparison.
>>>>>>
>>>>>> -Justin
>>>>>>
>>>>>> I have already done a considerable amount of work using md
>> integration,
>>>>>> and
>>>>>>
>>>>>>> I want to make sure that the free energy values I calculate are
>>>>>>> consistent
>>>>>>> with my previous work.
>>>>>>>
>>>>>>> If using the md integrator is not sound, is there an alternative way
>>> of
>>>>>>> calculating solvation energy that will be consistent?
>>>>>>>
>>>>>>> Best Regards,
>>>>>>>
>>>>>>> Dan
>>>>>>>
>>>>>>>
>>>>>>> --
>>>>>> ==================================================
>>>>>>
>>>>>> 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
>>>>>>
>>>>>> ==================================================
>>>>>> --
>>>>>> 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.
>>>>>>
>>>>>>
>>>> --
>>>> ==================================================
>>>>
>>>> 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
>>>>
>>>> ==================================================
>>>> --
>>>> 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.
>>
--
==================================================
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
==================================================
More information about the gromacs.org_gmx-users
mailing list