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.

Answers


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.')

Disclaimer

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.

PS

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.


Need Your Help

How can I bind different property types with Bond framework?

ios swift

I have a label (string) and my object has a property with NSNumber.

Long Polling in Django- Return the response when there occurs an activity in the DB

jquery python django long-polling

I've read on SO regarding Long Polling in Django, but my problem is not a complex one which requires using tornado or building a chat application.