Nonuniform Random Number Generators IV

Io has a normally distributed random number generator:

bash $ io
Io 20090105
Io> Random gaussian
==> 1.0989883665373978
Io> Random gaussian(100,10)
==> 88.2743995098800269

The underlying source code uses the Box-Muller transform:

double RandomGen_gaussian(RandomGen *self, double m, double s)
{
		// http://www.taygeta.com/random/gaussian.html
	double x1, x2, w, y1, y2;

	do {
		x1 = 2.0 * RandomGen_randomDouble(self) - 1.0;
		x2 = 2.0 * RandomGen_randomDouble(self) - 1.0;
		w = x1 * x1 + x2 * x2;
	} while ( w >= 1.0 );

	w = sqrt( (-2.0 * log( w ) ) / w );
	y1 = x1 * w;
	y2 = x2 * w;

		// The following code resulted in a lot of nans being returned. The
		// following code *should* also be slower.
		/*
	double x1 = 2.0 * RandomGen_randomDouble(self) - 1.0;
	double x2 = 2.0 * RandomGen_randomDouble(self) - 1.0;
	double y1 = sqrt( - 2.0 * log(x1) ) * cos( M_PI_2 * x2 );
		*/

	return ( m + y1 * s );
}

Factor also has a normally distributed generator:

bash $ factor
( scratchpad ) USE: random
( scratchpad ) 0 1 normal-random-float

--- Data stack:
-1.651647700885478

The interesting story here is the Factor source is algorithmically identical to what is commented out in the Io source:

: normal-random-float ( mean sigma -- n )
    0.0 1.0 uniform-random-float
    0.0 1.0 uniform-random-float
    [ 2 pi * * cos ]
    [ 1.0 swap - log -2.0 * sqrt ]
    bi* * * + ;