[
https://issues.apache.org/jira/browse/RNG-154?page=com.atlassian.jira.plugin.system.issuetabpanels:comment-tabpanel&focusedCommentId=17377416#comment-17377416
]
Alex Herbert commented on RNG-154:
----------------------------------
h2. Other normalised Gaussian samplers in the library
Infinite samples can also be returned by the BoxMullerNormalizedGaussianSampler
due to the use of Math.log(rng.nextDouble()).
The MarsagliaNormalizedGaussianSampler calls Math.log(r2) on a squared radius
of a sample (x,y) inside the unit circle. This only occurs when r2 is above
zero, thus the sampler eliminates the infinite extreme.
The ZigguratSampler.NormalisedGaussian does not use Math.log. It generates the
tail via an exponential distribution that uses recursive calls that add the
term 7.56 to repeat samples. It will create a stack overflow before ever
reaching infinity.
The Box-Muller method could also be fixed to avoid infinite samples. In this
case I would favour not doing so as the sampler is a reference implementation
and unlikely to be used in practice.
> Avoid infinite samples in the ZigguratNormalisedGaussianSampler
> ---------------------------------------------------------------
>
> Key: RNG-154
> URL: https://issues.apache.org/jira/browse/RNG-154
> Project: Commons RNG
> Issue Type: Improvement
> Components: sampling
> Affects Versions: 1.3
> Reporter: Alex Herbert
> Priority: Minor
>
> The ZigguratNormalisedGaussianSampler uses the following method to sample the
> tail of the Gaussian:
> {code:java}
> // 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;
> {code}
> 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:
> {code:java}
> Math.log(0x1.0p-53) == -36.7368005696771
> {code}
> Here the returned tail would be 3.44 + 36.7368005696771 / 3.44 = 14.11. 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.
>
--
This message was sent by Atlassian Jira
(v8.3.4#803005)
