1

C++/c++11の環境で疑似乱数生成のメルセンヌ・ツイスターを並列化して動作させる方法をご教示下さい。
boost/random.hppとopenmpの組合せで実装出来るのかとも考えいろいろ試したのですが、
簡単には並列化できないようです。

#include <boost/random.hpp>
#include <omp.h>

using namespace std;
using namespace boost;
mt19937 gen;
uniform_real<> uniform_dst(0,1);//min,max
variate_generator< mt19937,uniform_real<> > uniform(gen,uniform_dst);

void test1(){
  uniform.engine().seed(123);
  int NN=1600000;
  int i;
  double x;
#pragma omp parallel for private(i,x,gen)
  for(i=0;i<NN;++i){
    x=uniform();
    printf("%f in thread:%d\n",x,omp_get_thread_num());
  }
}
int main(){
    test1();
    return 0;
}

乱数生成器のmt19937をベースにしているのですが、openmpのループの中からこの乱数生成器を呼び出しても、乱数生成器自体は並列化されているわけではないようです。uniform()は個別のスレッドから呼び出されているようですが、乱数生成器の部分は共通して動作しており、並列化はされていないようです。この乱数生成器自体を並列化したいと考えています。

よろしくお願い致します。

4 件の回答 4

2

本質的でないuniform_realvariate_generatorを取り除いて説明します。

void test1() {
    printf("%s\n", __func__);
    std::mt19937 gen{ 123 };
    #pragma omp parallel for private(gen)
    for (int i = 0; i < 16; ++i) {
        unsigned int x = gen();
        printf("%u in thread:%d\n", x, omp_get_thread_num());
    }
}

メルセンヌ・ツイスターに限らず全ての疑似乱数生成器は

  • ある状態で得られる乱数は決まっている
  • ある状態から遷移する次の状態は決まっている
  • 最初の状態を決定する値をseedと呼ぶ

という特徴があります。
ですのでprivate(gen)を指定したとしても、得られる乱数はスレッド間で同じになりますし、遷移する次の状態も同じになります。結果的にスレッドごとに得られる乱数列は同一のものとなってしまいます。これを解消するにはまずはスレッド間排他を行うことです。

void test2() {
    printf("%s\n", __func__);
    std::mt19937 gen{ 123 };
#pragma omp parallel for
    for (int i = 0; i < 16; ++i) {
        unsigned int x;
#pragma omp critical
        {
            x = gen();
        }
        printf("%u in thread:%d\n", x, omp_get_thread_num());
    }
}

ただし、排他のコストが高すぎるため並列化している意義が薄くなります。そこで排他するよりもスレッド毎に疑似乱数生成器を用意することです。その際、個々の疑似乱数生成器のseedは異なる値に設定する必要があります。

void test3() {
    printf("%s\n", __func__);
#pragma omp parallel
    {
        std::mt19937 gen{ std::random_device{}() };
#pragma omp for
        for (int i = 0; i < 16; ++i) {
            unsigned int x = gen();
            printf("%u in thread:%d\n", x, omp_get_thread_num());
        }
    }
}

以上を踏まえて質問のコードは次のように書けます。

#include <random>
#include <omp.h>

void test1() {
    int const NN = 1600000;
#pragma omp parallel
    {
        std::mt19937 gen{ std::random_device{}() };
        std::uniform_real_distribution<> dst{ 0, 1 };
        auto uniform = [&] { return dst(gen); };
#pragma omp for
        for (int i = 0; i < NN; ++i) {
            auto x = uniform();
            printf("%f in thread:%d\n", x, omp_get_thread_num());
        }
    }
}
int main() {
    test1();
    return 0;
}
6
  • 基本的にはこの方針で良いと思うのですが、たとえばモンテカルロ法に使いたいなどの意図がある場合、MTの良い性質を保ったまま複数系列に分ける必要があるわけで、となると Fumu 7 さんの回答で触れられているページで紹介されている並列対応版MTや、この投稿で触れられている方法を使う必要があるのかもしれません。残念ながら詳しくないのできちんとは書けませんが……。
    – nekketsuuu
    Commented 2017年10月6日 1:35
  • 動作確認しました。これで乱数生成が並列化されています。
    – Abocado
    Commented 2017年10月7日 22:18
  • しかし@nekketsuuu さんが指摘しているようにモンテカルロ法での使用を考えると、さらに改善の余地があります。上のtest1()のような、ある程度の長さの乱数列の生成を、多数回呼び出すような場合があります。そのときに、このtest1()のように呼び出されるたびにseedをrandom_deviveで与え直すと、MTの良い性質が生かせないようになります。例えばtest1()の中のNNが小さな値で、test1()の呼び出し回数が多くなれば、得られる乱数列はrandom_deviceによるseedの生成に強く依存していることになります。(悪い影響が出るとしても、それを見極めるのは難しい状況になると思いますが)また乱数列の再現性を確保するために、seedを自分で設定したい場合もあります。
    – Abocado
    Commented 2017年10月7日 23:04
  • @Abocado さん、質問者さんに応えたのではなく、質問文に答えている点をご理解ください。つまりNNがそこそこ大きくtest1()が1回しか呼ばれない質問に答えました。後だしでコメントされるのではなく、質問文に含めていただけますでしょうか。とりあえずC++11のthread_localを使う方法もあるようですが、環境によってはOpenMPと共存できるかどうか分かれるようです。
    – sayuri
    Commented 2017年10月8日 0:53
  • @sayuriさん、失礼しました。 このtest1()の呼び出しを多数繰り返すという想定で、議論を続けたいと思います。質問文に追加しようと思います。あるいは新しいスレッドを立てたほうが良いですかね?
    – Abocado
    Commented 2017年10月8日 2:15
0

OpenMP並列処理でワーカスレッド毎にオブジェクトを複製したい場合、firstprivate節にて実現できます。このケースでは、それぞれのスレッドにて同一擬似乱数系列が得られることに注意ください。

#pragma omp parallel for private(i,x) firstprivate(uniform)
for(i=0;i<NN;++i){
  ...
}
-1

乱数呼び出し部分を排他制御するか、
あるいはスレッドセーフなメルセンヌ・ツイスターの実装を利用してはどうでしょうか。

例えば、mt-liteはスレッドセーフでありマルチスレッドで使えるとあります。
http://mt-lite.sourceforge.net/index.html.ja

2
  • 1
    mt-liteが謳うスレッドセーフ = グローバル変数や静的ローカル変数に状態を保持しない という意味ですから、Boost.Randomと同条件ですね。
    – yohjp
    Commented 2017年10月5日 5:47
  • マルチスレッドで使いたいのではなくて、スレッド毎に独立した乱数生成器を持ちたいという質問だったんですね。
    – shimitei
    Commented 2017年10月5日 8:56
-3

メルセンヌ・ツイスターについては、疑似乱数発生法に関する業績で"日本IBM科学賞"を受賞された松本 眞氏のMersenne Twister Home Pageが詳しいです。参考になる情報を得られると思います。

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

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