# Brew a Parallel RNG?

**URL:** <https://discourse.julialang.org/t/brew-a-parallel-rng/30286>\
**Category:** Performance\
**Created:** [October 25, 2019, 5:10am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286 "2019-10-25T05:10:11Z")\
**Posts on this page:** 20\
**Page:** 1

<div class="post-metadata">

**Author:** ![HaoLi111](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/haoli111/32/10982_2.png) [@HaoLi111](https://discourse.julialang.org/u/HaoLi111)\
**Post date:** [October 25, 2019, 5:10am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/1 "2019-10-25T05:10:11Z")

</div>

RNGs cannot be parallelized without modification, as the effectiveness of the randomness representation ceases if the sequences are related.  
Previous discussions have been found from

> [@Parallel Mersenne Twister](https://discourse.julialang.org/t/parallel-mersenne-twister/27567):
>
> Like many people, I’m pretty excited about the new multithreading capabilities in 1.3. But there’s one thing in [the announcement](https://julialang.org/blog/2019/07/multithreading) that I’m not sure about, so it seems worth some discussion. It says, The approach we’ve taken with Julia’s default global random number generator ( rand() and friends) is to make it thread-specific. On first use, each thread will create an independent instance of the default RNG type (currently MersenneTwister ) seeded from system entropy. I haven’t worked much wi…

But can we use some alternative solutions?

Consider some ancient algorithms that are probably not in use nowadays(most written from the book Numerical Analysis by T. Sauer)

```julia
# rng

#random number generators
#Computers are not capable of generating true random numbers
#but can generate sequences with statistical randomness

# 1 Julia's implemented solution
rand(100)

# LCP psuedo solutions

function rng_Park_Miller1998(n;x1=3)
    x=zeros(n)
    x[1] = x1
    u = zeros(n)
    for i in 2:n
        x[i] =mod(16807*x[i-1],2147483647)#7^5 \ 2^31-1 31th mason prime
    end
    u = x./2147483647
    return u
end
# implemented in Matlab 4 1990

rng_Park_Miller1998(100;x1=3)

function rng_randu(n;x1=3)
    x = zeros(n)
    x[1]=x1
    for i in 2:n
        x[i] = mod(65539*x[i-1],2147483648)
    end
    u = x./2147483648
    return u
end

rng_randu(100)

# iterative psudo solutions

function rng_LMap(n;x1=.4)
    x = zeros(n)
    x[1] = x1
    for i in 2:n
        x[i]=1-2*x[i-1]^2
    end
    return((x./2) .+.5)
end

rng_LMap(100)

# self - avoiding solutio

function rng_Halton(n;x1 = 3)
    b = zeros(Int(ceil(log(n)/log(x1))))
    u = zeros(n)
    for j in 1:n
        i = 1
        b[1] = b[1]+1
        while b[i]>(x1-1+eps())
            b[i]=0
            i=i+1
            b[i]=b[i]+1
        end

        u[j]=0
        for k=1:length(b)
            u[j] = u[j]+b[k]*Float64(x1)^(-k);
        end
    end
    return u
end

rng_Halton(100)

```

It seems not hard if a random seq can be generated with a “root”, so that if we pre-calculate some of the terms as starting values for each thread, then they may end up independent. Is this mistaken？

BTW the code is from my repository

> <https://github.com/HaoLi111/Julia_Numerical_Recipe/blob/master/RNG.jl>

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [October 25, 2019, 5:16am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/2 "2019-10-25T05:16:11Z")

</div>

I may be missing something, but I was under the impression that if you either

1. initialize all the parallel RNGs with a truly random seed (eg from hardware entropy), or

2. use `randjump` with a conservatively large step,

these issues are mitigated. See eg

[https://github.com/JuliaLang/julia/issues/94#issuecomment-515806131](https://github.com/JuliaLang/julia/issues/94#issuecomment-515806131)

---

<div class="post-metadata">

**Author:** ![HaoLi111](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/haoli111/32/10982_2.png) [@HaoLi111](https://discourse.julialang.org/u/HaoLi111)\
**Post date:** [October 25, 2019, 5:23am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/3 "2019-10-25T05:23:09Z")

</div>

Thank you! That’s a very neat solution.

---

<div class="post-metadata">

**Author:** ![cscherrer](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/cscherrer/32/7631_2.png) [@cscherrer](https://discourse.julialang.org/u/cscherrer)\
**Post date:** [October 25, 2019, 1:10pm UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/4 "2019-10-25T13:10:04Z")

</div>

Not sure why I didn’t think of this before, but there’s also [SPRNG](http://www.sprng.org/), which is designed specifically for this use case.

---

<div class="post-metadata">

**Author:** ![tkf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkf/32/17635_2.png) [@tkf](https://discourse.julialang.org/u/tkf)\
**Post date:** [October 26, 2019, 7:46am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/5 "2019-10-26T07:46:25Z")

</div>

I fond what Guy Steele mentioned in the last half of this talk interesting:

[![](https://global.discourse-cdn.com/julialang/original/3X/c/e/ce0b1e884ff882ad64dbd60f1729f65ccc3ffd4e.jpeg "Keynote with Guy Steele") ](https://www.youtube.com/watch?v=0hlBkQ5DjaY)

IIRC there is a type of SRNG that can be split in such a way that the distribution of resulting streams combined is equivalent to the original stream (or something like that). It would be nice to look into (probably some follow up studies of):

[http://supertech.lcs.mit.edu/papers/dprng.pdf](http://supertech.lcs.mit.edu/papers/dprng.pdf)

> **[Parallel random numbers: As easy as 1, 2, 3](https://ieeexplore.ieee.org/abstract/document/6114424)**
>
> Most pseudorandom number generators (PRNGs) scale poorly to massively parallel high-performance computation because they are designed as sequentially dependent state transformations. We demonstrate that independent, keyed transformations of counters...

which were mentioned in the talk.

My understanding is that this would guarantee that parallel computation result with SRNG to be deterministic even when the scheduler is dynamic (as in Julia), ~~if you construct the computation as some combination of map/filter/reduce/etc.~~ (edit: hmm… wait, does it work for any computation?)

---

<div class="post-metadata">

**Author:** ![jling](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jling/32/212909_2.png) [@jling](https://discourse.julialang.org/u/jling)\
**Post date:** [October 26, 2019, 8:06am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/6 "2019-10-26T08:06:00Z")

</div>

also in the coming up 1.3

> All operations that affect the random number state ( `rand` , `randn` , `seed!` , etc.) will then operate on only the current thread’s RNG state. This way, multiple independent code sequences that seed and then use random numbers will individually work as expected.

---

<div class="post-metadata">

**Author:** ![Elrod](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/elrod/32/22461_2.png) [@Elrod](https://discourse.julialang.org/u/Elrod)\
**Post date:** [October 26, 2019, 7:32pm UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/7 "2019-10-26T19:32:40Z")

</div>

Using different multipliers in a linear congruential generator (or a PCG, which uses an LCG to update state) will create unique, non-overlapping streams. You can take advantage of this for SIMD as well.  
This is what I have done with [VectorizedRNG](https://github.com/chriselrod/VectorizedRNG.jl/blob/master/src/multipliers.jl), and [ProbabilityModels](https://github.com/chriselrod/ProbabilityModels.jl/blob/master/src/rng.jl#L3), which initializes one per thread.

You can use [RNGTest](https://github.com/andreasnoack/RNGTest.jl) for validating an RNG or LCG multiplier (note that an LCG itself will fail a lot of tests, a PCG is much better in this regard).

---

<div class="post-metadata">

**Author:** ![tkf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkf/32/17635_2.png) [@tkf](https://discourse.julialang.org/u/tkf)\
**Post date:** [October 27, 2019, 2:09am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/8 "2019-10-27T02:09:39Z")

</div>

My understanding is that making RNG deterministic _and_ parallel is hard when the task scheduler is dynamic as in Julia \>= 1.3 (hence my last comment). Is there a packaged solution for Julia? Or maybe people actually don’t need it?

---

<div class="post-metadata">

**Author:** ![HaoLi111](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/haoli111/32/10982_2.png) [@HaoLi111](https://discourse.julialang.org/u/HaoLi111)\
**Post date:** [October 27, 2019, 5:17am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/9 "2019-10-27T05:17:57Z")

</div>

Yes. I am thinking the same. We have to specify scenarios.  
If the n on each core is not known then jumpers would not work…

---

<div class="post-metadata">

**Author:** ![HaoLi111](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/haoli111/32/10982_2.png) [@HaoLi111](https://discourse.julialang.org/u/HaoLi111)\
**Post date:** [October 27, 2019, 5:22am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/10 "2019-10-27T05:22:48Z")

</div>

Or should we just simply calculate the length of random sequences applied in each threads and apply the jumpers afterwards 😂  
Well…😅 that’s a sacrifice of the performance

---

<div class="post-metadata">

**Author:** ![HaoLi111](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/haoli111/32/10982_2.png) [@HaoLi111](https://discourse.julialang.org/u/HaoLi111)\
**Post date:** [October 27, 2019, 5:30am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/11 "2019-10-27T05:30:24Z")

</div>

Thanks!  
So does this mean the ‘jumper’ is already implemented or does it use somewhat different techniques to ensure the between-threads independence? So that  
If I  
Draw same length of sequences from each thread and then join them together  
Or if I  
Draw sequences from random length from each sequence and then join them together.  
Will the results be a good representation of random sequences?  
I guess the previous one is easier to solve. Plz point out if I am wrong

---

<div class="post-metadata">

**Author:** ![jling](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jling/32/212909_2.png) [@jling](https://discourse.julialang.org/u/jling)\
**Post date:** [October 27, 2019, 5:35am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/12 "2019-10-27T05:35:15Z")

</div>

If you just want thread-independent random number (or, sequence), `1.3` should have solved the problem for you (i.e this would be the default behavior). So you can just generate random sequences in different threads and combine afterwards if that’s what you want.

> [@HaoLi111](#):
>
> Draw sequences from random length

btw, I don’t quite understand what this means.

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [October 27, 2019, 5:48am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/13 "2019-10-27T05:48:56Z")

</div>

> [@tkf](#):
>
> making RNG deterministic

I am not sure I understand why one would want to do that.

---

<div class="post-metadata">

**Author:** ![tkf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkf/32/17635_2.png) [@tkf](https://discourse.julialang.org/u/tkf)\
**Post date:** [October 27, 2019, 5:53am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/14 "2019-10-27T05:53:36Z")

</div>

I think it’s highly useful that `sum(rand(MersenneTwister(0), 10))` returns the same value as long as I’m using same `julia` version and setting.

---

<div class="post-metadata">

**Author:** ![Tamas\_Papp](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tamas_papp/32/25949_2.png) [@Tamas\_Papp](https://discourse.julialang.org/u/Tamas_Papp)\
**Post date:** [October 27, 2019, 6:04am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/15 "2019-10-27T06:04:37Z")

</div>

I understand what it implies, but now how that is it so useful.

For example, I have seen people rely on this for testing, but I don’t think it is the right approach (as it is quite fragile to details of the algorithm).

Because getting this is quite difficult, before we focus on this as a desirable property of an RNG in a parallel setting, it would be great to understand the rationale.

---

<div class="post-metadata">

**Author:** ![jling](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/jling/32/212909_2.png) [@jling](https://discourse.julialang.org/u/jling)\
**Post date:** [October 27, 2019, 6:29am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/16 "2019-10-27T06:29:41Z")

</div>

the demand(?) maybe more common than you think, for example, I think most CMS(physics) papers would use `seed=123456` for quite a few Monte Carole/Inference tools to make the result ‘reproducible’(for debugging and peer review purpose, I only could imagine)

---

<div class="post-metadata">

**Author:** ![HaoLi111](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/haoli111/32/10982_2.png) [@HaoLi111](https://discourse.julialang.org/u/HaoLi111)\
**Post date:** [October 27, 2019, 6:53am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/18 "2019-10-27T06:53:11Z")

</div>

Agree

---

<div class="post-metadata">

**Author:** ![tkf](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/tkf/32/17635_2.png) [@tkf](https://discourse.julialang.org/u/tkf)\
**Post date:** [October 27, 2019, 6:57am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/19 "2019-10-27T06:57:06Z")

</div>

How about debugging some input-dependent bugs? If you have a deterministic program, you can reproduce it quite easily and re-do same execution as many as you want. If you have parallel code that depends on the exact scheduling of the tasks, I can imagine that it is extremely hard (impossible?) to reproduce the same execution reliably.

> [@Tamas\_Papp](#):
>
> For example, I have seen people rely on this for testing, but I don’t think it is the right approach (as it is quite fragile to details of the algorithm).
> 
> Because getting this is quite difficult

I do agree that relying on the exact output of RNG-depending code is not the right approach in the testing. That’s because getting same output from code using RNG _across library versions_ is hard. But I don’t think it is directly relevant for getting the same result of the same program with the same input if it only does pure computation. Also, there seems to be a way to do it already. Implementing it probably is a very hard work but I don’t think it is as difficult as improving libraries using RNG without changing the exact output.

---

<div class="post-metadata">

**Author:** ![HaoLi111](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/haoli111/32/10982_2.png) [@HaoLi111](https://discourse.julialang.org/u/HaoLi111)\
**Post date:** [October 27, 2019, 7:11am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/20 "2019-10-27T07:11:33Z")

</div>

That means if I don’t know the exact lengths of each random sequences that is called within each thread.  
For example thread one has a process requiring 100 random numbers, thread 2 requiring 120. Both of them are good representations of randomness for themselves as separated event lists -stochastic processes. But if I combine these they might not practically be good representations of true random as a whole.

Is this property algorithm dependent?Maybe I should read about pcg or some other algorithms?

---

<div class="post-metadata">

**Author:** ![HaoLi111](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/haoli111/32/10982_2.png) [@HaoLi111](https://discourse.julialang.org/u/HaoLi111)\
**Post date:** [October 27, 2019, 7:20am UTC](https://discourse.julialang.org/t/brew-a-parallel-rng/30286/21 "2019-10-27T07:20:16Z")

</div>

And I think at some point if the jumpers enable us to code and run RNG on general graphic cards… well that’s a different story… maybe forget about this:joy:

[Next page](https://discourse.julialang.org/t/brew-a-parallel-rng/30286.md?page=2)