diff --git a/src/java.base/share/classes/java/util/random/RandomGenerator.java b/src/java.base/share/classes/java/util/random/RandomGenerator.java index a7c6bcec339..9574fd18594 100644 --- a/src/java.base/share/classes/java/util/random/RandomGenerator.java +++ b/src/java.base/share/classes/java/util/random/RandomGenerator.java @@ -1,5 +1,5 @@ /* - * Copyright (c) 2021, 2022, Oracle and/or its affiliates. All rights reserved. + * Copyright (c) 2021, 2023, Oracle and/or its affiliates. All rights reserved. * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. * * This code is free software; you can redistribute it and/or modify it @@ -30,13 +30,14 @@ import java.security.SecureRandom; import java.util.Objects; import java.util.concurrent.ThreadLocalRandom; import jdk.internal.util.random.RandomSupport; -import jdk.internal.util.random.RandomSupport.*; import java.util.stream.DoubleStream; import java.util.stream.IntStream; import java.util.stream.LongStream; import java.util.stream.Stream; +import static java.lang.Math.*; + /** * The {@link RandomGenerator} interface is designed to provide a common * protocol for objects that generate random or (more typically) pseudorandom @@ -252,6 +253,174 @@ public interface RandomGenerator { return doubles(randomNumberOrigin, randomNumberBound).limit(streamSize); } + /** + * Returns an effectively unlimited stream of pseudorandomly chosen + * {@code double} values, where each value is between the specified + * {@code left} boundary and the specified {@code right} boundary. + * The {@code left} boundary is included as indicated by + * {@code isLeftIncluded}; similarly, the {@code right} boundary is included + * as indicated by {@code isRightIncluded}. + * + *

The stream potentially produces all multiples k δ + * (k integer) lying in the interval specified by the parameters, + * where δ > 0 is the smallest number for which all these multiples + * are exact {@code double}s. + * They are therefore all equidistant. + * The uniformity of the distribution of the {@code double}s produced by + * the stream depends on the quality of the underlying {@link #nextLong(long)}. + * + * @implSpec The default implementation first determines the δ above. + * It then computes both the smallest integer kl + * such that kl δ lies inside + * the given interval, and the smallest integer n > 0 such that + * (kl + n) δ lies + * outside the interval. + * Finally, it returns a stream which generates the {@code double}s + * according to (kl + {@code nextLong(}n{@code )}) + * δ. + * The stream never produces {@code -0.0}, although it may produce + * {@code 0.0} if the specified interval contains 0. + * + * @param left the left boundary + * @param right the right boundary + * @param isLeftIncluded whether the {@code left} boundary is included + * @param isRightIncluded whether the {@code right} boundary is included + * + * @return a stream of pseudorandomly chosen {@code double} values, each + * between {@code left} and {@code right}, as specified above. + * + * @throws IllegalArgumentException if {@code left} is not finite, + * or {@code right} is not finite, or if the specified interval + * is empty. + * + * @since 21 + */ + default DoubleStream equiDoubles(double left, double right, + boolean isLeftIncluded, boolean isRightIncluded) { + if (!(Double.NEGATIVE_INFINITY < left + && right < Double.POSITIVE_INFINITY + && (isLeftIncluded ? left : nextUp(left)) + <= (isRightIncluded ? right : nextDown(right)))) { + throw new IllegalArgumentException( + "the boundaries must be finite and the interval must not be empty"); + } + + /* + * Inspired by + * Goualard, "Drawing random floating-point numbers from an + * interval", ACM TOMACS, 2022, 32 (3) + * (https://hal.science/hal-03282794v4) + * although implemented differently. + * + * It is assumed that left <= right. + * Whether the boundaries of the interval I = are included + * is indicated by isLeftIncluded and isRightIncluded. + * + * delta > 0 is the smallest double such that every product k delta + * (k integer) that lies in I is an exact double as well. + * It turns out that delta is always a power of 2. + * + * kl is the smallest k such that k delta is inside I. + * kr > kl is the smallest k such that k delta is outside I. + * n is kr - kl + */ + double delta; // captured + long kl; // captured + long kr; + long n; // captured + + if (left <= -right) { + /* + * Here, + * left <= 0, left <= right <= -left + * P = Double.PRECISION + * + * delta is the distance from left to the next double in the + * direction of positive infinity. + * Most of the time, this is equivalent to the ulp of left, but not + * always. + * For example, for left == -1.0, Math.ulp(left) == 2.220446049250313E-16, + * whereas delta as computed here is 1.1102230246251565E-16. + * + * Every product k delta lying in [left, -left] is an exact double. + * Thus, every product k delta lying in I is an exact double, too. + * Any other positive eps < delta does not meet this property: + * some product k eps lying in I is not an exact double. + * On the other hand, any other eps > delta would generate more + * sparse products k eps, that is, fewer doubles in I. + * delta is therefore the best value to ensure the largest number + * of equidistant doubles in the interval I. + * + * left / delta is an exact double and an exact integer with + * -2^P <= left / delta <= 0 + * Thus, kl is computed exactly. + * + * Mathematically, + * kr = ceil(right / delta), if !isRightIncluded + * kr = floor(right / delta) + 1, if isRightIncluded + * The double division rd = right / delta never overflows and is + * exact, except in the presence of underflows. But even underflows + * do not affect the outcomes of ceil() and floor(), except, + * in turn, when the result drops to 0, that is, rd = 0. + * + * crd is a corrected version of rd when rd is zero. It is simply + * right / delta, but rounded away from 0 to preserve information + * ensuring correct outcomes in ceil() and floor(). + * + * We know that -2^P <= kl, so + * -2^P <= kl + nextLong(n) + * Also, since right <= -left, we know that + * kr <= -kl + 1 + * so that + * 0 < n <= -2 kl + 1 + * This implies + * kl + nextLong(n) <= kl + (-2 kl) = -kl <= 2^P + * and thus + * -2^P <= kl + nextLong(n) <= 2^P + * which shows that kl + nextLong(n) can be cast exactly to double. + * + * Further, if isLeftIncluded then left = kl delta, so that we get + * left = kl * delta <= (kl + nextLong(n)) * delta + * For any other k < kl, when nextLong(n) = 0 we would have + * (k + nextLong(n)) * delta < left + * Otherwise, left = (kl - 1) delta, and therefore + * left = (kl - 1) * delta < (kl + nextLong(n)) * delta + * For any other k < kl, when nextLong(n) = 0 we would get + * (k + nextLong(n)) * delta <= left + * Either way, the lhs expression would not belong to I. + * That is, kl is the smallest integer such that kl delta always + * lies in I (it is an exact double). + * + * Similar considerations show that kr is the smallest integer such + * that kr delta lies to the right of I (it is an exact double). + * + * All the above means that (kl + nextLong(n)) * delta is an exact + * double lying in I and that kl and kr, thus n, are the best + * possible choices to ensure the largest number of equidistant + * doubles in I. Uniform distribution relies on the guarantee + * afforded by nextLong(). + */ + delta = nextUp(left) - left; + double rd = right / delta; + double crd = rd != 0 || right == 0 ? rd : copySign(Double.MIN_VALUE, right); + kr = isRightIncluded ? (long) floor(crd) + 1 : (long) ceil(crd); + kl = (long) (left / delta) + (isLeftIncluded ? 0 : 1); + } else { + /* Here, + * right > 0, -right < left <= right + * + * Considerations similar to the ones above apply here as well. + */ + delta = right - nextDown(right); + double ld = left / delta; + double cld = ld != 0 || left == 0 ? ld : copySign(Double.MIN_VALUE, left); + kl = isLeftIncluded ? (long) ceil(cld) : (long) floor(cld) + 1; + kr = (long) (right / delta) + (isRightIncluded ? 1 : 0); + } + n = kr - kl; + return DoubleStream.generate(() -> (kl + nextLong(n)) * delta).sequential(); + } + /** * Returns an effectively unlimited stream of pseudorandomly chosen * {@code int} values. diff --git a/test/jdk/java/util/Random/EquiDoublesTest.java b/test/jdk/java/util/Random/EquiDoublesTest.java new file mode 100644 index 00000000000..ddc2cdfe8de --- /dev/null +++ b/test/jdk/java/util/Random/EquiDoublesTest.java @@ -0,0 +1,349 @@ +/* + * Copyright (c) 2023, Oracle and/or its affiliates. All rights reserved. + * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER. + * + * This code is free software; you can redistribute it and/or modify it + * under the terms of the GNU General Public License version 2 only, as + * published by the Free Software Foundation. + * + * This code is distributed in the hope that it will be useful, but WITHOUT + * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or + * FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License + * version 2 for more details (a copy is included in the LICENSE file that + * accompanied this code). + * + * You should have received a copy of the GNU General Public License version + * 2 along with this work; if not, write to the Free Software Foundation, + * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA. + * + * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA + * or visit www.oracle.com if you need additional information or have any + * questions. + */ + +import jdk.test.lib.RandomFactory; +import org.junit.jupiter.params.ParameterizedTest; +import org.junit.jupiter.params.provider.Arguments; +import org.junit.jupiter.params.provider.MethodSource; + +import java.util.Iterator; +import java.util.TreeSet; +import java.util.random.RandomGenerator; +import java.util.stream.DoubleStream; + +import static org.junit.jupiter.api.Assertions.*; +import static org.junit.jupiter.params.provider.Arguments.arguments; + +/** + * @test + * @bug 8302987 + * @key randomness + * + * @summary Check consistency of RandomGenerator::equiDoubles + * @library /test/lib + * @run junit EquiDoublesTest + * + */ + +public class EquiDoublesTest { + + private static final int SAMPLES = 100_000; + + /* + * A factor to use in the tight*() tests to make sure that + * all equidistant doubles are generated. + */ + private static final long SAFETY_FACTOR = 100L; + private static final RandomGenerator rnd = RandomFactory.getRandom(); + + private static double nextUp(double d, int steps) { + for (int i = 0; i < steps; ++i) { + d = Math.nextUp(d); + } + return d; + } + + private static double nextDown(double d, int steps) { + for (int i = 0; i < steps; ++i) { + d = Math.nextDown(d); + } + return d; + } + + static Arguments[] equi() { + return new Arguments[] { + arguments(0.0, 1e-9), + arguments(1.0, 1.1), + arguments(1.0e23, 1.1e23), + arguments(1.0e300, 1.1e300), + arguments(-1.2, 1.1), + arguments(-1.2e-30, 1.1e6), + arguments(-Double.MIN_VALUE, Double.MIN_VALUE), + arguments(-Double.MAX_VALUE, Double.MAX_VALUE), + }; + } + + @ParameterizedTest + @MethodSource + void equi(double l, double r) { + double[] minmax = new double[2]; + + resetMinmax(minmax); + DoubleStream equi = rnd.equiDoubles(l, r, true, false); + equi.limit(SAMPLES).forEach(d -> updateMinmax(minmax, d)); + assertTrue(l <= minmax[0]); + assertTrue(minmax[1] < r); + + resetMinmax(minmax); + equi = rnd.equiDoubles(l, r, true, true); + equi.limit(SAMPLES).forEach(d -> updateMinmax(minmax, d)); + assertTrue(l <= minmax[0]); + assertTrue(minmax[1] <= r); + + resetMinmax(minmax); + equi = rnd.equiDoubles(l, r, false, true); + equi.limit(SAMPLES).forEach(d -> updateMinmax(minmax, d)); + assertTrue(l < minmax[0]); + assertTrue(minmax[1] <= r); + + resetMinmax(minmax); + equi = rnd.equiDoubles(l, r, false, false); + equi.limit(SAMPLES).forEach(d -> updateMinmax(minmax, d)); + assertTrue(l < minmax[0]); + assertTrue(minmax[1] < r); + + /* with negated intervals */ + resetMinmax(minmax); + equi = rnd.equiDoubles(-r, -l, true, false); + equi.limit(SAMPLES).forEach(d -> updateMinmax(minmax, d)); + assertTrue(-r <= minmax[0]); + assertTrue(minmax[1] < -l); + + resetMinmax(minmax); + equi = rnd.equiDoubles(-r, -l, true, true); + equi.limit(SAMPLES).forEach(d -> updateMinmax(minmax, d)); + assertTrue(-r <= minmax[0]); + assertTrue(minmax[1] <= -l); + + resetMinmax(minmax); + equi = rnd.equiDoubles(-r, -l, false, true); + equi.limit(SAMPLES).forEach(d -> updateMinmax(minmax, d)); + assertTrue(-r < minmax[0]); + assertTrue(minmax[1] <= -l); + + resetMinmax(minmax); + equi = rnd.equiDoubles(-r, -l, false, false); + equi.limit(SAMPLES).forEach(d -> updateMinmax(minmax, d)); + assertTrue(-r < minmax[0]); + assertTrue(minmax[1] < -l); + } + + private void resetMinmax(double[] minmax) { + minmax[0] = Double.POSITIVE_INFINITY; + minmax[1] = Double.NEGATIVE_INFINITY; + } + + private void updateMinmax(double[] minmax, double d) { + if (d < minmax[0]) { + minmax[0] = d; + } + if (d > minmax[1]) { + minmax[1] = d; + } + } + + static Arguments[] tight() { + return new Arguments[] { + arguments(0.0, (short) 100), + arguments(1.0, (short) 100), + arguments(1.1, (short) 100), + arguments(1.0e23, (short) 100), + arguments(1.0e300, (short) 100), + arguments(-1.2, (short) 100), + arguments(-1.2e-30, (short) 100), + arguments(-Double.MIN_VALUE, (short) 100), + + arguments(-Double.MIN_VALUE, (short) 2), + arguments(-Double.MAX_VALUE, (short) 2), + }; + } + + /* + * All equidistant doubles in a tight range are expected to be generated. + * The arguments must be chosen as to not overlap a value with irregular + * spacing around it. + */ + @ParameterizedTest + @MethodSource + void tight(double l, short steps) { + double r = nextUp(l, steps); + + TreeSet set = new TreeSet<>(); + DoubleStream equi = rnd.equiDoubles(l, r, true, false); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(l, r, true, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps + 1, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(l, r, false, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(l, r, false, false); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps - 1, set.size()); + checkEquidistance(set); + + /* with negated intervals */ + set.clear(); + equi = rnd.equiDoubles(-r, -l, true, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps + 1, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(-r, -l, true, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps + 1, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(-r, -l, false, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(-r, -l, false, false); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps - 1, set.size()); + checkEquidistance(set); + } + + static Arguments[] tightWithIrregularSpacing() { + return new Arguments[] { + arguments(0x1p-1, (short) 15, (short) 23), + arguments(0x1p0, (short) 17, (short) 5), + arguments(0x1p1, (short) 7, (short) 8), + arguments(0x1p-600, (short) 28, (short) 33), + arguments(0x1p600, (short) 9, (short) 19), + }; + } + + /* + * m must be a power of 2 greater than Double.MIN_NORMAL + */ + @ParameterizedTest + @MethodSource + void tightWithIrregularSpacing(double m, short lSteps, short rSteps) { + double l = nextDown(m, 2 * lSteps); + double r = nextUp(m, rSteps); + int steps = lSteps + rSteps; + + TreeSet set = new TreeSet<>(); + DoubleStream equi = rnd.equiDoubles(l, r, true, false); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(l, r, true, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps + 1, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(l, r, false, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(l, r, false, false); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps - 1, set.size()); + checkEquidistance(set); + + /* with negated intervals */ + set.clear(); + equi = rnd.equiDoubles(-r, -l, true, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps + 1, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(-r, -l, true, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps + 1, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(-r, -l, false, true); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps, set.size()); + checkEquidistance(set); + + set.clear(); + equi = rnd.equiDoubles(-r, -l, false, false); + equi.limit(SAFETY_FACTOR * steps).forEach(set::add); + assertEquals(steps - 1, set.size()); + checkEquidistance(set); + } + + private void checkEquidistance(TreeSet set) { + if (set.size() < 3) { + return; + } + Iterator iter = set.iterator(); + double prev = iter.next(); + double curr = iter.next(); + double delta = curr - prev; + while (iter.hasNext()) { + prev = curr; + curr = iter.next(); + assertEquals(delta, curr - prev); + } + } + + static Arguments[] empty() { + return new Arguments[] { + arguments(1.0), + arguments(-1.0), + arguments(0.0), + arguments(nextDown(Double.MAX_VALUE, 1)), + arguments(nextUp(-Double.MAX_VALUE, 1)), + }; + } + + @ParameterizedTest + @MethodSource + void empty(double l) { + assertThrows(IllegalArgumentException.class, + () -> rnd.equiDoubles(l, l, true, false) + ); + assertThrows(IllegalArgumentException.class, + () -> rnd.equiDoubles(l, nextUp(l, 1), false, false) + ); + assertThrows(IllegalArgumentException.class, + () -> rnd.equiDoubles(nextDown(l, 1), l, false, false) + ); + assertThrows(IllegalArgumentException.class, + () -> rnd.equiDoubles(l, l, false, true) + ); + assertThrows(IllegalArgumentException.class, + () -> rnd.equiDoubles(l, nextDown(l, 1), true, true) + ); + assertThrows(IllegalArgumentException.class, + () -> rnd.equiDoubles(nextUp(l, 1), l, true, true) + ); + } + +}