Advice for improving Monte-Carlo code

New to Julia
, ,
18

Multi-threading version, which uses random number generation is more involved, it was discussed here

In a simplified version it may look like this

function meandist3(n, rng)
    d = zeros(Threads.nthreads())
    Threads.@threads for i in 1:n
        k = Threads.threadid()
        x = (rand(rng[k]),rand(rng[k]))
        y = (rand(rng[k]),rand(rng[k]))
        d[k] += sqrt((x[1]-y[1])^2+(x[2]-y[2])^2)
    end
    return sum(d)/n
end

And benchmarks

@btime meandist(100_000_000)
  # 1.690 s (4 allocations: 2.98 GiB)

rngs = Tuple(StableRNG.(2020 .+ collect(1:Threads.nthreads())))
@btime meandist3(100_000_000, $rngs)
  # 1.060 s (43 allocations: 6.13 KiB)

It is important, because if you have less then 30Gb it is the only way to run

julia> @time meandist3(1_000_000_000, rngs)
 13.683740 seconds (44 allocations: 6.141 KiB)
0.5214066471370079

EDIT: ha, interestingly enough, MersenneTwister suddenly is the best here

julia> rngs = Tuple(Random.MersenneTwister.(2020 .+ collect(1:Threads.nthreads())))
julia> @time meandist3(1_000_000_000, rngs)
  8.425097 seconds (48 allocations: 6.266 KiB)
0.521398839144987
2 Likes