Improve random numbers functionality and documentation

[navdeeprana]

While it's a joy to use CUDA.jl for GPU computing with Julia, I find the documentation for working with random numbers a bit confusing and lacking, in particular, on how to draw reproducible random numbers. Here's the list of specific points which can be improved. I am willing to submit a PR.

• rand!(array) or randn!(array) will use CURAND by default and not CUDA.default_rng() but it is not mentioned in the docs in Array programming.
• While launching a kernel, a unique seed is passed from host to the CUDA.default_rng() to ensure kernels draw different random numbers on multiple invocations. This is documented, but it uses Random.default_rng() to make seed, which means that one way to have reproducible results is to seed the TaskLocalRNG(), which I found amusing.
• Manual seeding is possible, seed!() but it seeds both CUDA.default_rng() and CURAND.default_rng() at the same time?
• I can produce reproducible random number arrays by manually seeding CUDA.default_rng() as

rng = CUDA.default_rng()
seed!(rng, 12345)
randn(rng, 10)

but inside kernels it does not work (because it is re-seeeded?). So

function _kernel(p)
    i = threadIdx().x
    p[i] = randn()
    return nothing
end
Random.seed!(rng, 12345)
@cuda threads=length(p) _kernel(p)

doesn't give the same random numbers. However, if I seed the TaskLocalRNG(), it does.

Random.seed!(Random.default_rng(), 12345)
@cuda threads=length(p) _kernel(p)

• Similarly this should give same random numbers using CURAND, but it does not?

CUDA.seed!(12345)
randn(4)

• It's unclear if it is possible to use CURAND backed rngs inside kernels.

I went through the following discourse post:
https://discourse.julialang.org/t/kernel-random-numbers-generation-entropy-randomness-issues/105637/4

[maleadt]

Yes, those are fair points that probably warrant some additional documentation. Feel free to open a PR with improvements you see fit. Some comments:

rand!(array) or randn!(array) will use CURAND by default and not CUDA.default_rng() but it is not mentioned in the docs in Array programming.

This is for performance reasons; CURAND is faster than the CUDA.jl native RNG last I checked.

While launching a kernel, a unique seed is passed from host to the CUDA.default_rng() to ensure kernels draw different random numbers on multiple invocations. This is documented, but it uses Random.default_rng() to make seed, which means that one way to have reproducible results is to seed the TaskLocalRNG(), which I found amusing.

Is there an alternative? I contemplated using a specific RNG inside of CUDA.jl, but people commented that it would be even more unexpected. Not seeding at all would results in kernels returning the same numbers when invoked close together temporally (as they would get seeded by the clock() function on the device), while now at least there's some level of determinism for kernels launched in sequence.

Manual seeding is possible, seed!() but it seeds both CUDA.default_rng() and CURAND.default_rng() at the same time?

Right, because users typically expect seeding to affect whatever RNG would end up using. Note that you listed this right after a point of the device-side RNG, so it needs to be said that you can separately seed that RNG by calling seed! device-side, which does not affect multiple RNGs. That also applies to your next point.

[navdeeprana]

Thanks @maleadt. I will soon open a PR.

Yes, I was thinking on similar lines for an alternative, i.e., having a specific host-side RNG that can be controlled instead of the global RNG. I guess both options are okay if they are documented?

The behavior of seed!() is also okay when documented, I think.

Can you also shed some light on my last two points? On the reproducibility and the possibility of using CURAND generators inside kernels?