Example usage for org.apache.commons.math3.analysis.integration RombergIntegrator RombergIntegrator

List of usage examples for org.apache.commons.math3.analysis.integration RombergIntegrator RombergIntegrator

Introduction

In this page you can find the example usage for org.apache.commons.math3.analysis.integration RombergIntegrator RombergIntegrator.

Prototype

public RombergIntegrator(final double relativeAccuracy, final double absoluteAccuracy,
        final int minimalIterationCount, final int maximalIterationCount)
        throws NotStrictlyPositiveException, NumberIsTooSmallException, NumberIsTooLargeException 

Source Link

Document

Build a Romberg integrator with given accuracies and iterations counts.

Usage

From source file:beast.math.distribution.GammaDistributionTest.java

@Test
public void testPdf() {

    final int numberOfTests = 100;
    double totErr = 0;
    double ptotErr = 0;
    int np = 0;/*from  w w w  .j  a va 2s  . co m*/
    double qtotErr = 0;

    Random random = new Random(37);

    for (int i = 0; i < numberOfTests; i++) {
        final double mean = .01 + (3 - 0.01) * random.nextDouble();
        final double var = .01 + (3 - 0.01) * random.nextDouble();

        final double scale = var / mean;
        final double shape = mean / scale;

        final GammaDistribution gamma = new GammaDistribution(shape, scale);

        final double value = gamma.nextGamma();

        final double mypdf = mypdf(value, shape, scale);
        final double pdf = gamma.pdf(value);
        if (Double.isInfinite(mypdf) && Double.isInfinite(pdf)) {
            continue;
        }

        assertFalse(Double.isNaN(mypdf));
        assertFalse(Double.isNaN(pdf));

        totErr += mypdf != 0 ? Math.abs((pdf - mypdf) / mypdf) : pdf;

        assertFalse("nan", Double.isNaN(totErr));
        //assertEquals("" + shape + "," + scale + "," + value, mypdf,gamma.pdf(value),1e-10);

        final double cdf = gamma.cdf(value);
        UnivariateFunction f = new UnivariateFunction() {
            public double value(double v) {
                return mypdf(v, shape, scale);
            }
        };
        final UnivariateIntegrator integrator = new RombergIntegrator(MachineAccuracy.SQRT_EPSILON, 1e-14, 1,
                16);

        double x;
        try {
            x = integrator.integrate(16, f, 0.0, value);
            ptotErr += cdf != 0.0 ? Math.abs(x - cdf) / cdf : x;
            np += 1;
            //assertTrue("" + shape + "," + scale + "," + value + " " + Math.abs(x-cdf)/x + "> 1e-6", Math.abs(1-cdf/x) < 1e-6);

            //System.out.println(shape + ","  + scale + " " + value);
        } catch (MaxCountExceededException e) {
            // can't integrate , skip test
            //  System.out.println(shape + ","  + scale + " skipped");
        }

        final double q = gamma.quantile(cdf);
        qtotErr += q != 0 ? Math.abs(q - value) / q : value;
        // assertEquals("" + shape + "," + scale + "," + value + " " + Math.abs(q-value)/value, q, value, 1e-6);
    }
    //System.out.println( !Double.isNaN(totErr) );
    // System.out.println(totErr);
    // bad test, but I can't find a good threshold that works for all individual cases 
    assertTrue("failed " + totErr / numberOfTests, totErr / numberOfTests < 1e-7);
    assertTrue("failed " + ptotErr / np, np > 0 ? (ptotErr / np < 1e-5) : true);
    assertTrue("failed " + qtotErr / numberOfTests, qtotErr / numberOfTests < 1e-7);
}