[gmx-users] potential energy of one molecule using GAFF: confused regarding Ryckaert-Bellemans and 1-4 potentials

Vedat Durmaz durmaz at zib.de
Thu Oct 20 15:52:58 CEST 2011

Am 20.10.2011 15:01, schrieb Mark Abraham:
> On 20/10/2011 7:25 PM, Vedat Durmaz wrote:
>> thanks mark. i guess it took some time answering all these questions. 
>> and i think you are right. trying to interpret each computed energy 
>> term seperately in a physical manner is senseless. especially since 
>> i'm not really deep inside the force field & implementation stuff. 
>> however, with one of your statements i cannot agree at all:
>> imho, free energy diff's on the basis of two conformations and their 
>> steady state distribution are independent from both the energy and 
>> type of the intermediate state. it's only the conversion's rate that 
>> does depend on it.
> I'm not sure with which statement of mine you are disagreeing. 

actually with this one: "Measuring the free energy difference with 
simulations is hard because you cannot model the intermediate stages of 
bond breaking and forming."

> The free energy is a state function, and so the free energy difference 
> between E and Z isomers does not depend on the path taken for the 
> conversion. The difficulty is just that there is no path for a cyclic 
> alkane from E to Z in a normal MD force field.
and that's the reason for the disagreement. also imho, the path is 
irrelevant. for a molecule given with a cis and a trans configuration, 
respectively (same atoms, quasi-same topology), the inner energy and 
entropy of both conformations (derived from classical force field 
simulations) is enough in order to get the free energy difference.
>> but the reason, why i refrained from extracting the potential energy 
>> for a subset of my system was simply due to the difficulties to find 
>> "detailed information" about how to play it out although many gromacs 
>> users have been asking for the issue over the last years. i needed 
>> nearly one full day to figure it out even though it's a simple series 
>> of about 4 gromacs commands only, given some md input & result files 
>> like md.xtc, md.gro, complex.top, index.ndx:
>> |grompp -f mdSubset.mdp -c md.gro -p complex.top -o mdSubsetTemp.tpr|
>> |echo ||"1"| || tpbconv -s mdSubsetTemp.tpr -nsteps ||0| |-n 
>> index.ndx -o mdSubset.tpr|
>> |echo ||"1"| || trjconv -s mdSubset.tpr -f md.xtc -o mdSubset.xtc|
>> |mdrun -s mdSubset.tpr -rerun mdSubset.xtc -v -deffnm mdSubset|
>> where "1" stands for the group to be extracted (from the list of 
>> groups in the index file) and "-nsteps 0" causes recomputation of for 
>> the extracted subset at the given time steps only. i mention it here 
>> so that other users looking for the details might find them a little 
>> bit faster ... it's especially this kind of "tiny tutorials" that i 
>> miss e. g. on the gromacs website.
> There's an art to trying to find detailed help. I know I've given the 
> advice to make subsets and use mdrun -rerun several times, but I don't 
> know how a newcomer is ever going to figure that out on their own!
Your somehow right, but unfortunately, i'm no artist at all. However, 
what's the reason for the availability of all the tutorials and how-tos 
on the gromacs documentation site? i'll tell you: it helps other 
researchers to reduce times of searching/finding/figuring out. and 
besides, the information retrieved from the tools' manpage or the manual 
is often far away from being helpful!

please, don't misunderstand me. i'm just saying: if there are how-tos 
about how to plot data, perform calculations of pKa values or diffusion 
constants, then the description of getting the potEnergy of one of the 
systems components has an a fortiori right to exist there or somewhere 
else at a central place within the gromacs world!

> BTW -nsteps 0 is not necessary. mdrun -rerun only computes on the 
> frames present in the input trajectory file. Also, if the index group 
> of interest is already present in the .tpr that ran the simulation, 
> then you can do the job with variations on only the latter three commands.

thanks for the hints!

> Mark


ps: don't kill me now ...
>> thanks for listening and kind regards,
>> vedat
>>>> Q8 does anyone have an idea, how to perform the simulation and on 
>>>> which energy terms to concentrate in order to get reliable results?
>>> You seem to be performing it OK, given that you've said very little 
>>> about any details...
>>> I think the problem is poorly constructed. You have some 
>>> experimental data that gives a general understanding of the size of 
>>> the free energy difference between the isomers. You can't 
>>> necessarily expect to reproduce that from (average) potential energy 
>>> differences between conformations of those isomers. Measuring the 
>>> free energy difference with simulations is hard because you cannot 
>>> model the intermediate stages of bond breaking and forming. There 
>>> are "alchemical" free energy methods that could in principle treat 
>>> this problem effectively, but there will be some significant issues 
>>> and you are best doing your own homework there.
>>> Mark
-------------- next part --------------
An HTML attachment was scrubbed...
URL: <http://maillist.sys.kth.se/pipermail/gromacs.org_gmx-users/attachments/20111020/1f4b22fc/attachment.html>

More information about the gromacs.org_gmx-users mailing list