/*
* 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.math.analysis.polynomials;
import java.util.ArrayList;
import org.apache.commons.math.fraction.BigFraction;
import org.apache.commons.math.util.FastMath;
/**
* A collection of static methods that operate on or return polynomials.
*
* @version $Revision: 990655 $ $Date: 2010-08-29 23:49:40 +0200 (dim. 29 août 2010) $
* @since 2.0
*/
public class PolynomialsUtils {
/** Coefficients for Chebyshev polynomials. */
private static final ArrayList<BigFraction> CHEBYSHEV_COEFFICIENTS;
/** Coefficients for Hermite polynomials. */
private static final ArrayList<BigFraction> HERMITE_COEFFICIENTS;
/** Coefficients for Laguerre polynomials. */
private static final ArrayList<BigFraction> LAGUERRE_COEFFICIENTS;
/** Coefficients for Legendre polynomials. */
private static final ArrayList<BigFraction> LEGENDRE_COEFFICIENTS;
static {
// initialize recurrence for Chebyshev polynomials
// T0(X) = 1, T1(X) = 0 + 1 * X
CHEBYSHEV_COEFFICIENTS = new ArrayList<BigFraction>();
CHEBYSHEV_COEFFICIENTS.add(BigFraction.ONE);
CHEBYSHEV_COEFFICIENTS.add(BigFraction.ZERO);
CHEBYSHEV_COEFFICIENTS.add(BigFraction.ONE);
// initialize recurrence for Hermite polynomials
// H0(X) = 1, H1(X) = 0 + 2 * X
HERMITE_COEFFICIENTS = new ArrayList<BigFraction>();
HERMITE_COEFFICIENTS.add(BigFraction.ONE);
HERMITE_COEFFICIENTS.add(BigFraction.ZERO);
HERMITE_COEFFICIENTS.add(BigFraction.TWO);
// initialize recurrence for Laguerre polynomials
// L0(X) = 1, L1(X) = 1 - 1 * X
LAGUERRE_COEFFICIENTS = new ArrayList<BigFraction>();
LAGUERRE_COEFFICIENTS.add(BigFraction.ONE);
LAGUERRE_COEFFICIENTS.add(BigFraction.ONE);
LAGUERRE_COEFFICIENTS.add(BigFraction.MINUS_ONE);
// initialize recurrence for Legendre polynomials
// P0(X) = 1, P1(X) = 0 + 1 * X
LEGENDRE_COEFFICIENTS = new ArrayList<BigFraction>();
LEGENDRE_COEFFICIENTS.add(BigFraction.ONE);
LEGENDRE_COEFFICIENTS.add(BigFraction.ZERO);
LEGENDRE_COEFFICIENTS.add(BigFraction.ONE);
}
/**
* Private constructor, to prevent instantiation.
*/
private PolynomialsUtils() {
}
/**
* Create a Chebyshev polynomial of the first kind.
* <p><a href="http://mathworld.wolfram.com/ChebyshevPolynomialoftheFirstKind.html">Chebyshev
* polynomials of the first kind</a> are orthogonal polynomials.
* They can be defined by the following recurrence relations:
* <pre>
* T<sub>0</sub>(X) = 1
* T<sub>1</sub>(X) = X
* T<sub>k+1</sub>(X) = 2X T<sub>k</sub>(X) - T<sub>k-1</sub>(X)
* </pre></p>
* @param degree degree of the polynomial
* @return Chebyshev polynomial of specified degree
*/
public static PolynomialFunction createChebyshevPolynomial(final int degree) {
return buildPolynomial(degree, CHEBYSHEV_COEFFICIENTS,
new RecurrenceCoefficientsGenerator() {
private final BigFraction[] coeffs = { BigFraction.ZERO, BigFraction.TWO, BigFraction.ONE };
/** {@inheritDoc} */
public BigFraction[] generate(int k) {
return coeffs;
}
});
}
/**
* Create a Hermite polynomial.
* <p><a href="http://mathworld.wolfram.com/HermitePolynomial.html">Hermite
* polynomials</a> are orthogonal polynomials.
* They can be defined by the following recurrence relations:
* <pre>
* H<sub>0</sub>(X) = 1
* H<sub>1</sub>(X) = 2X
* H<sub>k+1</sub>(X) = 2X H<sub>k</sub>(X) - 2k H<sub>k-1</sub>(X)
* </pre></p>
* @param degree degree of the polynomial
* @return Hermite polynomial of specified degree
*/
public static PolynomialFunction createHermitePolynomial(final int degree) {
return buildPolynomial(degree, HERMITE_COEFFICIENTS,
new RecurrenceCoefficientsGenerator() {
/** {@inheritDoc} */
public BigFraction[] generate(int k) {
return new BigFraction[] {
BigFraction.ZERO,
BigFraction.TWO,
new BigFraction(2 * k)};
}
});
}
/**
* Create a Laguerre polynomial.
* <p><a href="http://mathworld.wolfram.com/LaguerrePolynomial.html">Laguerre
* polynomials</a> are orthogonal polynomials.
* They can be defined by the following recurrence relations:
* <pre>
* L<sub>0</sub>(X) = 1
* L<sub>1</sub>(X) = 1 - X
* (k+1) L<sub>k+1</sub>(X) = (2k + 1 - X) L<sub>k</sub>(X) - k L<sub>k-1</sub>(X)
* </pre></p>
* @param degree degree of the polynomial
* @return Laguerre polynomial of specified degree
*/
public static PolynomialFunction createLaguerrePolynomial(final int degree) {
return buildPolynomial(degree, LAGUERRE_COEFFICIENTS,
new RecurrenceCoefficientsGenerator() {
/** {@inheritDoc} */
public BigFraction[] generate(int k) {
final int kP1 = k + 1;
return new BigFraction[] {
new BigFraction(2 * k + 1, kP1),
new BigFraction(-1, kP1),
new BigFraction(k, kP1)};
}
});
}
/**
* Create a Legendre polynomial.
* <p><a href="http://mathworld.wolfram.com/LegendrePolynomial.html">Legendre
* polynomials</a> are orthogonal polynomials.
* They can be defined by the following recurrence relations:
* <pre>
* P<sub>0</sub>(X) = 1
* P<sub>1</sub>(X) = X
* (k+1) P<sub>k+1</sub>(X) = (2k+1) X P<sub>k</sub>(X) - k P<sub>k-1</sub>(X)
* </pre></p>
* @param degree degree of the polynomial
* @return Legendre polynomial of specified degree
*/
public static PolynomialFunction createLegendrePolynomial(final int degree) {
return buildPolynomial(degree, LEGENDRE_COEFFICIENTS,
new RecurrenceCoefficientsGenerator() {
/** {@inheritDoc} */
public BigFraction[] generate(int k) {
final int kP1 = k + 1;
return new BigFraction[] {
BigFraction.ZERO,
new BigFraction(k + kP1, kP1),
new BigFraction(k, kP1)};
}
});
}
/** Get the coefficients array for a given degree.
* @param degree degree of the polynomial
* @param coefficients list where the computed coefficients are stored
* @param generator recurrence coefficients generator
* @return coefficients array
*/
private static PolynomialFunction buildPolynomial(final int degree,
final ArrayList<BigFraction> coefficients,
final RecurrenceCoefficientsGenerator generator) {
final int maxDegree = (int) FastMath.floor(FastMath.sqrt(2 * coefficients.size())) - 1;
synchronized (PolynomialsUtils.class) {
if (degree > maxDegree) {
computeUpToDegree(degree, maxDegree, generator, coefficients);
}
}
// coefficient for polynomial 0 is l [0]
// coefficients for polynomial 1 are l [1] ... l [2] (degrees 0 ... 1)
// coefficients for polynomial 2 are l [3] ... l [5] (degrees 0 ... 2)
// coefficients for polynomial 3 are l [6] ... l [9] (degrees 0 ... 3)
// coefficients for polynomial 4 are l[10] ... l[14] (degrees 0 ... 4)
// coefficients for polynomial 5 are l[15] ... l[20] (degrees 0 ... 5)
// coefficients for polynomial 6 are l[21] ... l[27] (degrees 0 ... 6)
// ...
final int start = degree * (degree + 1) / 2;
final double[] a = new double[degree + 1];
for (int i = 0; i <= degree; ++i) {
a[i] = coefficients.get(start + i).doubleValue();
}
// build the polynomial
return new PolynomialFunction(a);
}
/** Compute polynomial coefficients up to a given degree.
* @param degree maximal degree
* @param maxDegree current maximal degree
* @param generator recurrence coefficients generator
* @param coefficients list where the computed coefficients should be appended
*/
private static void computeUpToDegree(final int degree, final int maxDegree,
final RecurrenceCoefficientsGenerator generator,
final ArrayList<BigFraction> coefficients) {
int startK = (maxDegree - 1) * maxDegree / 2;
for (int k = maxDegree; k < degree; ++k) {
// start indices of two previous polynomials Pk(X) and Pk-1(X)
int startKm1 = startK;
startK += k;
// Pk+1(X) = (a[0] + a[1] X) Pk(X) - a[2] Pk-1(X)
BigFraction[] ai = generator.generate(k);
BigFraction ck = coefficients.get(startK);
BigFraction ckm1 = coefficients.get(startKm1);
// degree 0 coefficient
coefficients.add(ck.multiply(ai[0]).subtract(ckm1.multiply(ai[2])));
// degree 1 to degree k-1 coefficients
for (int i = 1; i < k; ++i) {
final BigFraction ckPrev = ck;
ck = coefficients.get(startK + i);
ckm1 = coefficients.get(startKm1 + i);
coefficients.add(ck.multiply(ai[0]).add(ckPrev.multiply(ai[1])).subtract(ckm1.multiply(ai[2])));
}
// degree k coefficient
final BigFraction ckPrev = ck;
ck = coefficients.get(startK + k);
coefficients.add(ck.multiply(ai[0]).add(ckPrev.multiply(ai[1])));
// degree k+1 coefficient
coefficients.add(ck.multiply(ai[1]));
}
}
/** Interface for recurrence coefficients generation. */
private static interface RecurrenceCoefficientsGenerator {
/**
* Generate recurrence coefficients.
* @param k highest degree of the polynomials used in the recurrence
* @return an array of three coefficients such that
* P<sub>k+1</sub>(X) = (a[0] + a[1] X) P<sub>k</sub>(X) - a[2] P<sub>k-1</sub>(X)
*/
BigFraction[] generate(int k);
}
}