001package squidpony.squidmath;
002
003import squidpony.StringKit;
004
005import java.io.Serializable;
006
007/**
008 * A mid-high-quality StatefulRandomness that is the second-fastest 64-bit generator in this library that is
009 * 1-dimensionally equidistributed across its 64-bit outputs. Has a period of 2 to the 64, and permits all states.
010 * Passes all statistical tests in PractRand up to 32TB of data. {@link DiverRNG} is a variant that is a little faster,
011 * while keeping the same other qualities, so it is currently recommended over LinnormRNG (however, the same structure
012 * is shared between LinnormRNG and {@link MizuchiRNG}, and Mizuchi has some extra features that Diver lacks). Has 64
013 * bits of state and natively outputs 64 bits at a time, changing the state with a basic linear congruential generator
014 * (it is simply {@code state = state * 3935559000370003845L + 1L}, with the large multiplier found by L'Ecuyer in a
015 * 1999 paper). Starting with that LCG's output, it xorshifts that output twice (using
016 * {@code z ^= (z >>> 23) ^ (z >>> 47);}), multiplies by a very large negative long, then returns another xorshift.
017 * Considering that some seeds don't have any anomalies in 8TB with Linnorm, the quality is probably fine except for the
018 * known potential issue that it can't return duplicate outputs until its period has cycled.
019 * {@link ThrustAltRNG} and {@link MiniMover64RNG} are faster (tied for first place), but unlike those, Linnorm can
020 * produce all long values as output; ThrustAltRNG bunches some outputs and makes producing them more likely while
021 * others can't be produced at all, while MiniMover64RNG cycles at some point before 2 to the 64 but after 2 to the 42
022 * (it doesn't produce any duplicates until then, but it also can't produce all values). Notably, this generator is
023 * faster than {@link LightRNG} with similar quality, and also faster than {@link XoRoRNG} while passing tests that
024 * XoRoRNG always or frequently fails (and fails early), such as binary matrix rank tests. It is slower than
025 * {@link DiverRNG}, which is a variant on the structure of LinnormRNG.
026 * <br>
027 * This generator is a StatefulRandomness but not a SkippingRandomness, so it can't (efficiently) have the skip() method
028 * that LightRNG has. A method could be written to run the generator's state backwards, though, as well as to get the
029 * state from an output of {@link #nextLong()}. {@link MizuchiRNG} uses the same algorithm except for the number added
030 * in the LCG state update; here this number is always 1, but in MizuchiRNG it can be any odd long. This means that any
031 * given MizuchiRNG object has two long values stored in it instead of the one here, but it allows two MizuchiRNG
032 * objects with different streams to produce different, probably-not-correlated sequences of results, even with the same
033 * seed. This property may be useful for cases where an adversary is trying to predict results in some way, though using
034 * different streams for this purpose isn't enough and should be coupled with truncation of a large part of output (see
035 * PCG-Random's techniques for this).
036 * <br>
037 * The name comes from LINear congruential generator this uses to change it state, while the rest is a NORMal
038 * SplitMix64-like generator. "Linnorm" is a Norwegian name for a kind of dragon, as well. 
039 * <br>
040 * Written June 29, 2019 by Tommy Ettinger. Thanks to M.E. O'Neill for her insights into the family of generators both
041 * this and her PCG-Random fall into, and to the team that worked on SplitMix64 for SplittableRandom in JDK 8. Chris
042 * Doty-Humphrey's work on PractRand has been invaluable. The LCG state multiplier is listed in a paper by L'Ecuyer from
043 * 1999, Tables of Linear Congruential Generators of Different Sizes and Good Lattice Structure. The other
044 * multiplier is from PCG-Random, and that's both the nothing-up-my-sleeve numbers used here. Thanks also to Sebastiano
045 * Vigna and David Blackwell for creating the incredibly fast xoroshiro128+ generator and also very fast
046 * <a href="http://xoshiro.di.unimi.it/hwd.php">HWD tool</a>; the former inspired me to make my code even faster and the
047 * latter tool seems useful so far in proving the quality of the generator (LinnormRNG passes over 100TB of HWD, and
048 * probably would pass much more if I gave it more days to run).
049 * @author Tommy Ettinger
050 */
051public final class LinnormRNG implements RandomnessSource, StatefulRandomness, Serializable {
052
053    private static final long serialVersionUID = 153186732328748834L;
054
055    private long state; /* The state can be seeded with any value. */
056
057    /**
058     * Constructor that seeds the generator with two calls to Math.random.
059     */
060    public LinnormRNG() {
061        this((long) ((Math.random() - 0.5) * 0x10000000000000L)
062                ^ (long) (((Math.random() - 0.5) * 2.0) * 0x8000000000000000L));
063    }
064
065    /**
066     * Constructor that uses the given seed as the state without changes; all long seeds are acceptable.
067     * @param seed any long, will be used exactly
068     */
069    public LinnormRNG(final long seed) {
070        state = seed;
071    }
072
073    /**
074     * Constructor that hashes seed with {@link CrossHash#hash64(CharSequence)} and uses the result as the state.
075     * @param seed any CharSequence, such as a String or StringBuilder; should probably not be null (it might work?)
076     */
077    public LinnormRNG(final CharSequence seed) {
078        state = CrossHash.hash64(seed);
079    }
080
081    @Override
082    public final int next(int bits)
083    {
084        long z = (state = state * 0x369DEA0F31A53F85L + 1L);
085        z = (z ^ z >>> 23 ^ z >>> 47) * 0xAEF17502108EF2D9L;
086        return (int)(z ^ z >>> 25) >>> (32 - bits);
087    }
088
089    /**
090     * Can return any long, positive or negative, of any size permissible in a 64-bit signed integer.
091     *
092     * @return any long, all 64 bits are random
093     */
094    @Override
095    public final long nextLong() {
096        long z = (state = state * 0x369DEA0F31A53F85L + 1L);
097        z = (z ^ z >>> 23 ^ z >>> 47) * 0xAEF17502108EF2D9L;
098        return (z ^ z >>> 25);
099    }
100
101    /**
102     * Produces a copy of this LinnormRNG that, if next() and/or nextLong() are called on this object and the
103     * copy, both will generate the same sequence of random numbers from the point copy() was called. This just need to
104     * copy the state so it isn't shared, usually, and produce a new value with the same exact state.
105     *
106     * @return a copy of this LinnormRNG
107     */
108    @Override
109    public LinnormRNG copy() {
110        return new LinnormRNG(state);
111    }
112
113    /**
114     * Can return any int, positive or negative, of any size permissible in a 32-bit signed integer.
115     *
116     * @return any int, all 32 bits are random
117     */
118    public final int nextInt() {
119        long z = (state = state * 0x369DEA0F31A53F85L + 1L);
120        z = (z ^ z >>> 23 ^ z >>> 47) * 0xAEF17502108EF2D9L;
121        return (int)(z ^ z >>> 25);
122    }
123
124    /**
125     * Exclusive on the outer bound.  The inner bound is 0.
126     * The bound can be negative, which makes this produce either a negative int or 0.
127     *
128     * @param bound the upper bound; should be positive
129     * @return a random int between 0 (inclusive) and bound (exclusive)
130     */
131    public final int nextInt(final int bound) {
132        long z = (state = state * 0x369DEA0F31A53F85L + 1L);
133        z = (z ^ z >>> 23 ^ z >>> 47) * 0xAEF17502108EF2D9L;
134        return (int)((bound * ((z ^ z >>> 25) & 0xFFFFFFFFL)) >> 32);
135    }
136
137    /**
138     * Inclusive inner, exclusive outer.
139     *
140     * @param inner the inner bound, inclusive, can be positive or negative
141     * @param outer the outer bound, exclusive, can be positive or negative, usually greater than inner
142     * @return a random int between inner (inclusive) and outer (exclusive)
143     */
144    public final int nextInt(final int inner, final int outer) {
145        return inner + nextInt(outer - inner);
146    }
147
148    /**
149     * Exclusive on bound (which may be positive or negative), with an inner bound of 0.
150     * If bound is negative this returns a negative long; if bound is positive this returns a positive long. The bound
151     * can even be 0, which will cause this to return 0L every time. This uses a biased technique to get numbers from
152     * large ranges, but the amount of bias is incredibly small (expected to be under 1/1000 if enough random ranged
153     * numbers are requested, which is about the same as an unbiased method that was also considered). It may have
154     * noticeable bias if the LinnormRNG's period is exhausted by only calls to this method, which would take months on
155     * 2018-era consumer hardware. Unlike all unbiased methods, this advances the state by an equivalent to exactly one
156     * call to {@link #nextLong()}, where rejection sampling would sometimes advance by one call, but other times by
157     * arbitrarily many more.
158     * <br>
159     * Credit for this method goes to <a href="https://oroboro.com/large-random-in-range/">Rafael Baptista's blog</a>
160     * for the original idea, and the JDK10 Math class' usage of Hacker's Delight code for the current algorithm. 
161     * This method is drastically faster than the previous implementation when the bound varies often (roughly 4x
162     * faster, possibly more). It also always gets exactly one random long, so by default it advances the state as much
163     * as {@link #nextLong()}.
164     *
165     * @param bound the outer exclusive bound; can be positive or negative
166     * @return a random long between 0 (inclusive) and bound (exclusive)
167     */
168    public long nextLong(long bound) {
169        long rand = (state = state * 0x369DEA0F31A53F85L + 1L);
170        rand = (rand ^ rand >>> 23 ^ rand >>> 47) * 0xAEF17502108EF2D9L;
171        rand ^= rand >>> 25;
172        final long randLow = rand & 0xFFFFFFFFL;
173        final long boundLow = bound & 0xFFFFFFFFL;
174        rand >>>= 32;
175        bound >>= 32;
176        final long t = rand * boundLow + (randLow * boundLow >>> 32);
177        return rand * bound + (t >> 32) + (randLow * bound + (t & 0xFFFFFFFFL) >> 32);
178    }
179    /**
180     * Inclusive inner, exclusive outer; lower and upper can be positive or negative and there's no requirement for one
181     * to be greater than or less than the other.
182     *
183     * @param lower the lower bound, inclusive, can be positive or negative
184     * @param upper the upper bound, exclusive, can be positive or negative
185     * @return a random long that may be equal to lower and will otherwise be between lower and upper
186     */
187    public final long nextLong(final long lower, final long upper) {
188        return lower + nextLong(upper - lower);
189    }
190
191    /**
192     * Gets a uniform random double in the range [0.0,1.0)
193     *
194     * @return a random double at least equal to 0.0 and less than 1.0
195     */
196    public final double nextDouble() {
197        long z = (state = state * 0x369DEA0F31A53F85L + 1L);
198        z = (z ^ z >>> 23 ^ z >>> 47) * 0xAEF17502108EF2D9L;
199        return ((z ^ z >>> 25) & 0x1FFFFFFFFFFFFFL) * 0x1p-53;
200
201    }
202
203    /**
204     * Gets a uniform random double in the range [0.0,outer) given a positive parameter outer. If outer
205     * is negative, it will be the (exclusive) lower bound and 0.0 will be the (inclusive) upper bound.
206     *
207     * @param outer the exclusive outer bound, can be negative
208     * @return a random double between 0.0 (inclusive) and outer (exclusive)
209     */
210    public final double nextDouble(final double outer) {
211        long z = (state = state * 0x369DEA0F31A53F85L + 1L);
212        z = (z ^ z >>> 23 ^ z >>> 47) * 0xAEF17502108EF2D9L;
213        return ((z ^ z >>> 25) & 0x1FFFFFFFFFFFFFL) * 0x1p-53 * outer;
214    }
215
216    /**
217     * Gets a uniform random float in the range [0.0,1.0)
218     *
219     * @return a random float at least equal to 0.0 and less than 1.0
220     */
221    public final float nextFloat() {
222        long z = (state = state * 0x369DEA0F31A53F85L + 1L);
223        return ((z ^ z >>> 23 ^ z >>> 47) * 0xAEF17502108EF2D9L >>> 40) * 0x1p-24f;
224    }
225
226    /**
227     * Gets a random value, true or false.
228     * Calls nextLong() once.
229     *
230     * @return a random true or false value.
231     */
232    public final boolean nextBoolean() {
233        long z = (state = state * 0x369DEA0F31A53F85L + 1L);
234        return ((z ^ z >>> 23 ^ z >>> 47) * 0xAEF17502108EF2D9L) < 0;
235    }
236
237    /**
238     * Given a byte array as a parameter, this will fill the array with random bytes (modifying it
239     * in-place). Calls nextLong() {@code Math.ceil(bytes.length / 8.0)} times.
240     *
241     * @param bytes a byte array that will have its contents overwritten with random bytes.
242     */
243    public final void nextBytes(final byte[] bytes) {
244        int i = bytes.length, n;
245        while (i != 0) {
246            n = Math.min(i, 8);
247            for (long bits = nextLong(); n-- != 0; bits >>>= 8) bytes[--i] = (byte) bits;
248        }
249    }
250
251    /**
252     * Sets the seed (also the current state) of this generator.
253     *
254     * @param seed the seed to use for this LinnormRNG, as if it was constructed with this seed.
255     */
256    @Override
257    public final void setState(final long seed) {
258        state = seed;
259    }
260
261    /**
262     * Gets the current state of this generator.
263     *
264     * @return the current seed of this LinnormRNG, changed once per call to nextLong()
265     */
266    @Override
267    public final long getState() {
268        return state;
269    }
270
271    @Override
272    public String toString() {
273        return "LinnormRNG with state 0x" + StringKit.hex(state) + 'L';
274    }
275
276    @Override
277    public boolean equals(Object o) {
278        if (this == o) return true;
279        if (o == null || getClass() != o.getClass()) return false;
280        return state == ((LinnormRNG) o).state;
281    }
282
283    @Override
284    public int hashCode() {
285        return (int) (state ^ (state >>> 32));
286    }
287
288    /**
289     * Static randomizing method that takes its state as a parameter; state is expected to change between calls to this.
290     * It is recommended that you use {@code LinnormRNG.determine(++state)} or {@code LinnormRNG.determine(--state)} to
291     * produce a sequence of different numbers, but you can also use {@code LinnormRNG.determine(state += 12345L)} or
292     * any odd-number increment. All longs are accepted by this method, and all longs can be produced; unlike several
293     * other classes' determine() methods, passing 0 here does not return 0.
294     * @param state any long; subsequent calls should change by an odd number, such as with {@code ++state}
295     * @return any long
296     */
297    public static long determine(long state)
298    {
299        return (state = ((state = ((state * 0x632BE59BD9B4E019L) ^ 0x9E3779B97F4A7C15L) * 0xC6BC279692B5CC83L) ^ state >>> 27) * 0xAEF17502108EF2D9L) ^ state >>> 25;
300    }
301
302    /**
303     * Static randomizing method that takes its state as a parameter and limits output to an int between 0 (inclusive)
304     * and bound (exclusive); state is expected to change between calls to this. It is recommended that you use
305     * {@code LinnormRNG.determineBounded(++state, bound)} or {@code LinnormRNG.determineBounded(--state, bound)} to
306     * produce a sequence of different numbers, but you can also use
307     * {@code LinnormRNG.determineBounded(state += 12345L, bound)} or any odd-number increment. All longs are accepted
308     * by this method, but not all ints between 0 and bound are guaranteed to be produced with equal likelihood (for any
309     * odd-number values for bound, this isn't possible for most generators). The bound can be negative.
310     * @param state any long; subsequent calls should change by an odd number, such as with {@code ++state}
311     * @param bound the outer exclusive bound, as an int
312     * @return an int between 0 (inclusive) and bound (exclusive)
313     */
314    public static int determineBounded(long state, final int bound)
315    {
316        return (int)((bound * (((state = ((state = (((state * 0x632BE59BD9B4E019L) ^ 0x9E3779B97F4A7C15L) * 0xC6BC279692B5CC83L)) ^ state >>> 27) * 0xAEF17502108EF2D9L) ^ state >>> 25) & 0x7FFFFFFFL)) >> 31);
317    }
318
319    /**
320     * Returns a random float that is deterministic based on state; if state is the same on two calls to this, this will
321     * return the same float. This is expected to be called with a changing variable, e.g. {@code determine(++state)},
322     * where the increment for state should be odd but otherwise doesn't really matter. This multiplies state by
323     * {@code 0x632BE59BD9B4E019L} within this method, so using a small increment won't be much different from using a
324     * very large one, as long as it is odd. The period is 2 to the 64 if you increment or decrement by 1, but there are
325     * less than 2 to the 30 possible floats between 0 and 1.
326     * @param state a variable that should be different every time you want a different random result;
327     *              using {@code determine(++state)} is recommended to go forwards or {@code determine(--state)} to
328     *              generate numbers in reverse order
329     * @return a pseudo-random float between 0f (inclusive) and 1f (exclusive), determined by {@code state}
330     */
331    public static float determineFloat(long state) { return ((((state = (((state * 0x632BE59BD9B4E019L) ^ 0x9E3779B97F4A7C15L) * 0xC6BC279692B5CC83L)) ^ state >>> 27) * 0xAEF17502108EF2D9L) >>> 40) * 0x1p-24f; }
332
333    /**
334     * Returns a random double that is deterministic based on state; if state is the same on two calls to this, this
335     * will return the same float. This is expected to be called with a changing variable, e.g.
336     * {@code determine(++state)}, where the increment for state should be odd but otherwise doesn't really matter. This
337     * multiplies state by {@code 0x632BE59BD9B4E019L} within this method, so using a small increment won't be much
338     * different from using a very large one, as long as it is odd. The period is 2 to the 64 if you increment or
339     * decrement by 1, but there are less than 2 to the 62 possible doubles between 0 and 1.
340     * @param state a variable that should be different every time you want a different random result;
341     *              using {@code determine(++state)} is recommended to go forwards or {@code determine(--state)} to
342     *              generate numbers in reverse order
343     * @return a pseudo-random double between 0.0 (inclusive) and 1.0 (exclusive), determined by {@code state}
344     */
345    public static double determineDouble(long state) { return (((state = ((state = (((state * 0x632BE59BD9B4E019L) ^ 0x9E3779B97F4A7C15L) * 0xC6BC279692B5CC83L)) ^ state >>> 27) * 0xAEF17502108EF2D9L) ^ state >>> 25) & 0x1FFFFFFFFFFFFFL) * 0x1p-53; }
346
347}