1

お世話になります。
現在モンテカルロ法を高速化しようとOpen MPを使って並列化処理を行なっているのですが、Open MPを使ったときと使わなかった時を比べると、実行速度は若干速くなり、逆に結果の精度は悪くなりました。
この原因は何でしょうか?わかる方教えてください。よろしくお願いいたします。

下記プログラムはモンテカルロ法を使って円周率を求めるものです。

【実行結果】
・openMPなしの場合
3.141556216 120948

・openMPありの場合(同じコードで複数回実行)
3.102521972 103864
3.104310832 108247
3.106061276 109388

【実行環境】
g++ (clangではなくGNUの方です)
CPUは2コア4スレッドです。

#include <iostream>
#include <random>
#include <chrono>
#include <iomanip>
#include <omp.h> //openMPを使わない時はコメントアウトします

using namespace std;
const long int N=1'000'000'000;

random_device seed_gen;
std::mt19937 engine(seed_gen());
std::uniform_real_distribution<double> dist(0.0,1.0);


int main(void)
{
    double x,y,val=0.0;
    long int count=0;
    chrono::system_clock::time_point start, end;

    start = chrono::system_clock::now();

    #pragma omp parallel for reduction(+:count) //同様
    for (long int loop=0; loop < N; ++loop)
    {
      x=dist(engine);
      y=dist(engine);
      if(x*x+y*y<=1.) count++;
    }
     end = std::chrono::system_clock::now();
     double elapsed = std::chrono::duration_cast<std::chrono::milliseconds>(end-start).count();
    val=4*double(count)/N;

  cout<<setprecision(10)<<val<<" "<<elapsed<<endl;

  return 0;
}

1 件の回答 1

3

疑似乱数生成(メルセンヌ・ツイスター)の並列化とほぼ同じ問題を抱えています。

OpenMPの並列ループ内で乱数engineを参照しています。マルチスレッド安全ではない気がします。最善でもスレッド間で同じ乱数値を取得することになりその時点で精度は1スレッド分程度に落ちます。つまり8並列であればN/8回の施行しか意味を持たなくなる。最悪な場合、std::mt19937の内部状態がスレッド間で不適切に更新され、乱数として機能しないもしくは乱数の精度が極端に落ちることも考えられます。

ところで、C言語及びC++言語ではBasic風に変数宣言を関数先頭に列挙すべきではありません。特にC++言語では宣言個所でコンストラクターが動作するため不適切です(更にはデストラクターの実行順が不要になった順ではなく、宣言と逆順に固定される点も不適切)。更に今回の場合、変数xyはループ内で並列に更新するわけですから、ループ内で宣言すべきです。

#pragma omp parallel
{
    std::mt19937 engine{ std::random_device{}() };
    std::uniform_real_distribution<double> dist{ 0.0, 1.0 };
#pragma omp for reduction(+:count)
    for (long int loop=0; loop < N; ++loop) {
       auto x=dist(engine);
       auto y=dist(engine);
       if(x*x+y*y<=1.) count++;
    }
}

この質問に回答するには、ログインする必要があります。

求めていた回答ではありませんか? のタグが付いた他の質問を参照する。