Looking through my 3 volumes of TAOCP I can not find the source of the short random sequence:
// Seed always in 1..0x10000
Seed = (Seed * const) % 0x10001
There was also a algorithm and possibly a MIX program to validate the const such that all 2^16 values will be returned. At least that’s what I remember. Also in the same general area was the fact that the above recursion works because (2^16)+1 is prime but alas, neither (2^32)+1 nor (2^64)-1 are primes.
FWIW, replacing const with iconst = 1/const (mod 0x10001) produces the sequence in reverse order. I.e. const*iconst%0x10001 = 1
The algorithm is a variant of linear congruential generators(LCGs) known as a Prime Modulus Multiplicative Linear Congruential Generator (PMMLCG - see Law's "Simulation Modeling and Analysis, 5e" p.400), also sometimes called a Lehmer generator. Both links I've provided give common parameterizations, the first link gives rules for choosing parameters which will give a full cycle. However, having full cycle is not a sufficient condition to have a good generator, as IBM devastatingly found out with the now infamous RANDU.
I strongly urge you to use only algorithms which have been well tested, it's amazingly hard to get this right. There are well-researched and readily available algorithms that are widely available now which are much better than any LCG.