RNG-50: PoissonSampler speed improvements Fixed checkstyle errors
Project: http://git-wip-us.apache.org/repos/asf/commons-rng/repo Commit: http://git-wip-us.apache.org/repos/asf/commons-rng/commit/9704699b Tree: http://git-wip-us.apache.org/repos/asf/commons-rng/tree/9704699b Diff: http://git-wip-us.apache.org/repos/asf/commons-rng/diff/9704699b Branch: refs/heads/1.1 Commit: 9704699b403f23939529e14b0560cabf00530d98 Parents: 7b4a414 Author: aherbert <a.herb...@sussex.ac.uk> Authored: Wed Aug 1 12:40:04 2018 +0100 Committer: aherbert <a.herb...@sussex.ac.uk> Committed: Wed Aug 1 12:40:04 2018 +0100 ---------------------------------------------------------------------- .../distribution/LargeMeanPoissonSampler.java | 45 +++++++++++++++----- .../distribution/SmallMeanPoissonSampler.java | 10 ++--- 2 files changed, 39 insertions(+), 16 deletions(-) ---------------------------------------------------------------------- http://git-wip-us.apache.org/repos/asf/commons-rng/blob/9704699b/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/LargeMeanPoissonSampler.java ---------------------------------------------------------------------- diff --git a/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/LargeMeanPoissonSampler.java b/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/LargeMeanPoissonSampler.java index c861802..a729686 100644 --- a/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/LargeMeanPoissonSampler.java +++ b/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/LargeMeanPoissonSampler.java @@ -31,7 +31,7 @@ import org.apache.commons.rng.sampling.distribution.InternalUtils.FactorialLog; * </blockquote> * </li> * </ul> - * + * * This sampler is suitable for {@code mean >= 40}. */ public class LargeMeanPoissonSampler @@ -39,7 +39,7 @@ public class LargeMeanPoissonSampler implements DiscreteSampler { /** Class to compute {@code log(n!)}. This has no cached values. */ - static private final InternalUtils.FactorialLog NO_CACHE_FACTORIAL_LOG; + private static final InternalUtils.FactorialLog NO_CACHE_FACTORIAL_LOG; static { // Create without a cache. @@ -52,17 +52,40 @@ public class LargeMeanPoissonSampler private final ContinuousSampler gaussian; /** Local class to compute {@code log(n!)}. This may have cached values. */ private final InternalUtils.FactorialLog factorialLog; - + // Working values + + /** Algorithm constant: {@code Math.floor(mean)}. */ private final double lambda; + /** Algorithm constant: {@code mean - lambda}. */ private final double lambdaFractional; + /** Algorithm constant: {@code Math.log(lambda)}. */ private final double logLambda; + /** Algorithm constant: {@code factorialLog((int) lambda)}. */ private final double logLambdaFactorial; + /** Algorithm constant: {@code Math.sqrt(lambda * Math.log(32 * lambda / Math.PI + 1))}. */ private final double delta; + /** Algorithm constant: {@code delta / 2}. */ private final double halfDelta; + /** Algorithm constant: {@code 2 * lambda + delta}. */ private final double twolpd; + /** + * Algorithm constant: {@code a1 / aSum} with + * <ul> + * <li>{@code a1 = Math.sqrt(Math.PI * twolpd) * Math.exp(c1)}</li> + * <li>{@code aSum = a1 + a2 + 1}</li> + * </ul> + */ private final double p1; + /** + * Algorithm constant: {@code a1 / aSum} with + * <ul> + * <li>{@code a2 = (twolpd / delta) * Math.exp(-delta * (1 + delta) / twolpd)}</li> + * <li>{@code aSum = a1 + a2 + 1}</li> + * </ul> + */ private final double p2; + /** Algorithm constant: {@code 1 / (8 * lambda)}. */ private final double c1; /** The internal Poisson sampler for the lambda fraction. */ @@ -73,13 +96,13 @@ public class LargeMeanPoissonSampler * @param mean Mean. * @throws IllegalArgumentException if {@code mean <= 0}. */ - public LargeMeanPoissonSampler(UniformRandomProvider rng, - double mean) { + public LargeMeanPoissonSampler(UniformRandomProvider rng, + double mean) { super(rng); if (mean <= 0) { throw new IllegalArgumentException(mean + " <= " + 0); } - + gaussian = new ZigguratNormalizedGaussianSampler(rng); exponential = new AhrensDieterExponentialSampler(rng, 1); // Plain constructor uses the uncached function. @@ -105,15 +128,15 @@ public class LargeMeanPoissonSampler null : // Not used. new SmallMeanPoissonSampler(rng, lambdaFractional); } - + /** {@inheritDoc} */ @Override public int sample() { - final int y2 = (smallMeanPoissonSampler == null) ? + final int y2 = (smallMeanPoissonSampler == null) ? 0 : // No lambda fraction smallMeanPoissonSampler.sample(); - + double x = 0; double y = 0; double v = 0; @@ -161,7 +184,7 @@ public class LargeMeanPoissonSampler break; } } - + return (int) Math.min(y2 + (long) y, Integer.MAX_VALUE); } @@ -172,7 +195,7 @@ public class LargeMeanPoissonSampler * @return {@code log(n!)} * @throws IllegalArgumentException if {@code n < 0}. */ - private final double factorialLog(int n) { + private double factorialLog(int n) { return factorialLog.value(n); } http://git-wip-us.apache.org/repos/asf/commons-rng/blob/9704699b/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/SmallMeanPoissonSampler.java ---------------------------------------------------------------------- diff --git a/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/SmallMeanPoissonSampler.java b/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/SmallMeanPoissonSampler.java index 9f47c76..283594e 100644 --- a/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/SmallMeanPoissonSampler.java +++ b/commons-rng-sampling/src/main/java/org/apache/commons/rng/sampling/distribution/SmallMeanPoissonSampler.java @@ -28,15 +28,15 @@ import org.apache.commons.rng.UniformRandomProvider; * The Poisson process (and hence, the returned value) is bounded by 1000 * mean. * </li> * </ul> - * + * * This sampler is suitable for {@code mean < 40}. */ public class SmallMeanPoissonSampler extends SamplerBase implements DiscreteSampler { - /** - * Pre-compute {@code Math.exp(-mean)}. + /** + * Pre-compute {@code Math.exp(-mean)}. * Note: This is the probability of the Poisson sample {@code P(n=0)}. */ private final double p0; @@ -54,7 +54,7 @@ public class SmallMeanPoissonSampler if (mean <= 0) { throw new IllegalArgumentException(mean + " <= " + 0); } - + p0 = Math.exp(-mean); // The returned sample is bounded by 1000 * mean or Integer.MAX_VALUE limit = (int) Math.ceil(Math.min(1000 * mean, Integer.MAX_VALUE)); @@ -76,7 +76,7 @@ public class SmallMeanPoissonSampler } return n; } - + /** {@inheritDoc} */ @Override public String toString() {