[PYTHON] Speed comparison of each language by Monte Carlo method

Introduction

[Monte Carlo method](https://ja.wikipedia.org/wiki/%E3%83%A2%E3%83%B3%E3%83%86%E3%82%AB%E3%83%AB%E3% 83% AD% E6% B3% 95) The speed comparison of each language is performed by calculating the pi in each language using the Monte Carlo method. The Monte Carlo method is to drop an arbitrary (random number) point at 1.0> = x> = 0, 1.0> = y> = 0, and the dropped point falls within the circle from the origin (0,0). The ratio between the thing and the thing that fell outside the circle is calculated, and that ratio becomes the pi. Of course, although it approximates the pi obtained by calculation, it does not give an accurate value, but I think it makes sense to implement it in a programming language. The reason is -The speed of floating point calculation can be compared. ・ Loop speed can be compared. -Memory usage is low and paging does not (probably) occur. -Since DiskIO does not occur, it is not easily affected by external storage devices. -The correct pi is known information and is easy to debug. And so on.

June 28, 2020 ruby was added. June 28, 2020 Changed the calculation method of python. June 28, 2020 Corrected an error in the rust loop count specification. June 28, 2020 Additional measurement of Xoroshiro version source provided by @ tatsuya6502 for rust

Below, each language source perl

Since perl had a very high-speed image I have tried so far, I prepared it as a benchmark (≒ standard). Since the source is simple, it should work at high speed without any overhead due to object orientation.

montecarlo.pl


use Time::HiRes qw( gettimeofday tv_interval);

local $ave = 0;

local $time_start, $time_end, $time_split = 0;

local $max_times    = 10000000;
local $test_times   = 100;

local $cnt = 0;

for ( $i = 1 ; $i <= $test_times ; $i++ ) {

    print($i."\n");

    $time_start = [gettimeofday];
    $cnt = 0;
    for ( $j = 0 ; $j < $max_times ; $j++) {
        $x = rand();
        $y = rand();

        if ( $x * $x + $y * $y <= 1) {
            $cnt++;
        }
    }

    $time_split = tv_interval($time_start);
    print("pi ".(4 * $cnt / $max_times)."\n");
    print("split:".$time_split."\n");

    $ave += $time_split;
}

print("AVE:".($ave / $test_times)."\n");
~                                                                                                                    ~                                                                                                                    

c language

For the time being, I prepared it because "C is fast, isn't it?" It may be faster if you tune it tightly, I tried to suppress it at this level because of the ease of coding.

montecarlo.c


#include <stdio.h>
#include <stdlib.h>
#include <time.h>

float montecarlo();

void main() {

    float ave = 0;
    int times = 100;

    for (int i = 1 ; i <= times ; i++){
        printf("%03d\n", i );
        ave += montecarlo();
    }

    printf("AVE %f\n\n", ave / times);

}

float montecarlo() {

    clock_t time_start, time_end;
    double  time_split;
    int times = 10000000;
    float x, y;

    time_start = clock();

    int cnt = 0;
    for ( int i = 0; i < times ; i++) {
        x = (float)rand() / (float)RAND_MAX;
        y = (float)rand() / (float)RAND_MAX;

        if ( (x * x + y * y) <= 1) {
            cnt++;
        }
    }

    time_end = clock();

    time_split = (double)(time_end - time_start) / 1000000;
    printf("%f\n", 4.0 * cnt / times);
    printf("time_split : %lf\n", time_split);

    return(time_split);

}

go language

Since go had a good reputation before, I prepared it as a benchmark different from perl. Since it is a compilation language, how thin can it be in C language? I am also interested in that.

montecarlo.go


package main

import (
        "fmt"
        "time"
        "math/rand"
)

var i, j, cnt int
var x, y float64
var time_start, time_end int64
var ave, time_split float64

const MAX_TIMES = 10000000
const TEST_TIMES = 100

func main() {
        //fmt.Printf("%03d\n", 10)

        //fmt.Printf("%d\n", time.Now())
        //fmt.Printf("%f\n", float64(time.Now().UnixNano()) / 1000000000.0)
        //fmt.Printf("%f\n", float64(time.Now().UnixNano()) / 1000000000.0)
        //fmt.Printf("%f\n", float64(time.Now().UnixNano()) / 1000000000.0)

        rand.Seed(time.Now().UnixNano())

        ave := 0.0
        for i := 1 ; i <= TEST_TIMES ; i++ {

                fmt.Printf("%03d times\n", i)

                time_start      := time.Now().UnixNano()
                cnt := 0
                for j := 0 ; j <  MAX_TIMES ; j++ {

                        x := rand.Float64()
                        y := rand.Float64()

                        if x * x + y * y <= 1 {
                                cnt++
                        }

                }
                fmt.Printf("pi : %f\n", 4 * float64(cnt) / float64(MAX_TIMES))

                time_end        := time.Now().UnixNano()
                time_split      := float64(time_end - time_start) / 1000000000
                fmt.Printf("time : %f\n", time_split)
                ave += time_split
                //ave := time_split

        }

        fmt.Printf("AVE : %f\n", ave / TEST_TIMES)

}

Java

As a compilation language, I think it is currently positioned as an iron plate. (Although each has their own opinions ...) It's expected to be slower than C with VMs in between, but past experience gives the impression that it's not fatally slow.

montecarlo.java


import java.util.Random;
import java.util.*;


class Montecarlo {
    public static void main(String[] args) {

        //System.out.println( getNowTime() );

                //Specify the radius of the circle
                double iR = 1.0;

                //Specify the number of draws
                int MAX_POINT = 10000000;

        //Specify the number of tests
        int TEST_TIMES = 100;

                //Random object definition
                Random rand = new Random();

        //Start and end times
        long time_start, time_end = 0;

        double time_split;
        double ave = 0.0;

                //Definition of shunt variable for Random value
                double ranValue = 0.0;
                double POS_X = 0.0;
                double POS_Y = 0.0;

        for ( int j = 1 ; j <= TEST_TIMES ; j ++) {

            //System.out.println("%03d", j);
            System.out.printf("%03d\n", j);

            //Start time
            time_start = System.nanoTime();
                //System.out.println( "start : " + time_start );

            //Count the number of cases that fell in an arc
            int HIT_COUNT = 0;

            //Initialization of work variables for output
            for ( int i = 0 ; i < MAX_POINT ; i++ ) {

                POS_X = rand.nextDouble();
                POS_Y = rand.nextDouble();

                if ( iR >= POS_X * POS_X + POS_Y * POS_Y )  {
                    HIT_COUNT++;
                }

            }
            // System.out.println(HIT_COUNT );
            System.out.println( (  4 * HIT_COUNT  / (double)MAX_POINT ) );

            //End time
            time_end = System.nanoTime();

                //System.out.println( "end   : " + time_end );

            time_split = (double)(time_end - time_start) / 1000000000.0;

                //System.out.println( "split : " + (time_end - time_start) );
                System.out.println( "split : " + time_split );

            ave += time_split;

        }

        //System.out.println( getNowTime() );
        System.out.println( "AVE : " + ave / TEST_TIMES );

    }

        ///Current time yyyy/mm/dd hh:mm:ss(jpn)Get in format
        public static String getNowTime() {

                return(String.format("%1$tF %1$tT" ,Calendar.getInstance()));

        }

}

lua It has a reputation for being light as a script language, so I added it. As a language specification, it is an image that it is fast because it uses instructions close to the OS.

montecarlo.lua


local socket = require("socket")

math.randomseed(os.time())

MAX_LOOP = 10000000
TEST_TIMES = 100

ave = 0.0
for i = 1, TEST_TIMES do
    print(string.format("%03d times", i))
    time_start  = socket.gettime()

    cnt = 0
    for j = 1 , MAX_LOOP do
        x = math.random()
        y = math.random()

        if x * x + y * y <= 1 then
            cnt = cnt + 1
        end

    end

    time_end    = socket.gettime()

    time_split  = time_end - time_start

    print(string.format("pi    : %f", 4 * cnt / MAX_LOOP))
    print(string.format("split : %f", time_split))

    ave = ave + time_split

end

print(string.format("ave : %f", ave / TEST_TIMES))

python Needless to say, it's python. This time, we will compare with 2nd and 3rd series, and also with pypy.

montecarlo.py


import random
import time
import math

_times = 10000000

def pi():
    _time_start = time.time()

    _cnt = 0
    for i in range(_times):

        _x = random.random()
        _y = random.random()

        #if  _x ** 2 + _y **2 <= 1:
        if  _x * _x + _y * _y <= 1:

            _cnt += 1

    print( 4.0 * _cnt / _times )

    _time_end = time.time()

    _time_split = _time_end - _time_start
    print(_time_split)

    return _time_split

_test = 100

_ave = 0.0
for _i in range(_test):
    print('{:03d} times'.format(_i+1))
    _ave += pi()

print('ave: {} s'.format(_ave/_test))

Rust Rust is a language that is attracting attention in AI and machine learning. It's faster than python and easier to develop than C.

Cargo.toml


[package]
name = "montecarlo"
version = "0.1.0"
authors = ["ubuntu"]
edition = "2018"

# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html

[dependencies]
rand = "0.6"
floating-duration="0.1.2"

src/main.rs


fn main() {

    use rand::Rng;                  //For random numbers
    use std::time::Instant;         //For time acquisition
    use floating_duration::TimeAsFloat; //For conversion from time to float

    println!("Hello, world!");

    let loop_max = 10000000;        //Number of RBIs
    let test_max = 100;             //Loop to average

    let mut x: f32;
    let mut y: f32;

    let mut rng = rand::thread_rng();//randseed settings

    let mut split_time: f32 = 0.0;
    //for _test_count in 1..test_max {
    for _test_count in 1..=test_max {
        println!("times:{:?}", _test_count );
        let start_time = Instant::now();

        let mut count = 0;
        //for _loot_count in 1..loop_max {
        for _loot_count in 1..=loop_max {
            x = rng.gen_range(0., 1.);  //Random number 0.0~1.0
            y = rng.gen_range(0., 1.);

            if x * x + y * y <= 1.0 {
                count = count + 1;
            }

        }
        println!("pi:{}", 4.0 * count as f32 / loop_max as f32); //as f32 is a type cast
        let end_time = start_time.elapsed().as_fractional_secs();//Convert time to seconds
        println!("time:{:?}", end_time );

        split_time = split_time + end_time as f32;

    }

    println!("AVE:{}", split_time / test_max as f32);

}

Xoroshiro version source provided by @ tatsuya6502

Cargo.toml


[package]
name = "montecarlo"
version = "0.1.0"
authors = ["ubuntu"]
edition = "2018"

# See more keys and their definitions at https://doc.rust-lang.org/cargo/reference/manifest.html

[dependencies]
rand = "0.7"
rand_xoshiro = "0.4"

src/main.rs


//For random numbers
use rand::Rng;
use rand_xoshiro::{rand_core::SeedableRng, Xoroshiro128Plus};

//For time acquisition
use std::time::{Instant, SystemTime};

const LOOP_MAX: usize = 10_000_000; //Number of RBIs
const TEST_MAX: usize = 100; //Loop to average

fn main() -> Result<(), Box<dyn std::error::Error>> {
    //Set random seed
    let now = SystemTime::now().duration_since(SystemTime::UNIX_EPOCH)?;
    let mut rng = Xoroshiro128Plus::seed_from_u64(now.as_millis() as u64);

    let mut split_time: f32 = 0.0;
    for _test_count in 1..=TEST_MAX {
        println!("times:{}", _test_count);
        let start_time = Instant::now();

        let mut count = 0;
        for _loot_count in 0..LOOP_MAX {
            //Generate random numbers[0.0, 1.0)・ ・ 0.0 or more 1.Less than 0
            let x: f32 = rng.gen_range(0.0, 1.0);
            let y: f32 = rng.gen_range(0.0, 1.0);

            if x * x + y * y <= 1.0 {
                count += 1;
            }
        }
        println!("pi:{}", 4.0 * count as f32 / LOOP_MAX as f32); //as f32 is a type cast
        let end_time = start_time.elapsed().as_secs_f32(); //Convert time to seconds
        println!("time:{}", end_time);

        split_time += end_time as f32;
    }

    println!("AVE:{}", split_time / TEST_MAX as f32);
    Ok(())
}

Nim For comparison ...

montecarlo.nim


import random
import times

randomize()

const TEST_MAX = 10000000
const LOOP_MAX = 100

var x, y = 0.0

var ave = 0.0

for loop in 1..LOOP_MAX:
    echo $loop & " times"

    var time_start = cpuTime()
    var count = 0
    for test in 1..TEST_MAX:

        x = rand(0.0..1.0)
        y = rand(0.0..1.0)

        # echo $x
        # echo $y

        if ( x * x + y * y <= 1):
            count += 1

    echo "pi:" &  repr(float(4.0 * float(count) / float(TEST_MAX)))

    var time_split = cpuTime() - time_start

    echo "split" & repr(time_split)
    ave += time_split

echo "AVE:" & repr(ave / float(LOOP_MAX))
                                                                                                                    ~                                                                                                                    

Ruby I added ruby because @scivola provided the source.

montecarlo.rb


TIMES = 10000000

def pi
  time_start = Time.now

  cnt = 0
  TIMES.times do
    x = rand
    y = rand

    if  x * x + y * y <= 1
      cnt += 1
    end
  end

  puts 4.0 * cnt / TIMES

  time_end = Time.now

  time_split = time_end - time_start
  puts time_split

   time_split
end

test = 100

ave = 0.0

test.times do |i|
    puts "%03d times" % (i + 1)
    ave += pi
end

puts "ave: #{ ave / test } s"

result

perl 5.30
	AVE:8.93247789s

java 14.0.1
	AVE: 1.4161036380199996s

c
	gcc 9.3.0

	gcc montecarlo.c -o pi

no gcc optimization
	AVE:1.501472s

optimization of gcc execution speed-O3
	AVE:1.425400s

python 2.7.18
	AVE:14.216013s

pypy 7.3.1 with GCC 9.3.0
	ave: 0.983056025505

python 3.8.2
	AVE:12.485190s
pypy3 7.3.1 with GCC 9.3.0
	ave: 1.1475282835960388 s

go 1.13.8
	AVE:2.132493

lua 5.3
	AVE:6.715329

rust (cargo 1.44.0)
	cargo run --release
	AVE 0.507196783
	
	cargo run
	AVE 16.449156

Nim 1.0.6
    nim c -r montecarlo.nim
    AVE 7.2063023393
    
    nim c -r -d:release montecarlo.nim
    AVE:0.3236947073200001

Added on June 28, 2020

Added ruby.
 Ruby 2.7.0dev (2019-12-25 master e1e1d92277) [aarch64-linux]
    ave: 7.5495061204099985 s

Fixed python source.
 if  _x ** 2 + _y **2 <= 1:    ①
 ↓
 if  _x * _x  + _y * _y <= 1:  ②

 ※python 3.8.Measurement of 2
① Re-acquisition
    ave: 11.132939925193787 s
 ②
    ave: 9.153075976371765 s
* Measurement of pypy3
① Re-acquisition
    ave: 1.0760090994834899 s
 ②
    ave: 1.3655683732032775 s
that?

Fixed a rust loop.
    for _test_count in 1..test_max {
        for _loot_count in 1..loop_max { ①
     ↓
    for _test_count in 1..=test_max {
        for _loot_count in 1..=loop_max { ②

 cargo run --Release measurement
① Re-acquisition
    AVE:0.5022585
 ②
    AVE:0.51400733

With rust@We measured the Xoroshiro version source provided by tatsuya6502.

 AVE:0.21163946

Conclusion

At this level, ~~ Nim ~~ rust was rated as the fastest. In terms of speed

~~ Nim> rust> pypy> C = java> pypy3> go> lua> ruby> perl> python2 series> python3 series ~~ rust (Xoroshiro version)> Nim> rust> pypy> C = java> = pypy3> go> lua> ruby> perl> python2 series> python3 series

have become. It is properly compiled language> scripting language. Of course, my coding may be bad and "this way it will be faster". However, it is unlikely that there will be a significant change in rankings, even if some ingenuity is used to change the ranking by one step.

Of particular note is python, and by using pypy, I was able to confirm (possible) that it works at a speed close to C language (of course, the source as it is may not work with pypy, and the library becomes pypy. It may not be supported). However, I think it's worth a try if you can speed it up with python's flexible description.

Added on June 28, 2020 ruby is faster than expected and I can't hide my surprise. When I tried it before, it was about c: ruby = 1: 20, so I expected it this time as well, but I was betrayed brilliantly. I ignored @ scivola's suggestions if I didn't. Thank you very much.

With this source fix, python is faster when started with python3, but slower when started with pypy3, which is a mysterious result. (Request for comments from experts) By the way, in the source of this article, I commented to show the difference, but since I replaced it in the measured source, I have not been late in the evaluation of the comment. (I think there was such a thing in the old Basic ...)

The difference between rust and the pseudo-random number generator is apparent. The Xoroshiro version source provided by @ tatsuya6502 revealed that the pseudo-random number generator is also in the right place. Thank you for pointing out.

This article will be updated occasionally. If there is a language that looks good, I will add it.

Recommended Posts

Speed comparison of each language by Monte Carlo method
Stack processing speed comparison by language
Estimating π by Monte Carlo method
Monte Carlo method
Find the ratio of the area of Lake Biwa by the Monte Carlo method
Dominion compression play analyzed by Monte Carlo method
The first Markov chain Monte Carlo method by PyStan
Calculation of the shortest path using the Monte Carlo method
Introduction to Monte Carlo Method
[Scientific / technical calculation by Python] Monte Carlo simulation of thermodynamics of 2D Ising spin system by Metropolis method
Simulate Monte Carlo method in Python
Comparison of calculation speed by implementation of python mpmath (arbitrary precision calculation) (Note)
Speed comparison of Python XML parsing
Speed improvement by self-implementation of numpy.random.multivariate_normal
I couldn't install pypy3.6-7.3.1 on macOS + pyenv, but I could install pypy3.6-7.3.0. I felt the wind of pypy by the Monte Carlo method.
Speed comparison when shifting by group by pandas
Test of uniqueness in paired comparison method
Summary of library hosting pages by language
Try speed comparison of BigQuery Storage API
Summary of SQLAlchemy connection method by DB
Gomoku AI by Monte Carlo tree search
Speed comparison of murmurhash3, md5 and sha1
[Statistics] Let's explain sampling by Markov chain Monte Carlo method (MCMC) with animation.
#Monte Carlo method to find pi using Python
[Python3] Coarse graining of numpy.ndarray Speed comparison etc.
A verification of AWS SDK performance by language
Try implementing the Monte Carlo method in Python