Undefined Behaviour in volk_32f_invsqrt_32f

[argilo]

UBSAN shows the following Undefined Behaviour in volk_32f_invsqrt_32f:

/home/argilo/git/volk/kernels/volk/volk_32f_invsqrt_32f.h:71:22: runtime error: signed integer overflow: 1597463007 - -569061536 cannot be represented in type 'int'

This is the problematic line:

volk/kernels/volk/volk_32f_invsqrt_32f.h

Line 71 in e853e9b

u.i = 0x5f3759df - (u.i >> 1); // what the fuck?

[jdemel]

This is: Wikipedia fast inverse square root
It is sometimes called John Carmack reverse because of its use in Q3A.

Fast inverse square root, sometimes referred to as Fast InvSqrt() or by the hexadecimal constant 0x5F3759DF

It is an approximation. It works. It is a hack.

Since this is actually a floating point operation, I suggest to ignore this specific error.

[argilo]

I suspect the signed integer overflow occurs because the input is negative. So this could just be a case where the test input (uniformly distributed floats in the range -1 .. +1) doesn't make sense.

[jdemel]

As long as the output signature is 32f, roots of negative values do not exist.

The C reference for sqrt says:

1. +-0 should return the input value
2. Values smaller +-0 are NaN. Since NaN == NaN is false, such input should make the test fail.

[michael-roe]

According to
https://en.wikipedia.org/wiki/Fast_inverse_square_root

.. the way to avoid undefined behaviour is to use std::bit_cast

Modern compilers are much more likely than older ones to generate crazy code for undefined behavior; unclear if this is actually a problem here.

[jdemel]

According to https://en.wikipedia.org/wiki/Fast_inverse_square_root

.. the way to avoid undefined behaviour is to use std::bit_cast

Modern compilers are much more likely than older ones to generate crazy code for undefined behavior; unclear if this is actually a problem here.

Unfortunately, the kernels as well as the library itself are written in C, unlike all the test code. Thus, we can't use std::bit_cast. There seem to be hacks according to SO. If available, we probably want to use CPU FPUs for the task anyways, or their vector equivalent.

This function is probably a good candidate to implement a test with the new googletest based CI, where we implement a naive invsqrt and compare against this result. As long as the result is approximately equal, we know the compiler didn't start to mess with our hack.