疑似乱数生成 メルセンヌ･ツイスター

良くないというのは、この乱数生成アルゴリズムが線形合同法という物を使用しているからです。実際調べたわけでないので、「らしい」としておきます。

基本的にコンピュータにおいて乱数は、乱数ではなく乱数に見せた数列です。(このことから、疑似乱数と呼ばれています。)
この数列のランダム性が、線形合同法だとあまり良くないです。それと、数列の周期が短く、すぐに繰り返しが起きます。
これだと、モンテカルロ法などの乱数を必要とする物に使ってしまうと、正しい結果が出にくくなります。(モンテカルロ法の結果に、正しいはありませんがね)

新たな疑似乱数生成アルゴリズムが必要です。
たくさん種類がありますが、ここで紹介するのは、メルセンヌ･ツイスターです。

メルセンヌ･ツイスターは、周期が非常に長くきれいな乱数を生成してくれます。
私は難しいことがわからないので、詳しくは、メルセンヌ･ツイスターのホームページを見てください。

使い方

事前準備

メルセンヌ･ツイスターのホームページからC言語のソースを入手して、関数をコンパイルし、それを使うだけです。

ソースは、mt19937ar.cを使用します。ここのページから入手してください。
現在このソースは、実行速度が遅いとされていますが、私の持つ環境では、高速化された物の方が動作が遅かったです。

下の方に、mt19937ar.sep.tgzがあり、これはすでにヘッダファイルが含まれているので、これを使うとかなり楽です。mt19937ar.cとmt19937ar.hをコピーします。

関数の使い方

まず、init_genrand関数で種を初期化します。種の初期化方法は、srand関数と同じように行います。

#include <time.h>

init_genrand((unsigned long)time(NULL));

これでも良いですが、1秒以内で実行すると同じ種になります。なので、以下のようにすると、それを防げます。

#include <windows.h>
#pragma comment(lib, "winmm.lib")

init_genrand((unsigned long)timeGetTime());

timeGetTime関数は、Windows APIに含まれています。Windowsが起動してからの時間を、ミリ秒単位で取得できます。このとき、winmm.libが必要となるので、リンカでリンクを追加するか、#pragmaのように記述を追加してください。(リンカに追加した場合、#pragmaの行は不要)

後は、genrand_*関数で乱数を生成して使うだけです。

unsigned long ul;
long l;
double d;
 
// 32ビットの整数乱数 0x0以上 0xffffffff以下
ul = genrand_int32();

// 31ビットの整数乱数 0x0以上 0x7fffffff以下
l = genrand_int31();

// 0以上 1以下の実数乱数
d = genrand_real1();

// 0以上 1未満の実数乱数
d = genrand_real2();

// 0より大きく 1未満の実数乱数
d = genrand_real3();

// 0以上 1未満の実数乱数
d = genrand_res53();

オーバーフローに注意して使用してください。genrand_real2関数とgenrand_res53関数との違いは、実数に変換する時に使う整数乱数の精度の違いです。realが32ビットでres53が53ビットです。
このソースにはバグがあり、genrand_real3関数で、まれに1が出てくることがあるそうです。使用時に注意してください。

私は、これをよく使うので、これらをクラス化して、ライブラリ化をしています。なかなか便利で、良い物です。