pi calculation

References

パフォーマンス一覧にSASが無かったのでやってみた

http://prunus1350.hatenablog.com/entry/2015/01/04/195603

Juliaのコードを勝手に最適化してみた

http://bicycle1885.hatenablog.com/entry/2015/01/04/175916

Juliaのススメ

http://www.slideshare.net/beeEaMa/lt-20140515

Code

Reference naive Julia code

tic()
circle_in = 0.0
for i in 1:100000000
    l = (rand()^2 + rand()^2) ^ 0.5
    if l <= 1
        circle_in = circle_in + 1
    end
end
println(4 * circle_in / 100000000)
toc()

R code

cat("### Naive\n")

## ### Naive

system.time({
    n = 100000000
    mat <- matrix(runif(n = 2*n), ncol = 2)
    print(4*sum((mat[,1]^2 + mat[,2]^2)^(0.5) <= 1)/n)
})

## [1] 3.141521

##    user  system elapsed 
##  17.721   2.142  19.871

cat("### sqrt\n")

## ### sqrt

system.time({
    n = 100000000
    mat <- matrix(runif(n = 2*n), ncol = 2)
    print(4*sum(sqrt(mat[,1]^2 + mat[,2]^2) <= 1)/n)
})

## [1] 3.141514

##    user  system elapsed 
##  13.903   2.204  16.119

cat("### no sqrt\n")

## ### no sqrt

system.time({
    n = 100000000
    mat <- matrix(runif(n = 2*n), ncol = 2)
    print(4*sum(mat[,1]^2 + mat[,2]^2 <= 1)/n)
})

## [1] 3.141891

##    user  system elapsed 
##  13.106   2.194  15.305

cat("### no sqrt\n")

## ### no sqrt

system.time({
    n = 100000000
    x <- runif(n = n)
    y <- runif(n = n)
    print(4*sum(x^2 + y^2 <= 1)/n)
})

## [1] 3.141777

##    user  system elapsed 
##   8.833   1.158  10.153