julia - 蒙特卡洛 - 模拟布丰投针实验

文章和代码参考：MATLAB模拟布丰投针实验

本文适合萌新学习蒙特卡洛入门

using Distributions
#布丰投针实验
function fn(n::Int64)
    """
    n:每次实验投针的次数（n次投针）
    """
    l = 0.6      #针的长度
    a = 1        #平行线的宽度    
    k = 0        #记录针与平行线相交的次数
    y = rand(Uniform(0.0, a/2), n)   #在[0, a/2]内服从均匀分布随机产生n个数
    theta = rand(Uniform(0.0,π), n) #在[0, pi]内服从均匀分布随机产生n个数
    for i=1:n
        if y[i] < (l/2)*sin(theta[i]) 
            k = k + 1
        end
    end
    f = k / n
    Pi = (2*l*n)/(a*k)
    return Pi
end

#测试代码
fn(1000) |> println
fn(10000000) |> println

3.076923076923077
3.139287847319689

#重复实验
function main(times::Int64,n::Int64)
    """
    times:实验重复的次数
    n:每次实验投针的次数（做n次投针）
    """
    res = Array{Float64}(undef,times)
    for j = 1:times      
        res[j] = fn(n)
    end
    return res
end

@time res = main(10000,10000)
res |> describe

2.048339 seconds (53.30 k allocations: 1.492 GiB, 6.01% gc time)
Summary Stats:
Length: 10000
Missing Count: 0
Mean: 3.141537 #所计算的10000个 Pi值求平均值
Minimum: 2.994012
1st Quartile: 3.114456
Median: 3.140539
3rd Quartile: 3.168108
Maximum: 3.298516
Type: Float64

总结：蒙特卡洛计算量比较大，本文如果进行main(10000,1000000)【一万次实验，每次投针一百万次】，耗时240秒，解决方案：多线程，多进程，Gpu计算