AMBER Archive (2005)

Subject: Re: AMBER: pmemd - bug report

From: David A. Case (case_at_scripps.edu)
Date: Tue Jul 05 2005 - 12:38:09 CDT

On Tue, Jul 05, 2005, Petr Kulhanek wrote:

>
> When I debug program I just found that some fractional coordinate is
> 1.0. I thought that round errors are responsible for this. OK, I
> analysed the problem in more details and I think that the calculation of
> fractions are not fully correct.
>
> I will assume only one-dimensional case in following consideration.
>
> * fractional coordinate has to be in range from 0.00 to 0.9999..
>
> * dnint(0.5) is 1.0 ( dnint(0.5) = dint(0.5 + 0.5) = dint(1.0) = 1.0 )
> * dnint(0.0) is 0.0 ( dnint(0.0) = dint(0.0 - 0.5) = dint(-0.5) = 0.0 )
> * dnint(-0.5) is -1.0 ( dint(-0.5) = dint(-0.5 - 0.5) = dint(-1.0) = -1.0 )
> (I checked result for dnint(-0.5) also by g77)
>
> * fractional coordinate is calculated by fract = f1 - dnint(f1) + 0.5d0
>
> * if f1 = 0.5
> fract = f1 - dnint(f1) + 0.5d0 = 0.5 - dnint(0.5) + 0.5 = 0.5 - 1.0 +
> 0.5 = 0.0
> this is correct
>
> * if f1 = 0.0
> fract = f1 - dnint(f1) + 0.5d0 = 0.0 - dnint(0.0) + 0.5 = 0.0 - 0.0 +
> 0.5 = 0.5
> this is correct
>
> * if f1 = -0.5
> fract = f1 - dnint(f1) + 0.5d0 = -0.5 - dnint(-0.5) + 0.5 = -0.5 -
> (-1.0) + 0.5 = 1.0
> this is NOT correct, fract has to be 0.0

The only problem I see here is a possible round-off problem: if the the
initial coordinate is exactly -0.5, then the fractional coordinate will be
exactly 1.0, and it may be that this is not allowed by some data structure.
But otherwise, I don't see anything wrong. [Note that a fractional coordiante
of 0.0 is logically the same as a fractional coordinate of 1.0.]

>
> It is my guess that fract could be calculated correctly only with
> floor() or ceiling() functions...

I don't see any easy way to use floor() or ceiling(); the "nearest integer"
(nint) is indeed the function that is required. Clearly, if a fractional
coordinate of exactly -0.5 causes problems, this needs to be trapped in the
code.

...regards...dac

-----------------------------------------------------------------------
The AMBER Mail Reflector
To post, send mail to amber_at_scripps.edu
To unsubscribe, send "unsubscribe amber" to majordomo_at_scripps.edu