Mersenne Twister のプチ並列化


Mersenne Twister は広く用いられる疑似乱数生成器です。疑似乱数生成器は乱数生成のために状態を持つので,並列化にはあまり向かないのですが, Mersenne Twister の内部状態の大きさを逆手にとって簡易的な並列化が可能です。

通常の 32 bit 版の Mersenne Twister の内部状態は,サイズ 624 の 32 bit 整数配列と,その配列のインデックスです。

type State = { Index : int; Vector : uint32 [] }

この状態から,次のような感じで乱数を返しています。

let bitop y =
   let mutable y = y
   y <- y ^^^ (y >>> 11)
   y <- y ^^^ ((y <<< 7) &&& 0x9D2C5680u)
   y <- y ^^^ ((y <<< 15) &&& 0xEFC60000u)
   y <- y ^^^ (y >>> 18)
   y
let nextRandom state =
   // 内部ベクトルを使い果たしたら状態の更新
   let state = if state.Index >= 624 then refresh state else state
   let index = state.Index
   let random = bitop state.Vector.[index]
   let nextState = { state with Index = index + 1 }
   random, nextState

最後の y が 32 bit 精度の乱数, nextState が次の内部状態です。

注目するのは,乱数に関係するのは内部ベクトルのたった 1 つの要素だけという点です。 Index が 0 である場合,以下の様に 624 個の乱数をまとめて生成することができます。

let random = Array.Parallel.map bitop state.Vector

さらに乱数 1 つから何かを行うシミュレーションであれば,次のようなことも可能です。

let simulationResult = Array.Parallel.map (bitop >> simulate) state.Vector

普通は倍精度浮動小数点数を利用すると思いますが,この方法では単精度になります。その場合は 64 bit 版の Mersenne Twister を利用すれば良いです。ちなみに 64 bit 版の場合は内部ベクトルのサイズは 312 と半減しますが,並列化を考慮するには十分なサイズだと思います。

乱数生成器の細かいチューニングで並列化するよりも,素直な並列化の方が良いとは思いますが,一応メモとして。