piruty.hatenablog.com

今を一所懸命生きる。

【ruby】線形合同法を用いた乱数生成

アルゴリズム事典よりC言語で書かれていた線形合同法を用いた乱数生成をrubyで書き直したのでまとめ。

線形合同法についてはwikipediaが詳しい。

線形合同法 - Wikipedia

これは、以下のような漸化式を解くことで乱数を生成するアルゴリズムである。

X_{n+1} = (A * X_n + C) mod M

ここで、A, C, Mは定数であり、0 < A < M かつ 0 ≤ C < Mである。
また、初期値X0はSeed値として与えられる。

ただし、アルゴリズム事典では更に式変形したものが紹介されており、それをC言語で実装している。

一般に0以上M未満の整数の一様乱数xから0以上L未満(L≤M)の整数の一様乱数yを作るには、 1≤D≤floor(M/L)の範囲のなるべく大きい整数Dを選んでおき、 y <- x/Dとし、y≥Lになったらやり直す。

従って、今回はこれを元に実装されたものをrubyで書き直してみた。

$seed = 1
MAX_INT = 680564733841876926926749214863536422911

def init_rnd(x)
  $seed = x
end

def irnd
  $seed = $seed * 1566083941 + 1
end

def rnd
  (1.0 / (MAX_INT + 1.0)) * irnd
end

デフォルトの初期値は1だが、init_rnd関数に値を渡すことで初期値を変更することができる。
元の実装ではULONG_MAXを利用していたが、rubyでは整数の最大値を取得するメソッドや定数が無いようなので、別途定数としてMAX_INTを定義している。
また、irnd関数内で使用されている定数については、「線形合同法の多次元分布の悪さ」を目立たないようにするために選ばれた値で、これは元の漸化式のAに相当する。

本の中で説明されているように、線形合同法を用いると値の組み合わせによっては生成される乱数が一定の周期をもってしまう。
それを改善するための文献も紹介されているので、そちらも読んでみたり、また関連した項目も確認して理解を深めたほうが良さそう。