[gmx-users] Free energy of transfer with accuracy

David Mobley dmobley at gmail.com
Wed Dec 17 19:19:09 CET 2008 

Eudes,

> But first I would want to know two things.
> 1) The parameter sc-sigma of 0.5 is reasonable?
> How this parameter affects the shape of dvdl plot?
> You performed tests about it? I found little about
> it in the forum.

Check out the manual on what this parameter does. IIRC it only applies
to atoms which have sigma parameters of zero to begin with (i.e.
hydrogens) and so it really doesn't do much.

> 2) In my tests, I realized that, simply changing the electrostatic
> interactions
> treatment (cut-off or PME), a large difference it was found at hydration
> free energy.

Of course. Long-range interactions make big differences in this sort
of thing. You should be using lattice-sum (i.e. PME) electrostatics if
you care about accuracy.

> With cutoff I have found 12.2kJ/mol while with PME, 24.0kJ/mol.
> The accepted value must be around 20kJ/mol).
> This to me was somewhat surprising, because, as mentioned above,
> C60 does not have charge, so that effect should be minimized or zero.
> In any case I assign it to the electrostatic of the solvent.
> I suppose the PME result is more reliable, it's ok?
> Do you know some work that explore the effect of the electrostatic of the
> solvent
> on the solvation free energy?

Treatment of long range electrostatics effects the solvent generally,
including the effective hydrogen bond strength and structure, solvent
density, etc. So it's not that surprising it should affect the
hydration free energy. Also, if you think about it, a hypothetical
water molecule at the center of a C60 will have very different
interactions with its environment depend on whether or not it sees
water far away. I would think it would cost more to create a cavity
with PME than without PME because of this effect, which I think is
consistent with the numbers you're giving (though you didn't specify
your sign convention).

The other thing that actually could affect this significantly is
treatment of long-range dispersion interactions. I have a separate
publication with Michael SHirts on that if you are interested (see
http://www.dillgroup.ucsf.edu/~dmobley).

In terms of effects of electrostatics on solvation free energy, I am
not sure whether anyone has studied this explicitly, but certainly the
effects of electrostatics on water structure generally are well known
and studied.

David


> Thank you
> eef
> _______________________________________
> Eudes Eterno Fileti
> Centro de Ciências Naturais e Humanas
> Universidade Federal do ABC
> Rua Santa Adélia, 166 - Bloco B, Sala 1048
> 09210-170  Santo André - SP Brasil
> +55.11.4437-8408
> skype: eefileti
> http://cromo.ufabc.edu.br/~fileti/
>>
>> Subject: Re: [gmx-users] Free energy of transfer with accuracy
>> To: "Discussion list for GROMACS users" <gmx-users at gromacs.org>
>> Message-ID:
>>        <bc2c99750812110724tc388a77p1577e9c999a6dd2c at mail.gmail.com>
>> Content-Type: text/plain; charset=ISO-8859-1
>>
>> Eudes,
>>
>> On Thu, Dec 11, 2008 at 8:27 AM, Eudes Fileti <fileti at ufabc.edu.br> wrote:
>> > Hello David, thanks for the reply!
>> > I am simulating the free energy of transfer of the C60, from ethanol to
>> > benzene.
>> > The softcore parameters, I used the are the same values of Tieleman
>> > et.al.
>> > (nature nanotechnology, 3, 363).
>> > At figure, each curve refers to a different solvent; one of them for
>> > ethanol
>> > and another for benzene.
>> > For lambda=0 I have the Coulomb and LJ parameters turn on and for
>> > lambda=1,
>> > turn off.
>> > I have found a large difference when simulating my system using
>> > sc_power=2.
>> > See in this link (
>> > http://cromo.ufabc.edu.br/~fileti/web/dgdl-scpower.jpg )
>> > the comparison between
>> > plots from sc_power=1 and sc_power=2 for solute in benzene.
>> > Near 0 and 1, I believe that problem cancels out at difference, as can
>> > seen
>> > at
>> > ( http://cromo.ufabc.edu.br/~fileti/web/dgdl-transfer.jpg ) . Then the
>> > main
>> > problem, in my view,
>> > is in the intermediate region, where the plot presents a deep (-500
>> > units)
>> > and noise minimum.
>> > Observe that in lambda=0, sc_power=2 reduce largely the derivative, but
>> > the
>> > critical region,
>> > one that makes the difference in calculation of the transfer free
>> > energy,
>> > remains the same,
>> > ie, the derivative still oscilates a lot and is not smooth.
>> > Are you believe that changing the value of alpha from 0.47 to 0.50 can
>> > resolve this?
>>
>> OK, there are several issues here.
>> 1) The "noisy" minimum in dV/dlambda
>> 2) The overall shape of dV/dlambda
>> 3) Other potential problems
>>
>> On (1), the noisiness of the minimum you refer to really means there
>> are convergence problems for some reason -- probably correlation times
>> are particularly long there. This may be helped by using the soft core
>> parameters I mention, but I am not sure. See item #3.
>>
>> On (2), in general, the overall shape of dV/dlambda seems to be
>> smoothest and best behaved with the particular soft core parameters I
>> mention. We did a lot of testing on this (i.e., many values of alpha
>> at each of several sc-power values).
>>
>> On (3), there are other thing that can go wrong when using soft core.
>> See http://www.alchemistry.org/wiki/index.php/Best_Practices for some
>> discussion, and see rule #2 in particular. If you have any partial
>> charges on atoms in your C60 this could be a potential source of error
>> and noise. Also, if your C60 is hollow, I can imagine that it might
>> exclude solvent atoms from the inside, in which case you are going to
>> have a big convergence issue as solvent begins to interpenetrate with
>> the C60, which may lead to "noise" problems when this starts to
>> happen, whatever soft core functional form you use. This just means
>> you will need to run longer simulations in this regime.
>>
>> Anyway, long story short, I suggest two things: (a) Make sure you turn
>> off partial charges on C60 (without using soft core) before turning
>> off the LJ using soft core, and (b) use the soft core functional form
>> I suggest.
>>
>> David
>>
>> > Please, ask me more details if you need.
>> > Thanks again.
>> > eef
>> >
>> > _______________________________________
>> > Eudes Eterno Fileti
>> > Centro de Ciências Naturais e Humanas
>> > Universidade Federal do ABC
>> > Rua Santa Adélia, 166 - Bloco B, Sala 1048
>> > 09210-170  Santo André - SP Brasil
>> > +55.11.4437-8408
>> > skype: eefileti
>> > http://cromo.ufabc.edu.br/~fileti/
>> >>
>> >>
>> >> ------------------------------
>> >>
>> >> Message: 4
>> >> Date: Thu, 11 Dec 2008 07:24:26 -0600
>> >> From: "David Mobley" <dmobley at gmail.com>
>> >> Subject: Re: [gmx-users] Free energy of transfer with accuracy
>> >> To: "Discussion list for GROMACS users" <gmx-users at gromacs.org>
>> >> Message-ID:
>> >>        <bc2c99750812110524r1aee080av4ac07eed79b4b0d2 at mail.gmail.com>
>> >> Content-Type: text/plain; charset=ISO-8859-1
>> >>
>> >> Hi,
>> >>
>> >> > Hi gmx-users,
>> >> > I have tried to calculate the free energy of transfer of
>> >> > a solute between two solvents by thermodynamic integration.
>> >> > However, the dgdl plot I have obtained is not smooth enough to ensure
>> >> > accuracy.
>> >> > I have used a non homogeneous spacing for lambda. The first is from
>> >> > 0.0
>> >> > to
>> >> > 0.4 (dl=0.04),
>> >> > the second is from 0.405 to 0.48 (dl=0.005) and the third is from
>> >> > 0.48
>> >> > to 1
>> >> > (dl=0.04).
>> >> > Overall, I was performed 40 simulations.
>> >> > As can be seen in the figure (
>> >> > http://cromo.ufabc.edu.br/~fileti/web/dgdl-transfer.jpg ),
>> >> > there is a large fluctuation in the region that goes from 0.40 to
>> >> > 0.48.
>> >> > Each simulation was carried out by 0.5ns preceded by 0.15ns of
>> >> > equilibration.
>> >> > I am using soft-core whose parameters are:
>> >> > ; Free energy control stuff
>> >> > free_energy              = yes
>> >> > init_lambda              = 0.00
>> >> > delta_lambda             = 0
>> >> > sc_alpha                 = 1.3
>> >> > sc-power                 = 1.0 (for this I also used 2).
>> >> > sc-sigma                 = 0.47
>> >> > In the example of the link sc-power=1, but I re-did the calculations
>> >> > with
>> >> > sc-power=1
>> >> > and I found the same behavior (although the format of the plot was a
>> >> > little
>> >> > different).
>> >>
>> >> > Can you suggest me a way to calculate this free energy difference
>> >> > with
>> >> > accuracy?
>> >> > The perturbation free energy method could give better results?
>> >> > Thank you
>> >> > EEF
>> >>
>> >> This does look alarming.
>> >>
>> >> I would make a couple comments here:
>> >> 1) Soft core parameters can be fairly sensitive. Michael Shirts and I
>> >> have both looked at these independently and we find that sc-power =
>> >> 1.0 with an alpha of 0.5 is typically substantially better than most
>> >> other options. Even alpha = 0.47 or 0.53 can be much worse than alpha
>> >> = 0.5. Your dV/dlambda curve looks to me like your soft core values
>> >> are far from ideal and that you are seeing a huge peak near lambda = 0
>> >> for this reason.
>> >> 2) Am I interpreting your e-mail correctly that lambda = 0 corresponds
>> >> to one solvent, and lambda = 1 corresponds to another? If so, I am
>> >> unclear on how exactly you're doing the transformation and could use
>> >> more detail on how you've set this up. A more straightforward (and
>> >> easier?) way to do this would be to transfer the solute to vacuum from
>> >> each solvent and then take the difference in the vacuum transfer free
>> >> energies. Please clarify which you're doing so we can help better.
>> >>
>> >>
>> >> David Mobley, Ph.D.
>> >> Assistant Professor of Chemistry
>> >> University of New Orleans
>> >> New Orleans, LA 70148
>> >> dlmobley at uno.edu
>> >> Office 504-280-6445
>> >> Fax 504-280-6860
>> >
>
> _______________________________________________
> gmx-users mailing list    gmx-users at gromacs.org
> http://www.gromacs.org/mailman/listinfo/gmx-users
> Please search the archive at http://www.gromacs.org/search before posting!
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-users-request at gromacs.org.
> Can't post? Read http://www.gromacs.org/mailing_lists/users.php
>

More information about the gromacs.org_gmx-users mailing list