[gmx-users] non zero free energy values in mass transformation

deepak bapat dubapat at gmail.com
Sat Aug 12 06:42:31 CEST 2017


Thank you Andre and Dan.

I'll read more about.

Deepak

On Fri, Aug 11, 2017 at 7:35 PM, Dan Gil <dan.gil9973 at gmail.com> wrote:

> Hi Deepak,
>
> I am battling the same problem too. I am doing free energy calculation for
> heptane -> perfluoro-heptane in water. Also SPC/E water, Gromacs 5.1.
>
> Like Andre said, there is expected to be a free energy change when you
> change the mass of a particle. But if you do the same transformation from
> Li+ to Cl- in the gas phase, hopefully the free energy change will be
> opposite and (approximately) equal. Have you tried doing this second
> calculation [Li+ -> Cl-] in the gas phase? I think this website:
> http://www.alchemistry.org/wiki/Thermodynamic_Cycle explains why this
> second step might be needed better.
>
> Dan
>
> On Fri, Aug 11, 2017 at 9:37 AM, André Farias de Moura <moura at ufscar.br>
> wrote:
>
> > Deepak,
> >
> > I never did a mass perturbation myself, but I would expect some free
> energy
> > change, since the translational partition function depends on the mass of
> > the molecules.
> >
> > Also, having the same kinetic energy is not meaningful at all, the
> correct
> > quantity you should always use in statistical mechanics reasoning is the
> > momentum of the particle, p=mv, which clearly changes during your
> > alchemical transformation. Put in another way: the phase space of the
> > system is a 6N dimensional space of positions and momenta, and you have
> > just perturbed the momenta part.
> >
> > And by the way: that is a typical case in which Monte Carlo methods would
> > fail to reproduce an MD result, since only the configurational phase
> space
> > is sampled in MC (masses effects would have to be added afterwards,
> > assuming that some analytical result could be computed)
> >
> > I hope it helps
> >
> > Andre
> >
> > On Fri, Aug 11, 2017 at 6:36 AM, deepak bapat <dubapat at gmail.com> wrote:
> >
> > > Dear Gmx Users
> > >
> > > I am trying to find out free energy of transforming one neutral LI to
> one
> > > CL in aqueous solution.
> > >
> > > Water model is SPC/E. Gromacs version: 5.1.2
> > >
> > > I am doing something like this
> > >
> > > Insert 1 Li  olecule in water. (Li parameters were taken from Kasper P.
> > > Jensen and William L. Jorgensen 2006.)
> > >
> > > step 1: increase mass of atom from 6.94100 to 35.45300.
> > >
> > > Step 2: change LJ parameters of LI to LJ of CL
> > >
> > >
> > > Now my doubt is in step 1:  i get non zero free energy values.
> > >
> > > *required nonbonded.itp file content is*
> > >
> > > ;LI+
> > > opls_406   LI  3   6.94100     0.000       A    2.8700e-01 2.0929e-03
> > > ; Li+ parameters from jensen and jorgensen JCTC 2006 2 1499-1509
> > >
> > > ;LI+ grown
> > > opls_407   LI  17   35.45300     0.000       A    2.8700e-01 2.0929e-03
> > > ; Li+ parameters from jensen and jorgensen JCTC 2006 2 1499-1509 ;;
> Mass
> > of
> > > chlorine
> > >
> > >
> > > *LI.itp file content of alchemical transformation for mass is*
> > > [ atoms ]
> > > ; id    at type     res nr  residu name at name  cg nr  charge   mass
> > > 1       opls_406    1       LI          LI       1      0.0
> 6.941
> > > opls_407 0.0 35.45300
> > >
> > >
> > > *.mdp file used is as follows*
> > >
> > >
> > > integrator             = sd
> > > tinit                      = 0
> > > dt                         = 0.001
> > > nsteps                  = 1000000
> > > nstcomm              = 100
> > > nstxout                 = 0
> > > nstvout                 = 0
> > > nstfout                  = 0
> > > nstlog                   = 5000
> > > nstenergy             = 5000
> > > nstxtcout              = 500
> > > cutoff-scheme            = verlet
> > > nstlist                  = 10
> > > ns_type                  = grid
> > > pbc                      = xyz
> > > rlist                    = 1.0
> > > coulombtype              = PME
> > > rcoulomb                 = 1.0
> > > vdwtype                  = cutoff
> > > rvdw                     = 1.0
> > > DispCorr                  = EnerPres
> > > fourierspacing           = 0.12
> > > pme_order                = 6
> > > tc_grps                  = system
> > > tau_t                    = 1.0
> > > ref_t                    = 300
> > > Pcoupl                   = Parrinello-Rahman
> > > tau_p                    = 1.0
> > > compressibility          = 4.5e-05
> > > ref_p                    = 1.0
> > > free_energy              = yes
> > > init_lambda  = 0.05 ;;; depending on state, lambda value to be changed
> > from
> > > 0 to 1 ;;;;;
> > > delta_lambda             = 0
> > > nstdhdl                  = 10
> > > gen_vel                  = no
> > > constraints              = h-bonds
> > > constraint-algorithm     = lincs
> > > continuation             = yes
> > > lincs-order              = 12
> > >
> > >
> > > Free enegy of step 1: i.e. mass transfomation is around +6.2 kJ/mol.
> > >
> > > In Gromacs 5.1.2 manual section 4.5 p.no 99 dependance of free energy
> on
> > > lambda equation is given.
> > >
> > > Now my doubt is if I am using the NPT ensemble then the total kinetic
> > > energy of the system should be constant i.e. summation(p^2/2m) which is
> > > being controlled by the thermostat.
> > >
> > > But i am getting non zero value in this case? Is it normal to get
> > something
> > > like this or i have blundered somewhere?
> > >
> > > Will be happy to share input files if it helps.
> > >
> > > Thanks and Regards
> > >
> > > --
> > > Deepak U. Bapat
> > > --
> > > 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.
> > >
> >
> >
> >
> > --
> > _____________
> >
> > Prof. Dr. André Farias de Moura
> > Department of Chemistry
> > Federal University of São Carlos
> > São Carlos - Brazil
> > phone: +55-16-3351-8090
> > --
> > 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.
>



-- 
Deepak U. Bapat


More information about the gromacs.org_gmx-users mailing list