とあるきっかけで自分が使っている乱数ライブラリを一部更新することにしました。

Xorshiftが高速で周期もある程度あるのでゲームなどで使いやすいらしい

と言うことをいろいろな所から聞いている内にやっている作業に絡む問題が見つかったので実装しようと思い実装してみました。

名前通りXorとShift演算だけを使った乱数の生成法で、ほぼどのCPUでも高速に計算できるのが売りでしょうか。

代表的な実装はXorshiftの128bitバージョンですが、それ以外のビット数でも乱数列として使用できるような物があるので気になる人はさらに調べてみるといいと思います。

いままで使っていた乱数ライブラリはさすがに10年くらい前に確認した組んだMersenneTwister法だけだったので拡張するためにいろいろと細工をするのにとまどりました。

Xorshiftでの乱数は初期状態に注意

どうせいろいろな場所に実装法やら(Wikipediaにコード本体が記述されているという珍しい状態)乱数の性質に関する考察やらがありますのでそちらに譲りましが、初期状態には注意です。

特に実装をそのまま使って乱数初期化を考えたときに何も考えないと初期化後しばらくの間は乱数としては使い物にならなくなります。

たとえば、Xorshift(128bit)の乱数初期化を以下のように組んだとします。

typedef unsigned int uint32_t;
static uint32_t x,y,z,w;
void srandxorshift(uint32_t seed)
{
	x = 123456789;
	y = 362436069;
	z = 521288629;
	w = seed;
}
uint32_t randxorshift(void)
{
	uint32_t t = x ^ (x << 11);
	x = y; y = z; z = w;
	return w = (w ^ (w >> 19)) ^ (t ^ (t >> 8));
}

このとき、乱数を初期化した直後だとseedの変更分はx,y,zの各値に影響がないため、影響が波及するまではどのようなseedで初期化されようが

seed値が似通っていれば乱数値も似通った値がしばらく計算されれしまうのが分かると思います。

そのため、Xorshift法でseedを与えて乱数列を初期化させるときはx,y,z,wのすべてに対してseedの変更を利かせるなり、先頭の乱数列をしばらく強制的に破棄してしまうのが良いと思われます。

先頭を破棄する方法の場合は「何回破棄すればよいか？」と言う問題が出てきますが・・・。定量的に説明しているページに譲ります。

たとえば、今回のXorshift法のルーチンの場合、実際に乱数をとってみると

というように、乱数の初期値による差が8回目までほぼ現れていないことが分かります。1回目なんかはこの表の値だと乱数初期値の差がそのまま出てきていますね。

もちろん、見て分かるとおり8回程度ではぜんぜん足りていませんね。あとは計算回数によるビットの変更状態(エントロピーやらハミング距離やら)を考えたり、実装者の気分次第と言ったところでしょうか。

SFMT法もついでに実装

MersenneTwister法の改良系として「SIMD-oriented Fast Mersenne Twister」という方法もあると言うことが分かりました。

MersenneTwister法で乱数の実装としてあまり良くない要素だった物を修正し、乱数列の周期を変更しやすくしたり、SIMD(この場合はSSE2)で高速に計算できるようにした物らしいです。

もちろん、逆に乱数列の周期を短くしたパターンもあるので、計算速度およびメモリ空間の関係からあまりにも長すぎる周期の乱数列を使いたくない、と言うときには下げるのも手だったりします。

演算は基本的に128bit単位(というよりは32bitx4単位がほとんど)で行われます。そのため、SSE2ととても親和性がよい計算方法です。

元のコードはエンディアンの関係で微妙に見づらくなっていますので、リトルエンディアンでしか使用しない、と決めうちすれば見やすくできます。

なお、元のソースコードのライセンスは修正BSDなので使いやすいと思います。

アルゴリズムも新しい物がいろいろと出てくる

特にゲームで使用するアルゴリズムは年々新しい物が出てきますので使える物はちゃんと取り込んだ方がいいです。

バグがあるライブラリを使い続けるというのもまた問題が出てきますし・・・。