1 /*

     2  * Copyright (c) 2013, 2019, Oracle and/or its affiliates. All rights reserved.

     3  * DO NOT ALTER OR REMOVE COPYRIGHT NOTICES OR THIS FILE HEADER.

     4  *

     5  * This code is free software; you can redistribute it and/or modify it

     6  * under the terms of the GNU General Public License version 2 only, as

     7  * published by the Free Software Foundation.  Oracle designates this

     8  * particular file as subject to the "Classpath" exception as provided

     9  * by Oracle in the LICENSE file that accompanied this code.

    10  *

    11  * This code is distributed in the hope that it will be useful, but WITHOUT

    12  * ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or

    13  * FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License

    14  * version 2 for more details (a copy is included in the LICENSE file that

    15  * accompanied this code).

    16  *

    17  * You should have received a copy of the GNU General Public License version

    18  * 2 along with this work; if not, write to the Free Software Foundation,

    19  * Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA.

    20  *

    21  * Please contact Oracle, 500 Oracle Parkway, Redwood Shores, CA 94065 USA

    22  * or visit www.oracle.com if you need additional information or have any

    23  * questions.

    24  */

    25 

    26 package java.util.random;

    27 

    28 import java.math.BigInteger;

    29 import java.util.concurrent.atomic.AtomicLong;

    30 import java.util.random.RandomGenerator.SplittableGenerator;

    31 import java.util.random.RandomSupport.AbstractSplittableWithBrineGenerator;

    32 

    33 /**

    34  * A generator of uniform pseudorandom values applicable for use in

    35  * (among other contexts) isolated parallel computations that may

    36  * generate subtasks.  Class {@link L64X1024MixRandom} implements

    37  * interfaces {@link RandomGenerator} and {@link SplittableGenerator},

    38  * and therefore supports methods for producing pseudorandomly chosen

    39  * numbers of type {@code int}, {@code long}, {@code float}, and {@code double}

    40  * as well as creating new split-off {@link L64X1024MixRandom} objects,

    41  * with similar usages as for class {@link java.util.SplittableRandom}.

    42  * <p>

    43  * Series of generated values pass the TestU01 BigCrush and PractRand test suites

    44  * that measure independence and uniformity properties of random number generators.

    45  * (Most recently validated with

    46  * <a href="http://simul.iro.umontreal.ca/testu01/tu01.html">version 1.2.3 of TestU01</a>

    47  * and <a href="http://pracrand.sourceforge.net">version 0.90 of PractRand</a>.

    48  * Note that TestU01 BigCrush was used to test not only values produced by the {@code nextLong()}

    49  * method but also the result of bit-reversing each value produced by {@code nextLong()}.)

    50  * These tests validate only the methods for certain

    51  * types and ranges, but similar properties are expected to hold, at

    52  * least approximately, for others as well.

    53  * <p>

    54  * {@link L64X1024MixRandom} is a specific member of the LXM family of algorithms

    55  * for pseudorandom number generators.  Every LXM generator consists of two

    56  * subgenerators; one is an LCG (Linear Congruential Generator) and the other is

    57  * an Xorshift generator.  Each output of an LXM generator is the result of

    58  * combining state from the LCG with state from the Xorshift generator by

    59  * using a Mixing function (and then the state of the LCG and the state of the

    60  * Xorshift generator are advanced).

    61  * <p>

    62  * The LCG subgenerator for {@link L64X1024MixRandom} has an update step of the

    63  * form {@code s = m * s + a}, where {@code s}, {@code m}, and {@code a} are all

    64  * of type {@code long}; {@code s} is the mutable state, the multiplier {@code m}

    65  * is fixed (the same for all instances of {@link L64X1024MixRandom}) and the addend

    66  * {@code a} is a parameter (a final field of the instance).  The parameter

    67  * {@code a} is required to be odd (this allows the LCG to have the maximal

    68  * period, namely 2<sup>64</sup>); therefore there are 2<sup>63</sup> distinct choices

    69  * of parameter.

    70  * <p>

    71  * The Xorshift subgenerator for {@link L64X1024MixRandom} is the {@code xoroshiro1024}

    72  * algorithm (parameters 25, 27, and 36), without any final scrambler such as "+" or "**".

    73  * Its state consists of an array {@code x} of sixteen {@code long} values,

    74  * which can take on any values provided that they are not all zero.

    75  * The period of this subgenerator is 2<sup>1024</sup>-1.

    76  * <p>

    77  * The mixing function for {@link L64X1024MixRandom} is {@link RandomSupport.mixLea64}

    78  * applied to the argument {@code (s + s0)}, where {@code s0} is the most recently computed

    79  * element of {@code x}.

    80  * <p>

    81  * Because the periods 2<sup>64</sup> and 2<sup>1024</sup>-1 of the two subgenerators

    82  * are relatively prime, the <em>period</em> of any single {@link L64X1024MixRandom} object

    83  * (the length of the series of generated 64-bit values before it repeats) is the product

    84  * of the periods of the subgenerators, that is, 2<sup>64</sup>(2<sup>1024</sup>-1),

    85  * which is just slightly smaller than 2<sup>1088</sup>.  Moreover, if two distinct

    86  * {@link L64X1024MixRandom} objects have different {@code a} parameters, then their

    87  * cycles of produced values will be different.

    88  * <p>

    89  * The 64-bit values produced by the {@code nextLong()} method are exactly equidistributed.

    90  * For any specific instance of {@link L64X1024MixRandom}, over the course of its cycle each

    91  * of the 2<sup>64</sup> possible {@code long} values will be produced 2<sup>1024</sup>-1 times.

    92  * The values produced by the {@code nextInt()}, {@code nextFloat()}, and {@code nextDouble()}

    93  * methods are likewise exactly equidistributed.

    94  * <p>

    95  * In fact, the 64-bit values produced by the {@code nextLong()} method are 16-equidistributed.

    96  * To be precise: for any specific instance of {@link L64X1024MixRandom}, consider

    97  * the (overlapping) length-16 subsequences of the cycle of 64-bit values produced by

    98  * {@code nextLong()} (assuming no other methods are called that would affect the state).

    99  * There are 2<sup>64</sup>(2<sup>1024</sup>-1) such subsequences, and each subsequence,

   100  * which consists of 16 64-bit values, can have one of 2<sup>1024</sup> values. Of those

   101  * 2<sup>1024</sup> subsequence values, nearly all of them (2<sup>1024</sup>-2<sup>64</sup>)

   102  * occur 2<sup>64</sup> times over the course of the entire cycle, and the other

   103  * 2<sup>64</sup> subsequence values occur only 2<sup>64</sup>-1 times.  So the ratio

   104  * of the probability of getting any specific one of the less common subsequence values and the

   105  * probability of getting any specific one of the more common subsequence values is 1-2<sup>-64</sup>.

   106  * (Note that the set of 2<sup>64</sup> less-common subsequence values will differ from

   107  * one instance of {@link L64X1024MixRandom} to another, as a function of the additive

   108  * parameter of the LCG.)  The values produced by the {@code nextInt()}, {@code nextFloat()},

   109  * and {@code nextDouble()} methods are likewise 16-equidistributed.

   110  * <p>

   111  * Method {@link #split} constructs and returns a new {@link L64X1024MixRandom}

   112  * instance that shares no mutable state with the current instance. However, with

   113  * very high probability, the values collectively generated by the two objects

   114  * have the same statistical properties as if the same quantity of values were

   115  * generated by a single thread using a single {@link L64X1024MixRandom} object.

   116  * This is because, with high probability, distinct {@link L64X1024MixRandom} objects

   117  * have distinct {@code a} parameters and therefore use distinct members of the

   118  * algorithmic family; and even if their {@code a} parameters are the same, with

   119  * very high probability they will traverse different parts of their common state

   120  * cycle.

   121  * <p>

   122  * As with {@link java.util.SplittableRandom}, instances of

   123  * {@link L64X1024MixRandom} are <em>not</em> thread-safe.

   124  * They are designed to be split, not shared, across threads. For

   125  * example, a {@link java.util.concurrent.ForkJoinTask} fork/join-style

   126  * computation using random numbers might include a construction

   127  * of the form {@code new Subtask(someL64X1024MixRandom.split()).fork()}.

   128  * <p>

   129  * This class provides additional methods for generating random

   130  * streams, that employ the above techniques when used in

   131  * {@code stream.parallel()} mode.

   132  * <p>

   133  * Instances of {@link L64X1024MixRandom} are not cryptographically

   134  * secure.  Consider instead using {@link java.security.SecureRandom}

   135  * in security-sensitive applications. Additionally,

   136  * default-constructed instances do not use a cryptographically random

   137  * seed unless the {@linkplain System#getProperty system property}

   138  * {@code java.util.secureRandomSeed} is set to {@code true}.

   139  *

   140  * @since 14

   141  */

   142 public final class L64X1024MixRandom extends AbstractSplittableWithBrineGenerator {

   143 

   144     /*

   145      * Implementation Overview.

   146      *

   147      * The split() operation uses the current generator to choose 18 new 64-bit

   148      * long values that are then used to initialize the parameter `a`, the

   149      * state variable `s`, and the array `x` for a newly constructed generator.

   150      *

   151      * With extremely high probability, no two generators so chosen

   152      * will have the same `a` parameter, and testing has indicated

   153      * that the values generated by two instances of {@link L64X1024MixRandom}

   154      * will be (approximately) independent if have different values for `a`.

   155      *

   156      * The default (no-argument) constructor, in essence, uses

   157      * "defaultGen" to generate 18 new 64-bit values for the same

   158      * purpose.  Multiple generators created in this way will certainly

   159      * differ in their `a` parameters.  The defaultGen state must be accessed

   160      * in a thread-safe manner, so we use an AtomicLong to represent

   161      * this state.  To bootstrap the defaultGen, we start off using a

   162      * seed based on current time unless the

   163      * java.util.secureRandomSeed property is set. This serves as a

   164      * slimmed-down (and insecure) variant of SecureRandom that also

   165      * avoids stalls that may occur when using /dev/random.

   166      *

   167      * File organization: First static fields, then instance

   168      * fields, then constructors, then instance methods.

   169      */

   170 

   171     /* ---------------- static fields ---------------- */

   172 

   173     /*

   174      * The length of the array x.

   175      */

   176 

   177     private static final int N = 16;

   178 

   179     /**

   180      * The seed generator for default constructors.

   181      */

   182     private static final AtomicLong defaultGen = new AtomicLong(RandomSupport.initialSeed());

   183 

   184     /*

   185      * The period of this generator, which is (2**1024 - 1) * 2**64.

   186      */

   187     private static final BigInteger PERIOD =

   188         BigInteger.ONE.shiftLeft(N*64).subtract(BigInteger.ONE).shiftLeft(64);

   189 

   190     /*

   191      * Multiplier used in the LCG portion of the algorithm.

   192      * Chosen based on research by Sebastiano Vigna and Guy Steele (2019).

   193      * The spectral scores for dimensions 2 through 8 for the multiplier 0xd1342543de82ef95

   194      * are [0.958602, 0.937479, 0.870757, 0.822326, 0.820405, 0.813065, 0.760215].

   195      */

   196 

   197     private static final long M = 0xd1342543de82ef95L;

   198 

   199     /* ---------------- instance fields ---------------- */

   200 

   201     /**

   202      * The parameter that is used as an additive constant for the LCG.

   203      * Must be odd.

   204      */

   205     private final long a;

   206 

   207     /**

   208      * The per-instance state: s for the LCG; the array x for the xorshift;

   209      * p is the rotating pointer into the array x.

   210      * At least one of the 16 elements of the array x must be nonzero.

   211      */

   212     private long s;

   213     private final long[] x;

   214     private int p = N - 1;

   215 

   216     /* ---------------- constructors ---------------- */

   217 

   218     /**

   219      * Basic constructor that initializes all fields from parameters.

   220      * It then adjusts the field values if necessary to ensure that

   221      * all constraints on the values of fields are met.

   222      *

   223      * @param a additive parameter for the LCG

   224      * @param s initial state for the LCG

   225      * @param x0 first word of the initial state for the xorshift generator

   226      * @param x1 second word of the initial state for the xorshift generator

   227      * @param x2 third word of the initial state for the xorshift generator

   228      * @param x3 fourth word of the initial state for the xorshift generator

   229      * @param x4 fifth word of the initial state for the xorshift generator

   230      * @param x5 sixth word of the initial state for the xorshift generator

   231      * @param x6 seventh word of the initial state for the xorshift generator

   232      * @param x7 eight word of the initial state for the xorshift generator

   233      * @param x8 ninth word of the initial state for the xorshift generator

   234      * @param x9 tenth word of the initial state for the xorshift generator

   235      * @param x10 eleventh word of the initial state for the xorshift generator

   236      * @param x11 twelfth word of the initial state for the xorshift generator

   237      * @param x12 thirteenth word of the initial state for the xorshift generator

   238      * @param x13 fourteenth word of the initial state for the xorshift generator

   239      * @param x14 fifteenth word of the initial state for the xorshift generator

   240      * @param x15 sixteenth word of the initial state for the xorshift generator

   241      */

   242     public L64X1024MixRandom(long a, long s,

   243                              long x0, long x1, long x2, long x3,

   244                              long x4, long x5, long x6, long x7,

   245                              long x8, long x9, long x10, long x11,

   246                              long x12, long x13, long x14, long x15) {

   247         // Force a to be odd.

   248         this.a = a | 1;

   249         this.s = s;

   250         this.x = new long[N];

   251         this.x[0] = x0;

   252         this.x[1] = x1;

   253         this.x[2] = x2;

   254         this.x[3] = x3;

   255         this.x[4] = x4;

   256         this.x[5] = x5;

   257         this.x[6] = x6;

   258         this.x[7] = x7;

   259         this.x[8] = x8;

   260         this.x[9] = x9;

   261         this.x[10] = x10;

   262         this.x[11] = x11;

   263         this.x[12] = x12;

   264         this.x[13] = x13;

   265         this.x[14] = x14;

   266         this.x[15] = x15;

   267         // If x0, x1, ..., x15 are all zero (very unlikely), we must choose nonzero values.

   268         if ((x0 | x1 | x2 | x3 | x4 | x5 | x6 | x7 | x8 | x9 | x10 | x11 | x12 | x13 | x14 | x15) == 0) {

   269 	    long v = s;

   270             // At least fifteen of the sixteen values generated here will be nonzero.

   271             for (int j = 0; j < N; j++) {

   272                 this.x[j] = RandomSupport.mixStafford13(v += RandomSupport.GOLDEN_RATIO_64);

   273             }

   274         }

   275     }

   276 

   277     /**

   278      * Creates a new instance of {@link L64X1024MixRandom} using the

   279      * specified {@code long} value as the initial seed. Instances of

   280      * {@link L64X1024MixRandom} created with the same seed in the same

   281      * program execution generate identical sequences of values.

   282      *

   283      * @param seed the initial seed

   284      */

   285     public L64X1024MixRandom(long seed) {

   286         // Using a value with irregularly spaced 1-bits to xor the seed

   287         // argument tends to improve "pedestrian" seeds such as 0 or

   288         // other small integers.  We may as well use SILVER_RATIO_64.

   289         //

   290         // The seed is hashed by mixMurmur64 to produce the `a` parameter.

   291         // The seed is hashed by mixStafford13 to produce the initial `x[0]`,

   292         // which will then be used to produce the first generated value.

   293         // The other x values are filled in as if by a SplitMix PRNG with

   294         // GOLDEN_RATIO_64 as the gamma value and mixStafford13 as the mixer.

   295         this(RandomSupport.mixMurmur64(seed ^= RandomSupport.SILVER_RATIO_64),

   296              1,

   297              RandomSupport.mixStafford13(seed),

   298              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   299              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   300              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   301              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   302              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   303              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   304              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   305              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   306              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   307              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   308              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   309              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   310              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   311              RandomSupport.mixStafford13(seed += RandomSupport.GOLDEN_RATIO_64),

   312              RandomSupport.mixStafford13(seed + RandomSupport.GOLDEN_RATIO_64));

   313     }

   314 

   315     /**

   316      * Creates a new instance of {@link L64X1024MixRandom} that is likely to

   317      * generate sequences of values that are statistically independent

   318      * of those of any other instances in the current program execution,

   319      * but may, and typically does, vary across program invocations.

   320      */

   321     public L64X1024MixRandom() {

   322         // Using GOLDEN_RATIO_64 here gives us a good Weyl sequence of values.

   323         this(defaultGen.getAndAdd(RandomSupport.GOLDEN_RATIO_64));

   324     }

   325 

   326     /**

   327      * Creates a new instance of {@link L64X1024MixRandom} using the specified array of

   328      * initial seed bytes. Instances of {@link L64X1024MixRandom} created with the same

   329      * seed array in the same program execution generate identical sequences of values.

   330      *

   331      * @param seed the initial seed

   332      */

   333     public L64X1024MixRandom(byte[] seed) {

   334         // Convert the seed to 18 long values, of which the last 16 are not all zero.

   335         long[] data = RandomSupport.convertSeedBytesToLongs(seed, 18, 16);

   336         long a = data[0], s = data[1];

   337         // Force a to be odd.

   338         this.a = a | 1;

   339         this.s = s;

   340         this.x = new long[N];

   341         for (int j = 0; j < N; j++) {

   342             this.x[j] = data[2+j];

   343         }

   344     }

   345 

   346     /* ---------------- public methods ---------------- */

   347     /**

   348      * Given 63 bits of "brine", constructs and returns a new instance of

   349      * {@code L64X1024MixRandom} that shares no mutable state with this instance.

   350      * However, with very high probability, the set of values collectively

   351      * generated by the two objects has the same statistical properties as if

   352      * same the quantity of values were generated by a single thread using

   353      * a single {@code L64X1024MixRandom} object.  Either or both of the two

   354      * objects may be further split using the {@code split} method,

   355      * and the same expected statistical properties apply to the

   356      * entire set of generators constructed by such recursive splitting.

   357      *

   358      * @param source a {@code SplittableGenerator} instance to be used instead

   359      *               of this one as a source of pseudorandom bits used to

   360      *               initialize the state of the new ones.

   361      * @param brine a long value, of which the low 63 bits are used to choose

   362      *              the {@code a} parameter for the new instance.

   363      * @return a new instance of {@code L64X1024MixRandom}

   364      */

   365     public SplittableGenerator split(SplittableGenerator source, long brine) {

   366 	// Pick a new instance "at random", but use the brine for `a`.

   367         return new L64X1024MixRandom(brine << 1, source.nextLong(),

   368 				     source.nextLong(), source.nextLong(),

   369                                      source.nextLong(), source.nextLong(),

   370                                      source.nextLong(), source.nextLong(),

   371                                      source.nextLong(), source.nextLong(),

   372                                      source.nextLong(), source.nextLong(),

   373                                      source.nextLong(), source.nextLong(),

   374                                      source.nextLong(), source.nextLong(),

   375                                      source.nextLong(), source.nextLong());

   376     }

   377 

   378     /**

   379      * Returns a pseudorandom {@code long} value.

   380      *

   381      * @return a pseudorandom {@code long} value

   382      */

   383     public long nextLong() {

   384         // First part of xoroshiro1024: fetch array data

   385         final int q = p;

   386         final long s0 = x[p = (p + 1) & (N - 1)];

   387         long s15 = x[q];

   388 

   389 	// Compute the result based on current state information

   390 	// (this allows the computation to be overlapped with state update).

   391 

   392 	final long result = RandomSupport.mixLea64(s + s0);

   393 	

   394 	// Update the LCG subgenerator

   395         s = M * s + a;  // LCG

   396 

   397         // Second part of xoroshiro1024: update array data

   398         s15 ^= s0;

   399         x[q] = Long.rotateLeft(s0, 25) ^ s15 ^ (s15 << 27);

   400         x[p] = Long.rotateLeft(s15, 36);

   401 

   402         return result;

   403     }

   404 

   405     /**

   406      * Returns the period of this random generator.

   407      *

   408      * @return a {@link BigInteger} whose value is the number of distinct possible states of this

   409      *         {@link RandomGenerator} object (2<sup>64</sup>(2<sup>1024</sup>-1)).

   410      */

   411     public BigInteger period() {

   412         return PERIOD;

   413     }

   414 }