[gmx-developers] Use of invsqrt with tables

[David van der Spoel]

On 23/03/16 14:07, Erik Lindahl wrote:
> Hi,
>
>
>> On 23 Mar 2016, at 14:00, Berk Hess <hess at kth.se> wrote:
>>>>
>>>> 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.
>
> I have planned infrastructure for the new kernels that can handle much of this with masked (SIMD) operations instead of conditionals, but my priority for the 2016 release is to get rid of group kernels, so this might not make it…
>

The issue is the interaction between two Gaussian-smeared charges. Both 
the energy and the force go to zero smoothly. For distances close to 
zero we have to multiply the scalar force by the normalized distance 
vector. If we do the normalization first that should not give numerical 
problems, right?

A further alternative (tables only) would be to store F/r in the tables 
if that is well behaved (which it is for the Gaussian charges).
> Cheers,
>
> Erik
>


-- 
David van der Spoel, Ph.D., Professor of Biology
Dept. of Cell & Molec. Biol., Uppsala University.
Box 596, 75124 Uppsala, Sweden. Phone:	+46184714205.
spoel at xray.bmc.uu.se    http://folding.bmc.uu.se