[gmx-users] Question about documentation of -a switch in g_hbond
Mark Abraham
Mark.Abraham at anu.edu.au
Mon Jun 11 08:31:11 CEST 2012
On 11/06/2012 12:41 AM, Justin A. Lemkul wrote:
>
>
> On 6/10/12 10:36 AM, Andrew DeYoung wrote:
>> Greetings,
>>
>> This might end up being a silly/embarrassing question, and if so, I
>> apologize. I feel like I may be making a conceptual mistake, but I'm
>> not
>> sure.
>>
>> Is it true that a hydrogen bond is of the following form?
>>
>> Donor---Hydrogen ... Acceptor
>>
>> Is this the correct order? I think so. For example, in my system, I
>> have:
>>
>> Oxygen---Hydrogen ... Fluorine
>>
>> In the documentation, the -a option has this description:
>>
>> "Cutoff angle (degrees, Acceptor - Donor - Hydrogen)"
>>
>> Why is the order Acceptor - Donor - Hydrogen (or, equivalently,
>> Hydrogen -
>> Donor - Acceptor)? Why isn't "Hydrogen" instead between "Donor" and
>> "Acceptor"? Which angle does -a specify?
>>
>
> Because it's easy to define from an algorithmic standpoint. A
> hydrogen bond is stronger depending on how close to linearity it is.
> Measuring an acute angle given by -a says that if the A-D-H angle is
> less than 35 degrees, you can vary between 0 and 35 easily to show the
> same thing as being within a certain range of being 180 degrees in
> terms of D-H-A.
The component of H-bond detection based on the cut-off angle is not
implemented as a cut-off in angle space, but rather in the space of the
cosine of the angle. Cosines are easily computed using the inner product
of the vectors that form them. There's no real advantage to computing
cos(A-D-H) vs cos(D-H-A) because in either case you just compare the
value in the appropriate sense with the appropriate value of the
(signed) cos(acutoff). Defining the H-bond in terms of A-D-H has the
advantage that if you later want to consider the angle distribution, the
relevant (signed) angle range is centered around zero degrees. D-H-A
would have a cluster below 180 and a cluster above -180, which is less
convenient to deal with.
For a given r(D-H), you can define the set of legal H-bond positions
equivalently from either A-D-H or D-H-A, of course. Since r(D-H) is
often constant for MD simulations, the conceptual difference between the
definitions becomes moot.
In neither case do you want to compute the angle with acos later,
because the numerical accuracy of acos is poor around 0 and +/-180 from
its near-linearity. Unfortunately g_hbond does do this, so I think I'll
patch it.
Mark
More information about the gromacs.org_gmx-users
mailing list