[gmx-users] non zero free energy values in mass transformation
Dan Gil
dan.gil9973 at gmail.com
Fri Aug 11 16:05:55 CEST 2017
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.
>
More information about the gromacs.org_gmx-users
mailing list