PCGと乱数生成について調べた

はじめに

テストでランダムな入力値を生成するのにどういうアルゴリズムを使うのが良いのかと調べていたのですが、 今回いろいろ知ったのでメモです。Goだとmath/randパッケージを使っておけば良いのですが、C言語だと下記のリンク先からコードをコピペ改変して使うのが良さそうです。

なお、暗号用の乱数はこの記事のスコープ外です。

pcg32が良さそう

今回調べるまではXorshift - Wikipediaが処理が軽くて周期も長くて良いかなと思ってました。

が、PCG, A Family of Better Random Number Generators | PCG, A Better Random Number Generatorの比較表を見ると、PCGのほうが良さそうです(ただ、PCGの作者によるサイトなので他の専門家の意見も聞いてみたいところではあります)。

Download the PCG Library | PCG, A Better Random Number Generatorに開発者によるCとC++の実装があり、最小限のC実装は5行と非常にコンパクトです。

Permuted congruential generator - Wikipediapcg32_fastというのも紹介されていますが、通常はpcg32の通常版を使うほうが良いらしいです。

乱数のテスト用にTestU01というライブラリがある

良い xorshift、悪い xorshiftで生成した乱数をプロットして規則的な模様が出るケースは良くないというのを見てなるほどと思いました。一方で良いほうは感覚的にしか分からないよなと思ったのですが、Xorshift - WikipediaPermuted congruential generator - Wikipediaで紹介されていたTestU01 - Empirical Testing of Random Number Generatorsが良いみたいです。

TestU01 - Wikipediaにも説明がありました。TestU01 - Empirical Testing of Random Number Generatorsgithub.com/umontreal-simul/TestU01-2009へのリンクがありました。論文のリンクもありました(が私は読んでないです)。

Ubuntu 22.04 LTSではtestu01 ソースパッケージがあったのですが、testu01-docパッケージ内にexamplesのソースが入っているけどgithub.com/umontreal-simul/TestU01-2009のexamples/にあるtestxoshiro128plusplus.cなど一部のファイルは含まれていませんでした。

というわけでdebパッケージは使わずにgithub.com/umontreal-simul/TestU01-2009のほうを取ってきて試してみました。

なぜかconfigureに実行パーミッションが付いていないので

sh configure

で実行して

make -j
sudo make install

でビルド・インストールします。

あとはcd examplesしてtestxoshiro128plusplus.ctestpcg32.cの先頭のコメントに書いてあるコマンドでコンパイルして実行します。

bbattery_SmallCrushだとすぐ終わるのですが、ソース内のコメントにbbattery_BigCrushに変えてもよいとあり、BigCrushXorshift - Wikipediaでも言及されていたなということで試してみました。どちらも2時間ぐらいかかってテストはパスしていました。実行中htopで見てみると、CPUの1つのコアが100%近く使われて、行程が変わると別の1つのコアが100%近く使われるという挙動になっていました。複数コアを並列では使ってくれないようです。

PCGのblog記事一覧でもTestU01が何度か出てきてました(が、今回は読んでないです)。

シード値について

横道: Xorshift用に良いシード値を生成する方法 SplitMix64

Xorshift - WikipediaのInitializationの項でSplitMix64 generatorというのを知りました。 良い xorshift、悪い xorshiftでも「レジスタ内に0のbitが固まって多く存在すると、しばらくの間0が多く含まれる値が出力されます」と説明されていました。 ということでSplitMix64でシード値を加工して使うのがお勧めらしいです。

pcg32のシード値について

https://www.pcg-random.org/using-pcg-c-basic.html#pcg32-srandom-r-rngptr-initstate-initseq によると/dev/randomが使えればそれを使うか、quick and dirtyにしたい場合は現在時刻とRNGの状態変数のアドレスを使うと良いらしいです。

指定した範囲の整数の乱数を効率よく生成するアルゴリズム

https://github.com/golang/go/blob/go1.20/src/math/rand/rand.go#L151-L173 のコードを見て https://lemire.me/blog/2016/06/27/a-fast-alternative-to-the-modulo-reduction/https://lemire.me/blog/2016/06/30/fast-random-shuffling/ の記事を読みました。 さらにEfficiently Generating a Number in a Range | PCG, A Better Random Number Generatorでも同じ題材について詳しく書かれていました。

[0, n)の整数の乱数を生成する際、バイアスありでも良い場合はpcg32[0, 2^32)の整数xを生成してx % nが一番シンプルです。が、剰余は遅いので、代わりに乗算とシフト演算で算出する方法が紹介されています。さらにRejection sampling - Wikipediaという手法でバイアス無しにする方法も紹介されていました。

記事からリンクされていた https://github.com/lemire/Code-used-on-Daniel-Lemire-s-blog/tree/master/2016/06/25 のコードを試してみた結果は以下の通りです。

~/ghq/github.com/lemire/Code-used-on-Daniel-Lemire-s-blog/2016/06/25$ ./fastrange
N = 31
modsum(z,N,accesses,nmbr):  30.01 cycles per operation
fastsum(z,N,accesses,nmbr):  2.36 cycles per operation
N = 1500
modsum(z,N,accesses,nmbr):  26.60 cycles per operation
fastsum(z,N,accesses,nmbr):  2.67 cycles per operation
N = 15000
modsum(z,N,accesses,nmbr):  21.27 cycles per operation
fastsum(z,N,accesses,nmbr):  2.05 cycles per operation
N = 32
modsum(z,N,accesses,nmbr):  23.44 cycles per operation
fastsum(z,N,accesses,nmbr):  1.86 cycles per operation
maskedsum(z,N,accesses,nmbr):  1.49 cycles per operation
N = 4096
modsum(z,N,accesses,nmbr):  19.65 cycles per operation
fastsum(z,N,accesses,nmbr):  1.74 cycles per operation
maskedsum(z,N,accesses,nmbr):  1.36 cycles per operation
N = 65536
modsum(z,N,accesses,nmbr):  17.98 cycles per operation
fastsum(z,N,accesses,nmbr):  1.98 cycles per operation
maskedsum(z,N,accesses,nmbr):  1.61 cycles per operation

あと、Efficiently Generating a Number in a Range | PCG, A Better Random Number Generatorの最後に書かれていたのですが、nが2^32未満の場合、uint64_tの乱数を生成するよりuint32_tの乱数を生成してそれを使うほうが速いそうです。Goのmath/randのIntnもそういう実装になっていました。 https://github.com/golang/go/blob/go1.20/src/math/rand/rand.go#L175-L185

自分でテストで使うケースを考えると、ほとんどの場合はnは2^32未満で使うと思うので、pcg32だけ使っておけば良いかもと思いました。

効率良くシャッフルするアルゴリズム Fisher–Yates shuffle

math/rand.Rand.ShuffleのコメントでFisher–Yates shuffle - Wikipediaを知りました。Fast random shuffling – Daniel Lemire’s blogの冒頭でも紹介されていました。

記事からリンクされていた https://github.com/lemire/Code-used-on-Daniel-Lemire-s-blog/tree/master/2016/06/29 のコードを試してみた結果は以下の通りです。

~/ghq/github.com/lemire/Code-used-on-Daniel-Lemire-s-blog/2016/06/29$ ./shuffle
Shuffling arrays of size 10000
Time reported in number of cycles per array element.
Tests assume that array is in cache as much as possible.
shuffle_pcg(testvalues,size)                                :  33.20 cycles per input key
shuffle_pcg_go(testvalues,size)                             :  33.15 cycles per input key
shuffle_pcg_java(testvalues,size)                           :  16.42 cycles per input key
shuffle_pcg_divisionless(testvalues,size)                   :  3.75 cycles per input key
shuffle_pcg_divisionless_with_slight_bias(testvalues,size)  :  3.43 cycles per input key

横道: Daniel Lemireさんはsimdjsonの作者でもある

https://lemire.me/blog/ から https://github.com/lemire を見て気づいたのですが、Daniel Lemireさんはsimdjsonの作者でもあったんですね。それと、他のブログ記事もいくつかチラ見したのですが、いろいろ面白そうな記事があったので、いつか読んで試してみたいところです。