/*
 * Licensed to the Apache Software Foundation (ASF) under one or more
 * contributor license agreements.  See the NOTICE file distributed with
 * this work for additional information regarding copyright ownership.
 * The ASF licenses this file to You under the Apache License, Version 2.0
 * (the "License"); you may not use this file except in compliance with
 * the License.  You may obtain a copy of the License at
 *
 *      http://www.apache.org/licenses/LICENSE-2.0
 *
 * Unless required by applicable law or agreed to in writing, software
 * distributed under the License is distributed on an "AS IS" BASIS,
 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
 * See the License for the specific language governing permissions and
 * limitations under the License.
 */
package org.apache.commons.math3.util;

import org.apache.commons.math3.exception.ConvergenceException;
import org.apache.commons.math3.exception.MaxCountExceededException;
import org.apache.commons.math3.exception.util.LocalizedFormats;

Provides a generic means to evaluate continued fractions. Subclasses simply provided the a and b coefficients to evaluate the continued fraction.

References:

/** * Provides a generic means to evaluate continued fractions. Subclasses simply * provided the a and b coefficients to evaluate the continued fraction. * * <p> * References: * <ul> * <li><a href="http://mathworld.wolfram.com/ContinuedFraction.html"> * Continued Fraction</a></li> * </ul> * </p> * */
public abstract class ContinuedFraction {
Maximum allowed numerical error.
/** Maximum allowed numerical error. */
private static final double DEFAULT_EPSILON = 10e-9;
Default constructor.
/** * Default constructor. */
protected ContinuedFraction() { super(); }
Access the n-th a coefficient of the continued fraction. Since a can be a function of the evaluation point, x, that is passed in as well.
Params:
  • n – the coefficient index to retrieve.
  • x – the evaluation point.
Returns:the n-th a coefficient.
/** * Access the n-th a coefficient of the continued fraction. Since a can be * a function of the evaluation point, x, that is passed in as well. * @param n the coefficient index to retrieve. * @param x the evaluation point. * @return the n-th a coefficient. */
protected abstract double getA(int n, double x);
Access the n-th b coefficient of the continued fraction. Since b can be a function of the evaluation point, x, that is passed in as well.
Params:
  • n – the coefficient index to retrieve.
  • x – the evaluation point.
Returns:the n-th b coefficient.
/** * Access the n-th b coefficient of the continued fraction. Since b can be * a function of the evaluation point, x, that is passed in as well. * @param n the coefficient index to retrieve. * @param x the evaluation point. * @return the n-th b coefficient. */
protected abstract double getB(int n, double x);
Evaluates the continued fraction at the value x.
Params:
  • x – the evaluation point.
Throws:
Returns:the value of the continued fraction evaluated at x.
/** * Evaluates the continued fraction at the value x. * @param x the evaluation point. * @return the value of the continued fraction evaluated at x. * @throws ConvergenceException if the algorithm fails to converge. */
public double evaluate(double x) throws ConvergenceException { return evaluate(x, DEFAULT_EPSILON, Integer.MAX_VALUE); }
Evaluates the continued fraction at the value x.
Params:
  • x – the evaluation point.
  • epsilon – maximum error allowed.
Throws:
Returns:the value of the continued fraction evaluated at x.
/** * Evaluates the continued fraction at the value x. * @param x the evaluation point. * @param epsilon maximum error allowed. * @return the value of the continued fraction evaluated at x. * @throws ConvergenceException if the algorithm fails to converge. */
public double evaluate(double x, double epsilon) throws ConvergenceException { return evaluate(x, epsilon, Integer.MAX_VALUE); }
Evaluates the continued fraction at the value x.
Params:
  • x – the evaluation point.
  • maxIterations – maximum number of convergents
Throws:
Returns:the value of the continued fraction evaluated at x.
/** * Evaluates the continued fraction at the value x. * @param x the evaluation point. * @param maxIterations maximum number of convergents * @return the value of the continued fraction evaluated at x. * @throws ConvergenceException if the algorithm fails to converge. * @throws MaxCountExceededException if maximal number of iterations is reached */
public double evaluate(double x, int maxIterations) throws ConvergenceException, MaxCountExceededException { return evaluate(x, DEFAULT_EPSILON, maxIterations); }
Evaluates the continued fraction at the value x.

The implementation of this method is based on the modified Lentz algorithm as described on page 18 ff. in:

Note: the implementation uses the terms ai and bi as defined in Continued Fraction @ MathWorld.

Params:
  • x – the evaluation point.
  • epsilon – maximum error allowed.
  • maxIterations – maximum number of convergents
Throws:
Returns:the value of the continued fraction evaluated at x.
/** * Evaluates the continued fraction at the value x. * <p> * The implementation of this method is based on the modified Lentz algorithm as described * on page 18 ff. in: * <ul> * <li> * I. J. Thompson, A. R. Barnett. "Coulomb and Bessel Functions of Complex Arguments and Order." * <a target="_blank" href="http://www.fresco.org.uk/papers/Thompson-JCP64p490.pdf"> * http://www.fresco.org.uk/papers/Thompson-JCP64p490.pdf</a> * </li> * </ul> * <b>Note:</b> the implementation uses the terms a<sub>i</sub> and b<sub>i</sub> as defined in * <a href="http://mathworld.wolfram.com/ContinuedFraction.html">Continued Fraction @ MathWorld</a>. * </p> * * @param x the evaluation point. * @param epsilon maximum error allowed. * @param maxIterations maximum number of convergents * @return the value of the continued fraction evaluated at x. * @throws ConvergenceException if the algorithm fails to converge. * @throws MaxCountExceededException if maximal number of iterations is reached */
public double evaluate(double x, double epsilon, int maxIterations) throws ConvergenceException, MaxCountExceededException { final double small = 1e-50; double hPrev = getA(0, x); // use the value of small as epsilon criteria for zero checks if (Precision.equals(hPrev, 0.0, small)) { hPrev = small; } int n = 1; double dPrev = 0.0; double cPrev = hPrev; double hN = hPrev; while (n < maxIterations) { final double a = getA(n, x); final double b = getB(n, x); double dN = a + b * dPrev; if (Precision.equals(dN, 0.0, small)) { dN = small; } double cN = a + b / cPrev; if (Precision.equals(cN, 0.0, small)) { cN = small; } dN = 1 / dN; final double deltaN = cN * dN; hN = hPrev * deltaN; if (Double.isInfinite(hN)) { throw new ConvergenceException(LocalizedFormats.CONTINUED_FRACTION_INFINITY_DIVERGENCE, x); } if (Double.isNaN(hN)) { throw new ConvergenceException(LocalizedFormats.CONTINUED_FRACTION_NAN_DIVERGENCE, x); } if (FastMath.abs(deltaN - 1.0) < epsilon) { break; } dPrev = dN; cPrev = cN; hPrev = hN; n++; } if (n >= maxIterations) { throw new MaxCountExceededException(LocalizedFormats.NON_CONVERGENT_CONTINUED_FRACTION, maxIterations, x); } return hN; } }