Вычисление объема n-шара методом Монте-Карло в Джулии

Хорошо, есть несколько бросающихся в глаза ошибок. У вас на один end больше, чем нужно, и строка

randoms(j) = rand(radius * 2) - radius;

имеет пару очевидных ошибок, потому что (1) вы индексируете массивы с [] в julia, а не (), и (2) rand(n) создает n случайные числа, а не случайные числа между 0 и 0 и rand(n) n>. Итак, если мы исправим эти проблемы

function MonteCarloHypersphereVolume(radius, dimension, number_of_generations)
    
    number_within_sphere = 0
    for i = 1 : number_of_generations
        randoms = zeros( 1, dimension )
        for j = 1 : dimension 
            randoms[j] = 2*radius*rand() - radius
        end

        if sum( randoms .^ 2 ) <= radius^2
            number_within_sphere = number_within_sphere + 1
        end
    end

    approximate_volume = (number_within_sphere/number_of_generations) * (2*radius)^dimension;

end

Теперь у нас есть подход, который, кажется, дает нам правильные числа, если мы тестируем простой трехмерный случай

julia> MonteCarloHypersphereVolume(1,3,100000)
4.18328

julia> 4/3*pi*1
4.1887902047863905

Теперь мы можем сделать несколько оптимизаций. В настоящее время это занимает около 13 мс и огромное количество распределений

julia> using BenchmarkTools

julia> @benchmark MonteCarloHypersphereVolume(1,3,100000)
BechmarkTools.Trial: 375 samples with 1 evaluations.
 Range (min … max):  11.016 ms … 23.764 ms  ┊ GC (min … max): 0.00% … 12.36%
 Time  (median):     12.893 ms              ┊ GC (median):    0.00%
 Time  (mean ± σ):   13.339 ms ±  1.647 ms  ┊ GC (mean ± σ):  4.24% ±  4.51%

      ▄▃▂█  ▆▇▅▇▁▁ ▂▅ ▁                                        
  ▃▄█▄█████▇███████████▆▇▆▅▇▅█▅▅▄▃▆▃▄▄▄▃▁▃▃▃▁▁▃▄▃▁▃▁▁▁▃▁▁▃▁▃▃ ▄
  11 ms           Histogram: frequency by time        18.9 ms <

 Memory estimate: 21.36 MiB, allocs estimate: 200000.

Мы можем сделать намного лучше, переместив randoms = zeros( 1, dimension ) за пределы цикла и используя sum с функцией, которая не должна выделять промежуточный массив для квадратов. Мы также можем использовать @inbounds, чтобы немного ускорить циклы for, учитывая, что здесь не следует ожидать индексов за пределами границ.

function montecarlohyperspherevolume(radius, dimension, number_of_generations)
    
    number_within_sphere = 0
    randoms = zeros( 1, dimension )
    @inbounds for i = 1 : number_of_generations
        for j = 1 : dimension 
            randoms[j] = 2*radius*rand() - radius
        end

        if sum(abs2, randoms) <= radius^2 # could also use anonymous function x->x^2 if you didn't wnat to use abs2
            number_within_sphere += 1
        end
    end

    return approximate_volume = (number_within_sphere/number_of_generations) * (2*radius)^dimension

end

Я также внес несколько стилистических изменений. Среди прочего, функции в Julia обычно пишутся строчными буквами по соглашению (в то время как типы могут быть CamelCase).

В любом случае это дает нам примерно 6-кратное ускорение и значительно снижает использование и выделение памяти

julia> @benchmark montecarlohyperspherevolume(1,3,100000)
BechmarkTools.Trial: 2249 samples with 1 evaluations.
 Range (min … max):  1.846 ms …   4.053 ms  ┊ GC (min … max): 0.00% … 0.00%
 Time  (median):     2.080 ms               ┊ GC (median):    0.00%
 Time  (mean ± σ):   2.217 ms ± 377.123 μs  ┊ GC (mean ± σ):  0.00% ± 0.00%

   █▆▅▂▁▂                                                      
  ▅████████▇▆▅▅▅▃▃▃▃▃▃▄▃▃▂▃▂▃▂▂▃▂▂▂▂▂▁▁▂▁▂▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁▁ ▂
  1.85 ms         Histogram: frequency by time        3.61 ms <

 Memory estimate: 112 bytes, allocs estimate: 1.