random: Fix γ-section implementation for Randomizer::getFloat()

[TimWolla]

The reference implementation of the "Drawing Random Floating-Point Numbers from an Interval" paper contains two mistakes that will result in erroneous values being returned under certain circumstances:

• For large values of g the multiplication of k * g might overflow to infinity.
• The value of ceilint() might exceed 2^53, possibly leading to a rounding error when promoting k to double within the multiplication of k * g.

This commit updates the implementation based on Prof. Goualard suggestions after reaching out to him. It will correctly handle inputs larger than 2^-1020 in absolute values. This limitation will be documented and those inputs possibly be rejected in a follow-up commit depending on performance concerns.

[Girgias]

This looks like black magic to me, but I trust that the implementation is sensible.

[zeriyoshi]

I have fully utilized my inadequate brain to check, and I think it's probably fine.
Although it would be better to obtain a review from a knowledgeable person in mathematics, Tim is probably already knowledgeable.

[TimWolla]

The changes are actually pretty simple: Both issues are caused by k being too large, leading to floating point imprecision, thus the multiplication is effectively performed in two steps with smaller k:

max - k * g

is turned into:

4 * ((max / 4) - floor(k / 4) * g) - (k % 4) * g
= max - 4 * floor(k / 4) * g - (k % 4) * g
= max - (4 * floor(k / 4) + (k % 4)) * g
= max - k * g

This transformation is correct if max is at least 2^-1020, any smaller max / 4 rounds down to 0 and thus the reverse direction of multiplying with 4 doesn't give back the original value.