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