While keeping the MT generator as the default, add new alternative random number generators as drop-in replacements.  Since MT was first introduced, PRNG technology has continued to advance.

I'm opening this feature request to be a centralized place to discuss which alternatives should be offered.

Since we already have a mostly-good-enough(tm) prng, any new generators need to be of top quality, be well researched, and offer significantly different performance characteristics than we have now (i.e. speed, cryptographic strength, simplicity, smaller state vectors, etc).

At least one of the new generators should be cryptographically strong (both to the left and to the right) while keeping reasonable speed for simulation, sampling, and gaming apps.  (The speed requirement precludes the likes of Blum Blum Shub for example.)  I believe there are several good candidates based on stream ciphers or that use block ciphers in a feedback mode.

George Marsaglia's latest random number generator KISS4691 is worth considering, though I am not sure the performance is that different from MT19937. 

Here is a link to Marsaglia's post on comp.lang.c. Marasglia passed away shortly after (Feb. 2011), and to my knowledge a paper on KISS4691 was never published:

http://www.rhinocerus.net/forum/lang-c/620168-kiss4691-potentially-top-ranked-rng.html

On my laptop, KISS4691 could produce about 110 million random numbers per second (148 millon if inlined), whereas MT19937 produced 118 million random numbers per second. Another user on comp.lang.c reported that (with my benchmark) KISS4691 was about twice as fast as MT19937 on his computer. As for quality, I have been told that MT19937 only failes a couple of obscure tests for randomness, whereas KISS4691 failes no single-seed test.

The source code I used for this test is available here:

http://folk.uio.no/sturlamo/prngtest.zip

(Requires Windows because that's what I use, sorry, might work with winelib on Linux though.)

Marsaglia has previously recommended several PRNGs that are considerably simpler and faster than MT19937. These are the ones used in the 3rd edition of "Numerical Receipes" (yes I know that not a sign of good quality). We can look at them too, with Marsaglia's comments:

https://groups.google.com/group/sci.stat.math/msg/edcb117233979602?hl=en&pli=1

https://groups.google.com/group/sci.math.num-analysis/msg/eb4ddde782b17051?hl=en&pli=1

There are also SIMD-oriented versions of MT19937, though for licensing and portability reasons they might not be suitable for Python's standard library.

High-performance PRNGs are also present in the Intel MKL and AMD ACML libraries. These could be used if Python was linked against these libraries at build-time.

Regards,
Sturla Molden

I'm posting the code for comparison of KISS4691 and MT19937. I do realize KISS4691 might not be sufficiently different from MT19937 in characteristics for Raymond Hettinger to consider it. But at least here it is for reference should it be of value.

Another (bug fix) post by Marsaglia on KISS4691:

http://www.phwinfo.com/forum/comp-lang-c/460292-ensuring-long-period-kiss4691-rng.html

Further suggestions to improve the random module:

** Object-oriented PRNG: Let it be an object which stores the random state internally, so we can create independent PRNG objects. I.e. not just one global generator.

** Generator for quasi-random Sobol sequences.  These fail most statistical tests for randomness. But for many practical uses of random numbers, e.g. simulations and numerical integration, convergence can be orders of magnitude faster (often 1000x faster convergence) but still numerically correct. That is, they are not designed to "look random", but fill the sample space as uniformly as possible.

** For cryptographically strong random numbers, os.urandom is a good back-end. It will use /dev/urandom on Linux and CryptGenRandom on Windows.  We still e.g. need a method to convert random bits from os.urandom to floats with given distributions. os.urandom is already in the standard library, so there is no reason the random module should not use it.

** Ziggurat method for generating normal, exponential and gamma deviates. This avoids most call to trancendental functions (sin, cos, log) in transforming from uniform random deviates, and is thus much faster.

** Option to return a buffer (bytearray?) with random numbers. Some times we don't need Python ints or floats, but rather the raw bytes.

** Support for more statistical distributions. Provide a much broader coverage than today. 

** Markov Chain Monte Carlo (MCMC) generator. Provide simple plug-in factories for the Gibbs sampler and the Metropolis-Hastings alorithm.

Please focus your thoughts.

I don't know if it would help you, but I wrote a C library offering a simple API and supporting various RNG (cryptographic, hardware, pseudo, ...). It reuses existing libraries like GSL, OpenSSL, glib, gcrypt, etc. It supports UNIX/BSD /dev/*random devices and the Windows CryptGen API. It has many tests and reuses external tools to test the RNG.

I added recently RAND_bytes to the Python ssl module, you may reuse it.

"On my laptop, KISS4691 could produce about 110 million random numbers per second (148 millon if inlined), whereas MT19937 produced 118 million random numbers per second."

The problem is that the Python API can only produce one number per call and a function call in Python is really slow (it creates a Python frame). If you want to speed it Python, it would be better to add methods to generate arrays to limit the overhead of Python function calls.

"The problem is that the Python API can only produce one number per call and a function call in Python is really slow (it creates a Python frame). If you want to speed it Python, it would be better to add methods to generate arrays to limit the overhead of Python function calls."


That is why I suggested an option to return a raw buffer (e.g. bytearray or memory_view). 

My experience with NumPy and MATLAB tells me that random numbers are nearly always needed in arrays with scripting languages.

But I cannot focus my thoughts...

Earlier this year I wrote Python wrappers for a number of generators:

https://github.com/douglasbagnall/riffle

They are mostly cryptographic stream ciphers from the ESTREAM[1] project, but I was also interested in dSFMT[2], which is a SIMD optimised descendant of MT19937 which runs several times faster and directly produces doubles using cunning bit tricks.

[1]http://www.ecrypt.eu.org/stream/
[2]http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/#dSFMT

Wrapped in Python, the stream ciphers ran about as fast as MT19937 on my laptop, while dSFMT took about two thirds the time to run a "random();random();random();..." test.  For a slightly more realistic test ("sum(random() for x in range(N))"), the performance levelled right off.  As expected.

The stream cipher generators have some good properties.  They generally generate random bytes using something analogous to hash('%s%s' % seed, counter), which means different seeds produce well separated streams, and to skip forward or back in the stream, you just adjust the counter.  This would allow the reinstatement of Random()'s stream-skipping function, which some people (e.g. L'Ecuyer) think is important. (incidentally, the MT people have come up with a jump-ahead algorithm for MT http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/JUMP/index.html).

Of the ciphers I tried, the chacha/salsa family and sosemanuk had the best combination of good testing and portable, reasonably fast, openly licensed C implementations.  HC128 and snow2 also perform well.

The chacha code is shorter than sosemanuk, so I would choose that.  It is used as a primitive in the BLAKE SHA3 candidate, which is a vote of confidence and an attractor of testing for the algorithm.

A bit more on the state size and period of the stream ciphers.

Chacha and Salsa use 64 bytes (512 bits) of state (vs ~2.5kB for MT19937).

Its counter is 64 bits, and its seed can be 320 bits (in cipher-speak, the seed is split between a 256 bit key and a 64 bit IV).

Each counter iteration produces 64 random bytes, or 8 doubles, so for any seed, you get a cycle of 2 ** 67, which would last in the order of 100 thousand years on current PCs.

Some of the other ciphers I looked at have smaller seeds and states, and some produce fewer bytes per iteration, but I don't think any of them will result in a cycle of smaller than 2 ** 64.

PS: Regarding the discussion of something like Random.getrandbytes(n): +1

Thanks Douglas.   Can you say what the cryptographic guarantees are for Chacha and Salsa (seeing a stream of randoms doesn't allow you to do deduce internal state, previous randoms, or future randoms)?  Is it suitably strong for gaming (dealing poker hands, lottery numbers, etc)?

I'm not sure I follow the notes on state size.  Is it 320 bits + 64 bits or is it 512 bits?  Also, I'm not sure that the smaller state is an advantage that users care about (unless they are pickling many instances of the prngs).

It's okay for jumpahead() to reappear in generators that support it, but   that method can't be a mandatory part of the Random API because it doesn't make sense for many PRNGs where a jumpahead function isn't known.

With respect to the SIMD optimizations and longlong to double operations, I'm curious to take a look at how it was done yet wonder if there is a provable, portable implementation and also wonder if it is worth it (the speed of generating a random() tends to be dwarfed by surrounding code that actually uses the result -- allocating the python object, etc).

I am no kind of crypto expert, but from what I read, there are no known attacks on chacha8 or salsa20/12 better than brute-forcing the key, and distinguishing the stream from random or deducing state would be considered an attack.  There's a summary of the ESTREAM cipher's security here:

http://cr.yp.to/streamciphers/attacks.html

-- be aware it was written by the chacha/salsa author, so may be biased.

> I'm not sure I follow the notes on state size.  Is it 320 bits + 64 bits or is it 512 bits?

Yeah. The state is contained u32[16], so the 512 is sizeof(that).  320 + 64 is the number of states I can see it getting into from the seeds and cycles.  I imagine the discrepancy is a convenience, just as the mt19937 struct uses a few more than 19937 bits.

> With respect to the SIMD optimizations and longlong to double operations, I'm curious to take a look at how it was done yet wonder if there is a provable, portable implementation and also wonder if it is worth it (the speed of generating a random() tends to be dwarfed by surrounding code that actually uses the result -- allocating the python object, etc).

I agree that it is not worth it.  However the dSFMT generator does seem quite portable and fall back to non-SIMD code (which is allegedly still faster), and its distribution is supposedly a bit better -- though not as good as WELL.

The bit magic is quite simple: if you set the top 12 bits to 0x7ff and randomise the other 52, you get a double in the range [1, 2).  So you subtract 1.  It costs one bit relative to the current method, which is equivalent to 53 bit fixed point.  They explain it reasonably well in these slides:

http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/dSFMT-slide-e.pdf

Before trying to find the best (CS)PRNG, can't we start with ssl.RAND_bytes() and ssl.RAND_pseudo_bytes()? I would be nice to use ssl.RAND_pseudo_bytes() to generate crypt.mksalt(): see issue #12858

Since this issue was opened, a few relevant changes have been made:

* firstly the random.SystemRandom API and subsequently the "secrets" module were added to provide ready access to the operating system's cryptographically secure PRNG. This addresses the "At least one of the new generators should be cryptographically strong" consideration in Raymond's original proposal
* secondly, the ensurepip module has been added and various other improvements to Python's 3rd party software distribution infrastructure have been made, lowering the barriers to retrieving alternate modelling and simulation RNGs from PyPI for folks that need or want them

Given those changes, would it make sense to close this enhancement proposal as out of date?

I agree with what Nick wrote. It became easy to install a third-party module, and we made progress on APIs to get random bytes from the system.

I'm not convinced neither that Mersenne Twister limitations are important enough to replace it.

If you want to see enhancements in the domain of RNG, I suggest to reopen a more specific issue. This issue looks like a general discussion around RNG, not really a concrete change.