メルセンヌ・ツイスタ乱数発生プログラムをJavaに移植してみる

会社の飲み会で、後輩と乱数の話になりました。「java.util.Randomはあまり良くない。もっと良いメルセンヌ・ツイスタ乱数という乱数があって、、、」。少し調べてみるとこれの事みたいです->Mersenne Twister: A random number generator (since 1997/10)。

Javaで実装したものは、ここにありました->yosaku: Java Mersenne Twister:。でも、少し動かしてみると、原典で紹介されているプログラムと3番目から結果が違う、、、(以下は、13579で初期化し、32bit整数の議事乱数を発生させてみた結果です)。

No	オリジナル(mt19937ar.c)	yosakuさん版(MersenneTwister.java)
1	720353472	720353472
2	341455330	341455330
3	1294704524	2969203518
4	1431604248	2827313322
...	...	...

違ってても良い気もしますが、、、。Mersenne Twister: A random number generator (since 1997/10)にある資料をみても、良し悪しが分かりません。というわけで、原典にあるmt19937ar.cと結果が一致するJava実装を作ってみました。

目次

• テストデータ作成
• テストコード作成
• プログラム本体を作る
• java.util.Randomに組み込む

テストデータ作成

まず、テストデータを作ります。Mersenne Twister with improved initializationにある、プロトタイプ宣言の入った版 mt19937ar.sep.tgz をダウンロードして、テストデータを作るプログラムを作ります。テストケースは7ケース。以下で初期化して、乱数を発生させています。

• 13579
• 0xffffffff
• 0x00000000
• 0xffff0000
• {0x123, 0x234, 0x345, 0x456}
• {0xffffffff, 0x00000000, 0xffffffff, 0x00000000}
• {0xffffffff, 0xffffffff, 0xffffffff, 0xffffffff}

#include 
#include "mt19937ar.h"
void long_random(void);
void array_random(void);
void long2_random(void);
void long3_random(void);
void long4_random(void);
void array2_random(void);
void array3_random(void);

int main(void)
{
    long_random();
    long2_random();
    long3_random();
    long4_random();
    
    array_random();
    array2_random();
    array3_random();
}

void long_random(void)
{
    FILE *fp;
    
    if ((fp = fopen("D:\\long.txt", "w")) == NULL) {
        return;
    }
    
    init_genrand(13579UL);
    for (int i = 0; i < 1000; i++) {
        fprintf(fp, "%d\n", genrand_int32());
    }
    
    fclose(fp);
}

void long2_random(void)
{
    FILE *fp;
    
    if ((fp = fopen("D:\\long2.txt", "w")) == NULL) {
        return;
    }
    
    init_genrand(0xffffffffUL);
    for (int i = 0; i < 1000; i++) {
        fprintf(fp, "%d\n", genrand_int32());
    }
    
    fclose(fp);
}

void long3_random(void)
{
    FILE *fp;
    
    if ((fp = fopen("D:\\long3.txt", "w")) == NULL) {
        return;
    }
    
    init_genrand(0x00000000UL);
    for (int i = 0; i < 1000; i++) {
        fprintf(fp, "%d\n", genrand_int32());
    }
    
    fclose(fp);
}

void long4_random(void)
{
    FILE *fp;
    
    if ((fp = fopen("D:\\long4.txt", "w")) == NULL) {
        return;
    }
    
    init_genrand(0xffff0000UL);
    for (int i = 0; i < 2000; i++) {
        fprintf(fp, "%d\n", genrand_int32());
    }
    
    fclose(fp);
}

void array_random(void)
{
    FILE *fp;
    
    if ((fp = fopen("D:\\array.txt", "w")) == NULL) {
        return;
    }
    
    unsigned long init[4] = {
          0x123
        , 0x234
        , 0x345
        , 0x456
    };
    int length = 4;
    init_by_array(init, length);
    for (int i = 0; i < 1000; i++) {
        fprintf(fp, "%d\n", genrand_int32());
    }
    
    fclose(fp);
}

void array2_random(void)
{
    FILE *fp;
    
    if ((fp = fopen("D:\\array2.txt", "w")) == NULL) {
        return;
    }
    
    unsigned long init[4] = {
          0xffffffffUL
        , 0x00000000UL
        , 0xffffffffUL
        , 0x00000000UL
    };
    int length = 4;
    init_by_array(init, length);
    for (int i = 0; i < 1000; i++) {
        fprintf(fp, "%d\n", genrand_int32());
    }
    
    fclose(fp);
}

void array3_random(void)
{
    FILE *fp;
    
    if ((fp = fopen("D:\\array3.txt", "w")) == NULL) {
        return;
    }
    
    unsigned long init[4] = {
          0xffffffffUL
        , 0xffffffffUL
        , 0xffffffffUL
        , 0xffffffffUL
    };
    int length = 4;
    init_by_array(init, length);
    for (int i = 0; i < 1000; i++) {
        fprintf(fp, "%d\n", genrand_int32());
    }
    
    fclose(fp);
}

テストコード作成

JUnit4で、テストコードを作ります。どうしてもApache Commonsを使いたかったので、Integer.parseInt(int)を使うところをNumberUtils.toInt(String)と書いています。

package mathpara.random.mt;

import static org.junit.Assert.*;

import java.io.BufferedReader;
import java.io.FileReader;
import java.io.IOException;

import org.apache.commons.lang.math.NumberUtils;
import org.junit.Test;

public class RandomIntGeneratorTest {
    @Test
    public void testNextInt_InitByInt() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(13579);
        BufferedReader reader = new BufferedReader(new FileReader("./testdata/long.txt"));

        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("Int型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();
    }

    @Test
    public void testNextInt_InitByInt2() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(0xffffffff);
        BufferedReader reader = new BufferedReader(new FileReader("./testdata/long2.txt"));

        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("Int型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();
    }

    @Test
    public void testNextInt_InitByInt3() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(0x00000000);
        BufferedReader reader = new BufferedReader(new FileReader("./testdata/long3.txt"));

        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("Int型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();
    }
    
    @Test
    public void testNextInt_IntByInt4() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(0xffff0000);
        BufferedReader reader = new BufferedReader(new FileReader("./testdata/long4.txt"));

        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("Int型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();        
    }

    @Test
    public void testNextInt_InitByArray() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(
            new int[] {
                  0x123
                , 0x234
                , 0x345
                , 0x456
            }
        );
        BufferedReader reader = new BufferedReader(new FileReader("./testdata/array.txt"));
        
        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("array型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();
    }
    
    @Test
    public void testNextInt_InitByArray2() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(
            new int[] {
                  0xffffffff
                , 0x00000000
                , 0xffffffff
                , 0x00000000
            }
        );
        BufferedReader reader = new BufferedReader(new FileReader("./testdata/array2.txt"));
        
        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("array型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();
    }
    
    @Test
    public void testNextInt_InitByArray3() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(
            new int[] {
                  0xffffffff
                , 0xffffffff
                , 0xffffffff
                , 0xffffffff
            }
        );
        BufferedReader reader = new BufferedReader(new FileReader("./testdata/array3.txt"));
        
        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("array型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();
    }

    @Test
    public void testInit_Int() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(Integer.MAX_VALUE);
        // aging
        for (int i = 0; i < 5000; i++) {
            generator.nextInt();
        }
        generator.init(13579);

        BufferedReader reader = new BufferedReader(new FileReader("./testdata/long.txt"));

        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("Int型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();
    }

    @Test
    public void testInit_Array() throws IOException {
        RandomIntGenerator generator = new RandomIntGenerator(Integer.MIN_VALUE);
        // aging
        for (int i = 0; i < 5000; i++) {
            generator.nextInt();
        }
        generator.init(
            new int[] {
                    0x123
                  , 0x234
                  , 0x345
                  , 0x456
              }
          );
        BufferedReader reader = new BufferedReader(new FileReader("./testdata/array.txt"));
        
        String line = null;
        while ((line = reader.readLine()) != null) {
            assertEquals("array型の初期化:乱数不一致", NumberUtils.toInt(line), generator.nextInt());
        }
        
        reader.close();
    }
}

プログラム本体を作る

プログラム本体を作ります。久しぶりにビット演算を使ったので、大変でした。synchronizedとか使って、マルチスレッド対策っぽい事をしてますが、自身がありません。ので、パクる人は注意して下さい。

package mathpara.random.mt;

public class RandomIntGenerator {
    private static final int N = 624;
    private static final int M = 397;
    
    private final int mt[] = new int[N];
    
    private int mtIndex = N + 1;
    
    /**
     * seedでRandomIntGeneratorのインスタンスを作る
     * 
     * @param seed シード
     */
    public RandomIntGenerator(int seed) {
        init(seed);
    }
    
    /**
     * 配列seedArrayでRandomIntGeneratorのインスタンスを作る
     * 
     * @param seed シード配列
     */
    public RandomIntGenerator(int[] seedArray) {
        init(seedArray);
    }
    
    /**
     * 一様分布の int 型の擬似乱数を返す。
     * 
     * @return 一様分布の int 型の次の擬似乱数値
     */
    public int nextInt() {
        int nextInt = generateInt();
        
        // 精度を上げる為の調律
        nextInt = temper(nextInt);
        
        return nextInt;
    }
    
    /**
     * seedでMersenne Twister配列を初期化する
     * 
     * @param seed シード
     */
    public void init(int seed) {
        
        long[] longMt = createLongMt(seed);
        
        setMt(longMt);
    }
    
    /**
     * 配列seedArrayでMersenne Twister配列を初期化する
     * 
     * @param seedArray シード配列
     */
    public void init(int[] seedArray) {
        long[] longMt = createLongMt(19650218);
        
        final int max = (N > seedArray.length ? N : seedArray.length);
        for (int i = 1, j = 0, counter = 0; counter < max; i++, j++, counter++) {
            if (N <= i) {
                longMt[0] = longMt[N - 1];
                i = 1;
            }
            if (seedArray.length <= j) {
                j = 0;
            }
            longMt[i] ^= ((longMt[i-1] ^ (longMt[i-1] >>> 30)) * 1664525);
            longMt[i] += seedArray[j] + j;
            longMt[i] &= 0xffffffffL;
        }
        
        final int initial = max % (N -1) + 1;
        for (int i = initial,counter = 0; counter < N - 1; i++, counter++) {
            if (N <= i) {
                longMt[0] = longMt[N - 1];
                i = 1;
            }
            longMt[i] ^= ((longMt[i-1] ^ (longMt[i-1] >>> 30)) * 1566083941);
            longMt[i] -= i;
            longMt[i] &= 0xffffffffL;
        }
        
        longMt[0] = 0x80000000L;
        
        setMt(longMt);
    }
    
    /**
     * Mersenne Twisterアルゴリズムで、一様分布の int 型の擬似乱数を返す。
     * 
     * @return 
     */
    private synchronized int generateInt() {
        twist();
        
        int ret = mt[mtIndex];
        mtIndex++;
        
        return ret;
    }
    
    /**
     * 配列をtwistする
     *
     */
    private void twist() {
        final int[] BIT_MATRIX = new int[] {0x0, 0x9908b0df};
        final int UPPER_MASK = 0x80000000;
        final int LOWER_MASK = 0x7fffffff;
        
        if (mtIndex < N) {
            return;
        }
        
        if(mtIndex > N) {
            init(5489);
        }
        
        for (int i = 0; i < N; i++) {
            int x = (mt[i] & UPPER_MASK) | (mt[(i + 1) % N] & LOWER_MASK);
            mt[i] = mt[(i + M) % N] ^ (x >>> 1) ^ BIT_MATRIX[x & 0x1];
        }
        
        mtIndex = 0;
    }
    
    /**
     * 調律する
     * 
     * @param num 調律するint
     * @return 調律されたint
     */
    private static int temper(int num) {
        
        num ^= (num >>> 11);
        num ^= (num << 7) & 0x9d2c5680;
        num ^= (num << 15) & 0xefc60000;
        num ^= (num >>> 18);
        
        return num;
    }
    
    /**
     * 2^32 - 1以下のlongを、右32ビット列だけintに変換する
     * 
     * @param num 2^32 - 1以下のlong
     * @return 変換したint
     */
    private static int toInt(long num) {
        return (int) (num > Integer.MAX_VALUE ? num - 0x100000000L : num);
    }
    
    /**
     * seedから、初期化したmtをlongで作る
     * 
     * @param seed シード
     * @return 初期化したmt
     */
    private static long[] createLongMt(int seed) {
        long[] longMt = new long[N];
        
        longMt[0] = seed & 0xffffffffL;
        
        for (int i = 1; i < N; i++) {
            longMt[i] = longMt[i - 1];
            
            longMt[i] >>>= 30;
            longMt[i] ^= longMt[i - 1];
            longMt[i] *= 0x6C078965L;
            longMt[i] += i;
            longMt[i] &= 0xffffffffL;
        }
        
        return longMt;
    }
    
    /**
     * long型の配列を、ビット配列としてintに変換して、インスタンス変数mtにセットする
     * 
     * @param longMt 配列
     */
    synchronized private void setMt(long[] longMt) {
        for (int i = 0; i < N; i++) {
            mt[i] = toInt(longMt[i]);
        }
     
        mtIndex = N;
    }
}

java.util.Randomに組み込む

　java.util.Randomを継承して、nextDouble()、nextFloat()とか利用できるようにしました。

package mathpara.random.mt;

public class Random extends java.util.Random {
    
    private static final long serialVersionUID = 2479573230251354962L;
    
    private RandomIntGenerator generator = null;
    
    public Random() {
        this(System.currentTimeMillis());
    }
    
    public Random(long seed) {
        super(seed);
        setSeed(seed);
    }
    
    public Random(int[] seedArray) {
        init(seedArray);
    }
    
    @Override
    synchronized public void setSeed(long seed) {
        int[] seedArray = new int[2];
        seedArray[0] = (int) (seed & 0xffffffff);
        seedArray[1] = (int) (seed >>> 32);
        
        if (seedArray[1] == 0) {
            init(seedArray[0]);
        } else {
            init(seedArray);
        }
    }
    
    @Override
    protected int next(int bits) {
        return generator.nextInt() >>> (32 - bits);
    }
    
    synchronized public void setSeed(int[] seedArray) {
        init(seedArray);
    }
    
    private void init(int seed) {
        if (generator == null) {
            generator = new RandomIntGenerator(seed);
        } else {
            generator.init(seed);
        }
    }
    
    private void init(int[] seedArray) {
        if (generator == null) {
            generator = new RandomIntGenerator(seedArray);
        } else {
            generator.init(seedArray);
        }
    }
}

ちょっと気にならないのが、2点。

　その１。コンストラクタRandom(long)で、setSeed(long)が2重に呼び出されてしまう。継承している以上、javaの仕組みから

1. スーパークラスがsetSeed(long)を呼び出す。
2. インスタンス変数を初期化する。
3. 自身のコンストラクタが実行される。

　その２。next(int)のオーバーライドが微妙。ドキュメントには

next の一般規約は、int 値を返し、引数 bits が 1 〜 32 (これを含む) の場合は戻り値のそれだけの個数の下位ビットがほぼ独立に選択されたビット値になる、というものです。

---

Up:目次へ戻る

Last Updated 2003/07/27

ご意見、ご感想はvya07536@nifty.comまで