C: Random number generation - what (if anything) is wrong with this

Question

For simple modeling in C, I need to create exponential random variables. I remember reading somewhere (but I can’t find it now, and I don’t remember why) that using the rand () function to generate random integers in a fixed range would create unevenly distributed integers. Because of this, I am wondering if this code has a similar problem:

//generate u ~ U[0,1] u = ( (double)rand() / ((double)(RAND_MAX)); //inverse of exponential CDF to get exponential random variable expon = -log(1-u) * mean; 

Thanks!

Answer 1

The problem with random numbers in a fixed range is that many people do this for numbers from 100 to 200, for example:

 100 + rand() % 100 

This is not uniform. But by doing this, it (or close enough to uniform, at least):

 u = 100 + 100 * ((double)rand() / ((double)(RAND_MAX)); 

Ever since you do, you must be safe.

Answer 2

In theory, at least rand () should give you a discrete uniform distribution from 0 to RAND_MAX ... in practice, it has some undesirable properties, such as a short period, so whether this is useful depends on using it.

Answer 3

RAND_MAX is usually 32k, while LCG rand () generates pseudo-random 32-bit numbers. Thus, the lack of uniformity, as well as low frequency, as a rule, go unnoticed.

If you need high-quality pseudo-random numbers, you can try George Marsaglia CMWC4096 ("Extra Carry Multiplication"). This is probably the best pseudo-random number generator, with extreme frequency and uniform distribution (you just need to choose good seeds for this). In addition, it flashes quickly (not as fast as the LCG, but about twice as fast as the Mersenne Twister).

Answer 4

Yes and no. The problem you are thinking is when you pin the output from rand() to a range that is less than RAND_MAX (i.e. there are fewer possible outputs than inputs).

In your case, you (usually) change this: you take a fairly small number of bits created by a random number generator and distribute them among what will usually be a large number of bits in your double mantissa. This means that there are usually several bit patterns in double (and therefore specific double values) that can never happen. For most people, this is not a problem.

As for "normal", it is always possible that you have a 64-bit random number generator, where the double usually has a 53-bit mantissa. In this case, you may have the same problem as when compressing the range with integers.

Answer 5

No, your algorithm will work; it uses a module function that makes things imperfect.
One of the problems is that since it is quantized, from time to time it will generate exactly RAND_MAX, and you will query log(1-1) . I would recommend at least (rand() + 0.5)/(RAND_MAX+1) , if not the best source, e.g. drand48() .

There is a much faster way to calculate the required numbers, for example. Ziggurat's algorithm .