Hello, I am re-running some of our gromacs 3 simulations using gromacs
4, and as far as I can tell the gromacs 4 pull code, while very nicely
enhanced from gromacs 3, has also lost some functionality.
I am calculating the PMF of a peptide across a bilayer and, to
simplify the issue, what I can't figure out how to do with gromacs 4
is to pull group1 to a position 1 nm more negative along the z-axis
than group 0 (the bilayer in this case) for some replicas and to 1 nm
more positive along the z-axis than group 0 in others. I'll focus on
the negative displacement here as it is the one that is giving me
problems.
This used to be possible with the following gromacs 3 pull.ppa options:
runtype = umbrella
reftype = com
pulldim = N N Y
reference_group = group0
group_1 = group1
K1 = 500
Pos1 = 0 0 -1.0
Sure, I could start group1 in the correct negative-z displacement from
group0 and use pull_init1=+1.0, but this will not work when, for
instance, I want to restrain it to -0.1 nm (using pull_init1=+0.1),
where the sampling will infrequently jump back and forth about z=0.
Just to be sure, I tried for gromacs 4 are the following pull code
.mdp options:
pull = umbrella
pull_geometry = distance
pull_dim = N N Y
pull_start = no
pull_ngroups = 1
pull_group0 = group0
pull_group1 = group1
pull_init1 = -1.0
pull_rate1 = 0
pull_k1 = 500.0
where mdrun complains:
"Pull reference distance for group 1 is negative (-1.000000)"
and it is pretty obvious why this doesn't work since it is asking for
a negative displacement. Nevertheless, I tried it and pull_init1
appears to get set to zero.
I also attempted the following.
pull = umbrella
pull_geometry = distance
pull_dim = N N Y
pull_start = no
pull_ngroups = 1
pull_group0 = group0
pull_group1 = group1
pull_init1 = 1.0
pull_rate1 = 0
pull_k1 = 500.0
pull_vec1 = 0 0 -1
where I would then use
pull_init1 = 1.0
pull_vec1 = 0 0 1
for the positive side of the bilayer.
However, when I look at the forces, I am getting the negative of what
I should get when pull_vec1 = 0 0 -1.
coord.xvg:
610.0000 4.9343 -1.02019
660.0001 4.91454 -0.949747
force.xvg:
610.0000 -10.0932
660.0001 25.1265
Although I am getting exactly what I should get when pull_vec1 = 0 0 1
(for intended positive displacements).
coord.xvg:
660.0001 4.9014 1.16304
force.xvg:
660.0001 -81.5201
Any ideas are greatly appreciated. I can probably mod the code for my
needs, but a standard gromacs binary is always preferable.
Thank you,
Chris.
--
gmx-users mailing list gmx-users@gromacs.org
http://lists.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-requ...@gromacs.org.
Can't post? Read http://www.gromacs.org/mailing_lists/users.php
