Regenerate vector of randoms in Python
I am trying to generate and regenerate vector of random normals
I want to be able to achieve the following by generating random normal matrix of size 100x3 and random normals having mean 0 and sd 1:
seed1 = '123' seed2 = 'asd' randMatrixrows = 100 randMatrixcols = 3 mu = 0 sd = 1 normRand1 = rekina_normRandomGenerator( seed1, randMatrixrows, randMatrixcols, mu, sd ) #normRand1 is of size 100x3 normRand2 = rekina_normRandomGenerator( seed2, randMatrixrows, randMatrixcols, mu, sd ) normRand3 = rekina_normRandomGenerator( seed1, randMatrixrows, randMatrixcols, mu, sd ) normRand4 = rekina_normRandomGenerator( seed2, randMatrixrows, randMatrixcols, mu, sd ) err1 = normRand1 - normRand3 err2 = normRand2 - normRand4
Each element of err1 and err2 should be 0
I have tried to read also but being new to Python I am totally lost with the implementation. I was expecting a simple implementation to show how to use the CBRNG.
Your initial assumption that you have to use reikna.cbrng.CBRNG is correct. It seems it offers a multitude of sources of pseudo randomness – based on counter-based RNGs – and distributions, which can be confusing. It also offers shortcuts for creating random number generators with given distribution.
from reikna.core import Type from reikna.cbrng import CBRNG from reikna.cluda import any_api import numpy as np # Get any supported API for this system, or an exception _api = any_api() # The global CLUDA Thread, wrapping context, shared by all # reikna_norm_rng's _thr = _api.Thread.create() def reikna_norm_rng(seed, rows, cols, mean, std, dtype=np.float32, generators_dim=1): """ A do-all generator function for creating a new Computation returning a stream of pseudorandom number arrays. """ # The Type of the output array randoms_arr = Type(dtype, (rows, cols)) # Shortcut for creating a Sampler for normally distributed # random numbers rng = CBRNG.normal_bm(randoms_arr=randoms_arr, generators_dim=generators_dim, sampler_kwds=dict(mean=mean, std=std), # Reikna expects a positive integer. This is a # quick and *dirty* solution. seed=abs(hash(seed))) compiled_comp = rng.compile(_thr) # RNG "state" counters = _thr.to_device(rng.create_counters()) # Output array randoms = _thr.empty_like(compiled_comp.parameter.randoms) while True: compiled_comp(counters, randoms) yield randoms.get()
To see it in action add:
if __name__ == '__main__': seed1 = '123' seed2 = 'asd' rows, cols = 100, 3 mu, sd = 0, 1 # This will create a new RNG, take 1 result and throw the RNG away r1 = next(reikna_norm_rng(seed1, rows, cols, mu, sd)) r2 = next(reikna_norm_rng(seed2, rows, cols, mu, sd)) r3 = next(reikna_norm_rng(seed1, rows, cols, mu, sd)) r4 = next(reikna_norm_rng(seed2, rows, cols, mu, sd)) err1 = r1 - r3 err2 = r2 - r4 print("all(err1 == 0):", np.all(err1 == 0)) print("all(err2 == 0):", np.all(err2 == 0))
I also want to make sure that randoms generated using the two different seeds have no correlation.
This depends on the quality of the implementation. Here's how 2 sets of numbers with seeds 0 and 1 plotted out:
rng1 = reikna_norm_rng(0, 100, 10000, 0, 1) rng2 = reikna_norm_rng(1, 100, 10000, 0, 1) A = next(rng1) B = next(rng2) A_r = A.ravel() B_r = B.ravel() for i in range(0, A_r.size, 1000): plot(A_r[i:i+1000], B_r[i:i+1000], 'b.')
This is my first run in with reikna. The above code may not release resources timely and/or leak like a sieve. It uses a global Thread, which may not be what you want in a larger application.
np.random.seed(seed) np.random.normal(0, 1, (100, 3))
also produces arrays of normally distributed random numbers with shape (100, 3), though it does not employ the GPU.