Can Fisher-Yates shuffle produce all playing card permutations?

Question

I'm using the standard Fisher-Yates algorithm to randomly shuffle a deck of cards in an array. However, I'm unsure if this will actually produce a true distribution of all possible permutations of a real-world shuffled deck of cards.

V8's Math.random only has 128-bits of internal state. Since there are 52 cards in a deck, 52 factorial would require 226-bits of internal state to generate all possible permutations.

However, I'm unsure if this applies when using Fisher-Yates since you aren't actually generating each possible but just getting one position randomly out of 52.

function shuffle(array) {
  var m = array.length, t, i;

  while (m) {
    i = Math.floor(Math.random() * m--);
    t = array[m];
    array[m] = array[i];
    array[i] = t;
  }

  return array;
}

Answer 1

In general, if a pseudorandom number generator admits fewer than 52 factorial different seeds, then there are some permutations that particular PRNG can't choose when it shuffles a 52-item list, and Fisher-Yates can't change that. (The set of permutations a particular PRNG can choose can be different from the set of permutations another PRNG can choose, even if both PRNGs are initialized with the same seed.) See also this question.

Note that although the Math.random algorithm used by V8 admits any of about 2^128 seeds at the time of this writing, no particular random number algorithm is mandated by the ECMAScript specification of Math.random, which states only that that method uses an "implementation-dependent algorithm or strategy" to generate random numbers (see ECMAScript sec. 20.2.2.27).

---

A PRNG's period can be extended with the Bays-Durham shuffle, which effectively increases that PRNG's state length (see Severin Pappadeux's answer). However, if you merely initialize the Bays-Durham table entries with outputs of the PRNG (rather than use the seed to initialize those entries), it will still be the case that that particular PRNG (which includes the way in which it initializes those entries and selects those table entries based on the random numbers it generates) can't choose more permutations than the number of possible seeds to initialize its original state, because there would be only one way to initialize the Bays-Durham entries for a given seed — unless, of course, the PRNG actually shuffles an exorbitant amount of lists, so many that it generates more random numbers without cycling than it otherwise would without the Bays-Durham shuffle.

For example, if the PRNG is 128 bits long, there are only 2^128 possible seeds, so there are only 2^128 ways to initialize the Bays-Durham shuffle, one for each seed, unless a seed longer than 128 bits extends to the Bays-Durham table entries and not just the PRNG's original state. (This is not to imply that the set of permutations that PRNG can choose is always the same no matter how it selects table entries in the Bays-Durham shuffle.)

EDIT (Aug. 7): Clarifications.

EDIT (Jan. 7, 2020): Edited.

Answer 2

You are right. With 128 bits of starting state, you can only generate at most 2¹²⁸ different permutations. It doesn't matter how often you use this state (call Math.random()), the PRNG is deterministic after all.

Where the number of calls to Math.random() actually matter is when

• each call would draw some more entropy (e.g. from hardware random) into the system, instead of relying on the internal state that is initialised only once
• the entropy of a single call result is so low that you don't use the entire internal state over the run of the algorithm

Answer 3

Well, you definitely need RNG with 226bits period for all permutation to be covered, @PeterO answer is correct in this regard. But you could extend period using Bays-Durham shuffle, paying by effectively extending state of RNG. There is an estimate of the period of the B-D shuffled RNG and it is

P = sqrt(Pi * N! / (2*O))

where Pi=3.1415..., N is B-D table size, O is period of the original generator. If you take log₂ of the whole expression, and use Stirling formula for factorial, and assume P=2²²⁶ and O=2¹²⁸, you could get estimate for N, size of the table in B-D algorithm. From back-of-the envelope calculation N=64 would be enough to get all your permutations.

UPDATE

Ok, here is an example implementation of RNG extended with B-D shuffle. First, I implemented in Javascript Xorshift128+, using BigInt, which is apparently default RNG in V8 engine as well. Compared with C++ one, they produced identical output for first couple of dozen calls. 128bits seed as two 64bits words. Windows 10 x64, NodeJS 12.7.

const WIDTH = 2n ** 64n;
const MASK  = WIDTH - 1n; // to keep things as 64bit values

class XorShift128Plus { // as described in https://v8.dev/blog/math-random
    _state0 = 0n;
    _state1 = 0n;

    constructor(seed0, seed1) { // 128bit seed as 2 64bit values 
        this._state0 = BigInt(seed0) & MASK;
        this._state1 = BigInt(seed1) & MASK;
        if (this._state0 <= 0n)
            throw new Error('seed 0 non-positive');
        if (this._state1 <= 0n)
            throw new Error('seed 1 non-positive');
    }

    next() {
        let s1 = this._state0;
        let s0 = this._state1;
        this._state0 = s0;
        s1  = ((s1 << 23n) ^ s1 ) & MASK;
        s1 ^= (s1 >> 17n);
        s1 ^= s0;
        s1 ^= (s0 >> 26n);
        this._state1 = s1;
        return (this._state0 + this._state1) & MASK; // modulo WIDTH
    }
}

Ok, then on top of XorShift128+ I've implemented B-D shuffle, with table of size 4. For your purpose you'll need table more than 84 entries, and power of two table is much easier to deal with, so let's say 128 entries table (7bit index) shall be good enough. Anyway, even with 4 entries table and 2bit index we need to know which bits to pick to form index. In original paper B-D discussed picking them from the back of rv as well as from front of rv etc. Here is where B-D shuffle needs another seed value - telling algorithm to pick say, bits from position 2 and 6.

class B_D_XSP {
    _xsprng;
    _seedBD = 0n;

    _pos0 = 0n;
    _pos1 = 0n;

    _t; // B-D table, 4 entries
    _Z = 0n;

    constructor(seed0, seed1, seed2) { // note third seed for the B-D shuffle
        this._xsprng = new XorShift128Plus(seed0, seed1);

        this._seedBD = BigInt(seed2) & MASK;
        if (this._seedBD <= 0n)
            throw new Error('B-D seed non-positive');

        this._pos0 = findPosition(this._seedBD);                         // first  non-zero bit position
        this._pos1 = findPosition(this._seedBD & (~(1n << this._pos0))); // second non-zero bit position

        // filling up table and B-D shuffler
        this._t = new Array(this._xsprng.next(), this._xsprng.next(), this._xsprng.next(), this._xsprng.next());
        this._Z = this._xsprng.next();
    }

    index(rv) { // bit at first position plus 2*bit at second position
        let idx = ((rv >> this._pos0) & 1n) + (((rv >> this._pos1) & 1n) << 1n);
        return idx;
    }

    next() {
        let retval = this._Z;
        let j = this.index(this._Z);
        this._Z = this._t[j];
        this._t[j] = this._xsprng.next();
        return retval;
    }
}

Use example is as follow.

let rng = new B_D_XSP(1, 2, 4+64); // bits at second and sixth position to make index

console.log(rng._pos0.toString(10));
console.log(rng._pos1.toString(10));

console.log(rng.next());
console.log(rng.next());
console.log(rng.next());

Obviously, third seed value of say 8+128 would produce different permutation from what is shown in the example, you could play with it.

Last step would be to make 226bit random value by calling several (3 of 4) times B-D shuffled rng and combine 64bit values (and potential carry over) to make 226 random bits and then convert them to the deck shuffle.