Examples with RandomEngine cern.jet.random.engine.RandomEngine used on opensource projects

Example 11 with RandomEngine

use of cern.jet.random.engine.RandomEngine in project tdq-studio-se by Talend.

the class Poisson method nextInt.

/**
 * Returns a random number from the distribution; bypasses the internal state.
 */
public int nextInt(double theMean) {
    /**
     ****************************************************************
     *                                                                *
     * Poisson Distribution - Patchwork Rejection/Inversion           *
     *                                                                *
     ******************************************************************
     *                                                                *
     * For parameter  my < 10  Tabulated Inversion is applied.        *
     * For my >= 10  Patchwork Rejection is employed:                 *
     * The area below the histogram function f(x) is rearranged in    *
     * its body by certain point reflections. Within a large center   *
     * interval variates are sampled efficiently by rejection from    *
     * uniform hats. Rectangular immediate acceptance regions speed   *
     * up the generation. The remaining tails are covered by          *
     * exponential functions.                                         *
     *                                                                *
     ****************************************************************
     */
    RandomEngine gen = this.randomGenerator;
    double my = theMean;
    double t, g, my_k;
    double gx, gy, px, py, e, x, xx, delta, v;
    int sign;
    // static double p,q,p0,pp[36];
    // static long ll,m;
    double u;
    int k, i;
    if (my < SWITCH_MEAN) {
        // CASE B: Inversion- start new table and calculate p0
        if (my != my_old) {
            my_old = my;
            llll = 0;
            p = Math.exp(-my);
            q = p;
            p0 = p;
        // for (k=pp.length; --k >=0; ) pp[k] = 0;
        }
        m = (my > 1.0) ? (int) my : 1;
        for (; ; ) {
            // Step U. Uniform sample
            u = gen.raw();
            k = 0;
            if (u <= p0)
                return (k);
            if (llll != 0) {
                // Step T. Table comparison
                i = (u > 0.458) ? Math.min(llll, m) : 1;
                for (k = i; k <= llll; k++) if (u <= pp[k])
                    return (k);
                if (llll == 35)
                    continue;
            }
            for (k = llll + 1; k <= 35; k++) {
                // Step C. Creation of new prob.
                p *= my / (double) k;
                q += p;
                pp[k] = q;
                if (u <= q) {
                    llll = k;
                    return (k);
                }
            }
            llll = 35;
        }
    } else // end my < SWITCH_MEAN
    if (my < MEAN_MAX) {
        // CASE A: acceptance complement
        // static double        my_last = -1.0;
        // static long int      m,  k2, k4, k1, k5;
        // static double        dl, dr, r1, r2, r4, r5, ll, lr, l_my, c_pm,
        // f1, f2, f4, f5, p1, p2, p3, p4, p5, p6;
        int Dk, X, Y;
        double Ds, U, V, W;
        m = (int) my;
        if (my != my_last) {
            // set-up
            my_last = my;
            // approximate deviation of reflection points k2, k4 from my - 1/2
            Ds = Math.sqrt(my + 0.25);
            // mode m, reflection points k2 and k4, and points k1 and k5, which
            // delimit the centre region of h(x)
            k2 = (int) Math.ceil(my - 0.5 - Ds);
            k4 = (int) (my - 0.5 + Ds);
            k1 = k2 + k2 - m + 1;
            k5 = k4 + k4 - m;
            // range width of the critical left and right centre region
            dl = (double) (k2 - k1);
            dr = (double) (k5 - k4);
            // recurrence constants r(k) = p(k)/p(k-1) at k = k1, k2, k4+1, k5+1
            r1 = my / (double) k1;
            r2 = my / (double) k2;
            r4 = my / (double) (k4 + 1);
            r5 = my / (double) (k5 + 1);
            // reciprocal values of the scale parameters of expon. tail envelopes
            // expon. tail left
            ll = Math.log(r1);
            // expon. tail right
            lr = -Math.log(r5);
            // Poisson constants, necessary for computing function values f(k)
            l_my = Math.log(my);
            c_pm = m * l_my - Arithmetic.logFactorial(m);
            // function values f(k) = p(k)/p(m) at k = k2, k4, k1, k5
            f2 = f(k2, l_my, c_pm);
            f4 = f(k4, l_my, c_pm);
            f1 = f(k1, l_my, c_pm);
            f5 = f(k5, l_my, c_pm);
            // area of the two centre and the two exponential tail regions
            // area of the two immediate acceptance regions between k2, k4
            // immed. left
            p1 = f2 * (dl + 1.0);
            // centre left
            p2 = f2 * dl + p1;
            // immed. right
            p3 = f4 * (dr + 1.0) + p2;
            // centre right
            p4 = f4 * dr + p3;
            // expon. tail left
            p5 = f1 / ll + p4;
            // expon. tail right
            p6 = f5 / lr + p5;
        }
        for (; ; ) {
            // case distinction corresponding to U
            if ((U = gen.raw() * p6) < p2) {
                // immediate acceptance region R2 = [k2, m) *[0, f2),  X = k2, ... m -1
                if ((V = U - p1) < 0.0)
                    return (k2 + (int) (U / f2));
                // immediate acceptance region R1 = [k1, k2)*[0, f1),  X = k1, ... k2-1
                if ((W = V / dl) < f1)
                    return (k1 + (int) (V / f1));
                // computation of candidate X < k2, and its counterpart Y > k2
                // either squeeze-acceptance of X or acceptance-rejection of Y
                Dk = (int) (dl * gen.raw()) + 1;
                if (W <= f2 - Dk * (f2 - f2 / r2)) {
                    // X = k2 - Dk
                    return (k2 - Dk);
                }
                if ((V = f2 + f2 - W) < 1.0) {
                    // quick reject of Y
                    Y = k2 + Dk;
                    if (V <= f2 + Dk * (1.0 - f2) / (dl + 1.0)) {
                        // Y = k2 + Dk
                        return (Y);
                    }
                    // final accept of Y
                    if (V <= f(Y, l_my, c_pm))
                        return (Y);
                }
                X = k2 - Dk;
            } else if (U < p4) {
                // immediate acceptance region R3 = [m, k4+1)*[0, f4), X = m, ... k4
                if ((V = U - p3) < 0.0)
                    return (k4 - (int) ((U - p2) / f4));
                // immediate acceptance region R4 = [k4+1, k5+1)*[0, f5)
                if ((W = V / dr) < f5)
                    return (k5 - (int) (V / f5));
                // computation of candidate X > k4, and its counterpart Y < k4
                // either squeeze-acceptance of X or acceptance-rejection of Y
                Dk = (int) (dr * gen.raw()) + 1;
                if (W <= f4 - Dk * (f4 - f4 * r4)) {
                    // X = k4 + Dk
                    return (k4 + Dk);
                }
                if ((V = f4 + f4 - W) < 1.0) {
                    // quick reject of Y
                    Y = k4 - Dk;
                    if (V <= f4 + Dk * (1.0 - f4) / dr) {
                        // Y = k4 - Dk
                        return (Y);
                    }
                    // final accept of Y
                    if (V <= f(Y, l_my, c_pm))
                        return (Y);
                }
                X = k4 + Dk;
            } else {
                W = gen.raw();
                if (U < p5) {
                    // expon. tail left
                    Dk = (int) (1.0 - Math.log(W) / ll);
                    // 0 <= X <= k1 - 1
                    if ((X = k1 - Dk) < 0)
                        continue;
                    // W -- U(0, h(x))
                    W *= (U - p4) * ll;
                    // quick accept of X
                    if (W <= f1 - Dk * (f1 - f1 / r1))
                        return (X);
                } else {
                    // expon. tail right
                    Dk = (int) (1.0 - Math.log(W) / lr);
                    // X >= k5 + 1
                    X = k5 + Dk;
                    // W -- U(0, h(x))
                    W *= (U - p5) * lr;
                    // quick accept of X
                    if (W <= f5 - Dk * (f5 - f5 * r5))
                        return (X);
                }
            }
            // log f(X) = (X - m)*log(my) - log X! + log m!
            if (Math.log(W) <= X * l_my - Arithmetic.logFactorial(X) - c_pm)
                return (X);
        }
    } else {
        // mean is too large
        return (int) my;
    }
}

Example 12 with RandomEngine

use of cern.jet.random.engine.RandomEngine in project tdq-studio-se by Talend.

the class Benchmark method random.

/**
 * Benchmarks all subclasses
 * @param size the number of random numbers to be generated per subclass.
 * @param print <tt>true</tt> prints each generated number, <tt>false</tt> does not print generated numbers (use this setting for benchmarking).
 * @param mean the mean for distributions that require a mean.
 */
public static void random(int size, boolean print, double mean, String generatorName) {
    System.out.println("Generating " + size + " random numbers per distribution...\n");
    // int large = 100000000;
    int largeVariance = 100;
    // = new MersenneTwister();
    RandomEngine gen;
    try {
        gen = (RandomEngine) Class.forName(generatorName).newInstance();
    } catch (Exception exc) {
        throw new InternalError(exc.getMessage());
    }
    /*
	randomInstance(size,print,new Zeta(10.0, 10.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Zeta(1.0, 1.0, (RandomEngine)gen.clone()));
	randomInstance(size,print,new Zeta(mean, mean, (RandomEngine)gen.clone()));
	randomInstance(size,print,new Zeta(mean, 1/mean, (RandomEngine)gen.clone()));
	//randomInstance(size,print,new Zeta(1/mean, mean, (RandomEngine)gen.clone()));
*/
    /*
	
	randomInstance(size,print,new Beta(10.0, 10.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Beta(1.0, 1.0, (RandomEngine)gen.clone()));
	randomInstance(size,print,new Beta(mean, mean, (RandomEngine)gen.clone()));
	randomInstance(size,print,new Beta(mean, 1/mean, (RandomEngine)gen.clone()));
	randomInstance(size,print,new Beta(1/mean, mean, (RandomEngine)gen.clone()));
	
	randomInstance(size,print,new Uniform((RandomEngine)gen.clone()));
	*/
    randomInstance(size, print, new Poisson(mean, (RandomEngine) gen.clone()));
/*
	randomInstance(size,print,new PoissonSlow(mean,(RandomEngine)gen.clone()));
	
	randomInstance(size,print,new Poisson(3.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new PoissonSlow(3.0,(RandomEngine)gen.clone()));
	
	randomInstance(size,print,new Binomial(1,0.5,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Binomial(5,0.3,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Binomial((int)mean,0.999999999,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Binomial((int)mean,1.0/mean,(RandomEngine)gen.clone()));
	
	randomInstance(size,print,new Exponential(1.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Exponential(3.0,(RandomEngine)gen.clone()));
	
	randomInstance(size,print,new Normal(0.0,1.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Normal(3.0,1.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Normal(mean,largeVariance,(RandomEngine)gen.clone()));
	
	randomInstance(size,print,new BreitWigner(1.0, 0.2, Double.NEGATIVE_INFINITY, (RandomEngine)gen.clone()));
	randomInstance(size,print,new BreitWigner(1.0, 0.2, 1.0, (RandomEngine)gen.clone()));
	
	randomInstance(size,print,new BreitWignerMeanSquare(1.0, 0.2, Double.NEGATIVE_INFINITY, (RandomEngine)gen.clone()));	
	randomInstance(size,print,new BreitWignerMeanSquare(1.0, 0.2, 1.0, (RandomEngine)gen.clone()));
	
	randomInstance(size,print,new ChiSquare(1.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new ChiSquare(5.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new ChiSquare(mean,(RandomEngine)gen.clone()));
	
	randomInstance(size,print,new Gamma(0.2,1.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Gamma(1.0,1.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Gamma(3.0,0.5,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Gamma(mean,0.5,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Gamma(mean,1.0/mean,(RandomEngine)gen.clone()));
	randomInstance(size,print,new Gamma(mean,mean,(RandomEngine)gen.clone()));
	
	randomInstance(size,print,new StudentT(1.0,(RandomEngine)gen.clone()));
	randomInstance(size,print,new StudentT(2.5,(RandomEngine)gen.clone()));
	randomInstance(size,print,new StudentT(mean,(RandomEngine)gen.clone()));
	randomInstance(size,print,new StudentT(1.0/mean,(RandomEngine)gen.clone()));

	int probs = 10000;
	double[] pdf = new double[probs];
	for (int i=0; i<probs; i++) pdf[i]=i*i; // prepare f(x)=x^2 distrib.
	randomInstance(size,print,new Empirical(pdf,Empirical.NO_INTERPOLATION, (RandomEngine)gen.clone()));
	randomInstance(size,print,new Empirical(pdf,Empirical.LINEAR_INTERPOLATION, (RandomEngine)gen.clone()));
	*/
}

Example 13 with RandomEngine

use of cern.jet.random.engine.RandomEngine in project stream-lib by addthis.

the class TestConcurrentStreamSummary method testGeometricDistribution.

@Test
public void testGeometricDistribution() {
    ConcurrentStreamSummary<Integer> vs = new ConcurrentStreamSummary<Integer>(10);
    RandomEngine re = RandomEngine.makeDefault();
    for (int i = 0; i < NUM_ITERATIONS; i++) {
        int z = Distributions.nextGeometric(0.25, re);
        vs.offer(z);
    }
    List<Integer> top = vs.peek(5);
    System.out.println("Geometric:");
    for (Integer e : top) {
        System.out.println(e);
    }
    int tippyTop = top.get(0);
    assertEquals(0, tippyTop);
    System.out.println(vs);
}

Example 14 with RandomEngine

use of cern.jet.random.engine.RandomEngine in project stream-lib by addthis.

the class TestStochasticTopper method testRandomEngine.

@Test
public void testRandomEngine() {
    int[] maxcounts = new int[10];
    int[] counts = new int[20];
    RandomEngine re = RandomEngine.makeDefault();
    for (int i = 0; i < NUM_ITERATIONS; i++) {
        //            int z = Distributions.nextZipfInt(1.2D, re);
        int z = Distributions.nextGeometric(0.25, re);
        if (z > Integer.MAX_VALUE - 9) {
            maxcounts[Integer.MAX_VALUE - z]++;
        }
        if (z < 20) {
            counts[z]++;
        }
    }
    for (int i = 0; i < 20; i++) {
        System.out.println(i + ": " + counts[i]);
    }
    for (int i = 9; i >= 0; i--) {
        System.out.println((Integer.MAX_VALUE - i) + ": " + maxcounts[i]);
    }
}

Example 15 with RandomEngine

use of cern.jet.random.engine.RandomEngine in project processdash by dtuma.

the class EVScheduleConfidenceIntervals method runSimulation.

private void runSimulation() {
    long start = System.currentTimeMillis();
    RandomEngine random = new MersenneTwister();
    int sampleCount = Settings.getInt("ev.simulationSize", BOOTSTRAP_SIZE);
    if (USE_RATIO && indivDates == null) {
        double factor = Math.exp(0.75 * Math.log(randomObjects.size()));
        sampleCount = (int) (sampleCount / factor);
        if (sampleCount < 100)
            sampleCount = 100;
    }
    for (int i = 0; i < sampleCount; i++) runOneTest(random);
    cost.samplesDone();
    date.samplesDone();
    if (optimizedDate != null)
        optimizedDate.samplesDone();
    if (indivDates != null)
        for (int i = 0; i < indivDates.length; i++) indivDates[i].samplesDone();
    long finish = System.currentTimeMillis();
    long elapsed = finish - start;
    System.out.println("schedule simulation took " + elapsed + " ms.");
    if (optimizedDate != null)
        date.debug = true;
}