Extending#
TheBitGenerators have been designed to be extendable using standard toolsfor high-performance Python – numba and Cython. TheGenerator object canalso be used with user-providedBitGenerators as long as these export asmall set of required functions.
Numba#
Numba can be used with eitherCTypesorCFFI.The current iteration of theBitGenerators all export a small set of functions through both interfaces.
This example shows how Numba can be used to produce Gaussian samples usinga pure Python implementation which is then compiled. The random numbers areprovided byctypes.next_double.
fromtimeitimporttimeitimportnumbaasnbimportnumpyasnpfromnumpy.randomimportPCG64bit_gen=PCG64()next_d=bit_gen.cffi.next_doublestate_addr=bit_gen.cffi.state_addressdefnormals(n,state):out=np.empty(n)foriinrange((n+1)//2):x1=2.0*next_d(state)-1.0x2=2.0*next_d(state)-1.0r2=x1*x1+x2*x2whiler2>=1.0orr2==0.0:x1=2.0*next_d(state)-1.0x2=2.0*next_d(state)-1.0r2=x1*x1+x2*x2f=np.sqrt(-2.0*np.log(r2)/r2)out[2*i]=f*x1if2*i+1<n:out[2*i+1]=f*x2returnout# Compile using Numbanormalsj=nb.jit(normals,nopython=True)# Must use state address not state with numban=10000defnumbacall():returnnormalsj(n,state_addr)rg=np.random.Generator(PCG64())defnumpycall():returnrg.normal(size=n)# Check that the functions workr1=numbacall()r2=numpycall()assertr1.shape==(n,)assertr1.shape==r2.shapet1=timeit(numbacall,number=1000)print(f'{t1:.2f} secs for{n} PCG64 (Numba/PCG64) gaussian randoms')t2=timeit(numpycall,number=1000)print(f'{t2:.2f} secs for{n} PCG64 (NumPy/PCG64) gaussian randoms')
Both CTypes and CFFI allow the more complicated distributions to be useddirectly in Numba after compiling the file distributions.c into aDLL orso. An example showing the use of a more complicated distribution is intheExamples section below.
Cython#
Cython can be used to unpack thePyCapsule provided by aBitGenerator.This example usesPCG64 and the example from above. The usual caveatsfor writing high-performance code using Cython – removing bounds checks andwrap around, providing array alignment information – still apply.
#cython: language_level=3"""This file shows how the to use a BitGenerator to create a distribution."""importnumpyasnpcimportnumpyasnpcimportcythonfromcpython.pycapsulecimportPyCapsule_IsValid,PyCapsule_GetPointerfromlibc.stdintcimportuint16_t,uint64_tfromnumpy.randomcimportbitgen_tfromnumpy.randomimportPCG64fromnumpy.random.c_distributionscimport(random_standard_uniform_fill,random_standard_uniform_fill_f)np.import_array()@cython.boundscheck(False)@cython.wraparound(False)defuniforms(Py_ssize_tn):""" Create an array of `n` uniformly distributed doubles. A 'real' distribution would want to process the values into some non-uniform distribution """cdefPy_ssize_ticdefbitgen_t *rngcdefconstchar *capsule_name="BitGenerator"cdefdouble[::1]random_valuesx=PCG64()capsule=x.capsule# Optional check that the capsule if from a BitGeneratorifnotPyCapsule_IsValid(capsule,capsule_name):raiseValueError("Invalid pointer to anon_func_state")# Cast the pointerrng=<bitgen_t*>PyCapsule_GetPointer(capsule,capsule_name)random_values=np.empty(n,dtype='float64')withx.lock,nogil:foriinrange(n):# Call the functionrandom_values[i]=rng.next_double(rng.state)randoms=np.asarray(random_values)returnrandoms
TheBitGenerator can also be directly accessed using the members of thebitgen_tstruct.
@cython.boundscheck(False)@cython.wraparound(False)defuint10_uniforms(Py_ssize_tn):"""Uniform 10 bit integers stored as 16-bit unsigned integers"""cdefPy_ssize_ticdefbitgen_t *rngcdefconstchar *capsule_name="BitGenerator"cdefuint16_t[::1]random_valuescdefintbits_remainingcdefintwidth=10cdefuint64_tbuff,mask=0x3FFx=PCG64()capsule=x.capsuleifnotPyCapsule_IsValid(capsule,capsule_name):raiseValueError("Invalid pointer to anon_func_state")rng=<bitgen_t*>PyCapsule_GetPointer(capsule,capsule_name)random_values=np.empty(n,dtype='uint16')# Best practice is to release GIL and acquire the lockbits_remaining=0withx.lock,nogil:foriinrange(n):ifbits_remaining<width:buff=rng.next_uint64(rng.state)random_values[i]=buff&maskbuff>>=widthrandoms=np.asarray(random_values)returnrandoms
Cython can be used to directly access the functions innumpy/random/c_distributions.pxd. This requires linking with thenpyrandom library located innumpy/random/lib.
defuniforms_ex(bit_generator,Py_ssize_tn,dtype=np.float64):""" Create an array of `n` uniformly distributed doubles via a "fill" function. A 'real' distribution would want to process the values into some non-uniform distribution Parameters ---------- bit_generator: BitGenerator instance n: int Output vector length dtype: {str, dtype}, optional Desired dtype, either 'd' (or 'float64') or 'f' (or 'float32'). The default dtype value is 'd' """cdefPy_ssize_ticdefbitgen_t *rngcdefconstchar *capsule_name="BitGenerator"cdefnp.ndarrayrandomscapsule=bit_generator.capsule# Optional check that the capsule if from a BitGeneratorifnotPyCapsule_IsValid(capsule,capsule_name):raiseValueError("Invalid pointer to anon_func_state")# Cast the pointerrng=<bitgen_t*>PyCapsule_GetPointer(capsule,capsule_name)_dtype=np.dtype(dtype)randoms=np.empty(n,dtype=_dtype)if_dtype==np.float32:withbit_generator.lock:random_standard_uniform_fill_f(rng,n,<float*>np.PyArray_DATA(randoms))elif_dtype==np.float64:withbit_generator.lock:random_standard_uniform_fill(rng,n,<double*>np.PyArray_DATA(randoms))else:raiseTypeError('Unsupported dtype%r for random'%_dtype)returnrandoms
SeeExtending numpy.random via Cython for the complete listings of these examplesand a minimalsetup.py to build the c-extension modules.
CFFI#
CFFI can be used to directly access the functions ininclude/numpy/random/distributions.h. Some “massaging” of the headerfile is required:
"""Use cffi to access any of the underlying C functions from distributions.h"""importosimportcffiimportnumpyasnpfrom.parseimportparse_distributions_hffi=cffi.FFI()inc_dir=os.path.join(np.get_include(),'numpy')# Basic numpy typesffi.cdef(''' typedef intptr_t npy_intp; typedef unsigned char npy_bool;''')parse_distributions_h(ffi,inc_dir)
Once the header is parsed byffi.cdef, the functions can be accesseddirectly from the_generator shared object, using theBitGenerator.cffi interface.
lib=ffi.dlopen(np.random._generator.__file__)# Compare the distributions.h random_standard_normal_fill to# Generator.standard_randombit_gen=np.random.PCG64()rng=np.random.Generator(bit_gen)state=bit_gen.stateinterface=rng.bit_generator.cffin=100vals_cffi=ffi.new('double[%d]'%n)lib.random_standard_normal_fill(interface.bit_generator,n,vals_cffi)# reset the statebit_gen.state=statevals=rng.standard_normal(n)foriinrange(n):assertvals[i]==vals_cffi[i]
New BitGenerators#
Generator can be used with user-providedBitGenerators. The simplestway to write a newBitGenerator is to examine the pyx file of one of theexistingBitGenerators. The key structure that must be provided is thecapsule which contains aPyCapsule to a struct pointer of typebitgen_t,
typedefstructbitgen{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);}bitgen_t;
which provides 5 pointers. The first is an opaque pointer to the data structureused by theBitGenerators. The next three are function pointers whichreturn the next 64- and 32-bit unsigned integers, the next random double andthe next raw value. This final function is used for testing and so can be setto the next 64-bit unsigned integer function if not needed. Functions insideGenerator use this structure as in
bitgen_state->next_uint64(bitgen_state->state)