モンテカルロ法で円周率を求める

Last update: 2024/12/30

動作確認環境

• OS：Alma Linux9
• コンパイラ：g++ (GCC) 11.5.0

モンテカルロ法で円周率を求めるc++コード

今後並列計算や可視化を勉強するにあたり何かしら題材があったほうがやりやすいので、モンテカルロ法で円周率を求めるC++のプログラムをつくる（モンテカルロ法についての説明は割愛する）。以下がソースコード。

#include <iostream>
#include <random>
#include <cstdlib>

int main(int argc, char *argv[]) {
        int n = std::stoi(argv[1]);
        double square = 0;
        double circle = 0;
        std::random_device seed;
        std::mt19937 engine(seed());
        std::uniform_real_distribution<> dist(0.0, 1.0);

        for(int i = 0; i < n; i++){
                double x = dist(engine);
                double y = dist(engine);
                square += 1;
                if(x*x + y*y <= 1) {
                        circle += 1;
                }
        }
        double pi = (circle / square) * 4;
        std::cout << "pi:" << pi << std::endl;

        return 0;
}

これをpi.cppという名前で保存し、g++でコンパイルする。

$ g++ pi.cpp

サンプル数と結果の比較

上記のコードはプログラムの実行時にサンプル数を引数として与えるようにしており、サンプル数を増やして行ったときに結果が以下となる。サンプル数が増えるほど円周率の誤差が少なくなっていることがわかる。

$ ./a.out 10000
pi:3.168
$ ./a.out 100000
pi:3.1354
$ ./a.out 1000000
pi:3.14255
$ ./a.out 10000000
pi:3.14176

参考文献

一週間でなれる！スパコンプログラマ Day 3 : 自明並列

Copyright (c) 2024 kd