Take a look here: https://bashtage.github.io/ng-numpy-randomstate/doc/index.html On Mon, Apr 3, 2017 at 9:45 AM Jaime Fernández del Río <jaime.frio@gmail.com> wrote:
On Mon, Apr 3, 2017 at 3:20 PM, Pierre Haessig <pierre.haessig@crans.org> wrote:
Hello, Le 30/03/2017 à 13:31, Pierre Haessig a écrit :
[....]
But how come Julia is 4-5x faster since Numpy uses C implementation for the entire process ? (Mersenne Twister -> uniform double -> Box-Muller transform to get a Gaussian https://github.com/numpy/numpy/blob/master/numpy/random/mtrand/randomkit.c). Also I noticed that Julia uses a different algorithm (Ziggurat Method from Marsaglia and Tsang , https://github.com/JuliaLang/julia/blob/master/base/random.jl#L700) but this doesn't explain the difference for uniform rng.
Any ideas?
This <https://github.com/JuliaLang/julia/blob/7fb758a275a0b4cf0e3f4cbf482c065cb32f...> says that Julia uses this library <http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/#dSFMT>, which is different from the home brewed version of the Mersenne twister in NumPy. The second link I posted claims their speed comes from generating double precision numbers directly, rather than generating random bytes that have to be converted to doubles, as is the case of NumPy through this magical incantation <https://github.com/numpy/numpy/blob/master/numpy/random/mtrand/randomkit.c#L...>. They also throw the SIMD acronym around, which likely means their random number generation is parallelized.
My guess is that most of the speed-up comes from the SIMD parallelization: the Mersenne algorithm does a lot of work <https://github.com/numpy/numpy/blob/master/numpy/random/mtrand/randomkit.c#L...> to produce 32 random bits, so that likely dominates over a couple of arithmetic operations, even if divisions are involved.
Jaime
Do you think Stackoverflow would be a better place for my question?
best,
Pierre
_______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@python.org https://mail.python.org/mailman/listinfo/numpy-discussion
-- (\__/) ( O.o) ( > <) Este es Conejo. Copia a Conejo en tu firma y ayúdale en sus planes de dominación mundial. _______________________________________________ NumPy-Discussion mailing list NumPy-Discussion@python.org https://mail.python.org/mailman/listinfo/numpy-discussion