# [NumPy-Tickets] [NumPy] #1450: Patch with Ziggurat method for Normal distribution

NumPy Trac numpy-tickets@scipy....
Fri Apr 9 20:18:40 CDT 2010

```#1450: Patch with Ziggurat method for Normal distribution
---------------------------------+------------------------------------------
Reporter:  ilan                 |       Owner:  somebody
Type:  enhancement          |      Status:  new
Priority:  normal               |   Milestone:
Component:  numpy.random         |     Version:  1.4.0
Keywords:  normal distribution  |
---------------------------------+------------------------------------------
I've written a patch which replaces the Box-Muller transform
in numpy/random/mtrand/randomkit.c with the faster Ziggurat method.
The patch also includes updates to doc-strings which include the
relevant references:

Doornik, J.A. (2005), "An Improved Ziggurat Method to Generate Normal
Random Samples", mimeo, Nuffield College, University of Oxford, and
www.doornik.com/research.

Marsaglia, G. and Tsang, W. W. (2000), "The Ziggurat Method for
Generating Random Variables", Journal of Statistcal Software 5,
http://www.jstatsoft.org/v05/i08/

Marsaglia, G. (1964) "Generating a variable from the tail of the
normal distribution", Technometrics 6, 101-102.

My patch is based on the support code which J. Doornik has made available
on his web page (first reference).  I've contacted him, and he has given
permission to include his code in numpy on the condition that his paper
is referenced (which the patch does).

The patch was written against http://svn.scipy.org/svn/numpy/tags/1.4.0,
and I've just checked to make sure it can still be applied to the
trunk (revision 8324).  I've tested the patch on Windows, MacOSX, Linux
and Solaris.

The speedup is a bit over 2 times, when creating a large number of
Standard
Normal distributed random numbers.  To verify that the method works and
produces good random numbers, I have:
* created a large array, using randn()
* transformed (using the cumulative distribution function) it back to
uniform distributed random numbers
* written the data to a large file
* performed the diehard tests (I used the dieharder program) on the data

Feel free to hammer the patch.  I hope this patch will make its way
into numpy soon.

--
Ticket URL: <http://projects.scipy.org/numpy/ticket/1450>
NumPy <http://projects.scipy.org/numpy>
My example project
```