Note that, in Julia 1.7, @threads’s scheduling is context-dependent and thus the task tree is not quite predictable. If you use it @threads in, e.g., at the “root task,” you get multiple tasks (on multiple worker threads):
However, inside of another @threads, it uses single task (on single worker thread):
julia> Threads.@threads for i in 1:Threads.nthreads()
if i == 1
Threads.@threads for j in 1:Threads.nthreads()
@show j => current_task()
end
end
end
j => current_task() = 1 => Task (runnable) @0x00007fb46e331880
j => current_task() = 2 => Task (runnable) @0x00007fb46e331880
j => current_task() = 3 => Task (runnable) @0x00007fb46e331880
j => current_task() = 4 => Task (runnable) @0x00007fb46e331880
Thus, it implies that Random.seed!(0); f() may not have the same result, depending on the context. This is fixed in Julia 1.8 thanks to PRs such as
As for manual task creation that @lmiq is mentioning, I think the most important benefit is that it lets you have multiple streams of RNG and also with arbitrary algorithms.
Thanks a lot, @tkf@lmiq@fredrikekre@stevengj I now understand the global rng and how I could use the “task pool pattern” mentioned by @lmiq . I have one more question. I used the pattern in the opening post
rngs = [make_rng(tid) for tid in 1:Threads.nthreads()]
quite a bit in old code that I don’t want to rewrite.
Is this pattern safe in Julia 1.7?
Is there a super simple way to make it safe? Say by annotating each loop like @threads :static for ... or something?
You could use task-local storage, e.g. something like:
rng = get!(task_local_storage(), _RNGS_KEY) do
MersenneTwister()
end::MersenneTwister
where const _RNGS_KEY = gensym(:rngs)
But since I think new tasks are spawned every time you do @threads, it’s not clear what the advantage of this would be over just rng = MersenneTwister() in the type of usage you seem like you’re talking about.
I actually recommend this, if you just want a very easy and conservative approach. There’s no change in Julia 1.7 about this but it could have been invalid already in 1.5 and later if you happen to start using old code in nested tasks. Furthermore, there are changes in the upcoming 1.8 which makes old programs invalid in more cases.
But the short answer is that this is because the following code may throw in 1.8 and later
Threads.@threads for _ in 1:Threads.nthreads()
i = Threads.threadid()
yield()
j = Threads.threadid()
@assert i == j
end
The problem is rather how this is used. If you access it using “task id” you issued yourself, using a pattern like Random numbers and threads - #20 by lmiq, there’s no problem. But if you use threadid(), there is a problem in 1.8 or later as I mentioned above. If you use rngs as rand(rngs[threadid()]) or something similar, it may be fine (but it depends on the RNG type implementation). But, if you cache rngs access in a variable, it’d be problematic:
Threads.@threads for _ in 1:Threads.nthreads()
rng = rngs[Threads.threadid()]
f() # some function that may contain a yield point (e.g., printing, locking, ...)
rand(rng) # not OK because we may have `rng !== rngs[Threads.threadid()]`
end
The rule of thumb for good high-level multi-tasking code in Julia is “don’t talk about threads” although this attitude is rather too idealistic ATM.
I would like to mention that you shouldn’t combine the random number sequence of Mersenne Twister RNGs initialized with different seeds. The theoretical guarantees of Mersenne Twisters don’t work in those settings. Instead, it is necessary to use RNGs specialized for parallel RNG.
Some typical choices for parallel RNG are the PCG and Random123 families. These RNGs work by having another layer of “seeds” known as keys. You start from a single seed, and generate multiple RNG sequences by passing a key for each parallel sequence.
While both are not “cryptographically secure”, they both generate good quality random numbers with much less computation compared to Mersenne Twisters. They are also quite simple to use. Personally, I prefer Random123, and there is a nice implementation in Julia.
No, what I’m saying is that you shouldn’t use Mersenne Twisters in parallel even if they are thread safe. It’s “mathematically” unsafe. This is a very common mistake that is even stated in Wikipedia:
Multiple Mersenne Twister instances that differ only in seed value (but not other parameters) are not generally appropriate for Monte-Carlo simulations that require independent random number generators, though there exists a method for choosing multiple sets of parameter values.
No, it’s a problem with any RNG that is not specifically designed to be used in parallel. To do proper parallel RNG generation, you have to use RNGs specialized for parallel generation (like Random123) in the correct way.
Are the arguments surrounding Mersenne Twister even applicable to Xoshiro? According to the authors, they do provide jump polynomials for this exact use case:
All generators, being based on linear recurrences, provide jump functions that make it possible to simulate any number of calls to the next-state function in constant time, once a suitable jump polynomial has been computed. We provide ready-made jump functions for a number of calls equal to the square root of the period, to make it easy generating non-overlapping sequences for parallel computations, and equal to the cube of the fourth root of the period, to make it possible to generate independent sequences on different parallel processors.
Though admittedly, I don’t think they are exposed in julia. Maybe they should be added to the Random stdlib.
That comes as a surprise (something that can be optimized, or at least faster in all real-word code?). Xoshiro supposed to be fast, it’s short code, with much smaller state than MersenneTwister.
Right, as of 1.7.
If you mention Xoshiro by name in your code, it’s no longer compatible with Julia 1.6 LTS. If you just use a global thread, I suppose you can have the same code working in both 1.7/1.8 and 1.6 LTS, but might there be a way for 1.7 code using threads to keep some (partial) compatibility with 1.6? In both cases the stream of numbers will be different.
[If you want to target 1.6 and later, then you can use MersenneTwister, that’s I believe why it was kept around, for compatibility, for those who mentioned it my name (thus will get same stream).]