[gmx-users] non zero free energy values in mass transformation
deepak bapat
dubapat at gmail.com
Mon Aug 14 05:09:06 CEST 2017
Hi
To study the efect of particle mass on excess chemical potential i did test
particle insertion (tpi) in SPC/E water. I created a lennard jones
particle with very small LJ parameters so that the probability of gettining
it inside the cavity increases.
In one case I kept mass of tpi particle equal to that of Lithium. In second
simulation the mass was equal to that of Chlorine atome.
LJ parameters for tpi particles in both cases were 2.0e-02 7.64793e-02.
In both cases the excess chemical potential was in the range of 0.0798
kJ/mol.
tpi only considers energy as a function of coordinates only so i assume
that it would depend only on LJ parameters as long as we are delaing with
neutral systems.
So looking at this I again come back to the mass perturbation which started
all this.
Since I was doing NPT calculations at 300 K I expected that average total
kinetic energy of the system should not be different. The themostat would
just adjust velocity distribution so that it brings the temperature of the
system to the set temperature. Any increase or decrease in mass should be
accounted for by the thermostat.
As I said, in free energy calculations I got non zero value. (~ +6.2
kJ/mol). So I am still worried about this non zero behavior.
Thanks and regards
Deepak
On Sat, Aug 12, 2017 at 10:11 AM, deepak bapat <dubapat at gmail.com> wrote:
> 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
>
>
--
Deepak U. Bapat
More information about the gromacs.org_gmx-users
mailing list