How good are current implementations of rand() in C?

Question

There is a lot of literature that indicates using rand() for simulations might be harmful. A couple of them are as follows:

• https://www.sciencedirect.com/science/article/pii/S0167947398800029
• https://channel9.msdn.com/Events/GoingNative/2013/rand-Considered-Harmful

A common recommendation has been to use the ran1/ran2/ran3 functions as defined in Numerical Recipes.

But these seem pretty old criticisms of the rand function. So my questions is, how do its current implementations (GCC/clang) measure up for the purpose of being a good RNG in, say Monte Carlo simulations?

This is the current implementation of the function used by rand() in glibc: (rand calls random, which calls random_r)

• https://sourceware.org/git/?p=glibc.git;a=blob;f=stdlib/rand.c
• https://sourceware.org/git/?p=glibc.git;a=blob;f=stdlib/random.c
• https://sourceware.org/git/?p=glibc.git;a=blob;f=stdlib/random_r.c

Answer 1

But these seem pretty old criticisms of the rand function

This may be a nitpick, but I want to point out what I think is a flaw in this logic. Compilers are often extremely conservative about changing program behaviour, even when that behaviour (foolishly) depends on implementation details.

This may or may not be true for the big compilers you’re familiar with, but it is conceivable that a compiler would keep a broken legacy RNG for compatibility reasons.

The following program on my macbook:

#include <cstdlib>
#include <iostream>
using namespace std;

int main() {
  srand(1);
  for (int i = 0; i < 20; ++i) {
    cout << (i > 0 ? ", " : "") << rand();
  }
  cout << endl;
  return 0;
}

produces the output

16807, 282475249, 1622650073, 984943658, 1144108930, 470211272, 101027544, 1457850878, 1458777923, 2007237709, 823564440, 1115438165, 1784484492, 74243042, 114807987, 1137522503, 1441282327, 16531729, 823378840, 143542612

which you can look up in OEIS: https://oeis.org/A096550, it is the sequence $$16807^n\bmod(2^{31}-1).$$

How this does under the SmallCrush RNG test (in the TestU01 library):

========= Summary results of SmallCrush =========

 Version:          TestU01 1.2.3
 Generator:        ulcg_CreateLCG
 Number of statistics:  15
 Total CPU time:   00:00:08.05
 The following tests gave p-values outside [0.001, 0.9990]:
 (eps  means a value < 1.0e-300):
 (eps1 means a value < 1.0e-15):

       Test                          p-value
 ----------------------------------------------
  1  BirthdaySpacings                 eps
  2  Collision                        eps
  6  MaxOft                           eps
 ----------------------------------------------
 All other tests were passed

As far as I know there is simply no good reason to assume that rand() might be okay. If you assume otherwise, you'll usually end up writing some extremely non-portable code that will break under different compilers on different systems.

Answer 2

Though the question is about C, I still consider the following information regarding C++11 relevant. In C++11, <random> was added that is both designed to fix many of C rand() implementation flaws, as well as to provide additional capabilities. So, if you have an option of using C++11, I would point your attention to <random> and its benefits.

In most cases, <random> function family would provide sufficient quality random numbers (as opposed to rand(), which is mentioned in the materials you linked in your question. The critique to rand() is still relevant (as ANSI standard did not change); however, certainly exaggerated.

My advice would be to first assess the need in higher quality random numbers and go step-by-step. There is a high chance, that for particular purposes <random> or even rand() will provide decent enough results. Otherwise (unlikely), you might have to go to special libraries or even hardware generators.

Answer 3

There is a finite number of correlation metrics relevant to your problem which you ideally should know before trying to implement MC!

Regarding rand(): I have written (parallel) MC methods which could easily draw $>10^9$ PRNs at which point you cycle the sequence of PRNs provided by rand(), whose period is $2^{32}-1\approx 4.3\cdot 10^9$ if I understand the comments in the code you link to correctly. Now you either have to analyse whether potential correlations will impact your problem, or you can choose a PRNG with a longer period, thereby avoiding the problem (or at least pushing it off to larger problem sizes).

Especially in parallel environments (the common way to run a MC simulation) the usage of rand() harbors the additional problem that multiple processes may end up initializing the PRNG to either start at the same position, or at positions at which their sequences will eventually overlap. The latter will result in two processes using the same PRNs for (supposedly) independent simulations. Other PRNGs, like the Mersenne Twister engine used in C++11 suffer the same problem, but to a lesser degree since their (limited) parameter sets produce different sequences, so as long as you do not initialize them with the same seed they will yield non-overlapping sequences.

A C code I have worked with used rand() for a parallel MC method used to approximate an inverse of a (very large) matrix. I have empirically determined that replacing it by a TRNG ( https://www.numbercrunch.de/trng/ ) - had to write a C-wrapper for it, though - measurably reduced the errors of the final result. Especially for large matrices (e.g., those that needed a lot of PRNs). Note here: Only because you are limited to work with C does not mean you can not use C++ libraries, you just have to write wrappers for the functions you need.