Java Code Examples for org.apache.commons.math3.util.FastMath#exp()

@Override
public double[][] computeJacobian(double[] variables) {
    double   x1 = variables[0];
    double   x2 = variables[1];
    double   x3 = variables[2];
    double   x4 = variables[3];
    double[][] jacobian = new double[m][];
    for (int i = 0; i < m; ++i) {
        double temp = (i + 1) / 5.0;
        double ti   = FastMath.sin(temp);
        double tmp1 = x1 + temp * x2 - FastMath.exp(temp);
        double tmp2 = x3 + ti   * x4 - FastMath.cos(temp);
        jacobian[i] = new double[] {
            2 * tmp1, 2 * temp * tmp1, 2 * tmp2, 2 * ti * tmp2
        };
    }
    return jacobian;
}
 

/**
 * Creates a Gamma distribution.
 *
 * @param rng Random number generator.
 * @param shape the shape parameter
 * @param scale the scale parameter
 * @param inverseCumAccuracy the maximum absolute error in inverse
 * cumulative probability estimates (defaults to
 * {@link #DEFAULT_INVERSE_ABSOLUTE_ACCURACY}).
 * @throws NotStrictlyPositiveException if {@code shape <= 0} or
 * {@code scale <= 0}.
 * @since 3.1
 */
public GammaDistribution(RandomGenerator rng,
                         double shape,
                         double scale,
                         double inverseCumAccuracy)
    throws NotStrictlyPositiveException {
    super(rng);

    if (shape <= 0) {
        throw new NotStrictlyPositiveException(LocalizedFormats.SHAPE, shape);
    }
    if (scale <= 0) {
        throw new NotStrictlyPositiveException(LocalizedFormats.SCALE, scale);
    }

    this.shape = shape;
    this.scale = scale;
    this.solverAbsoluteAccuracy = inverseCumAccuracy;
    this.shiftedShape = shape + Gamma.LANCZOS_G + 0.5;
    final double aux = FastMath.E / (2.0 * FastMath.PI * shiftedShape);
    this.densityPrefactor2 = shape * FastMath.sqrt(aux) / Gamma.lanczos(shape);
    this.densityPrefactor1 = this.densityPrefactor2 / scale *
            FastMath.pow(shiftedShape, -shape) *
            FastMath.exp(shape + Gamma.LANCZOS_G);
    this.minY = shape + Gamma.LANCZOS_G - FastMath.log(Double.MAX_VALUE);
    this.maxLogY = FastMath.log(Double.MAX_VALUE) / (shape - 1.0);
}
 

/**
 * Computing the log probability in a recursive way on the tree. We return a list
 * here and not the actual number so we can use the logSumExp trick (Link: ??) to
 * avoid underflow.
 * 
 * @param e	- Event to estimate the log probability for
 * @param ball - The region to start from.
 * @return	list of all log probabilities that were computed.
 */
private static List<Double> logPdfRecurse(Event e, Ball ball)
{
	List<Double> logValues = new ArrayList<Double>();
	if(ball.events != null)
	{
		// This is a leaf ball - compute pdf of all events
		logValues.addAll(computeLogKernel(e.getPoint(), ball.events));
		return logValues;
	}
	
	// Else, we need to check if we need to recurse further down or this
	// is a good stopping point
	double minPdf = ball.minPdf(e);
	double maxPdf = ball.maxPdf(e);
	
	if(FastMath.exp(maxPdf) - FastMath.exp(minPdf) < 0.001)
	{
		// No recursion needed 
		logValues.add(FastMath.log(ball.numPoints) + (maxPdf + minPdf)/2);
		
		numSkipped += ball.numPoints;
		
		return logValues;
	}

	// No recursion is definitely needed
	logValues.add(computeLogKernel(e.getPoint(), ball.event.getPoint(), ball.event.getH()));
	logValues.addAll(logPdfRecurse(e, ball.leftBall));
	logValues.addAll(logPdfRecurse(e, ball.rightBall));
	
	return logValues;
}
 

/** {@inheritDoc} */
public double probability(int x) {
    double ret;
    if (x < 0 || x == Integer.MAX_VALUE) {
        ret = 0.0;
    } else if (x == 0) {
        ret = FastMath.exp(-mean);
    } else {
        ret = FastMath.exp(-SaddlePointExpansion.getStirlingError(x) -
              SaddlePointExpansion.getDeviancePart(x, mean)) /
              FastMath.sqrt(MathUtils.TWO_PI * x);
    }
    return ret;
}
 

public static StatisticalReferenceDataset createMGH17()
    throws IOException {
    final BufferedReader in = createBufferedReaderFromResource("MGH17.dat");
    StatisticalReferenceDataset dataset = null;
    try {
        dataset = new StatisticalReferenceDataset(in) {

            @Override
            public double getModelValue(final double x, final double[] a) {
                return a[0] + a[1] * FastMath.exp(-a[3] * x) + a[2] *
                       FastMath.exp(-a[4] * x);
            }

            @Override
            public double[] getModelDerivatives(final double x,
                                                final double[] a) {
                final double[] dy = new double[5];
                dy[0] = 1.0;
                dy[1] = FastMath.exp(-x * a[3]);
                dy[2] = FastMath.exp(-x * a[4]);
                dy[3] = -x * a[1] * dy[1];
                dy[4] = -x * a[2] * dy[2];
                return dy;
            }
        };
    } finally {
        in.close();
    }
    return dataset;
}
 

/** {@inheritDoc} */
public double density(double x) {
    if (x < 0) {
        return 0;
    }
    return FastMath.exp(-x / mean) / mean;
}
 

@Test
public void testErfInv() {
    for (double x = -5.9; x < 5.9; x += 0.01) {
        final double y = Erf.erf(x);
        final double dydx = 2 * FastMath.exp(-x * x) / FastMath.sqrt(FastMath.PI);
        Assert.assertEquals(x, Erf.erfInv(y), 1.0e-15 / dydx);
    }
}
 

public static StatisticalReferenceDataset createMGH17()
    throws IOException {
    final BufferedReader in = createBufferedReaderFromResource("MGH17.dat");
    StatisticalReferenceDataset dataset = null;
    try {
        dataset = new StatisticalReferenceDataset(in) {

            @Override
            public double getModelValue(final double x, final double[] a) {
                return a[0] + a[1] * FastMath.exp(-a[3] * x) + a[2] *
                       FastMath.exp(-a[4] * x);
            }

            @Override
            public double[] getModelDerivatives(final double x,
                                                final double[] a) {
                final double[] dy = new double[5];
                dy[0] = 1.0;
                dy[1] = FastMath.exp(-x * a[3]);
                dy[2] = FastMath.exp(-x * a[4]);
                dy[3] = -x * a[1] * dy[1];
                dy[4] = -x * a[2] * dy[2];
                return dy;
            }
        };
    } finally {
        in.close();
    }
    return dataset;
}
 

/** {@inheritDoc} */
public double cumulativeProbability(double x) {
    double ret;
    if (x <= 0.0) {
        ret = 0.0;
    } else {
        ret = 1.0 - FastMath.exp(-FastMath.pow(x / scale, shape));
    }
    return ret;
}
 

/**
 * Computes the value of the gradient at {@code x}.
 * The components of the gradient vector are the partial
 * derivatives of the function with respect to each of the
 * <em>parameters</em>.
 *
 * @param x Value at which the gradient must be computed.
 * @param param Values for {@code k}, {@code m}, {@code b}, {@code q},
 * {@code a} and  {@code n}.
 * @return the gradient vector at {@code x}.
 * @throws NullArgumentException if {@code param} is {@code null}.
 * @throws DimensionMismatchException if the size of {@code param} is
 * not 6.
 */
public double[] gradient(double x, double ... param)
    throws NullArgumentException,
           DimensionMismatchException,
           NotStrictlyPositiveException {
    validateParameters(param);

    final double b = param[2];
    final double q = param[3];

    final double mMinusX = param[1] - x;
    final double oneOverN = 1 / param[5];
    final double exp = FastMath.exp(b * mMinusX);
    final double qExp = q * exp;
    final double qExp1 = qExp + 1;
    final double factor1 = (param[0] - param[4]) * oneOverN / FastMath.pow(qExp1, oneOverN);
    final double factor2 = -factor1 / qExp1;

    // Components of the gradient.
    final double gk = Logistic.value(mMinusX, 1, b, q, 0, oneOverN);
    final double gm = factor2 * b * qExp;
    final double gb = factor2 * mMinusX * qExp;
    final double gq = factor2 * exp;
    final double ga = Logistic.value(mMinusX, 0, b, q, 1, oneOverN);
    final double gn = factor1 * Math.log(qExp1) * oneOverN;

    return new double[] { gk, gm, gb, gq, ga, gn };
}
 

@Override
public double[][] computeJacobian(double[] variables) {
    double   x1 = variables[0];
    double   x2 = variables[1];
    double   x3 = variables[2];
    double[][] jacobian = new double[m][];
    for (int i = 0; i < m; ++i) {
        double temp = 5.0 * (i + 1) + 45.0 + x3;
        double tmp1 = x2 / temp;
        double tmp2 = FastMath.exp(tmp1);
        double tmp3 = x1 * tmp2 / temp;
        jacobian[i] = new double[] { tmp2, tmp3, -tmp1 * tmp3 };
    }
    return jacobian;
}
 

/** {@inheritDoc} */
public SparseGradient exp() {
    final double e = FastMath.exp(value);
    return new SparseGradient(e, e, derivatives);
}
 

/** {@inheritDoc} */
@Override
public double sample()  {
    final double n = random.nextGaussian();
    return FastMath.exp(scale + shape * n);
}
 

/** {@inheritDoc}
 * @since 3.1
 */
public DerivativeStructure value(final DerivativeStructure t)
    throws DimensionMismatchException {

    double[] f = new double[t.getOrder() + 1];
    final double exp = FastMath.exp(-t.getValue());
    if (Double.isInfinite(exp)) {

        // special handling near lower boundary, to avoid NaN
        f[0] = lo;
        Arrays.fill(f, 1, f.length, 0.0);

    } else {

        // the nth order derivative of sigmoid has the form:
        // dn(sigmoid(x)/dxn = P_n(exp(-x)) / (1+exp(-x))^(n+1)
        // where P_n(t) is a degree n polynomial with normalized higher term
        // P_0(t) = 1, P_1(t) = t, P_2(t) = t^2 - t, P_3(t) = t^3 - 4 t^2 + t...
        // the general recurrence relation for P_n is:
        // P_n(x) = n t P_(n-1)(t) - t (1 + t) P_(n-1)'(t)
        final double[] p = new double[f.length];

        final double inv   = 1 / (1 + exp);
        double coeff = hi - lo;
        for (int n = 0; n < f.length; ++n) {

            // update and evaluate polynomial P_n(t)
            double v = 0;
            p[n] = 1;
            for (int k = n; k >= 0; --k) {
                v = v * exp + p[k];
                if (k > 1) {
                    p[k - 1] = (n - k + 2) * p[k - 2] - (k - 1) * p[k - 1];
                } else {
                    p[0] = 0;
                }
            }

            coeff *= inv;
            f[n]   = coeff * v;

        }

        // fix function value
        f[0] += lo;

    }

    return t.compose(f);

}
 

public double value(double[] point) {
    final double x = point[0], y = point[1];
    final double twoS2 = 2.0 * std * std;
    return 1.0 / (twoS2 * FastMath.PI) * FastMath.exp(-(x * x + y * y) / twoS2);
}
 

/** {@inheritDoc}
* <p>
* From Wikipedia: The probability density function of the L&eacute;vy distribution
* over the domain is
* </p>
* <pre>
* f(x; &mu;, c) = &radic;(c / 2&pi;) * e<sup>-c / 2 (x - &mu;)</sup> / (x - &mu;)<sup>3/2</sup>
* </pre>
* <p>
* For this distribution, {@code X}, this method returns {@code P(X < x)}.
* If {@code x} is less than location parameter &mu;, {@code Double.NaN} is
* returned, as in these cases the distribution is not defined.
* </p>
*/
public double density(final double x) {
    if (x < mu) {
        return Double.NaN;
    }

    final double delta = x - mu;
    final double f     = halfC / delta;
    return FastMath.sqrt(f / FastMath.PI) * FastMath.exp(-f) /delta;
}
 

/**
 * Compute the
 * <a href="http://mathworld.wolfram.com/ExponentialFunction.html" TARGET="_top">
 * exponential function</a> of this complex number.
 * Implements the formula:
 * <pre>
 *  <code>
 *   exp(a + bi) = exp(a)cos(b) + exp(a)sin(b)i
 *  </code>
 * </pre>
 * where the (real) functions on the right-hand side are
 * {@link FastMath#exp}, {@link FastMath#cos}, and
 * {@link FastMath#sin}.
 * <br/>
 * Returns {@link Complex#NaN} if either real or imaginary part of the
 * input argument is {@code NaN}.
 * <br/>
 * Infinite values in real or imaginary parts of the input may result in
 * infinite or NaN values returned in parts of the result.
 * <pre>
 *  Examples:
 *  <code>
 *   exp(1 &plusmn; INFINITY i) = NaN + NaN i
 *   exp(INFINITY + i) = INFINITY + INFINITY i
 *   exp(-INFINITY + i) = 0 + 0i
 *   exp(&plusmn;INFINITY &plusmn; INFINITY i) = NaN + NaN i
 *  </code>
 * </pre>
 *
 * @return <code><i>e</i><sup>this</sup></code>.
 * @since 1.2
 */
public Complex exp() {
    if (isNaN) {
        return NaN;
    }

    double expReal = FastMath.exp(real);
    return createComplex(expReal *  FastMath.cos(imaginary),
                         expReal * FastMath.sin(imaginary));
}
 

/**
 * Compute the
 * <a href="http://mathworld.wolfram.com/ExponentialFunction.html" TARGET=
 * "_top"> exponential function</a> of this complex number. Implements the
 * formula:
 * 
 * <pre>
 *  <code>
 *   exp(a + bi) = exp(a)cos(b) + exp(a)sin(b)i
 *  </code>
 * </pre>
 * 
 * where the (real) functions on the right-hand side are
 * {@link java.lang.Math#exp}, {@link java.lang.Math#cos}, and
 * {@link java.lang.Math#sin}. <br/>
 * Returns {@link Complex#NaN} if either real or imaginary part of the input
 * argument is {@code NaN}. <br/>
 * Infinite values in real or imaginary parts of the input may result in
 * infinite or NaN values returned in parts of the result.
 * 
 * <pre>
 *  Examples:
 *  <code>
 *   exp(1 &plusmn; INFINITY i) = NaN + NaN i
 *   exp(INFINITY + i) = INFINITY + INFINITY i
 *   exp(-INFINITY + i) = 0 + 0i
 *   exp(&plusmn;INFINITY &plusmn; INFINITY i) = NaN + NaN i
 *  </code>
 * </pre>
 *
 * @return <code><i>e</i><sup>this</sup></code>.
 * @since 1.2
 */
public Complex exp() {
	if (isNaN) {
		return NaN;
	}

	double expReal = FastMath.exp(real);
	return createComplex(expReal * FastMath.cos(imaginary), expReal * FastMath.sin(imaginary));
}
 

/** {@inheritDoc}
* <p>
* From Wikipedia: The probability density function of the L&eacute;vy distribution
* over the domain is
* </p>
* <pre>
* f(x; &mu;, c) = &radic;(c / 2&pi;) * e<sup>-c / 2 (x - &mu;)</sup> / (x - &mu;)<sup>3/2</sup>
* </pre>
* <p>
* For this distribution, {@code X}, this method returns {@code P(X < x)}.
* If {@code x} is less than location parameter &mu;, {@code Double.NaN} is
* returned, as in these cases the distribution is not defined.
* </p>
*/
public double density(final double x) {
    if (x < mu) {
        return Double.NaN;
    }

    final double delta = x - mu;
    final double f     = halfC / delta;
    return FastMath.sqrt(f / FastMath.PI) * FastMath.exp(-f) /delta;
}
 

/**
 * Compute the
 * <a href="http://mathworld.wolfram.com/ExponentialFunction.html" TARGET="_top">
 * exponential function</a> of this complex number.
 * Implements the formula:
 * <pre>
 *  <code>
 *   exp(a + bi) = exp(a)cos(b) + exp(a)sin(b)i
 *  </code>
 * </pre>
 * where the (real) functions on the right-hand side are
 * {@link java.lang.Math#exp}, {@link java.lang.Math#cos}, and
 * {@link java.lang.Math#sin}.
 * <br/>
 * Returns {@link Complex#NaN} if either real or imaginary part of the
 * input argument is {@code NaN}.
 * <br/>
 * Infinite values in real or imaginary parts of the input may result in
 * infinite or NaN values returned in parts of the result.
 * <pre>
 *  Examples:
 *  <code>
 *   exp(1 &plusmn; INFINITY i) = NaN + NaN i
 *   exp(INFINITY + i) = INFINITY + INFINITY i
 *   exp(-INFINITY + i) = 0 + 0i
 *   exp(&plusmn;INFINITY &plusmn; INFINITY i) = NaN + NaN i
 *  </code>
 * </pre>
 *
 * @return <code><i>e</i><sup>this</sup></code>.
 * @since 1.2
 */
public Complex exp() {
    if (isNaN) {
        return NaN;
    }

    double expReal = FastMath.exp(real);
    return createComplex(expReal *  FastMath.cos(imaginary),
                         expReal * FastMath.sin(imaginary));
}