[gmx-developers] Shake question?

[hess at sbc.su.se]

> On 6/21/10 8:57 AM, Shirts, Michael (mrs5pt) wrote:
>> Hi, all-
>>
>> A couple questions about the shake implementation.
>>
>> Right now, it appears that shake tolerance (set as shake_tol) is the
>> tolerance on the relative distance from the constraint squared -- from
>> the
>> manual, one might guess that shake_tol is the maximum relative deviation
>> of
>> a constrained bond.  I think expressing an the parameter in terms of
>> distances instead of distances**2 would make more sense.
>>
>> Also, currently in the cshake code, the constraint is considered
>> satisfied
>> if iconv = 0, where:
>>
>> iconv   = fabs(diff)*tt[ll];
>>
>> Diff is the difference between the square bond length and the square
>> desired
>> bond length, and tt[ll] currently is 1/(2*tol) -- rather than the
>> tol*tol
>> which I think would be a bit more clear where tol is equal to shake_tol.
>>
>> If tt is too big, then its possible that there's a 1 in about 4 billion
>> chance that iconv will be zero even if the tolerance isn't satisfied.
>> Actually, this appears to be the case even with the current code - if
>> ir->shake_tol is less than 10^-13, then you'd start seeing the looping
>> around the integers in a long int.
>>
>> I'd propose changing the definition of 1/(2*tol) to 1/(tol*tol), and to
>> avoid the problem of integer wrapping by capping shake_tol internally in
>> readir.c , so that  the hardware limit for "long int" isn't violated.
>>
>> I can go ahead and make these changes, but I wanted to ask around and
>> verify
>> there wasn't something I was missing!
>>
>
> The routine that is being used is vec_shake which says:
>
>      tol2=2.0*tol;
>      tt[ll] = 1.0/(toler*tol2);
>
>

Or to clarify things a bit more:
the codes uses (1+a)^2 = 1+2a+a^2 approx 1+2a for small a.
I have been staring at this code a few years ago,
but it really does what the mdp option says.

Berk

>
>> Best,
>> ~~~~~~~~~~~~
>> Michael Shirts
>> Assistant Professor
>> Department of Chemical Engineering
>> University of Virginia
>> michael.shirts at virginia.edu
>> (434)-243-1821
>>
>
>
> --
> David van der Spoel, Ph.D., Professor of Biology
> Dept. of Cell & Molec. Biol., Uppsala University.
> Box 596, 75124 Uppsala, Sweden. Phone:	+46184714205. Fax: +4618511755.
> spoel at xray.bmc.uu.se	spoel at gromacs.org   http://folding.bmc.uu.se
> --
> gmx-developers mailing list
> gmx-developers at gromacs.org
> http://lists.gromacs.org/mailman/listinfo/gmx-developers
> Please don't post (un)subscribe requests to the list. Use the
> www interface or send it to gmx-developers-request at gromacs.org.
>