Xoroshiro128psp (plus star plus) PRNG - 2^64 De-correlated, Jump-less Streams Each With Period 2^128-1 !!!

[Chris Rutz]

I am unable to identify any significant fault with this code or the random numbers it produces.

Still in Beta while evaluating statistics and constants.

Let me know if you disagree and/or find it useful.

// Xoroshiro128psp Beta Test Code - Copyright Christopher Rutz - Free for all uses.
// Based on xoroshiro++ by S. Vigna / D. Blackman, and discussions of it with the Parallax development team.
// Designed specifically to address all common issues and criticisms of PRNGs, and to allow for most real usage scenarios without caveat.
// 1-dimensional equidistribution per stream, with one 64-bit output occurring 2^64-1 times and the rest 2^64 times.
// All 2^64 streams are non-overlapping, nearly-perfectly de-correlated from each other, and form a super-set with near-perfect 2-dimensional equidistribution.
// Though not intended for cryptographic uses, the basic design is reasonably secure when seeded properly (note: intentionally non-trivially invertible).
// Any subset of output bits may be used for floating-point conversion.
// Pipeline optimized for Intel CPUs to meet the approximate performance speed of xoroshiro128** (within +/- 10%, using MSVC '/Ox' and GCC '-O3').

#include <stdint.h>

// Current state (including two 64-bit state values and a 64-bit stream selector), seeded by calling xoroshiro128psp_seed

// Stream 0 is usable (i.e. base xoroshiro engine), but increases worst-case all-bits de-correlation by at least 1 cycle


uint64_t s[3] = { 0 , 0 , 1 }; // 1 is the default stream, but escape-from-zero is instantaneous with ~50% set bits


// 64-bit barrel-rotation, simplifies to a single ROL assembly language op-code by most compilers
uint64_t rotl(const uint64_t x, int k) {
    return (x << k) | (x >> (64 - k));
}

// Returns next 64-bit output value from the currently selected stream
uint64_t xoroshiro128psp() {
    uint64_t s0 = s[0];
    uint64_t s1 = s[1] ^ s[2];
    uint64_t result = s0 * 5 + rotl((s0 + s1) * 5, 33); // 33 is significantly better than {31,32,34} (using a TestU01 BigCrush meta-analysis)
    s1 ^= s0;
    s[0] = rotl(s0, 24) ^ s1 ^ (s1 << 16);
    s[1] = rotl(s1, 37);
    return result;
}

// Must use this function to pass seeds to state variables and stream selector
// Stream is not intended to be modified after an initial value has already been selected
// It is recommended to use SplitMix to generate the seeds and stream selector for passing to this function
// May also be used to return hash or encryption values via seed and / or stream in counter mode: experimental
uint64_t xoroshiro128psp_seed(uint64_t seed_0, uint64_t seed_1, uint64_t stream) {
    s[0] = seed_0; s[1] = seed_1; s[2] = stream;                      // Initialize state and stream
    uint64_t result = xoroshiro128psp();                              // Pump engine once to advance state and store result
    s[0] ^= ((s[0]==seed_0) & (s[1]==seed_1));                        // Flip state LSB only if state did not advance

//  xoroshiro128psp(); xoroshiro128psp(); result = xoroshiro128psp(); // Experimental: pump for hash/encryption/counter mode: single bit change in seeds/stream de-correlates all output bits


    return result;                                                    // Return result for inspection, or (by un-commenting line above) optional hash/encryption


}