On 2009-12-19 09:14 AM, Carl Johan Rehn wrote:
On Dec 19, 2:49 pm, sturlamolden<sturlamol...@yahoo.no> wrote:
On 19 Des, 11:05, Carl Johan Rehn<car...@gmail.com> wrote:
I plan to port a Monte Carlo engine from Matlab to Python. However,
when I timed randn(N1, N2) in Python and compared it with Matlab's
randn, Matlab came out as a clear winner with a speedup of 3-4 times.
This was truly disappointing. I ran tthis test on a Win32 machine and
without the Atlas library.
This is due to the algorithm. Matlab is using Marsaglia's ziggurat
method. Is is the fastest there is for normal and gamma random
variates. NumPy uses the Mersenne Twister to produce uniform random
deviates, and then applies trancendental functions to transform to the
normal distribution. Marsaglia's C code for ziggurat is freely
available, so you can compile it yourself and call from ctypes, Cython
or f2py.
The PRNG does not use BLAS/ATLAS.
Thank you, this was very informative. I know about the Mersenne
Twister, but had no idea about Marsaglia's ziggurat method. I will
definitely try f2py or cython.
It's worth noting that the ziggurat method can be implemented incorrectly, and
requires some testing before I will accept such in numpy. That's why we don't
use the ziggurat method currently. C.f.
http://www.doornik.com/research/ziggurat.pdf
Requests for the ziggurat method come up occasionally on the numpy list, but so
far no one has shown code or test results. Charles Harris and Bruce Carneal seem
to have gotten closest. You can search numpy-discussion's archive for their
email addresses and previous threads. Additionally, you should ask future numpy
questions on the numpy-discussion mailing list.
http://www.scipy.org/Mailing_Lists
--
Robert Kern
"I have come to believe that the whole world is an enigma, a harmless enigma
that is made terrible by our own mad attempt to interpret it as though it had
an underlying truth."
-- Umberto Eco
--
http://mail.python.org/mailman/listinfo/python-list