[gmx-users] restraints

[David Mobley]

Dear Gromacs users,

I need to apply a number of restraints (angles, distances, and
torsions) between some reference atoms in two molecules. I have been
looking at how to do this in the manual, and find a discussion of
distance and angle restraints (i.e. on page 99 of the Gromacs 3.2
manual). The distance restraints are straightforward; the angle
restraints take four atoms.

Upon further reading, I found on page 64 that the angle is computed by

cos theta=((r_j-r_i)/|r_j-r_i| -dot- (r_l-r_k)/|r_l-r_k|)

where i, j, k, and l are presumably the atoms I specify in the
restraints section.

Here's my question: Does this give the torsional angle one would
expect? It doesn't seem like it does: It's clear that in the case
where my four atoms are ABBC, the angle is the angle between the unit
vectors AB and BC, which is exactly what I want for an ordinary angle.

But for the ABCD case, this doesn't seem to be a torsion: One wants in
that case the angle between the planes ABC and BCD; this gives the
angle between unit vectors AB and CD, which is generally not the same
as the angle between the planes ABC and BCD, unless I'm missing
something obvious.

If this doesn't give a torsional angle, does anyone have other
suggestions on how best to impose torsional restraints between two
molecules in Gromacs? I had thought of trying to simply add another
torsional parameter, but it seems more desirable to treat restraints
as restraints.

Thanks in advance for your help.

David Mobley
Postdoc, University of California, San Francisco.