[gmx-users] Rugged dG/dlambda when turning off charges

[Maik Goette]

OK, forget my first reply and the mentioned topology splitting ;)
It MAY be a problem, that you switch off two charged molecules, but I 
really dont know, if this is the problem.
Try some totally artificial process, where you build a neutral base and 
use no counterion and morph the neutral base away.
 From my experience with DNA I can tell you, that these huge 
fluctuations are quite intrinsic due to the size of the perturbation. 
Actually, I don't have an idea, how to get rid of them.
Trying to sample much longer could help to converge the system. 
Furthermore, using more lambda-steps could help you in getting a better 
total error and stuff.
Again, from my experience I can tell you, that calculating free energies 
from letting nucleotides appear/disappear is really no fun.
If you come up with a solution yourself, I would be highly interested in it.

Regards

Maik Goette, Dipl. Biol.
Max Planck Institute for Biophysical Chemistry
Theoretical & computational biophysics department
Am Fassberg 11
37077 Goettingen
Germany
Tel.  : ++49 551 201 2310
Fax   : ++49 551 201 2302
Email : mgoette[at]mpi-bpc.mpg.de
         mgoette2[at]gwdg.de
WWW   : http://www.mpibpc.gwdg.de/groups/grubmueller/


Robert Johnson wrote:
> Thanks David,
> Since I'm turning off the VdW parameters as a second step, there
> really isn't any reason to use the soft core potentials when turning
> off the charges right? Thus, sc_alpha should be set to zero when
> turning off charges?
> Thanks,
> Bob
> 
> On Jan 25, 2008 5:38 PM, David Mobley <dmobley at gmail.com> wrote:
>> It will only be smooth if you don't use soft core (i.e. sc-alpha =0 ).
>>
>> Best,
>> David Mobley
>> http://www.dillgroup.ucsf.edu/~dmobley
>>
>>
>> On Jan 25, 2008 12:40 PM, Robert Johnson <bobjohnson1981 at gmail.com> wrote:
>>> Hello everyone,
>>> As I mentioned in a previous post, I'm trying to compute the free
>>> energy of binding of a DNA base on a carbon nanotube. My system
>>> consists of a single DNA nucleotide (base, sugar, phosphate group) on
>>> a rigid carbon nanotube in aqueous solution. A single Na+ counterion
>>> is also included for charge neutralization. I'm proceeding with this
>>> by computing the free energy associated with the following
>>> transformations:
>>>
>>> Water + Nanotube + DNA + Na+ -> Water + Nanotube + Nothing
>>> Water + DNA + Na+ -> Water + Nothing
>>>
>>> For each of these transformations I FIRST turn off the charges on the
>>> base and counterion over a set of 5 lambda values (0.0, 0.25, 0.5,
>>> 0.75, 1.0). From literature and other posts on the GMX mailing list,
>>> it seems like this small set of values should be adequate because
>>> dG/dlambda tends to be pretty smooth for discharging the molecule.
>>>
>>> Unfortunately, that doesn't seem to be the case in my system. Here is
>>> a table of the dG/dlambda values I get from running 5 ns trajectories
>>> of my system:
>>>
>>> lambda     <dG/dlambda>
>>> 0.0           6211.105
>>> 0.25         1055.254
>>> 0.5           1230.675
>>> 0.75         1128.359
>>> 1.0           756.2904
>>>
>>> If you plot this you will see a large drop from lambda=0 to 0.25. Then
>>> dG/dlambda increase slightly from lambda=0.25 to 0.5 and then smoothly
>>> decreases after that. This seems like a pretty strange result. Also,
>>> it seems that something weird is going on around lambda=0. If I
>>> compute dG/dlambda for lambda=0.05, I get something around 1200-1300,
>>> which is again a very sharp decrease from the value of 6211.105 I get
>>> at lambda=0. Thus, it doesn't seem like simply adding more lambda
>>> values is going to help me.
>>>
>>> Here are the values I'm using for the free energy calculation:
>>> free_energy         =  yes
>>> init_lambda         =  0.05
>>> delta_lambda        =  0
>>> sc_alpha            =  0.5
>>> sc_power            =  1
>>> sc_sigma            =  0.3
>>>
>>> Also, here is an excerpt from my topology file:
>>> [ atoms ]
>>> 1       P       1       DG      P       1       1.1659  30.9700
>>>  P       0.0000  30.9700
>>> 2       O       1       DG      O1P     1       -0.7761 16.0000
>>>  O       0.0000  16.0000
>>> 3       O       1       DG      O2P     1       -0.7761 16.0000
>>>  O       0.0000  16.0000
>>>
>>> I'm not perturbing any of the bonding interactions or VdW parameters -
>>> I'm just turning off all charges on the molecule.
>>>
>>> Does anyone have any ideas or suggestions about how this can be improved?
>>>
>>> Thanks,
>>> Bob Johnson
>>> _______________________________________________
>>> 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
>>>
>> _______________________________________________
>> 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
>>
> _______________________________________________
> 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
> 
> .
>