While it seems to work and the output looks sensible, I have this lingering suspicion that rand() is not threadsafe and that using it in a Thread.@threads loop is volatile.
Browsing through base/random.jl and base/dSFMT.jl I saw no obvious sign of thread protection, which makes me wonder why it seems to work so well in the first place. Is dSFMT itself threadsafe maybe? I did not find any info on that.
So my question is if it is safe to call rand() in a Threads.@threads loop, and if not what the encourage way of doing so is.
Something like this?
julia> const trand_mutex = Threads.Mutex()
Base.Threads.Mutex(0, Ptr{Void} @0x00007ff8949c3fc0)
julia> function trand(args...)
lock(trand_mutex)
result = rand(args...)
unlock(trand_mutex)
result
end
trand (generic function with 1 method)
Won’t this make every thread have identical random numbers?
@noinline function init_thread_rng()
# Allocate the random number generator on the thread's own heap lazily
# instead of the master thread heap to minimize memory conflict.
ThreadRNG[Threads.threadid()] = MersenneTwister(0)
end
I assume for separate random numbers per thread, it’d be fine to initialize with MersenneTwister(Threads.threadid())?
My understanding of seeding is that method (of using consecutive small integers) isn’t necessarily any random than MersenneTwister(0), but that you can jump one RNG ahead to effectively split it into two (or N) independent streams (https://github.com/JuliaLang/julia/pull/12498). I could be wrong though.
I understand that the default jump in julia is 10^20.
So if I need 36 random generators for 36 processes and on each thread I anticipate to generate around 10^13 random numbers. Am I safe by using the code below?
Will each MersennTwister be 10^20 ‘apart’ from the next one? (the documentation is not too easy to read)?
seed=2001
mta = MersenneTwister(seed)
mts = randjump(mta, 36)
# -> use mts[x] on process x (for x in 1:36), rand(mts[x]) respectively
Oh. My bad.
This is embarrassing, I must have confused 10^6001 (which should be close to 2^19937) with 10^61.
I guess this can only happen to someone with no intuition for power series.
You are correct.
Seems like I have plenty of ‘room’ available with 36*20=720 << 6001
That is to say: It is my understanding that the second element of the array would effectively start at the ‘random element’ 1+10^20 which the first MersenneTwister of the array produces.