add square with carry

[lemire]

Stefan Kanthak wrote:

While Lehmer’s multiplicative congruential generator is fast, it does NOT pass PractRand 0.94. In does it beat the minimal standard and too big to fail Melissa O’Neill tested these generators with PractRand 0.93, where they passed.

Consider an idea from George Marsaglia instead: save and add the “carry” of the multiplication, giving a multiply-with-carry generator:

uint64_t seed = 0, carry = 0x...;
uint64_t square_with_carry() {
__uint128_t mwc = (__uint128_t) seed * 0x... + carry;
seed = mwc;
carry = square >> 64;
return seed;
}

uint64_t seed = 0, carry = 0xb5ad4eceda1ce2a9;
uint64_t square_with_carry() {
__uint128_t square = (__uint128_t) seed * seed + carry;
seed = square;
carry = square >> 64;
return seed;
}

This 64-bit generator passes the PractRand test suite (version 0.94) at least up to 4 TB! See source and source for the code generated by GCC 8.2 and clang 7.0 respectively, plus VC source for an implementation using Microsofts Visual C compiler (which does not support a 128 bit integer type). I used the latter for the PractRand test.