Example usage for org.apache.commons.math.random MersenneTwister MersenneTwister

List of usage examples for org.apache.commons.math.random MersenneTwister MersenneTwister

Introduction

In this page you can find the example usage for org.apache.commons.math.random MersenneTwister MersenneTwister.

Prototype

public MersenneTwister() 

Source Link

Document

Creates a new random number generator.

Usage

From source file:net.sourceforge.jags.model.ModelTest.java

static double[] normalSample() {
    MersenneTwister rng = new MersenneTwister();
    double rval[] = new double[1000];
    for (int i = 0; i < rval.length; ++i) {
        rval[i] = rng.nextGaussian();/*from www. j  av a2s  .  c o m*/
    }
    return rval;
}

From source file:jmona.random.RandomUtilsTester.java

/**
 * Test method for {@link jmona.random.RandomUtils#setRandomData(RandomData)}.
 *//*from  w w w.  j ava 2  s . co  m*/
@Test
public void testSetRandomData() {
    final RandomData newRandomData = new RandomDataImpl(new MersenneTwister());
    RandomUtils.setRandomData(newRandomData);
    assertSame(newRandomData, RandomUtils.randomData());
}

From source file:cz.paulrz.montecarlo.random.Sobol.java

/**
 * @param dim dimension of the Sobol sequence.
 *//*from  w ww  . ja va 2 s.c  om*/
public Sobol(int dim) throws Exception {
    this.dim = dim;
    if (dim > 300) {
        throw new Exception("Sobol sequence only implemented for dimension at most 300." + "\nExiting.");
    }

    // array of polynomial degrees and coefficient arrays
    g = new int[dim];
    p = new int[dim][];
    // degree zero
    g[0] = 0;
    p[0] = new int[1];
    p[0][0] = 1;
    // positive degree polynomials read from pp
    int k = 1;
    for (int i = 0; i < pp.length; i++)
        for (int j = 0; j < pp[i].length; j++)
            if (k < dim) {
                g[k] = i + 1;
                read_prim_pol(g[k], pp[i][j], k);
                k++;
            } else
                break;

    /* DEBUG: print the primitive polynomials:
    System.out.println("\nPrimitive polynomials:\n");
    for(k=0;k<50;k++)
    {  
        for(int j=0;j<=g[k];j++)System.out.print(p[k][j]);
        System.out.println(""); 
    }
    */

    // initialize the array of direction integers
    v = new long[dim][bits];
    for (int j = 0; j < bits; j++)
        v[0][j] = (1L << (bits - j - 1));

    if (dim > 1)
        v[1][0] = (1L << bits - 1);

    if (dim > 2) {
        v[2][0] = (1L << bits - 1);
        v[2][1] = (1L << bits - 2);
    }

    if (dim > 3) {
        v[3][0] = (1L << bits - 1);
        v[3][1] = (3L << bits - 2);
        v[3][2] = (7L << bits - 3);
    }

    if (dim > 4) {
        v[4][0] = (1L << bits - 1);
        v[4][1] = (1L << bits - 2);
        v[4][2] = (5L << bits - 3);
    }

    if (dim > 5) {
        v[5][0] = (1L << bits - 1);
        v[5][1] = (3L << bits - 2);
        v[5][2] = (1L << bits - 3);
        v[5][3] = (1L << bits - 4);
    }

    if (dim > 6) {
        v[6][0] = (1L << bits - 1);
        v[6][1] = (1L << bits - 2);
        v[6][2] = (3L << bits - 3);
        v[6][3] = (7L << bits - 4);
    }

    if (dim > 7) {
        v[7][0] = (1L << bits - 1);
        v[7][1] = (3L << bits - 2);
        v[7][2] = (3L << bits - 3);
        v[7][3] = (9L << bits - 4);
        v[7][4] = (9L << bits - 5);
    }

    if (dim > 8) {
        v[8][0] = (1L << bits - 1);
        v[8][1] = (3L << bits - 2);
        v[8][2] = (7L << bits - 3);
        v[8][3] = (7L << bits - 4);
        v[8][4] = (21L << bits - 5);
    }

    if (dim > 9) {
        v[9][0] = (1L << bits - 1);
        v[9][1] = (1L << bits - 2);
        v[9][2] = (5L << bits - 3);
        v[9][3] = (11L << bits - 4);
        v[9][4] = (27L << bits - 5);
    }

    if (dim > 10) {
        v[10][0] = (1L << bits - 1);
        v[10][1] = (1L << bits - 2);
        v[10][2] = (7L << bits - 3);
        v[10][3] = (3L << bits - 4);
        v[10][4] = (29L << bits - 5);
    }

    if (dim > 11) {
        v[11][0] = (1L << bits - 1);
        v[11][1] = (3L << bits - 2);
        v[11][2] = (7L << bits - 3);
        v[11][3] = (13L << bits - 4);
        v[11][4] = (3L << bits - 5);
    }

    if (dim > 12) {
        v[12][0] = (1L << bits - 1);
        v[12][1] = (3L << bits - 2);
        v[12][2] = (5L << bits - 3);
        v[12][3] = (1L << bits - 4);
        v[12][4] = (15L << bits - 5);
    }

    if (dim > 13) {
        v[13][0] = (1L << bits - 1);
        v[13][1] = (1L << bits - 2);
        v[13][2] = (1L << bits - 3);
        v[13][3] = (9L << bits - 4);
        v[13][4] = (23L << bits - 5);
        v[13][5] = (37L << bits - 6);
    }

    if (dim > 14) {
        v[14][0] = (1L << bits - 1);
        v[14][1] = (1L << bits - 2);
        v[14][2] = (3L << bits - 3);
        v[14][3] = (13L << bits - 4);
        v[14][4] = (11L << bits - 5);
        v[14][5] = (7L << bits - 6);
    }

    if (dim > 15) {
        v[15][0] = (1L << bits - 1);
        v[15][1] = (3L << bits - 2);
        v[15][2] = (3L << bits - 3);
        v[15][3] = (5L << bits - 4);
        v[15][4] = (19L << bits - 5);
        v[15][5] = (33L << bits - 6);
    }

    if (dim > 16) {
        v[16][0] = (1L << bits - 1);
        v[16][1] = (1L << bits - 2);
        v[16][2] = (7L << bits - 3);
        v[16][3] = (13L << bits - 4);
        v[16][4] = (25L << bits - 5);
        v[16][5] = (5L << bits - 6);
    }

    if (dim > 17) {
        v[17][0] = (1L << bits - 1);
        v[17][1] = (1L << bits - 2);
        v[17][2] = (1L << bits - 3);
        v[17][3] = (13L << bits - 4);
        v[17][4] = (15L << bits - 5);
        v[17][5] = (39L << bits - 6);
    }

    if (dim > 18) {
        v[18][0] = (1L << bits - 1);
        v[18][1] = (3L << bits - 2);
        v[18][2] = (5L << bits - 3);
        v[18][3] = (11L << bits - 4);
        v[18][4] = (7L << bits - 5);
        v[18][5] = (11L << bits - 6);
    }

    if (dim > 19) {
        v[19][0] = (1L << bits - 1);
        v[19][1] = (3L << bits - 2);
        v[19][2] = (1L << bits - 3);
        v[19][3] = (7L << bits - 4);
        v[19][4] = (3L << bits - 5);
        v[19][5] = (23L << bits - 6);
        v[19][6] = (79L << bits - 7);
    }

    if (dim > 20) {
        v[20][0] = (1L << bits - 1);
        v[20][1] = (3L << bits - 2);
        v[20][2] = (1L << bits - 3);
        v[20][3] = (15L << bits - 4);
        v[20][4] = (17L << bits - 5);
        v[20][5] = (63L << bits - 6);
        v[20][6] = (13L << bits - 7);
    }

    if (dim > 21) {
        v[21][0] = (1L << bits - 1);
        v[21][1] = (3L << bits - 2);
        v[21][2] = (3L << bits - 3);
        v[21][3] = (3L << bits - 4);
        v[21][4] = (25L << bits - 5);
        v[21][5] = (17L << bits - 6);
        v[21][6] = (115L << bits - 7);
    }

    if (dim > 22) {
        v[22][0] = (1L << bits - 1);
        v[22][1] = (3L << bits - 2);
        v[22][2] = (7L << bits - 3);
        v[22][3] = (9L << bits - 4);
        v[22][4] = (31L << bits - 5);
        v[22][5] = (29L << bits - 6);
        v[22][6] = (17L << bits - 7);
    }

    if (dim > 23) {
        v[23][0] = (1L << bits - 1);
        v[23][1] = (1L << bits - 2);
        v[23][2] = (3L << bits - 3);
        v[23][3] = (15L << bits - 4);
        v[23][4] = (29L << bits - 5);
        v[23][5] = (15L << bits - 6);
        v[23][6] = (41L << bits - 7);
    }

    if (dim > 24) {
        v[24][0] = (1L << bits - 1);
        v[24][1] = (3L << bits - 2);
        v[24][2] = (1L << bits - 3);
        v[24][3] = (9L << bits - 4);
        v[24][4] = (5L << bits - 5);
        v[24][5] = (21L << bits - 6);
        v[24][6] = (119L << bits - 7);
    }

    if (dim > 25) {
        v[25][0] = (1L << bits - 1);
        v[25][1] = (1L << bits - 2);
        v[25][2] = (5L << bits - 3);
        v[25][3] = (5L << bits - 4);
        v[25][4] = (1L << bits - 5);
        v[25][5] = (27L << bits - 6);
        v[25][6] = (33L << bits - 7);
    }

    if (dim > 26) {
        v[26][0] = (1L << bits - 1);
        v[26][1] = (1L << bits - 2);
        v[26][2] = (3L << bits - 3);
        v[26][3] = (1L << bits - 4);
        v[26][4] = (23L << bits - 5);
        v[26][5] = (13L << bits - 6);
        v[26][6] = (75L << bits - 7);
    }

    if (dim > 27) {
        v[27][0] = (1L << bits - 1);
        v[27][1] = (1L << bits - 2);
        v[27][2] = (7L << bits - 3);
        v[27][3] = (7L << bits - 4);
        v[27][4] = (19L << bits - 5);
        v[27][5] = (25L << bits - 6);
        v[27][6] = (105L << bits - 7);
    }

    if (dim > 28) {
        v[28][0] = (1L << bits - 1);
        v[28][1] = (3L << bits - 2);
        v[28][2] = (5L << bits - 3);
        v[28][3] = (5L << bits - 4);
        v[28][4] = (21L << bits - 5);
        v[28][5] = (9L << bits - 6);
        v[28][6] = (7L << bits - 7);
    }

    if (dim > 29) {
        v[29][0] = (1L << bits - 1);
        v[29][1] = (1L << bits - 2);
        v[29][2] = (1L << bits - 3);
        v[29][3] = (15L << bits - 4);
        v[29][4] = (5L << bits - 5);
        v[29][5] = (49L << bits - 6);
        v[29][6] = (59L << bits - 7);
    }

    if (dim > 30) {
        v[30][0] = (1L << bits - 1);
        v[30][1] = (3L << bits - 2);
        v[30][2] = (5L << bits - 3);
        v[30][3] = (15L << bits - 4);
        v[30][4] = (17L << bits - 5);
        v[30][5] = (19L << bits - 6);
        v[30][6] = (21L << bits - 7);
    }

    if (dim > 31) {
        v[31][0] = (1L << bits - 1);
        v[31][1] = (1L << bits - 2);
        v[31][2] = (7L << bits - 3);
        v[31][3] = (11L << bits - 4);
        v[31][4] = (13L << bits - 5);
        v[31][5] = (29L << bits - 6);
        v[31][6] = (3L << bits - 7);
    }

    // random initialization in dimension bigger than 32
    MersenneTwister mt = new MersenneTwister();
    for (k = 32; k < dim; k++) {
        for (int l = 0; l < g[k]; l++) {
            double u = mt.nextDouble();
            long f = (1L << l + 1), n = (int) (f * u);
            while (n % 2 == 0) {
                u = mt.nextDouble();
                n = (int) (f * u);
            }

            v[k][l] = (n << (bits - l - 1));
        }
    } // end direction integer initialization

    // computation of direction integer v_kl for k>=degree[k]
    for (k = 1; k < dim; k++)
        for (int l = g[k]; l < bits; l++) {
            long n = (v[k][l - g[k]] >> g[k]);
            for (int j = 1; j <= g[k]; j++)
                if (p[k][j] != 0)
                    n = n ^ v[k][l - j];

            v[k][l] = n;

        }

    // initialize the vector of Sobol integers and Sobol points
    index = 1;
    x_int = new long[dim];
    for (k = 0; k < dim; k++)
        x_int[k] = v[k][0];

    x = new double[dim];
    for (k = 0; k < dim; k++)
        x[k] = ((double) x_int[k]) / N;

}