モンテカルロ法をＭＰＩ並列化させた時に乱数がの一様性を保証するには？
私はモンテカルロ法によるπの計算をＭＰＩにより並列化させるプログラムを作り、速度向上を図りました。各プロセッサの乱数発生の初期値seedはそれぞれ変えているので全く同じ乱数を発生させることはありません。ですが何回かは同じ点を取ってしまう可能性があります。私は乱数の検定として乱数に偏りがないかを確認するために度数分布を作ることで偏りがないことは保証されたのですがその問題の確認方法と解決方法が分かりません。ご教授宜しくお願い致します。
この辺りかな？。
>http://madia.world.coocan.jp/cgi-bin/DelphiBBS/wwwlng.cgi?search
こっちだった。
https://www.petitmonte.com/bbs/answers?question_id=141
https://www.petitmonte.com/bbs/answers?question_id=1828
https://www.petitmonte.com/bbs/answers?question_id=5319
https://www.petitmonte.com/bbs/answers?question_id=5207
まずは投稿する前に文章を読み返そうよ！
「偏りがないことは保証された」って自分で書いているジャン。
何が問題なのさ
第一 高度な数学問題に rand みたいな線形合同法を用いた乱数はとか使うべきじゃないって常識でしょ。
それは理解しています。
乱数の発生回数を複数台により分割することで並列化させています。プログラムを並列化させた時それぞれのプロセッサに違うseedを与え乱数を発生させているのですがそれぞれのプロセッサが発生させた乱数は時に全く同じ点を指してしまう可能性があると思うのです。それを解消、または検定によって発見する方法が分からなくて。
if文を使ってプロセッサごとに乱数のでる範囲を制限すれば計算時間がかかり並列化の意味がなくなってしまいます。
また配列を使って発生させた乱数全て格納しかぶっていないかを確認する方法も時間がかかりすぎてしまい良い方法がないかと・・・。
未だ質問者の意図が掴めませんが、
乱数発生器A,Bがあって、それぞれ完全に理想的にランダムな数値を発生させることが（仮に）可能な場合、AとBから、たまたま同じ値が出てしまうことはあるわけで、これがイヤといわれても。
次に、A,Bとも同じ線形合同法アルゴリズムを使っている場合、
運が悪いと何百万回か、あるいは何十万回か後に、かつてAが出していた乱数列にBが乗ってしまう、というケース。
これがまずい場合もあるでしょうが、種を適宜変更することで、便宜的に避けることは出来るでしょう。その為の種ですから。
そもそも線形合同法とはそういう危険性があるものですから。違う乱数発生器を使いましょう。
複数の乱数が同じタイミングで同じ値を出す事を回避したいだけではなく
自分もしくは他の乱数が過去に出した数を出す事も回避したいって事かな？
乱数の発生範囲分の配列を用意できれば・・・・時間はかかりませんよ
var
R1 : Integer;
H1 : array[0..99] of Boolean;
begin
R1 := Random(100);
while H1[R1] = True do //既に乱数が発生していれば再度乱数発生
begin
R1 := Random(100);
end;
H1[R1] := True;
end;
乱数の発生結果に意図的に介入するのは、長期的にみて一様性が乱れることに
なりませんか？ 質問者の意図が不明ですが、一様性が証明されているなら
結果を気にしないことです。あと、なるべく長周期の発生機構を採用されては
どうでしょうか。
えーとさんが書かれているように
乱数発生結果に介入するのはナンセンスですよね？
あと、過去に出した数値を禁止するっていうのはどうなの？
取り扱う問題によっては何億〜何十億個ものデータを確保しなければなりませんしそれを検索するのも恐ろしく大変な気がしますが。
それとも検索だけなら高速アルゴリズムが存在するのかしら？
質問の説明が稚拙になってすいません。
KHE00221のおっしゃられている通りのことです。
いま私はKHE00221のように配列を用意して
for(i=0;i<POINT;i++){
x=drand48();
y=drand48();
if(lengthx[x]=NULL && iength[y]=NULL){
lengthx[]=x:
lengthy[]=y;
}
というように生成し格納しようとしました。
しかしこれでは並列化する特にすべての乱数情報を１つに集めなければいけなく通信のコストが膨大にかかってしまう。ような気がします。
また無論さんがおっしゃている問題もでてくると思うんです。
私はdrand48を使って乱数を生成しています。メルセンヌ・ツイスタという非常に周期が長く高速な乱数も調べてはみたのですが使い方が理解できませんでした。
私の場合最高７台のプロセッサでの並列化ですのでrandは少し悪すぎるのでdrand48にしました
またネットで円周率をモンテカルロを求める問題ででこのように配列を使った方法がありました。
1:length=(double*)malloc(sizeof(double)*(POINT+1));
2:
3:cnt=0;
4:for(i=0;i<POINT;i++){
5:p.x=(double)rand()/(double)RAND_MAX;
6:p.y=(double)rand()/(double)RAND_MAX;
7:r=sqrt(p.x*p.x+p.y*p.y);
8:length[i]=r;
9:}
10:for(i=0;i<POINT;i++){
11:r=length[i];
というようなプログラムがありましたが第8行と第11行はそのようなことは意識されていなく、ひとまとめにするのと同じことなのでしょうか？
すみませんが、これはモンテカルロ法をつかって円周率を計算するとして
乱数で座標を発生し、半径を求めてるということは、円の面積を計算して
それから円周率を出す、ってことでいいでしょうか？
それならば、複数の乱数発生器があって、それらの一様性が保障されていれば
その発生結果についてなんら考慮する必要はないのでは？ 同じ座標が出ようが
出まいが、トータルとして円内に収まる確率だけが問題です。そもそもモンテカルロ法の
正当性は、使う乱数が一様であることと、確率的に大数の法則が成り立つことだけ
です。したがって、一様性が保障されている乱数発生の結果に意図的に介入するのは
まったく意味がありませんし、そんな計算結果は信用できません。
違っていれば、申し訳ありません。
使用メモリを減らす為にビットで使用済みを判定するようにしてみました
1024x1024(1M) 確保で 最大数 8388608 の乱数を確保できます
var
H1 : array[0..1023*1023-1] of Byte;
begin
//データ初期化
for I:=0 to High(H1)-1 do
begin
H1[I] := 0;
end;
//乱数発生
for I:=0 to (High(H1) *8) - 1 do
begin
R1 := Random(High(H1)*8);
R2 := R1 div 8;
R3 := 1 shl (R1 mod 8);
while (H1[R2] and R3) = R3 do //既に乱数が発生していれば再度乱数発生
begin
R1 := Random(High(H1)*8);
R2 := R1 div 8;
R3 := 1 shl (R1 mod 8);
end;
H1[R2] := H1[R2] or R3;
case I mod 7 of
0:// 1台目に送信
1:// 2台目に送信
2:// 3台目に送信
3:// 4台目に送信
4:// 5台目に送信
5:// 6台目に送信
6:// 7台目に送信
end;
end;
乱数発生は１台にして乱数の数値を各プロセッサーに送信して計算だけさせるといのはどうですか？
質問者さんは、典型的な「ギャンブラーの誤謬」をおかしています。
過去に同じ値が出ていたら・・・という心配をしてるんですが、逆を言えば同じ値が２度と出ない、あるいは出にくい事が当たり前と思っているんですね。
これは母集団が籤とか麻雀牌ならば当てはまりますが、繰り返し振るサイコロやルーレットでは当てはまりません。
赤が10回連続で出たら次は黒が出やすい・・・ワケではありません。
ルーレットは過去の結果に左右されません。
これを俗に「ギャンブラーの誤謬（ごびゅう）」というのです。
モンテカルロ法でπを求めるとき、使うべき母集団はどっちでしょうか？
引いたら捨てちゃう籤？それとも繰り返し振るサイコロ？
| ツイート |  |