diff --git a/core/src/main/scala/spire/random/Ziggurat.scala b/core/src/main/scala/spire/random/Ziggurat.scala index 42393ce7a..33a071bd9 100644 --- a/core/src/main/scala/spire/random/Ziggurat.scala +++ b/core/src/main/scala/spire/random/Ziggurat.scala @@ -76,8 +76,8 @@ object Ziggurat { x = -log(g.nextDouble()) * r1 y = -log(g.nextDouble()) y + y < x * x - }) - return if (hz > 0) r + x else -r - x + }) () + return if (hz > 0) r + x else -r - x } if (fn(iz) + g.nextDouble() * (fn(iz - 1) - fn(iz)) < exp(-.5 * x * x)) return x diff --git a/tests/shared/src/test/scala/spire/random/GaussianSuite.scala b/tests/shared/src/test/scala/spire/random/GaussianSuite.scala index 0895041ee..91bff27ad 100644 --- a/tests/shared/src/test/scala/spire/random/GaussianSuite.scala +++ b/tests/shared/src/test/scala/spire/random/GaussianSuite.scala @@ -44,6 +44,14 @@ class GaussianSuite extends munit.FunSuite { } } + test("spire.random.Gaussian does not crash Ziggurat (regression)") { + try { + spire.random.Gaussian(0.0, 1.0).histogram(975)(rng.Lcg64.fromTime(42L)) + } catch { + case e: Throwable => fail(f"crashed with $e") + } + } + test("MarsagliaGaussian[Float] is normal")(checkMarsagliaGaussian[Float]) test("MarsagliaGaussian[Double] is normal")(checkMarsagliaGaussian[Double]) test("MarsagliaGaussian[BigDecimal] is normal")(checkMarsagliaGaussian[BigDecimal])