[gmx-users] Free energy and and-to-end distances
wgallin at ualberta.ca
Thu May 13 23:33:16 CEST 2010
I am trying to calculate the relative free energy differences of conformations of a peptide as a function of end-to-end distance (linear distance between N-terminal nitrogen atom to C-terminal carbon atom).
According to the manual pull code should not be used because the two atoms being separated are joined by the constraints of the bonds (Section 6.3, Limitations). However, section 6.4 mentions two ways to get around this limitation, but I have not been able to get a clear sense of how to implement these ideas.
I am looking at a 15 amino acid peptide, I have run a simple simulation, and it appears that the closest the ends get are 0.5 nm and the farthest apart that they get is 3.5 nm. So I thought that the approach would be to collect a set of 11 conformations from the original simulation with the ends at different distances apart (in steps of about 0.3 nm), add a type 2 constraint to the topology of each with the A state being 0.5 and the B state being 3.5, in the run.mdp file for each run set:
free_energy = yes
init_lambda = <varied between 0 and 1 in 0.1 increments>
delta_lambda = 0
and run a 20 ns simulation on each of the 11 systems.
This kind of run, I gather will yield a dgdl.xvg file for each run.
However, in section 6.4 of the manual it says:
When the position of the minimum
or the constraint length is 1 nm more in state B than in state A, the restraint or constraint force is
given by dH/dlambda. The distance between the atoms can be changed as a function of and time by
setting delta-lambda in the .mdp file. The results should be identical (although not numerically
due to the different implementations) to the results of the pull code with umbrella sampling
and constraint pulling. Unlike the pull code, the free energy code can also handle atoms that are
connected by constraints.
But I was going to set my constraints at 0.5 and 3.5, so the difference is not 1 nm, so my numbers will not be useful.
I have been around in circles with this a number of times, and I was wondering if anyone could point me at an example of doing this kind of analysis, or else perhaps clarify how to generate the equivalent of pullx and pullf files for g_wham analysis without using the pull code.
All advice gratefully accepted.
More information about the gromacs.org_gmx-users