[gmx-developers] Use of invsqrt with tables

[Berk Hess]

On 2016-03-23 13:42, David van der Spoel wrote:
> On 23/03/16 13:31, Berk Hess wrote:
>> Hi,
>>
>> But avoiding the inversion is not possible, since you need to multiply
>> the scaling force by the normalized distance vector.
> The force should be zero at r=0 for symmetry reasons and hence the 
> code should not crash.
I know that. But that doesn't help in avoiding the division needed to 
get the normalized distance vector. You can avoid the division by 
exactly zero using a conditional. But you will still have issues for 
distances very close to zero. A solution I use for forces between 
excluded atoms in the non-bonded kernel, which can potentially overlap, 
is to add a very small value epsilon to r^2, such that invsqrt produces 
a large number that not close to max real. Multiplying this with r gives 
a value close to zero that's slightly off from 1, but the error will be 
less than real precision for potentials that are symmetric around 0.

Berk
>>
>> Berk
>>
>> On 2016-03-23 13:00, David van der Spoel wrote:
>>> Hi,
>>>
>>> In all versions of gromacs the distance between two nonbonded atoms is
>>> computed as
>>>             rinv00           = gmx_invsqrt(rsq00);
>>>             r00              = rsq00*rinv00;
>>>
>>> This is fine when we use Lennard Jones and Coulomb, but not when using
>>> table potentials which do not necessarily have singularities at r = 0.
>>> I would therefore propose to replace the statements by a normal sqrt
>>> in all innerloops using tables.
>>>
>>> Thoughts?
>>
>
>