Example usage for java.lang Math exp

List of usage examples for java.lang Math exp

Introduction

In this page you can find the example usage for java.lang Math exp.

Prototype

@HotSpotIntrinsicCandidate
public static double exp(double a) 

Source Link

Document

Returns Euler's number e raised to the power of a double value.

Usage

From source file:compilateur.test_independance.TestG2.java

private double pochisq(double x, int df) {
    double a, s;/*from ww w. j av  a  2s.c om*/
    double e, c, z;

    if (x <= 0.0 || df < 1) {
        return 1.0;
    }
    a = 0.5 * x;
    boolean even = (df & 1) == 0;
    double y = 0;
    if (df > 1) {
        y = Math.exp(-a);
    }
    s = (even ? y : (2.0 * poz(-Math.sqrt(x))));
    if (df > 2) {
        x = 0.5 * (df - 1.0);
        z = (even ? 1.0 : 0.5);
        if (a > MAX_X) {
            e = (even ? 0.0 : LOG_SQRT_PI);
            c = Math.log(a);
            while (z <= x) {
                e = Math.log(z) + e;
                s += Math.exp(c * z - a - e);
                z += 1.0;
            }
            return s;
        } else {
            e = (even ? 1.0 : (I_SQRT_PI / Math.sqrt(a)));
            c = 0.0;
            while (z <= x) {
                e = e * (a / z);
                c = c + e;
                z += 1.0;
            }
            return c * y + s;
        }
    } else {
        return s;
    }
}

From source file:com.opengamma.analytics.financial.model.volatility.BlackScholesFormulaRepository.java

/**
* The spot delta//from   w ww .j  a  v a 2s  .  c om
* @param spot The spot value of the underlying
* @param strike The Strike
* @param timeToExpiry The time-to-expiry
* @param lognormalVol The log-normal volatility
* @param interestRate The interest rate 
* @param costOfCarry The cost-of-carry rate
* @param isCall true for call
* @return The spot delta
*/
@ExternalFunction
public static double delta(final double spot, final double strike, final double timeToExpiry,
        final double lognormalVol, final double interestRate, final double costOfCarry, final boolean isCall) {
    ArgumentChecker.isTrue(spot >= 0.0, "negative/NaN spot; have {}", spot);
    ArgumentChecker.isTrue(strike >= 0.0, "negative/NaN strike; have {}", strike);
    ArgumentChecker.isTrue(timeToExpiry >= 0.0, "negative/NaN timeToExpiry; have {}", timeToExpiry);
    ArgumentChecker.isTrue(lognormalVol >= 0.0, "negative/NaN lognormalVol; have {}", lognormalVol);
    ArgumentChecker.isFalse(Double.isNaN(interestRate), "interestRate is NaN");
    ArgumentChecker.isFalse(Double.isNaN(costOfCarry), "costOfCarry is NaN");

    double coef = 0.;
    if ((interestRate > LARGE && costOfCarry > LARGE) || (-interestRate > LARGE && -costOfCarry > LARGE)
            || Math.abs(costOfCarry - interestRate) < SMALL) {
        coef = 1.; //ref value is returned
    } else {
        final double rate = costOfCarry - interestRate;
        if (rate > LARGE) {
            return isCall ? Double.POSITIVE_INFINITY : (costOfCarry > LARGE ? 0. : Double.NEGATIVE_INFINITY);
        }
        if (-rate > LARGE) {
            return 0.;
        }
        coef = Math.exp(rate * timeToExpiry);
    }

    if (spot > LARGE * strike) {
        return isCall ? coef : 0.;
    }
    if (spot < SMALL * strike) {
        return isCall ? 0. : -coef;
    }

    final int sign = isCall ? 1 : -1;
    final double rootT = Math.sqrt(timeToExpiry);
    double sigmaRootT = lognormalVol * rootT;
    if (Double.isNaN(sigmaRootT)) {
        sigmaRootT = 1.; //ref value is returned
    }

    double factor = Math.exp(costOfCarry * timeToExpiry);
    if (Double.isNaN(factor)) {
        factor = 1.; //ref value is returned
    }
    double rescaledSpot = spot * factor;

    double d1 = 0.;
    if (Math.abs(spot - strike) < SMALL || sigmaRootT > LARGE || (spot > LARGE && strike > LARGE)) {
        final double coefD1 = (costOfCarry / lognormalVol + 0.5 * lognormalVol);
        final double tmp = coefD1 * rootT;
        d1 = Double.isNaN(tmp) ? 0. : tmp;
    } else {
        if (sigmaRootT < SMALL) {
            return isCall ? (rescaledSpot > strike ? coef : 0.) : (rescaledSpot < strike ? -coef : 0.);
        }
        final double tmp = costOfCarry * rootT / lognormalVol;
        final double sig = (costOfCarry >= 0.) ? 1. : -1.;
        final double scnd = Double.isNaN(tmp)
                ? ((lognormalVol < LARGE && lognormalVol > SMALL) ? sig / lognormalVol : sig * rootT)
                : tmp;
        d1 = Math.log(spot / strike) / sigmaRootT + scnd + 0.5 * sigmaRootT;
    }
    //    if (Double.isNaN(d1)) {
    //      throw new IllegalArgumentException("NaN found");
    //    }
    final double norm = NORMAL.getCDF(sign * d1);

    return norm < SMALL ? 0. : sign * coef * norm;
}

From source file:etymology.util.EtyMath.java

public static double approxGamma(double x) {
    return Math.exp(approxLogGamma(x));
}

From source file:com.opengamma.analytics.financial.model.option.pricing.analytic.formula.BlackBarrierPriceFunction.java

/**
 * Computes the price of a barrier option in the Black world.
 * @param option The underlying European vanilla option.
 * @param barrier The barrier.//from   ww w .  j  av  a  2 s.  c o  m
 * @param rebate The rebate.
 * @param spot The spot price.
 * @param costOfCarry The cost of carry.
 * @param rate The interest rate.
 * @param sigma The Black volatility.
 * @param derivatives Array used to return the derivatives. Will be changed during the call. The derivatives are [0] spot, [1] strike, [2] rate, [3] cost-of-carry, [4] volatility.
 * @return The price.
 */
public double getPriceAdjoint(final EuropeanVanillaOption option, final Barrier barrier, final double rebate,
        final double spot, final double costOfCarry, final double rate, final double sigma,
        final double[] derivatives) {
    Validate.notNull(option, "option");
    Validate.notNull(barrier, "barrier");
    for (int loopder = 0; loopder < 5; loopder++) { // To clean the array.
        derivatives[loopder] = 0.0;
    }
    final boolean isKnockIn = (barrier.getKnockType() == KnockType.IN);
    final boolean isDown = (barrier.getBarrierType() == BarrierType.DOWN);
    final double h = barrier.getBarrierLevel();
    Validate.isTrue(!(barrier.getBarrierType() == BarrierType.DOWN && spot < barrier.getBarrierLevel()),
            "The Data is not consistent with an alive barrier (DOWN and spot<barrier).");
    Validate.isTrue(!(barrier.getBarrierType() == BarrierType.UP && spot > barrier.getBarrierLevel()),
            "The Data is not consistent with an alive barrier (UP and spot>barrier).");
    final boolean isCall = option.isCall();
    final double t = option.getTimeToExpiry();
    final double strike = option.getStrike();
    final int phi = isCall ? 1 : -1;
    final double eta = isDown ? 1 : -1;
    final double df1 = Math.exp(t * (costOfCarry - rate));
    final double df2 = Math.exp(-rate * t);
    if (CompareUtils.closeEquals(sigma, 0, 1e-16)) {
        final double priceBar = 1.0;
        final double df1Bar = rebate * priceBar;
        derivatives[2] = -t * df1 * df1Bar;
        derivatives[3] = t * df1 * df1Bar;
        return df1 * rebate;
    }
    final double sigmaSq = sigma * sigma;
    final double sigmaT = sigma * Math.sqrt(t);
    final double mu = (costOfCarry - 0.5 * sigmaSq) / sigmaSq;
    final double lambda = Math.sqrt(mu * mu + 2 * rate / sigmaSq);
    final double m1 = sigmaT * (1 + mu);
    final double x1 = Math.log(spot / strike) / sigmaT + m1;
    final double x2 = Math.log(spot / h) / sigmaT + m1;
    final double y1 = Math.log(h * h / spot / strike) / sigmaT + m1;
    final double y2 = Math.log(h / spot) / sigmaT + m1;
    final double z = Math.log(h / spot) / sigmaT + lambda * sigmaT;
    final double[] aDerivatives = new double[6];
    final double xA = getAAdjoint(spot, strike, df1, df2, x1, sigmaT, phi, aDerivatives);
    final double[] bDerivatives = new double[6];
    final double xB = getAAdjoint(spot, strike, df1, df2, x2, sigmaT, phi, bDerivatives);
    final double[] cDerivatives = new double[7];
    final double xC = getCAdjoint(spot, strike, df1, df2, y1, sigmaT, h, mu, phi, eta, cDerivatives);
    final double[] dDerivatives = new double[7];
    final double xD = getCAdjoint(spot, strike, df1, df2, y2, sigmaT, h, mu, phi, eta, dDerivatives);
    final double[] fDerivatives = new double[5];
    final double[] eDerivatives = new double[6];
    final double xE = isKnockIn ? getEAdjoint(spot, rebate, df2, x2, y2, sigmaT, h, mu, eta, eDerivatives)
            : getFAdjoint(spot, rebate, z, sigmaT, h, mu, lambda, eta, fDerivatives);
    double xEBar = 0.0;
    double xDBar = 0.0;
    double xCBar = 0.0;
    double xBBar = 0.0;
    double xABar = 0.0;
    double price;
    if (isKnockIn) { // IN start
        if (isDown) { // DOWN start
            if (isCall) { // Call start
                if (strike > h) {
                    xCBar = 1.0;
                    xEBar = 1.0;
                    price = xC + xE;
                } else {
                    xABar = 1.0;
                    xBBar = -1.0;
                    xDBar = 1.0;
                    xEBar = 1.0;
                    price = xA - xB + xD + xE;
                }
            } else { // Put start
                if (strike > h) {
                    xBBar = 1.0;
                    xCBar = -1.0;
                    xDBar = 1.0;
                    xEBar = 1.0;
                    price = xB - xC + xD + xE;
                } else {
                    xABar = 1.0;
                    xEBar = -1.0;
                    price = xA + xE;
                }
            } // DOWN end
        } else { // UP start
            if (isCall) {
                if (strike > h) {
                    xABar = 1.0;
                    xEBar = -1.0;
                    price = xA + xE;
                } else {
                    xBBar = 1.0;
                    xCBar = -1.0;
                    xDBar = 1.0;
                    xEBar = 1.0;
                    price = xB - xC + xD + xE;
                }
            } else {
                if (strike > h) {
                    xABar = 1.0;
                    xBBar = -1.0;
                    xDBar = 1.0;
                    xEBar = 1.0;
                    price = xA - xB + xD + xE;
                } else {
                    xCBar = 1.0;
                    xEBar = 1.0;
                    price = xC + xE;
                }
            } // UP end
        } // IN end
    } else { // OUT start
        if (isDown) { // DOWN start
            if (isCall) { // CALL start
                if (strike > h) {
                    xABar = 1.0;
                    xCBar = -1.0;
                    xEBar = 1.0;
                    price = xA - xC + xE;
                } else {
                    xBBar = 1.0;
                    xDBar = -1.0;
                    xEBar = 1.0;
                    price = xB - xD + xE;
                }
            } else { // PUT start
                if (strike > h) {
                    xABar = 1.0;
                    xBBar = -1.0;
                    xCBar = 1.0;
                    xDBar = -1.0;
                    xEBar = 1.0;
                    price = xA - xB + xC - xD + xE;
                } else {
                    xEBar = 1.0;
                    price = xE;
                } // PUT end
            } // DOWN end
        } else { // UP start
            if (isCall) {
                if (strike > h) {
                    xEBar = 1.0;
                    price = xE;
                } else {
                    xABar = 1.0;
                    xBBar = -1.0;
                    xCBar = 1.0;
                    xDBar = -1.0;
                    xEBar = 1.0;
                    price = xA - xB + xC - xD + xE;
                }
            } else {
                if (strike > h) {
                    xBBar = 1.0;
                    xDBar = -1.0;
                    xEBar = 1.0;
                    price = xB - xD + xE;
                } else {
                    xABar = 1.0;
                    xCBar = -1.0;
                    xEBar = 1.0;
                    price = xA - xC + xE;
                } // PUT end
            } // UP end
        } // OUT end
    }
    // Backward sweep (first step in the foward sweep)
    double zBar = 0.0;
    double y2Bar = 0.0;
    double x2Bar = 0.0;
    double lambdaBar = 0.0;
    double muBar = 0.0;
    double sigmaTBar = 0.0;
    double df2Bar = 0.0;
    if (isKnockIn) {
        y2Bar = eDerivatives[3] * xEBar;
        x2Bar = eDerivatives[2] * xEBar;
        muBar = eDerivatives[5] * xEBar;
        sigmaTBar = eDerivatives[4] * xEBar;
        df2Bar = eDerivatives[1] * xEBar;
        derivatives[0] = eDerivatives[0] * xEBar;
    } else {
        zBar = fDerivatives[1] * xEBar;
        lambdaBar = fDerivatives[4] * xEBar;
        muBar = fDerivatives[3] * xEBar;
        sigmaTBar = fDerivatives[2] * xEBar;
        derivatives[0] = fDerivatives[0] * xEBar;
    }
    y2Bar += dDerivatives[4] * xDBar;
    final double y1Bar = cDerivatives[4] * xCBar;
    x2Bar += bDerivatives[4] * xBBar;
    final double x1Bar = aDerivatives[4] * xABar;
    final double m1Bar = x1Bar + x2Bar + y1Bar + y2Bar;
    muBar += cDerivatives[6] * xCBar + dDerivatives[6] * xDBar + sigmaT * m1Bar + mu / lambda * lambdaBar;
    sigmaTBar += aDerivatives[5] * xABar + bDerivatives[5] * xBBar + cDerivatives[5] * xCBar
            + dDerivatives[5] * xDBar + (lambda - Math.log(h / spot) / (sigmaT * sigmaT)) * zBar
            - Math.log(h / spot) / (sigmaT * sigmaT) * y2Bar
            - Math.log(h * h / spot / strike) / (sigmaT * sigmaT) * y1Bar
            - Math.log(spot / h) / (sigmaT * sigmaT) * x2Bar
            - Math.log(spot / strike) / (sigmaT * sigmaT) * x1Bar + (1 + mu) * m1Bar;
    final double sigmaSqBar = -costOfCarry / (sigmaSq * sigmaSq) * muBar
            - rate / (sigmaSq * sigmaSq) / lambda * lambdaBar;
    df2Bar += aDerivatives[3] * xABar + bDerivatives[3] * xBBar + cDerivatives[3] * xCBar
            + dDerivatives[3] * xDBar;
    final double df1Bar = aDerivatives[2] * xABar + bDerivatives[2] * xBBar + cDerivatives[2] * xCBar
            + dDerivatives[2] * xDBar;
    derivatives[0] += aDerivatives[0] * xABar + bDerivatives[0] * xBBar + cDerivatives[0] * xCBar
            + dDerivatives[0] * xDBar + 1.0 / spot / sigmaT * x1Bar + 1.0 / spot / sigmaT * x2Bar
            + -1.0 / spot / sigmaT * y1Bar + -1.0 / spot / sigmaT * y2Bar - 1.0 / spot / sigmaT * zBar;
    derivatives[1] = aDerivatives[1] * xABar + bDerivatives[1] * xBBar + cDerivatives[1] * xCBar
            + dDerivatives[1] * xDBar + -1.0 / strike / sigmaT * x1Bar - 1 / strike / sigmaT * y1Bar;
    derivatives[2] = -t * df1 * df1Bar - t * df2 * df2Bar + 1.0 / lambda / sigmaSq * lambdaBar;
    derivatives[3] = t * df1 * df1Bar + 1.0 / sigmaSq * muBar;
    derivatives[4] = 2 * sigma * sigmaSqBar + Math.sqrt(t) * sigmaTBar;
    return price;
}

From source file:mastodon.algorithms.SALinearAlgorithm.java

protected void setNewBest() {
    if (currScore[0] > maxScore[0]) { //set new optimum
        maxScore = currScore; //might need a clone here
        maxScorePruning.clear();//from   ww w.ja v  a  2 s .com
        maxScorePruning.put((BitSet) currPruning.clone(), currScore.clone());
    } else if (currScore[0] == maxScore[0] && currScore[1] != 1) { //save variations with same score, but no need to if it produces no results
        maxScorePruning.put((BitSet) currPruning.clone(), currScore.clone());
    }

    if (Random.nextDouble() < Math.exp((currScore[0] - prevScore[0]) / currTemp)) {
        prevPruning = (BitSet) currPruning.clone();
        prevScore = currScore.clone();

        for (int a = currPruning.nextSetBit(0); a >= 0; a = currPruning.nextSetBit(a + 1)) {
            pruningFreq.put(a, pruningFreq.get(a) + 1);
        }
        totalPruningFreq++;

    } //try different pruning otherwise
}

From source file:org.fhcrc.cpl.viewer.quant.gui.LogRatioHistMouseListener.java

/**
 * Undraw the previous selected region (if it was drawn), calculate the new regions, draw again, save
 * the points, and draw the numeric ratio in its little box
 * @param e//from  www  .  j a v a 2 s .c  o m
 */
public void mouseDragged(MouseEvent e) {

    if (this.selectedRegionStart == null || e.getX() < this.selectedRegionStart.getX()) {
        return;
    }

    if (this.selectedRegion != null)
        drawOrUndrawRegion();

    // Erase the previous zoom rectangle (if any)...
    Rectangle2D scaledDataArea = _chartPanel.getScreenDataArea();

    this.selectedRegion = new Rectangle2D.Double(this.selectedRegionStart.getX(), scaledDataArea.getMinY(),
            Math.min(Math.abs(e.getX() - selectedRegionStart.getX()),
                    _chartPanel.getWidth() - this.selectedRegionStart.getX()),
            scaledDataArea.getHeight());
    transformAndSaveSelectedRegion();

    // Draw the new zoom rectangle...
    drawOrUndrawRegion();

    lastMousedRatio = Rounder.round(Math.exp(transformMouseXValue(e.getX())), 2);
    drawRatioInBox(getChartPanelGraphics());
}

From source file:com.opengamma.analytics.financial.model.finitedifference.ForwardPDETest.java

@Test(enabled = false)
public void debugTest() {

    final Function<Double, Double> lvFunc = new Function<Double, Double>() {
        @Override//w  w w.  jav a 2  s  .c o m
        public Double evaluate(final Double... tm) {
            Validate.isTrue(tm.length == 2);
            final double t = tm[0];
            final double m = tm[1];
            final double x = Math.log(m);
            if (Math.abs(x) > Math.sqrt(t) * 1.2) {
                return 0.0;
            }
            return 0.4;
        }
    };
    final LocalVolatilitySurfaceMoneyness lv = new LocalVolatilitySurfaceMoneyness(
            FunctionalDoublesSurface.from(lvFunc), new ForwardCurve(1.0));

    final Function1D<Double, Double> initCon = new Function1D<Double, Double>() {
        @Override
        public Double evaluate(final Double x) {
            return Math.max(0, 1 - Math.exp(x));
        }
    };

    final ConvectionDiffusionPDESolver solver = new ThetaMethodFiniteDifference(0.50, true);
    final ConvectionDiffusionPDE1DStandardCoefficients pde = getForwardLocalVol(lv);
    final double xMin = -2.5;
    final double xMax = 2.5;

    PDEUtilityTools.printSurface("lv", lv.getSurface(), 0, 2.0, Math.exp(xMin), Math.exp(xMax));

    final DirichletBoundaryCondition lower = new DirichletBoundaryCondition(initCon.evaluate(xMin), xMin);
    final DirichletBoundaryCondition upper = new DirichletBoundaryCondition(initCon.evaluate(xMax), xMax);
    final MeshingFunction timeMesh = new ExponentialMeshing(0, 2.0, 12, 0.0);
    final MeshingFunction spaceMesh = new ExponentialMeshing(xMin, xMax, 17, 0.0);
    final PDEGrid1D grid = new PDEGrid1D(timeMesh, spaceMesh);
    final PDE1DDataBundle<ConvectionDiffusionPDE1DCoefficients> db = new PDE1DDataBundle<ConvectionDiffusionPDE1DCoefficients>(
            pde, initCon, lower, upper, grid);

    final PDEFullResults1D prices = (PDEFullResults1D) solver.solve(db);
    PDEUtilityTools.printSurface("prices", prices);
}

From source file:com.opengamma.analytics.financial.model.option.pricing.fourier.FFTModelGreeks.java

/**
 * //from ww w .  ja va2  s.c  o m
 * @param forward The forward value of the underlying
 * @param discountFactor 
 * @param t Time to expiry
 * @param isCall true for call 
 * @param ce The Characteristic Exponent (log of characteristic function) of the returns of the underlying
 * @param nStrikesBelowATM maximum number of strikes below ATM to be returned 
 * @param nStrikesAboveATM maximum number of strikes above ATM to be returned 
 * @param alpha Regularization factor. Values of 0 or -1 are not allowed. -0.5 is recommended  
 * @param delta The spacing for sampling the function 
 * @param n The (zero padded) array of sample values. <b>Use a power of 2</b>
 * @param m The actual number of samples. Need n >= 2m-1
 * @return  an array of arrays where is first array is the strikes, the second the prices, the first the derivatives of price wrt the first parameter etc 
 */
public double[][] getGreeks(final double forward, final double discountFactor, final double t,
        final boolean isCall, final MartingaleCharacteristicExponent ce, final int nStrikesBelowATM,
        final int nStrikesAboveATM, final double alpha, final double delta, final int n, final int m) {

    Validate.notNull(ce, "characteristic exponent");
    Validate.isTrue(nStrikesBelowATM >= 0, "nStrikesBelowATM >= 0");
    Validate.isTrue(nStrikesAboveATM >= 0, "nStrikesAboveATM >= 0");
    Validate.isTrue(alpha != 0.0 && alpha != -1.0, "alpha cannot be -1 or 0");
    Validate.isTrue(delta > 0.0, "need delta > 0");
    Validate.isTrue(m > 0, "need m > 0");
    Validate.isTrue(n >= 2 * m - 1, "need n > 2m-1");

    final Function1D<ComplexNumber, ComplexNumber[]> func = ce.getAdjointFunction(t);
    final int halfN = n % 2 == 0 ? n / 2 : (n + 1) / 2;
    final double a = -(halfN - 1) * delta;
    final ComplexNumber[][] z = getPaddedArrays(alpha, delta, n, m, func, halfN);
    int size = z.length;
    final ComplexNumber[][] x = new ComplexNumber[size][];
    for (int i = 0; i < size; i++) {
        x[i] = JTransformsWrapper.transform1DComplex(z[i]);
    }
    final int nLowStrikes = Math.min(halfN, nStrikesBelowATM);
    final int nHighStrikes = Math.min(n - halfN, nStrikesAboveATM);
    final int p = 1 + nLowStrikes + nHighStrikes;
    final double[][] res = new double[size + 1][p];
    final double deltaK = 2 * Math.PI / delta / n;
    for (int i = 0; i < nLowStrikes; i++) {
        final double k = (i - nLowStrikes) * deltaK;
        res[0][i] = forward * Math.exp(k);
        res[1][i] = discountFactor * forward
                * getReducedPrice(x[0][i + n - nLowStrikes], alpha, delta, k, a, isCall);
        for (int j = 1; j < size; j++) {
            res[j + 1][i] = discountFactor * forward
                    * getReducedGreek(x[j][i + n - nLowStrikes], alpha, delta, k, a);
        }
    }
    for (int i = nLowStrikes; i < p; i++) {
        final double k = (i - nLowStrikes) * deltaK;
        res[0][i] = forward * Math.exp(k);
        res[1][i] = discountFactor * forward
                * getReducedPrice(x[0][i - nLowStrikes], alpha, delta, k, a, isCall);
        for (int j = 1; j < size; j++) {
            res[j + 1][i] = discountFactor * forward
                    * getReducedGreek(x[j][i - nLowStrikes], alpha, delta, k, a);
        }
    }
    return res;
}

From source file:beast.evomodel.speciation.BirthDeathModel.java

public double logNodeProbability(Tree tree, NodeRef node) {
    final double height = tree.getNodeHeight(node);
    final double r = getR();
    final double mrh = -r * height;
    final double a = getA();
    final double rho = getRho();

    if (conditionalOnRoot && tree.isRoot(node)) {
        return (tree.getTaxonCount() - 2) * logConditioningTerm(height);
    }//w ww .j av  a  2 s .c  om

    final double z = Math.log(rho + ((1 - rho) - a) * Math.exp(mrh));
    double l = -2 * z + mrh;

    if (!conditionOnOrigin && !conditionalOnRoot && tree.isRoot(node)) {
        l += mrh - z;
    }

    return l;
}

From source file:etomica.math.SpecialFunctions.java

public static void main(String[] args) {
    System.out.println(confluentHypergeometric1F1(-0.25, 0.5, 1.0));
    System.out.println();/*from  w w w .ja  v  a  2  s  .  c  om*/

    for (int i = -10; i < 2; i++) {
        double g = SpecialFunctions.gamma(i - 0.5);
        System.out.println((i - 0.5) + " " + g);
    }
    System.out.println(0 + " " + SpecialFunctions.gamma(0));
    System.out.println(
            0.25 + " " + SpecialFunctions.gamma(0.25) + " " + Math.exp(SpecialFunctions.lnGamma(0.25)));
    System.out.println(0.5 + " " + SpecialFunctions.gamma(0.5) + " " + Math.exp(SpecialFunctions.lnGamma(0.5)));
    System.out.println(1 + " " + SpecialFunctions.gamma(1));
    System.out.println();
    //       for (int i=2; i<1000; i++) {
    //           double lnfac = SpecialFunctions.lnFactorial(i);
    //           System.out.println(i+" "+lnfac+" "+(SpecialFunctions.lnGamma(i+1)-lnfac)/lnfac);
    //       }
    double w = SpecialFunctions.wigner3J(4, 2, 2, -2, 2, 0);
    System.out.println(w);
    //        System.out.println(besselI(0.32,0.25)+" "+besselI(-1.59, 0.25));
    //        double ap3 = SpecialFunctions.besselI(3, 200);
    //        System.out.println("ap3 = " + ap3);
    double a = 1E-24 / Mole.UNIT.fromSim(1);
    System.out.println(a * a * 3007.044183696613 + " " + a * a * 198.22529144211072);
    System.out.println(a * a * (-28532.603135414935 + 11560.282978358246) + " " + a * a
            * Math.sqrt(198.22529144211072 * 198.22529144211072 + 33.767924918154755 * 33.767924918154755));

}