diff options
Diffstat (limited to 'doc/source')
26 files changed, 1661 insertions, 6 deletions
diff --git a/doc/source/conf.py b/doc/source/conf.py index 072a3b44e..dec8fff05 100644 --- a/doc/source/conf.py +++ b/doc/source/conf.py @@ -19,11 +19,19 @@ needs_sphinx = '1.0' sys.path.insert(0, os.path.abspath('../sphinxext')) -extensions = ['sphinx.ext.autodoc', 'numpydoc', - 'sphinx.ext.intersphinx', 'sphinx.ext.coverage', - 'sphinx.ext.doctest', 'sphinx.ext.autosummary', - 'sphinx.ext.graphviz', 'sphinx.ext.ifconfig', - 'matplotlib.sphinxext.plot_directive'] +extensions = [ + 'sphinx.ext.autodoc', + 'numpydoc', + 'sphinx.ext.intersphinx', + 'sphinx.ext.coverage', + 'sphinx.ext.doctest', + 'sphinx.ext.autosummary', + 'sphinx.ext.graphviz', + 'sphinx.ext.ifconfig', + 'matplotlib.sphinxext.plot_directive', + 'IPython.sphinxext.ipython_console_highlighting', + 'IPython.sphinxext.ipython_directive', +] if sphinx.__version__ >= "1.4": extensions.append('sphinx.ext.imgmath') @@ -234,7 +242,7 @@ numpydoc_use_plots = True # ----------------------------------------------------------------------------- import glob -autosummary_generate = glob.glob("reference/*.rst") +autosummary_generate = True # ----------------------------------------------------------------------------- # Coverage checker @@ -355,3 +363,8 @@ def linkcode_resolve(domain, info): else: return "https://github.com/numpy/numpy/blob/v%s/numpy/%s%s" % ( numpy.__version__, fn, linespec) + +doctest_global_setup = ''' +import numpy as np +from numpy.random import randomgen +''' diff --git a/doc/source/reference/randomgen/brng/dsfmt.rst b/doc/source/reference/randomgen/brng/dsfmt.rst new file mode 100644 index 000000000..f9de48d61 --- /dev/null +++ b/doc/source/reference/randomgen/brng/dsfmt.rst @@ -0,0 +1,43 @@ +Double SIMD Mersenne Twister (dSFMT) +------------------------------------ + +.. module:: numpy.random.randomgen.dsfmt + +.. currentmodule:: numpy.random.randomgen.dsfmt + + +.. autoclass:: DSFMT + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~DSFMT.seed + ~DSFMT.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~DSFMT.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~DSFMT.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~DSFMT.cffi + ~DSFMT.ctypes + + diff --git a/doc/source/reference/randomgen/brng/index.rst b/doc/source/reference/randomgen/brng/index.rst new file mode 100644 index 000000000..aceecc792 --- /dev/null +++ b/doc/source/reference/randomgen/brng/index.rst @@ -0,0 +1,40 @@ +Basic Random Number Generators +------------------------------ + +The random values produced by :class:`~randomgen.generator.RandomGenerator` +are produced by a basic RNG. These basic RNGs do not directly provide +random numbers and only contains methods used for seeding, getting or +setting the state, jumping or advancing the state, and for accessing +low-level wrappers for consumption by code that can efficiently +access the functions provided, e.g., `numba <https://numba.pydata.org>`_. + +Stable RNGs +=========== +These RNGs will be included in future releases. + + +.. toctree:: + :maxdepth: 1 + + DSFMT <dsfmt> + MT19937 <mt19937> + PCG64 <pcg64> + Philox <philox> + ThreeFry <threefry> + XoroShiro128+ <xoroshiro128> + Xorshift1024*φ <xorshift1024> + Xoshiro256** <xoshiro256starstar> + Xoshiro512** <xoshiro512starstar> + + +Experimental RNGs +================= + +These RNGs are currently included for testing but are may not be +permanent. + +.. toctree:: + :maxdepth: 1 + + PCG32 <pcg32> + ThreeFry32 <threefry32> diff --git a/doc/source/reference/randomgen/brng/mt19937.rst b/doc/source/reference/randomgen/brng/mt19937.rst new file mode 100644 index 000000000..7739e16ce --- /dev/null +++ b/doc/source/reference/randomgen/brng/mt19937.rst @@ -0,0 +1,42 @@ +Mersenne Twister (MT19937) +-------------------------- + +.. module:: numpy.random.randomgen.mt19937 + +.. currentmodule:: numpy.random.randomgen.mt19937 + +.. autoclass:: MT19937 + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~MT19937.seed + ~MT19937.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~MT19937.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~MT19937.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~MT19937.cffi + ~MT19937.ctypes + + diff --git a/doc/source/reference/randomgen/brng/pcg32.rst b/doc/source/reference/randomgen/brng/pcg32.rst new file mode 100644 index 000000000..aaf3929e8 --- /dev/null +++ b/doc/source/reference/randomgen/brng/pcg32.rst @@ -0,0 +1,43 @@ +Parallel Congruent Generator (32-bit, PCG32) +-------------------------------------------- + +.. module:: numpy.random.randomgen.pcg32 + +.. currentmodule:: numpy.random.randomgen.pcg32 + +.. autoclass:: PCG32 + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~PCG32.seed + ~PCG32.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~PCG32.advance + ~PCG32.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~PCG32.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~PCG32.cffi + ~PCG32.ctypes + + diff --git a/doc/source/reference/randomgen/brng/pcg64.rst b/doc/source/reference/randomgen/brng/pcg64.rst new file mode 100644 index 000000000..94e73e491 --- /dev/null +++ b/doc/source/reference/randomgen/brng/pcg64.rst @@ -0,0 +1,43 @@ +Parallel Congruent Generator (64-bit, PCG64) +-------------------------------------------- + +.. module:: numpy.random.randomgen.pcg64 + +.. currentmodule:: numpy.random.randomgen.pcg64 + +.. autoclass:: PCG64 + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~PCG64.seed + ~PCG64.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~PCG64.advance + ~PCG64.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~PCG64.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~PCG64.cffi + ~PCG64.ctypes + + diff --git a/doc/source/reference/randomgen/brng/philox.rst b/doc/source/reference/randomgen/brng/philox.rst new file mode 100644 index 000000000..091c4d3e0 --- /dev/null +++ b/doc/source/reference/randomgen/brng/philox.rst @@ -0,0 +1,43 @@ +Philox Counter-based RNG +------------------------ + +.. module:: numpy.random.randomgen.philox + +.. currentmodule:: numpy.random.randomgen.philox + +.. autoclass:: Philox + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~Philox.seed + ~Philox.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~Philox.advance + ~Philox.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~Philox.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~Philox.cffi + ~Philox.ctypes + + diff --git a/doc/source/reference/randomgen/brng/threefry.rst b/doc/source/reference/randomgen/brng/threefry.rst new file mode 100644 index 000000000..4f5c56bae --- /dev/null +++ b/doc/source/reference/randomgen/brng/threefry.rst @@ -0,0 +1,43 @@ +ThreeFry Counter-based RNG +-------------------------- + +.. module:: numpy.random.randomgen.threefry + +.. currentmodule:: numpy.random.randomgen.threefry + +.. autoclass:: ThreeFry + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~ThreeFry.seed + ~ThreeFry.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~ThreeFry.advance + ~ThreeFry.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~ThreeFry.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~ThreeFry.cffi + ~ThreeFry.ctypes + + diff --git a/doc/source/reference/randomgen/brng/threefry32.rst b/doc/source/reference/randomgen/brng/threefry32.rst new file mode 100644 index 000000000..bd85db4a7 --- /dev/null +++ b/doc/source/reference/randomgen/brng/threefry32.rst @@ -0,0 +1,43 @@ +ThreeFry32 Counter-based RNG +---------------------------- + +.. module:: numpy.random.randomgen.threefry32 + +.. currentmodule:: numpy.random.randomgen.threefry32 + +.. autoclass:: ThreeFry32 + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~ThreeFry32.seed + ~ThreeFry32.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~ThreeFry32.advance + ~ThreeFry32.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~ThreeFry32.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~ThreeFry32.cffi + ~ThreeFry32.ctypes + + diff --git a/doc/source/reference/randomgen/brng/xoroshiro128.rst b/doc/source/reference/randomgen/brng/xoroshiro128.rst new file mode 100644 index 000000000..6796c4457 --- /dev/null +++ b/doc/source/reference/randomgen/brng/xoroshiro128.rst @@ -0,0 +1,42 @@ +Xoroshiro128+ +------------- + +.. module:: numpy.random.randomgen.xoroshiro128 + +.. currentmodule:: numpy.random.randomgen.xoroshiro128 + +.. autoclass:: Xoroshiro128 + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~Xoroshiro128.seed + ~Xoroshiro128.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~Xoroshiro128.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~Xoroshiro128.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~Xoroshiro128.cffi + ~Xoroshiro128.ctypes + + diff --git a/doc/source/reference/randomgen/brng/xorshift1024.rst b/doc/source/reference/randomgen/brng/xorshift1024.rst new file mode 100644 index 000000000..64df7e050 --- /dev/null +++ b/doc/source/reference/randomgen/brng/xorshift1024.rst @@ -0,0 +1,42 @@ +Xorshift1024*φ +-------------- + +.. module:: numpy.random.randomgen.xorshift1024 + +.. currentmodule:: numpy.random.randomgen.xorshift1024 + +.. autoclass:: Xorshift1024 + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~Xorshift1024.seed + ~Xorshift1024.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~Xorshift1024.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~Xorshift1024.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~Xorshift1024.cffi + ~Xorshift1024.ctypes + + diff --git a/doc/source/reference/randomgen/brng/xoshiro256starstar.rst b/doc/source/reference/randomgen/brng/xoshiro256starstar.rst new file mode 100644 index 000000000..7603e6f1b --- /dev/null +++ b/doc/source/reference/randomgen/brng/xoshiro256starstar.rst @@ -0,0 +1,42 @@ +Xoshiro256** +------------ + +.. module:: numpy.random.randomgen.xoshiro256starstar + +.. currentmodule:: numpy.random.randomgen.xoshiro256starstar + +.. autoclass:: Xoshiro256StarStar + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~Xoshiro256StarStar.seed + ~Xoshiro256StarStar.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~Xoshiro256StarStar.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~Xoshiro256StarStar.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~Xoshiro256StarStar.cffi + ~Xoshiro256StarStar.ctypes + + diff --git a/doc/source/reference/randomgen/brng/xoshiro512starstar.rst b/doc/source/reference/randomgen/brng/xoshiro512starstar.rst new file mode 100644 index 000000000..64f95f750 --- /dev/null +++ b/doc/source/reference/randomgen/brng/xoshiro512starstar.rst @@ -0,0 +1,42 @@ +Xoshiro512** +------------ + +.. module:: numpy.random.randomgen.xoshiro512starstar + +.. currentmodule:: numpy.random.randomgen.xoshiro512starstar + +.. autoclass:: Xoshiro512StarStar + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~Xoshiro512StarStar.seed + ~Xoshiro512StarStar.state + +Parallel generation +=================== +.. autosummary:: + :toctree: generated/ + + ~Xoshiro512StarStar.jump + +Random Generator +================ +.. autosummary:: + :toctree: generated/ + + ~Xoshiro512StarStar.generator + +Extending +========= +.. autosummary:: + :toctree: generated/ + + ~Xoshiro512StarStar.cffi + ~Xoshiro512StarStar.ctypes + + diff --git a/doc/source/reference/randomgen/change-log.rst b/doc/source/reference/randomgen/change-log.rst new file mode 100644 index 000000000..f791c8f54 --- /dev/null +++ b/doc/source/reference/randomgen/change-log.rst @@ -0,0 +1,39 @@ +Change Log +---------- +v1.16.1 +======= +- Synchronized with upstream changes. +- Fixed a bug in gamma generation if the shape parameters is 0.0. + +v1.16.0 +======= +- Fixed a bug that affected :class:`~randomgen.dsfmt.DSFMT` when calling + :func:`~randomgen.dsfmt.DSFMT.jump` or :func:`~randomgen.dsfmt.DSFMT.seed` + that failed to reset the buffer. This resulted in upto 381 values from the + previous state being used before the buffer was refilled at the new state. +- Fixed bugs in :class:`~randomgen.xoshiro512starstar.Xoshiro512StarStar` + and :class:`~randomgen.xorshift1024.Xorshift1024` where the fallback + entropy initialization used too few bytes. This bug is unlikely to be + encountered since this path is only encountered if the system random + number generator fails. +- Synchronized with upstream changes. + +v1.15.1 +======= +- Added Xoshiro256** and Xoshiro512**, the preferred generators of this class. +- Fixed bug in `jump` method of Random123 generators which did nto specify a default value. +- Added support for generating bounded uniform integers using Lemire's method. +- Synchronized with upstream changes, which requires moving the minimum supported NumPy to 1.13. + +v1.15 +===== +- Synced empty choice changes +- Synced upstream docstring changes +- Synced upstream changes in permutation +- Synced upstream doc fixes +- Added absolute_import to avoid import noise on Python 2.7 +- Add legacy generator which allows NumPy replication +- Improve type handling of integers +- Switch to array-fillers for 0 parameter distribution to improve performance +- Small changes to build on manylinux +- Build wheels using multibuild diff --git a/doc/source/reference/randomgen/entropy.rst b/doc/source/reference/randomgen/entropy.rst new file mode 100644 index 000000000..6814edfbe --- /dev/null +++ b/doc/source/reference/randomgen/entropy.rst @@ -0,0 +1,6 @@ +System Entropy +============== + +.. module:: numpy.random.randomgen.entropy + +.. autofunction:: random_entropy diff --git a/doc/source/reference/randomgen/extending.rst b/doc/source/reference/randomgen/extending.rst new file mode 100644 index 000000000..c9d987b59 --- /dev/null +++ b/doc/source/reference/randomgen/extending.rst @@ -0,0 +1,163 @@ +Extending +--------- +The basic RNGs have been designed to be extendable using standard tools for +high-performance Python -- numba and Cython. +The :class:`randomgen.generator.RandomGenerator` object can also be used with +user-provided basic RNGs as long as these export a small set of required +functions. + +Numba +===== +Numba can be used with either CTypes or CFFI. The current iteration of the +basic RNGs all export a small set of functions through both interfaces. + +This example shows how numba can be used to produce Box-Muller normals using +a pure Python implementation which is then compiled. The random numbers are +provided by ``ctypes.next_double``. + +.. code-block:: python + + from randomgen import Xoroshiro128 + import numpy as np + import numba as nb + + x = Xoroshiro128() + f = x.ctypes.next_double + s = x.ctypes.state + state_addr = x.ctypes.state_address + + def normals(n, state): + out = np.empty(n) + for i in range((n+1)//2): + x1 = 2.0*f(state) - 1.0 + x2 = 2.0*f(state) - 1.0 + r2 = x1*x1 + x2*x2 + while r2 >= 1.0 or r2 == 0.0: + x1 = 2.0*f(state) - 1.0 + x2 = 2.0*f(state) - 1.0 + r2 = x1*x1 + x2*x2 + g = np.sqrt(-2.0*np.log(r2)/r2) + out[2*i] = g*x1 + if 2*i+1 < n: + out[2*i+1] = g*x2 + return out + + # Compile using Numba + print(normals(10, s).var()) + # Warm up + normalsj = nb.jit(normals, nopython=True) + # Must use state address not state with numba + normalsj(1, state_addr) + %timeit normalsj(1000000, state_addr) + print('1,000,000 Box-Muller (numba/Xoroshiro128) randoms') + %timeit np.random.standard_normal(1000000) + print('1,000,000 Box-Muller (NumPy) randoms') + + +Both CTypes and CFFI allow the more complicated distributions to be used +directly in Numba after compiling the file distributions.c into a DLL or so. +An example showing the use of a more complicated distribution is in the +examples folder. + +Cython +====== + +Cython can be used to unpack the ``PyCapsule`` provided by a basic RNG. +This example uses :class:`~randomgen.xoroshiro128.Xoroshiro128` and +``random_gauss_zig``, the Ziggurat-based generator for normals, to fill an +array. The usual caveats for writing high-performance code using Cython -- +removing bounds checks and wrap around, providing array alignment information +-- still apply. + +.. code-block:: cython + + import numpy as np + cimport numpy as np + cimport cython + from cpython.pycapsule cimport PyCapsule_IsValid, PyCapsule_GetPointer + from randomgen.common cimport * + from randomgen.distributions cimport random_gauss_zig + from randomgen.xoroshiro128 import Xoroshiro128 + + + @cython.boundscheck(False) + @cython.wraparound(False) + def normals_zig(Py_ssize_t n): + cdef Py_ssize_t i + cdef brng_t *rng + cdef const char *capsule_name = "BasicRNG" + cdef double[::1] random_values + + x = Xoroshiro128() + capsule = x.capsule + # Optional check that the capsule if from a Basic RNG + if not PyCapsule_IsValid(capsule, capsule_name): + raise ValueError("Invalid pointer to anon_func_state") + # Cast the pointer + rng = <brng_t *> PyCapsule_GetPointer(capsule, capsule_name) + random_values = np.empty(n) + for i in range(n): + # Call the function + random_values[i] = random_gauss_zig(rng) + randoms = np.asarray(random_values) + return randoms + + +The basic RNG can also be directly accessed using the members of the basic +RNG structure. + +.. code-block:: cython + + @cython.boundscheck(False) + @cython.wraparound(False) + def uniforms(Py_ssize_t n): + cdef Py_ssize_t i + cdef brng_t *rng + cdef const char *capsule_name = "BasicRNG" + cdef double[::1] random_values + + x = Xoroshiro128() + capsule = x.capsule + # Optional check that the capsule if from a Basic RNG + if not PyCapsule_IsValid(capsule, capsule_name): + raise ValueError("Invalid pointer to anon_func_state") + # Cast the pointer + rng = <brng_t *> PyCapsule_GetPointer(capsule, capsule_name) + random_values = np.empty(n) + for i in range(n): + # Call the function + random_values[i] = rng.next_double(rng.state) + randoms = np.asarray(random_values) + return randoms + +These functions along with a minimal setup file are included in the +examples folder. + +New Basic RNGs +============== +:class:`~randomgen.generator.RandomGenerator` can be used with other +user-provided basic RNGs. The simplest way to write a new basic RNG is to +examine the pyx file of one of the existing basic RNGs. The key structure +that must be provided is the ``capsule`` which contains a ``PyCapsule`` to a +struct pointer of type ``brng_t``, + +.. code-block:: c + + typedef struct brng { + void *state; + uint64_t (*next_uint64)(void *st); + uint32_t (*next_uint32)(void *st); + double (*next_double)(void *st); + uint64_t (*next_raw)(void *st); + } brng_t; + +which provides 5 pointers. The first is an opaque pointer to the data structure +used by the basic RNG. The next three are function pointers which return the +next 64- and 32-bit unsigned integers, the next random double and the next +raw value. This final function is used for testing and so can be set to +the next 64-bit unsigned integer function if not needed. Functions inside +:class:`~randomgen.generator.RandomGenerator` use this structure as in + +.. code-block:: c + + brng_state->next_uint64(brng_state->state) diff --git a/doc/source/reference/randomgen/generator.rst b/doc/source/reference/randomgen/generator.rst new file mode 100644 index 000000000..d59efd68c --- /dev/null +++ b/doc/source/reference/randomgen/generator.rst @@ -0,0 +1,88 @@ +Random Generator +---------------- +The :class:`~randomgen.generator.RandomGenerator` provides access to +a wide range of distributions, and served as a replacement for +:class:`~numpy.random.RandomState`. The main difference between +the two is that :class:`~randomgen.generator.RandomGenerator` relies +on an additional basic RNG to manage state and generate the random +bits which are then transformed into random values from useful +distributions. The default basic RNG used by +:class:`~randomgen.generator.RandomGenerator` is +:class:`~randomgen.xoroshiro128.Xoroshiro128`. The basic RNG can be +changed by passing an instantized basic RNG to +:class:`~randomgen.generator.RandomGenerator`. + +.. currentmodule:: numpy.random.randomgen.generator + +.. autoclass:: RandomGenerator + :exclude-members: + +Seed and State Manipulation +=========================== +.. autosummary:: + :toctree: generated/ + + ~RandomGenerator.seed + ~RandomGenerator.state + +Simple random data +================== +.. autosummary:: + :toctree: generated/ + + ~RandomGenerator.rand + ~RandomGenerator.randn + ~RandomGenerator.randint + ~RandomGenerator.random_integers + ~RandomGenerator.random_sample + ~RandomGenerator.choice + ~RandomGenerator.bytes + +Permutations +============ +.. autosummary:: + :toctree: generated/ + + ~RandomGenerator.shuffle + ~RandomGenerator.permutation + +Distributions +============= +.. autosummary:: + :toctree: generated/ + + ~RandomGenerator.beta + ~RandomGenerator.binomial + ~RandomGenerator.chisquare + ~RandomGenerator.dirichlet + ~RandomGenerator.exponential + ~RandomGenerator.f + ~RandomGenerator.gamma + ~RandomGenerator.geometric + ~RandomGenerator.gumbel + ~RandomGenerator.hypergeometric + ~RandomGenerator.laplace + ~RandomGenerator.logistic + ~RandomGenerator.lognormal + ~RandomGenerator.logseries + ~RandomGenerator.multinomial + ~RandomGenerator.multivariate_normal + ~RandomGenerator.negative_binomial + ~RandomGenerator.noncentral_chisquare + ~RandomGenerator.noncentral_f + ~RandomGenerator.normal + ~RandomGenerator.pareto + ~RandomGenerator.poisson + ~RandomGenerator.power + ~RandomGenerator.rayleigh + ~RandomGenerator.standard_cauchy + ~RandomGenerator.standard_exponential + ~RandomGenerator.standard_gamma + ~RandomGenerator.standard_normal + ~RandomGenerator.standard_t + ~RandomGenerator.triangular + ~RandomGenerator.uniform + ~RandomGenerator.vonmises + ~RandomGenerator.wald + ~RandomGenerator.weibull + ~RandomGenerator.zipf diff --git a/doc/source/reference/randomgen/index.rst b/doc/source/reference/randomgen/index.rst new file mode 100644 index 000000000..67d0441a2 --- /dev/null +++ b/doc/source/reference/randomgen/index.rst @@ -0,0 +1,228 @@ +Randomgen.RandomGen +=================== +This package contains replacements for the NumPy +:class:`~numpy.random.RandomState` object that allows the core random number +generator be be changed. + +.. current_module numpy.random.randomgen + +Quick Start +----------- + +Like :mod:`numpy.random`, RandomGen can be used at the module level. +This uses the default :class:`~randomgen.generator.RandomGenerator` which +uses normals provided by :class:`~randomgen.xoroshiro128.Xoroshiro128`. + +.. code-block:: python + + # As replacement for numpy.random + import randomgen.generator as random + random.standard_normal() + +:class:`~randomgen.generator.RandomGenerator` can also be used as a +replacement for :class:`~numpy.random.RandomState`, although the random +values are generated by :class:`~randomgen.xoroshiro128.Xoroshiro128`. It +also isn't possible to directly seed a +:class:`~randomgen.generator.RandomGenerator`. + + +.. code-block:: python + + # As replacement for RandomState() + from randomgen import RandomGenerator + rg = RandomGenerator() + rg.standard_normal() + + +Seeds can be passed to any of the basic RNGs. Here :class:`~randomgen.mt19937.MT19937` +is used and the :class:`~randomgen.generator.RandomGenerator` is accessed via +the property :attr:`~randomgen.mt19937.MT19937.generator`. + +.. code-block:: python + + from randomgen import MT19937 + rg = MT19937(12345).generator + rg.standard_normal() + + +Introduction +------------ +RandomGen takes a different approach to producing random numbers from the +:class:`numpy.random.RandomState` object used in NumPy. Random number +generation is separated into two components, a basic RNG and a random +generator. + +The basic RNG has a limited set of responsibilities -- it manages the +underlying RNG state and provides functions to produce random doubles and +random unsigned 32- and 64-bit values. The basic random generator also handles +all seeding since this varies when using alternative basic RNGs. + +The random generator (:class:`~randomgen.generator.RandomGenerator`) takes the +basic RNG-provided functions and transforms them into more useful +distributions, e.g., simulated normal random values. This structure allows +alternative basic RNGs to be used without code duplication. + +The :class:`~randomgen.generator.RandomGenerator` is the user-facing object +that is nearly identical to :class:`~numpy.random.RandomState`. The canonical +method to initialize a generator passes a basic RNG -- +:class:`~randomgen.mt19937.MT19937`, the underlying RNG in NumPy -- as the +sole argument. Note that the basic RNG must be instantized. + +.. code-block:: python + + from randomgen import RandomGenerator, MT19937 + rg = RandomGenerator(MT19937()) + rg.random_sample() + +Seed information is directly passed to the basic RNG. + +.. code-block:: python + + rg = RandomGenerator(MT19937(12345)) + rg.random_sample() + +A shorthand method is also available which uses the +:meth:`~randomgen.mt19937.MT19937.generator` property from a basic RNG to +access an embedded random generator. + +.. code-block:: python + + rg = MT19937(12345).generator + rg.random_sample() + +What's New or Different +~~~~~~~~~~~~~~~~~~~~~~~ +.. warning:: + + The Box-Muller method used to produce NumPy's normals is no longer available + in :class:`~randomgen.generator.RandomGenerator`. It is not possible to + reproduce the random values using :class:`~randomgen.generator.RandomGenerator` + for the normal distribution or any other distribution that relies on the + normal such as the gamma or student's t. If you require backward compatibility, a + legacy generator, :class:`~randomgen.legacy.LegacyGenerator`, has been created + which can fully reproduce the sequence produced by NumPy. + +* The normal, exponential and gamma generators use 256-step Ziggurat + methods which are 2-10 times faster than NumPy's Box-Muller or inverse CDF + implementations. +* Optional ``dtype`` argument that accepts ``np.float32`` or ``np.float64`` + to produce either single or double prevision uniform random variables for + select distributions +* Optional ``out`` argument that allows existing arrays to be filled for + select distributions +* Simulate from the complex normal distribution + (:meth:`~randomgen.generator.RandomGenerator.complex_normal`) +* :func:`~randomgen.entropy.random_entropy` provides access to the system + source of randomness that is used in cryptographic applications (e.g., + ``/dev/urandom`` on Unix). +* All basic random generators functions to produce doubles, uint64s and + uint32s via CTypes (:meth:`~randomgen.xoroshiro128.Xoroshiro128.ctypes`) + and CFFI (:meth:`~randomgen.xoroshiro128.Xoroshiro128.cffi`). This allows + these basic RNGs to be used in numba. +* The basic random number generators can be used in downstream projects via + Cython. +* Support for Lemire’s method [Lemire]_ of generating uniform integers on an + arbitrary interval by setting ``use_masked=True`` in + (:meth:`~randomgen.generator.RandomGenerator.randint`). + + +See :ref:`new-or-different` for a complete list of improvements and +differences. + +Parallel Generation +~~~~~~~~~~~~~~~~~~~ + +The included generators can be used in parallel, distributed applications in +one of two ways: + +* :ref:`independent-streams` +* :ref:`jump-and-advance` + +Supported Generators +-------------------- +The main innovation is the inclusion of a number of alternative pseudo-random number +generators, 'in addition' to the standard PRNG in NumPy. The included PRNGs are: + +* MT19937 - The standard NumPy generator. Produces identical results to NumPy + using the same seed/state. Adds a jump function that advances the generator + as-if 2**128 draws have been made (:meth:`~randomgen.mt19937.MT19937.jump`). + See `NumPy's documentation`_. +* dSFMT - SSE2 enabled versions of the MT19937 generator. Theoretically + the same, but with a different state and so it is not possible to produce a + sequence identical to MT19937. Supports ``jump`` and so can + be used in parallel applications. See the `dSFMT authors' page`_. +* XoroShiro128+ - Improved version of XorShift128+ with better performance + and statistical quality. Like the XorShift generators, it can be jumped + to produce multiple streams in parallel applications. See + :meth:`~randomgen.xoroshiro128.Xoroshiro128.jump` for details. + More information about this PRNG is available at the + `xorshift, xoroshiro and xoshiro authors' page`_. +* XorShift1024*φ - Fast fast generator based on the XSadd + generator. Supports ``jump`` and so can be used in + parallel applications. See the documentation for + :meth:`~randomgen.xorshift1024.Xorshift1024.jump` for details. More information + about these PRNGs is available at the + `xorshift, xoroshiro and xoshiro authors' page`_. +* Xorshiro256** and Xorshiro512** - The most recently introduced XOR, + shift, and rotate generator. Supports ``jump`` and so can be used in + parallel applications. See the documentation for + :meth:`~randomgen.xoshiro256starstar.Xoshirt256StarStar.jump` for details. More + information about these PRNGs is available at the + `xorshift, xoroshiro and xoshiro authors' page`_. +* PCG-64 - Fast generator that support many parallel streams and + can be advanced by an arbitrary amount. See the documentation for + :meth:`~randomgen.pcg64.PCG64.advance`. PCG-64 has a period of + :math:`2^{128}`. See the `PCG author's page`_ for more details about + this class of PRNG. +* ThreeFry and Philox - counter-based generators capable of being advanced an + arbitrary number of steps or generating independent streams. See the + `Random123`_ page for more details about this class of PRNG. + +.. _`NumPy's documentation`: https://docs.scipy.org/doc/numpy/reference/routines.random.html +.. _`dSFMT authors' page`: http://www.math.sci.hiroshima-u.ac.jp/~m-mat/MT/SFMT/ +.. _`xorshift, xoroshiro and xoshiro authors' page`: http://xoroshiro.di.unimi.it/ +.. _`PCG author's page`: http://www.pcg-random.org/ +.. _`Random123`: https://www.deshawresearch.com/resources_random123.html + +Random Generator +---------------- +.. toctree:: + :maxdepth: 1 + + generator + legacy + +Basic Random Number Generators +------------------------------ + +.. toctree:: + :maxdepth: 1 + + Basic Random Number Generators <brng/index> + +New Features +------------ +.. toctree:: + :maxdepth: 2 + + Parallel Applications <parallel> + Multithreaded Generation <multithreading> + new-or-different + Comparing Performance <performance> + extending + Reading System Entropy <entropy> + references + +Changes +~~~~~~~ +.. toctree:: + :maxdepth: 2 + + Change Log <change-log> + +Indices and tables +~~~~~~~~~~~~~~~~~~ + +* :ref:`genindex` +* :ref:`modindex` +* :ref:`search` diff --git a/doc/source/reference/randomgen/legacy.rst b/doc/source/reference/randomgen/legacy.rst new file mode 100644 index 000000000..7e87f871c --- /dev/null +++ b/doc/source/reference/randomgen/legacy.rst @@ -0,0 +1,110 @@ +Legacy Random Generation +------------------------ +The :class:`~randomgen.legacy.LegacyGenerator` provides access to +some legacy generators. These all depend on Box-Muller normals or +inverse CDF exponentials or gammas. This class should only be used +if it is essential to have randoms that are identical to what +would have been produced by NumPy. + +:class:`~randomgen.legacy.LegacyGenerator` add additional information +to the state which is required when using Box-Muller normals since these +are produced in pairs. It is important to use +:attr:`~randomgen.legacy.LegacyGenerator.state` +when accessing the state so that these extra values are saved. + +.. warning:: + + :class:`~randomgen.legacy.LegacyGenerator` only contains functions + that have changed. Since it does not contain other functions, it + is not direclty possible to replace :class:`~numpy.random.RandomState`. + In order to full replace :class:`~numpy.random.RandomState`, it is + necessary to use both :class:`~randomgen.legacy.LegacyGenerator` + and :class:`~randomgen.generator.RandomGenerator` both driven + by the same basic RNG. Methods present in :class:`~randomgen.legacy.LegacyGenerator` + must be called from :class:`~randomgen.legacy.LegacyGenerator`. Other Methods + should be called from :class:`~randomgen.generator.RandomGenerator`. + + +.. code-block:: python + + from randomgen import RandomGenerator, MT19937 + from randomgen.legacy import LegacyGenerator + from numpy.random import RandomState + # Use same seed + rs = RandomState(12345) + mt19937 = MT19937(12345) + rg = RandomGenerator(mt19937) + lg = LegacyGenerator(mt19937) + + # Identical output + rs.standard_normal() + lg.standard_normal() + + rs.random_sample() + rg.random_sample() + + rs.standard_exponential() + lg.standard_exponential() + + +.. currentmodule:: numpy.random.randomgen.legacy + +.. autoclass:: LegacyGenerator + :exclude-members: + +Seeding and State +================= + +.. autosummary:: + :toctree: generated/ + + ~LegacyGenerator.get_state + ~LegacyGenerator.set_state + +Simple random data +================== +.. autosummary:: + :toctree: generated/ + + ~LegacyGenerator.randn + ~LegacyGenerator.randint + ~LegacyGenerator.random_integers + ~LegacyGenerator.random_sample + ~LegacyGenerator.choice + ~LegacyGenerator.bytes + +Permutations +============ +.. autosummary:: + :toctree: generated/ + + ~LegacyGenerator.shuffle + ~LegacyGenerator.permutation + +Distributions +============= +.. autosummary:: + :toctree: generated/ + + ~LegacyGenerator.beta + ~LegacyGenerator.chisquare + ~LegacyGenerator.dirichlet + ~LegacyGenerator.exponential + ~LegacyGenerator.f + ~LegacyGenerator.gamma + ~LegacyGenerator.lognormal + ~LegacyGenerator.multivariate_normal + ~LegacyGenerator.negative_binomial + ~LegacyGenerator.noncentral_chisquare + ~LegacyGenerator.noncentral_f + ~LegacyGenerator.normal + ~LegacyGenerator.pareto + ~LegacyGenerator.power + ~LegacyGenerator.standard_cauchy + ~LegacyGenerator.standard_exponential + ~LegacyGenerator.standard_gamma + ~LegacyGenerator.standard_normal + ~LegacyGenerator.standard_t + ~LegacyGenerator.wald + ~LegacyGenerator.weibull + ~LegacyGenerator.zipf diff --git a/doc/source/reference/randomgen/multithreading.rst b/doc/source/reference/randomgen/multithreading.rst new file mode 100644 index 000000000..6efbcdbe7 --- /dev/null +++ b/doc/source/reference/randomgen/multithreading.rst @@ -0,0 +1,106 @@ +Multithreaded Generation +======================== + +The four core distributions all allow existing arrays to be filled using the +``out`` keyword argument. Existing arrays need to be contiguous and +well-behaved (writable and aligned). Under normal circumstances, arrays +created using the common constructors such as :meth:`numpy.empty` will satisfy +these requirements. + +This example makes use of Python 3 :mod:`concurrent.futures` to fill an array +using multiple threads. Threads are long-lived so that repeated calls do not +require any additional overheads from thread creation. The underlying PRNG is +xorshift2014 which is fast, has a long period and supports using ``jump`` to +advance the state. The random numbers generated are reproducible in the sense +that the same seed will produce the same outputs. + +.. code-block:: ipython + + from randomgen import Xorshift1024 + import multiprocessing + import concurrent.futures + import numpy as np + + class MultithreadedRNG(object): + def __init__(self, n, seed=None, threads=None): + rg = Xorshift1024(seed) + if threads is None: + threads = multiprocessing.cpu_count() + self.threads = threads + + self._random_generators = [] + for _ in range(0, threads-1): + _rg = Xorshift1024() + _rg.state = rg.state + self._random_generators.append(_rg.generator) + rg.jump() + self._random_generators.append(rg.generator) + + self.n = n + self.executor = concurrent.futures.ThreadPoolExecutor(threads) + self.values = np.empty(n) + self.step = np.ceil(n / threads).astype(np.int) + + def fill(self): + def _fill(random_state, out, first, last): + random_state.standard_normal(out=out[first:last]) + + futures = {} + for i in range(self.threads): + args = (_fill, + self._random_generators[i], + self.values, + i * self.step, + (i + 1) * self.step) + futures[self.executor.submit(*args)] = i + concurrent.futures.wait(futures) + + def __del__(self): + self.executor.shutdown(False) + + +The multithreaded random number generator can be used to fill an array. +The ``values`` attributes shows the zero-value before the fill and the +random value after. + +.. code-block:: ipython + + In [2]: mrng = MultithreadedRNG(10000000, seed=0) + ...: print(mrng.values[-1]) + 0.0 + + In [3]: mrng.fill() + ...: print(mrng.values[-1]) + 3.296046120254392 + +The time required to produce using multiple threads can be compared to +the time required to generate using a single thread. + +.. code-block:: ipython + + In [4]: print(mrng.threads) + ...: %timeit mrng.fill() + + 4 + 32.8 ms ± 2.71 ms per loop (mean ± std. dev. of 7 runs, 10 loops each) + +The single threaded call directly uses the PRNG. + +.. code-block:: ipython + + In [5]: values = np.empty(10000000) + ...: rg = Xorshift1024().generator + ...: %timeit rg.standard_normal(out=values) + + 99.6 ms ± 222 µs per loop (mean ± std. dev. of 7 runs, 10 loops each) + +The gains are substantial and the scaling is reasonable even for large that +are only moderately large. The gains are even larger when compared to a call +that does not use an existing array due to array creation overhead. + +.. code-block:: ipython + + In [6]: rg = Xorshift1024().generator + ...: %timeit rg.standard_normal(10000000) + + 125 ms ± 309 µs per loop (mean ± std. dev. of 7 runs, 10 loops each) diff --git a/doc/source/reference/randomgen/new-or-different.rst b/doc/source/reference/randomgen/new-or-different.rst new file mode 100644 index 000000000..6598c13fe --- /dev/null +++ b/doc/source/reference/randomgen/new-or-different.rst @@ -0,0 +1,98 @@ +.. _new-or-different: + +What's New or Different +----------------------- + +.. warning:: + + The Box-Muller method used to produce NumPy's normals is no longer available + in :class:`~randomgen.generator.RandomGenerator`. It is not possible to + reproduce the random values using :class:`~randomgen.generator.RandomGenerator` + for the normal distribution or any other distribution that relies on the + normal such as the gamma or student's t. If you require backward compatibility, a + legacy generator, :class:`~randomgen.legacy.LegacyGenerator`, has been created + which can fully reproduce the sequence produced by NumPy. + + +* :func:`~randomgen.entropy.random_entropy` provides access to the system + source of randomness that is used in cryptographic applications (e.g., + ``/dev/urandom`` on Unix). +* Simulate from the complex normal distribution + (:meth:`~randomgen.generator.RandomGenerator.complex_normal`) +* The normal, exponential and gamma generators use 256-step Ziggurat + methods which are 2-10 times faster than NumPy's default implementation in + :meth:`~randomgen.generator.RandomGenerator.standard_normal`, + :meth:`~randomgen.generator.RandomGenerator.standard_exponential` or + :meth:`~randomgen.generator.RandomGenerator.standard_gamma`. +* The Box-Muller used to produce NumPy's normals is no longer available. +* All basic random generators functions to produce doubles, uint64s and + uint32s via CTypes (:meth:`~randomgen.xoroshiro128.Xoroshiro128.ctypes`) + and CFFI (:meth:`~randomgen.xoroshiro128.Xoroshiro128.cffi`). This allows + these basic RNGs to be used in numba. +* The basic random number generators can be used in downstream projects via + Cython. + + +.. ipython:: python + + from randomgen import Xoroshiro128 + import numpy.random + rg = Xoroshiro128().generator + %timeit rg.standard_normal(100000) + %timeit numpy.random.standard_normal(100000) + +.. ipython:: python + + %timeit rg.standard_exponential(100000) + %timeit numpy.random.standard_exponential(100000) + +.. ipython:: python + + %timeit rg.standard_gamma(3.0, 100000) + %timeit numpy.random.standard_gamma(3.0, 100000) + +* Optional ``dtype`` argument that accepts ``np.float32`` or ``np.float64`` + to produce either single or double prevision uniform random variables for + select distributions + + * Uniforms (:meth:`~randomgen.generator.RandomGenerator.random_sample` and + :meth:`~randomgen.generator.RandomGenerator.rand`) + * Normals (:meth:`~randomgen.generator.RandomGenerator.standard_normal` and + :meth:`~randomgen.generator.RandomGenerator.randn`) + * Standard Gammas (:meth:`~randomgen.generator.RandomGenerator.standard_gamma`) + * Standard Exponentials (:meth:`~randomgen.generator.RandomGenerator.standard_exponential`) + +.. ipython:: python + + rg.seed(0) + rg.random_sample(3, dtype='d') + rg.seed(0) + rg.random_sample(3, dtype='f') + +* Optional ``out`` argument that allows existing arrays to be filled for + select distributions + + * Uniforms (:meth:`~randomgen.generator.RandomGenerator.random_sample`) + * Normals (:meth:`~randomgen.generator.RandomGenerator.standard_normal`) + * Standard Gammas (:meth:`~randomgen.generator.RandomGenerator.standard_gamma`) + * Standard Exponentials (:meth:`~randomgen.generator.RandomGenerator.standard_exponential`) + + This allows multithreading to fill large arrays in chunks using suitable + PRNGs in parallel. + +.. ipython:: python + + existing = np.zeros(4) + rg.random_sample(out=existing[:2]) + print(existing) + +.. * For changes since the previous release, see the :ref:`change-log` + +* Support for Lemire’s method of generating uniform integers on an + arbitrary interval by setting ``use_masked=True`` in + (:meth:`~randomgen.generator.RandomGenerator.randint`). + +.. ipython:: python + + %timeit rg.randint(0, 1535, use_masked=False) + %timeit numpy.random.randint(0, 1535) diff --git a/doc/source/reference/randomgen/parallel.rst b/doc/source/reference/randomgen/parallel.rst new file mode 100644 index 000000000..df6f58d75 --- /dev/null +++ b/doc/source/reference/randomgen/parallel.rst @@ -0,0 +1,141 @@ +Parallel Random Number Generation +================================= + +There are three strategies implemented that can be used to produce +repeatable pseudo-random numbers across multiple processes (local +or distributed). + +.. _independent-streams: + +Independent Streams +------------------- + +:class:`~randomgen.pcg64.PCG64`, :class:`~randomgen.threefry.ThreeFry` +and :class:`~randomgen.philox.Philox` support independent streams. This +example shows how many streams can be created by passing in different index +values in the second input while using the same seed in the first. + +.. code-block:: python + + from randomgen.entropy import random_entropy + from randomgen import PCG64 + + entropy = random_entropy(4) + # 128-bit number as a seed + seed = sum([int(entropy[i]) * 2 ** (32 * i) for i in range(4)]) + streams = [PCG64(seed, stream) for stream in range(10)] + + +:class:`~randomgen.philox.Philox` and :class:`~randomgen.threefry.ThreeFry` are +counter-based RNGs which use a counter and key. Different keys can be used +to produce independent streams. + +.. code-block:: python + + import numpy as np + from randomgen import ThreeFry + + key = random_entropy(8) + key = key.view(np.uint64) + key[0] = 0 + step = np.zeros(4, dtype=np.uint64) + step[0] = 1 + streams = [ThreeFry(key=key + stream * step) for stream in range(10)] + +.. _jump-and-advance: + +Jump/Advance the PRNG state +--------------------------- + +Jump +**** + +``jump`` advances the state of the PRNG *as-if* a large number of random +numbers have been drawn. The specific number of draws varies by PRNG, and +ranges from :math:`2^{64}` to :math:`2^{512}`. Additionally, the *as-if* +draws also depend on the size of the default random number produced by the +specific PRNG. The PRNGs that support ``jump``, along with the period of +the PRNG, the size of the jump and the bits in the default unsigned random +are listed below. + ++-----------------+-------------------------+-------------------------+-------------------------+ +| PRNG | Period | Jump Size | Bits | ++=================+=========================+=========================+=========================+ +| DSFMT | :math:`2^{19937}` | :math:`2^{128}` | 53 | ++-----------------+-------------------------+-------------------------+-------------------------+ +| MT19937 | :math:`2^{19937}` | :math:`2^{128}` | 32 | ++-----------------+-------------------------+-------------------------+-------------------------+ +| PCG64 | :math:`2^{128}` | :math:`2^{64}` | 64 | ++-----------------+-------------------------+-------------------------+-------------------------+ +| Philox | :math:`2^{256}` | :math:`2^{128}` | 64 | ++-----------------+-------------------------+-------------------------+-------------------------+ +| ThreeFry | :math:`2^{256}` | :math:`2^{128}` | 64 | ++-----------------+-------------------------+-------------------------+-------------------------+ +| Xoroshiro128 | :math:`2^{128}` | :math:`2^{64}` | 64 | ++-----------------+-------------------------+-------------------------+-------------------------+ +| Xorshift1024 | :math:`2^{1024}` | :math:`2^{512}` | 64 | ++-----------------+-------------------------+-------------------------+-------------------------+ + +``jump`` can be used to produce long blocks which should be long enough to not +overlap. + +.. code-block:: python + + from randomgen.entropy import random_entropy + from randomgen import Xorshift1024 + + entropy = random_entropy(2).astype(np.uint64) + # 64-bit number as a seed + seed = entropy[0] * 2**32 + entropy[1] + blocked_rng = [] + for i in range(10): + rng = Xorshift1024(seed) + rng.jump(i) + blocked_rng.append(rng) + + +Advance +******* +``advance`` can be used to jump the state an arbitrary number of steps, and so +is a more general approach than ``jump``. :class:`~randomgen.pcg64.PCG64`, +:class:`~randomgen.threefry.ThreeFry` and :class:`~randomgen.philox.Philox` +support ``advance``, and since these also support independent +streams, it is not usually necessary to use ``advance``. + +Advancing a PRNG updates the underlying PRNG state as-if a given number of +calls to the underlying PRNG have been made. In general there is not a +one-to-one relationship between the number output random values from a +particular distribution and the number of draws from the core PRNG. +This occurs for two reasons: + +* The random values are simulated using a rejection-based method + and so, on average, more than one value from the underlying + PRNG is required to generate an single draw. +* The number of bits required to generate a simulated value + differs from the number of bits generated by the underlying + PRNG. For example, two 16-bit integer values can be simulated + from a single draw of a 32-bit PRNG. + +Advancing the PRNG state resets any pre-computed random numbers. This is +required to ensure exact reproducibility. + +This example uses ``advance`` to advance a :class:`~randomgen.pcg64.PCG64` +generator 2 ** 127 steps to set a sequence of random number generators. + +.. code-block:: python + + from randomgen import PCG64 + brng = PCG64() + brng_copy = PCG64() + brng_copy.state = brng.state + + advance = 2**127 + brngs = [brng] + for _ in range(9): + brng_copy.advance(advance) + brng = PCG64() + brng.state = brng_copy.state + brngs.append(brng) + +.. end block + diff --git a/doc/source/reference/randomgen/performance.py b/doc/source/reference/randomgen/performance.py new file mode 100644 index 000000000..12cbbc5d3 --- /dev/null +++ b/doc/source/reference/randomgen/performance.py @@ -0,0 +1,74 @@ +from collections import OrderedDict +from timeit import repeat + +import numpy as np +import pandas as pd + +from randomgen import MT19937, DSFMT, ThreeFry, PCG64, Xoroshiro128, \ + Xorshift1024, Philox, Xoshiro256StarStar, Xoshiro512StarStar + +PRNGS = [DSFMT, MT19937, Philox, PCG64, ThreeFry, Xoroshiro128, Xorshift1024, + Xoshiro256StarStar, Xoshiro512StarStar] + +funcs = {'32-bit Unsigned Ints': 'random_uintegers(size=1000000,bits=32)', + '64-bit Unsigned Ints': 'random_uintegers(size=1000000,bits=32)', + 'Uniforms': 'random_sample(size=1000000)', + 'Complex Normals': 'complex_normal(size=1000000)', + 'Normals': 'standard_normal(size=1000000)', + 'Exponentials': 'standard_exponential(size=1000000)', + 'Gammas': 'standard_gamma(3.0,size=1000000)', + 'Binomials': 'binomial(9, .1, size=1000000)', + 'Laplaces': 'laplace(size=1000000)', + 'Poissons': 'poisson(3.0, size=1000000)', } + +setup = """ +from randomgen import {prng} +rg = {prng}().generator +""" + +test = "rg.{func}" +table = OrderedDict() +for prng in PRNGS: + print(prng) + col = OrderedDict() + for key in funcs: + t = repeat(test.format(func=funcs[key]), + setup.format(prng=prng().__class__.__name__), + number=1, repeat=3) + col[key] = 1000 * min(t) + col = pd.Series(col) + table[prng().__class__.__name__] = col + +npfuncs = OrderedDict() +npfuncs.update(funcs) +npfuncs['32-bit Unsigned Ints'] = 'randint(2**32,dtype="uint32",size=1000000)' +npfuncs['64-bit Unsigned Ints'] = 'tomaxint(size=1000000)' +del npfuncs['Complex Normals'] +setup = """ +from numpy.random import RandomState +rg = RandomState() +""" +col = {} +for key in npfuncs: + t = repeat(test.format(func=npfuncs[key]), + setup.format(prng=prng().__class__.__name__), + number=1, repeat=3) + col[key] = 1000 * min(t) +table['NumPy'] = pd.Series(col) + +table = pd.DataFrame(table) +table = table.reindex(table.mean(1).sort_values().index) +order = np.log(table).mean().sort_values().index +table = table.T +table = table.reindex(order) +table = table.T +print(table.to_csv(float_format='%0.1f')) + +rel = table.loc[:, ['NumPy']].values @ np.ones((1, table.shape[1])) / table +rel.pop(rel.columns[0]) +rel = rel.T +rel['Overall'] = np.exp(np.log(rel).mean(1)) +rel *= 100 +rel = np.round(rel) +rel = rel.T +print(rel.to_csv(float_format='%0d')) diff --git a/doc/source/reference/randomgen/performance.rst b/doc/source/reference/randomgen/performance.rst new file mode 100644 index 000000000..2dfb32101 --- /dev/null +++ b/doc/source/reference/randomgen/performance.rst @@ -0,0 +1,75 @@ +Performance +----------- + +.. py:module:: randomgen + +Recommendation +************** +The recommended generator for single use is +:class:`~randomgen.xoroshiro128.Xoroshiro128`. The recommended generator +for use in large-scale parallel applications is +:class:`~randomgen.xorshift1024.Xorshift1024` +where the `jump` method is used to advance the state. For very large scale +applications -- requiring 1,000+ independent streams, +:class:`~randomgen.pcg64.PCG64` or :class:`~randomgen.threefry.ThreeFry` are +the best choices. + +Timings +******* + +The timings below are the time in ms to produce 1,000,000 random values from a +specific distribution. :class:`~randomgen.xoroshiro128.Xoroshiro128` is the +fastest, followed by :class:`~randomgen.xorshift1024.Xorshift1024` and +:class:`~randomgen.pcg64.PCG64`. The original :class:`~randomgen.mt19937.MT19937` +generator is much slower since it requires 2 32-bit values to equal the output +of the faster generators. + +Integer performance has a similar ordering although `dSFMT` is slower since +it generates 53-bit floating point values rather than integer values. On the +other hand, it is very fast for uniforms, although slower than `xoroshiro128+`. + +The pattern is similar for other, more complex generators. The normal +performance of NumPy's MT19937 is much lower than the other since it +uses the Box-Muller transformation rather than the Ziggurat generator. The +performance gap for Exponentials is also large due to the cost of computing +the log function to invert the CDF. + +.. csv-table:: + :header: ,Xoroshiro128,Xorshift1024,PCG64,DSFMT,MT19937,Philox,ThreeFry,NumPy + :widths: 14,14,14,14,14,14,14,14,14 + + 32-bit Unsigned Ints,3.0,3.0,3.0,3.5,3.7,6.8,6.6,3.3 + 64-bit Unsigned Ints,2.6,3.0,3.1,3.4,3.8,6.9,6.6,8.8 + Uniforms,3.2,3.8,4.4,5.0,7.4,8.9,9.9,8.8 + Normals,11.0,13.9,13.7,15.8,16.9,17.8,18.8,63.0 + Exponentials,7.0,8.4,9.0,11.2,12.5,14.1,15.0,102.2 + Binomials,20.9,22.6,22.0,21.2,26.7,27.7,29.2,26.5 + Complex Normals,23.2,28.7,29.1,33.2,35.4,37.6,38.6, + Gammas,35.3,38.6,39.2,41.3,46.7,49.4,51.2,98.8 + Laplaces,97.8,99.9,99.8,96.2,104.1,104.6,104.8,104.1 + Poissons,104.8,113.2,113.3,107.6,129.7,135.6,138.1,131.9 + + +The next table presents the performance relative to `xoroshiro128+` in +percentage. The overall performance was computed using a geometric mean. + +.. csv-table:: + :header: ,Xorshift1024,PCG64,DSFMT,MT19937,Philox,ThreeFry,NumPy + :widths: 14,14,14,14,14,14,14,14 + + 32-bit Unsigned Ints,102,99,118,125,229,221,111 + 64-bit Unsigned Ints,114,116,129,143,262,248,331 + Uniforms,116,137,156,231,275,306,274 + Normals,126,124,143,153,161,170,572 + Exponentials,121,130,161,179,203,215,1467 + Binomials,108,105,101,128,133,140,127 + Complex Normals,124,125,143,153,162,166, + Gammas,109,111,117,132,140,145,280 + Laplaces,102,102,98,106,107,107,106 + Poissons,108,108,103,124,129,132,126 + Overall,113,115,125,144,172,177,251 + + +.. note:: + + All timings were taken using Linux on a i5-3570 processor. diff --git a/doc/source/reference/randomgen/references.rst b/doc/source/reference/randomgen/references.rst new file mode 100644 index 000000000..0dc99868f --- /dev/null +++ b/doc/source/reference/randomgen/references.rst @@ -0,0 +1,5 @@ +References +---------- + +.. [Lemire] Daniel Lemire., "Fast Random Integer Generation in an Interval", + CoRR, Aug. 13, 2018, http://arxiv.org/abs/1805.10941. diff --git a/doc/source/reference/routines.rst b/doc/source/reference/routines.rst index a9e80480b..0ed99cbda 100644 --- a/doc/source/reference/routines.rst +++ b/doc/source/reference/routines.rst @@ -42,6 +42,7 @@ indentation. routines.padding routines.polynomials routines.random + randomgen/index routines.set routines.sort routines.statistics |