I ran equilibrium simulations, and then a longer simulations with the recommended sd integrator. Still, the mass_lambda dH/dL is nonzero and sometimes it is the dominant contribution at certain lambdas. I don't think it is a convergence issue because there aren't huge fluctuations in dH/dL values seen in dhdl.xvg.
I also ran simulations where only the mass changed as a function of lambda. Others (i.e. Electrostatic, Van der Waals, bonded, and restraints) were kept constant. The mass still has an effect on the Helmholtz free energy change. Is there possibly a mistake with how I've prepared the 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 ... 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 ... [ bonds ] ... [ angles ] ... etc... On Thu, May 11, 2017 at 11:33 AM, Justin Lemkul <jalem...@vt.edu> wrote: > > > On 5/11/17 9:49 AM, Dan Gil wrote: > >> Thank you Dr. Lemkul, >> >> The simulations ran, but I have some a question about the results. >> >> Is it natural to have a nonzero mass_lambda dH/dL? My intuition was that >> the potential energy does not depend on mass, and my thinking was that the >> mass contribution should be very small, if not zero. >> >> > It should. Maybe it's a convergence issue? > > > -Justin > > Best Regards, >> >> Dan >> >> On Wed, May 3, 2017 at 10:46 AM, Justin Lemkul <jalem...@vt.edu> wrote: >> >> >>> >>> 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 < >>>> tasneemkausa...@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.gil9...@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 <jalem...@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 <jalem...@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 >>>>>>>>> >>>>>>>>> jalem...@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-requ...@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 >>>>>>> >>>>>>> jalem...@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-requ...@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-requ...@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-requ...@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 >>> >>> jalem...@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-requ...@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 > > jalem...@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-requ...@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-requ...@gromacs.org.