Cannot reproduce claimed MRG32k3a failures

[vigna]

I happened to read your comment about MRG32k3a in the PRNG section. I'm not able to reproduce your failures.

rng=MRG32k3a, seed=0xdb3b0bab
length= 2 terabytes (2^41 bytes), time= 100581 seconds
  no anomalies in 2194 test result(s)

rng=MRG32k3a, seed=0xdb3b0bab
length= 4 terabytes (2^42 bytes), time= 194231 seconds
  no anomalies in 2255 test result(s)

(It is possible that faults may appear later—I'll let this go for a couple of weeks.)

I see two issues in your test code:

• The generator, by design, will never generate zero. It's a feature—not a bug, as it is desirable in a number of scientific environments. PractRand expects a generator that is uniform on 32 bits. A missing zero very 2³² elements will skew its counters and produce false positives.

• The generator, by design, emits a double in (0 . . 1) by dividing part of its internal state by 4294967088. You are testing the raw state. That was never intended to be used as a pseudorandom number.

If you wanna test MRG32k3a with PractRand, you must modify the difference normalization so that it can also generate zero, and test the double result multiplied by 2³² (that would be the part of the state you're testing now divided by 4294967088, and multiplied by 2³²). Otherwise you will get false positives from a wrong assumption of PractRand and a mistake in the use of the generator.

You might consider using my branchless version, which is about twice as fast. In that case the only modification to obtain also zero is to replace ((r - 1) >> 63) with (r >> 63).

[vigna]

It definitely looks like the problem is in the test code:

rng=MRG32k3a, seed=0xdb3b0bab
length= 8 terabytes (2^43 bytes), time= 380875 seconds
  no anomalies in 2313 test result(s)

rng=MRG32k3a, seed=0xdb3b0bab
length= 16 terabytes (2^44 bytes), time= 753934 seconds
  no anomalies in 2366 test result(s)

[alvoskov]

@vigna Yes, I agree that the problem was with the renormalization. After addition of such renormalization to the [0; 2**32 - 1] interval that PractRand failure disappeared. I've added this code to the release (version 0.48), see the https://github.com/alvoskov/SmokeRand/blob/main/generators/mrg32k3a.c file. PractRand stopped to fail at 4 TiB.

The conversion/renormalization is made without floating point numbers at all, the underlying theory together with Python code is given inside the mrg32k3a Doxygen comments.

static inline uint64_t get_bits_raw(Mrg32k3aState *obj)
{
    const uint32_t r = Mrg32k3aState_get_u32_raw(obj) - 1;
    const uint32_t x = (uint32_t) (r + ((1753219158U * (uint64_t) r) >> 55));
    return x;
}