comp.lang.c FAQ list · Question 13.20

Q: How can I generate random numbers with a normal or Gaussian distribution?

A: There are a number of ways of doing this.

1. Exploit the Central Limit Theorem (``law of large numbers'') and add up several uniformly-distributed random numbers:

   #include <stdlib.h>
   #include <math.h>
   
   #define NSUM 25
   
   double gaussrand()
   {
   	double x = 0;
   	int i;
   	for(i = 0; i < NSUM; i++)
   		x += (double)rand() / RAND_MAX;
   
   	x -= NSUM / 2.0;
   	x /= sqrt(NSUM / 12.0);
   
   	return x;
   }
   

   (Don't overlook the sqrt(NSUM / 12.) correction, though it's easy to do so accidentally, especially when NSUM is 12.)
2. Use a method described by Abramowitz and Stegun:

   #include <stdlib.h>
   #include <math.h>
   
   #define PI 3.141592654
   
   double gaussrand()
   {
   	static double U, V;
   	static int phase = 0;
   	double Z;
   
   	if(phase == 0) {
   		U = (rand() + 1.) / (RAND_MAX + 2.);
   		V = rand() / (RAND_MAX + 1.);
   		Z = sqrt(-2 * log(U)) * sin(2 * PI * V);
   	} else
   		Z = sqrt(-2 * log(U)) * cos(2 * PI * V);
   
   	phase = 1 - phase;
   
   	return Z;
   }
   

3. Use a method discussed in Knuth and due originally to Marsaglia:

   #include <stdlib.h>
   #include <math.h>
   
   double gaussrand()
   {
   	static double V1, V2, S;
   	static int phase = 0;
   	double X;
   
   	if(phase == 0) {
   		do {
   			double U1 = (double)rand() / RAND_MAX;
   			double U2 = (double)rand() / RAND_MAX;
   
   			V1 = 2 * U1 - 1;
   			V2 = 2 * U2 - 1;
   			S = V1 * V1 + V2 * V2;
   			} while(S >= 1 || S == 0);
   
   		X = V1 * sqrt(-2 * log(S) / S);
   	} else
   		X = V2 * sqrt(-2 * log(S) / S);
   
   	phase = 1 - phase;
   
   	return X;
   }
   

Additional links

References: Knuth Sec. 3.4.1 p. 117
Box and Muller, ``A Note on the Generation of Random Normal Deviates''
Marsaglia and Bray, ``A Convenient Method for Generating Normal Variables''
Abramowitz and Stegun, Handbook of Mathematical Functions
Press et al., Numerical Recipes in C Sec. 7.2 pp. 288-290

Hosted by