[RNG-154] Avoid infinite samples in the normalised Gaussian samplers

• Type: Improvement
• Status: Closed
• Priority: Minor
• Resolution: Implemented
• Affects Version/s: 1.3
• Fix Version/s: 1.4
• Component/s: sampling
• Assignee: Unassigned
• Reporter: Alex Herbert
• Created: 08/Jul/21 13:47
• Updated: 13/Sep/21 10:41
• Resolved: 09/Jul/21 09:45

Description

The ZigguratNormalisedGaussianSampler uses the following method to sample the tail of the Gaussian:

// Note
// R = 3.442619855899
// 1 / R = 0.2905...
double y;
double x;
do {
    y = -Math.log(rng.nextDouble());
    x = -Math.log(rng.nextDouble()) * ONE_OVER_R;
} while (y + y < x * x);

final double out = R + x;
return hz > 0 ? out : -out;

In the unlikely event the RNG produces 0 then Math.log(0) is infinity. Two zeros in a row can result in an infinite sample for the tail. This is very unlikely but would be unexpected for a user of the library since a sample should be roughly within +/-3.

Note that if the RNG is sampling from the 2^53 dyadic rationals in [0, 1) then the next value is:

Math.log(0x1.0p-53) == -36.7368005696771

The largest value x where 2y < x^2 is false is sqrt(2*36.74) = 8.571 and the returned tail would be +/- 12.01. This is very far from the extreme of infinity.

To avoid infinity this can be fixed by:
1. Assuming the RNG is returning a value in [0, 1) and using Math.log(1.0 - rng.nextDouble())
2. Generating the double u from a long to ensure the value is in [0, 1) and using 1.0 - u.

 

Issue Links

• is related to: RNG-148 UnitSphereSampler does not check for infinite length dimension