1/* 2 * Licensed to the Apache Software Foundation (ASF) under one or more 3 * contributor license agreements. See the NOTICE file distributed with 4 * this work for additional information regarding copyright ownership. 5 * The ASF licenses this file to You under the Apache License, Version 2.0 6 * (the "License"); you may not use this file except in compliance with 7 * the License. You may obtain a copy of the License at 8 * 9 * http://www.apache.org/licenses/LICENSE-2.0 10 * 11 * Unless required by applicable law or agreed to in writing, software 12 * distributed under the License is distributed on an "AS IS" BASIS, 13 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 14 * See the License for the specific language governing permissions and 15 * limitations under the License. 16 */ 17package org.apache.commons.math.analysis.polynomials; 18 19import java.util.ArrayList; 20 21import org.apache.commons.math.fraction.BigFraction; 22import org.apache.commons.math.util.FastMath; 23 24/** 25 * A collection of static methods that operate on or return polynomials. 26 * 27 * @version $Revision: 990655 $ $Date: 2010-08-29 23:49:40 +0200 (dim. 29 août 2010) $ 28 * @since 2.0 29 */ 30public class PolynomialsUtils { 31 32 /** Coefficients for Chebyshev polynomials. */ 33 private static final ArrayList<BigFraction> CHEBYSHEV_COEFFICIENTS; 34 35 /** Coefficients for Hermite polynomials. */ 36 private static final ArrayList<BigFraction> HERMITE_COEFFICIENTS; 37 38 /** Coefficients for Laguerre polynomials. */ 39 private static final ArrayList<BigFraction> LAGUERRE_COEFFICIENTS; 40 41 /** Coefficients for Legendre polynomials. */ 42 private static final ArrayList<BigFraction> LEGENDRE_COEFFICIENTS; 43 44 static { 45 46 // initialize recurrence for Chebyshev polynomials 47 // T0(X) = 1, T1(X) = 0 + 1 * X 48 CHEBYSHEV_COEFFICIENTS = new ArrayList<BigFraction>(); 49 CHEBYSHEV_COEFFICIENTS.add(BigFraction.ONE); 50 CHEBYSHEV_COEFFICIENTS.add(BigFraction.ZERO); 51 CHEBYSHEV_COEFFICIENTS.add(BigFraction.ONE); 52 53 // initialize recurrence for Hermite polynomials 54 // H0(X) = 1, H1(X) = 0 + 2 * X 55 HERMITE_COEFFICIENTS = new ArrayList<BigFraction>(); 56 HERMITE_COEFFICIENTS.add(BigFraction.ONE); 57 HERMITE_COEFFICIENTS.add(BigFraction.ZERO); 58 HERMITE_COEFFICIENTS.add(BigFraction.TWO); 59 60 // initialize recurrence for Laguerre polynomials 61 // L0(X) = 1, L1(X) = 1 - 1 * X 62 LAGUERRE_COEFFICIENTS = new ArrayList<BigFraction>(); 63 LAGUERRE_COEFFICIENTS.add(BigFraction.ONE); 64 LAGUERRE_COEFFICIENTS.add(BigFraction.ONE); 65 LAGUERRE_COEFFICIENTS.add(BigFraction.MINUS_ONE); 66 67 // initialize recurrence for Legendre polynomials 68 // P0(X) = 1, P1(X) = 0 + 1 * X 69 LEGENDRE_COEFFICIENTS = new ArrayList<BigFraction>(); 70 LEGENDRE_COEFFICIENTS.add(BigFraction.ONE); 71 LEGENDRE_COEFFICIENTS.add(BigFraction.ZERO); 72 LEGENDRE_COEFFICIENTS.add(BigFraction.ONE); 73 74 } 75 76 /** 77 * Private constructor, to prevent instantiation. 78 */ 79 private PolynomialsUtils() { 80 } 81 82 /** 83 * Create a Chebyshev polynomial of the first kind. 84 * <p><a href="http://mathworld.wolfram.com/ChebyshevPolynomialoftheFirstKind.html">Chebyshev 85 * polynomials of the first kind</a> are orthogonal polynomials. 86 * They can be defined by the following recurrence relations: 87 * <pre> 88 * T<sub>0</sub>(X) = 1 89 * T<sub>1</sub>(X) = X 90 * T<sub>k+1</sub>(X) = 2X T<sub>k</sub>(X) - T<sub>k-1</sub>(X) 91 * </pre></p> 92 * @param degree degree of the polynomial 93 * @return Chebyshev polynomial of specified degree 94 */ 95 public static PolynomialFunction createChebyshevPolynomial(final int degree) { 96 return buildPolynomial(degree, CHEBYSHEV_COEFFICIENTS, 97 new RecurrenceCoefficientsGenerator() { 98 private final BigFraction[] coeffs = { BigFraction.ZERO, BigFraction.TWO, BigFraction.ONE }; 99 /** {@inheritDoc} */ 100 public BigFraction[] generate(int k) { 101 return coeffs; 102 } 103 }); 104 } 105 106 /** 107 * Create a Hermite polynomial. 108 * <p><a href="http://mathworld.wolfram.com/HermitePolynomial.html">Hermite 109 * polynomials</a> are orthogonal polynomials. 110 * They can be defined by the following recurrence relations: 111 * <pre> 112 * H<sub>0</sub>(X) = 1 113 * H<sub>1</sub>(X) = 2X 114 * H<sub>k+1</sub>(X) = 2X H<sub>k</sub>(X) - 2k H<sub>k-1</sub>(X) 115 * </pre></p> 116 117 * @param degree degree of the polynomial 118 * @return Hermite polynomial of specified degree 119 */ 120 public static PolynomialFunction createHermitePolynomial(final int degree) { 121 return buildPolynomial(degree, HERMITE_COEFFICIENTS, 122 new RecurrenceCoefficientsGenerator() { 123 /** {@inheritDoc} */ 124 public BigFraction[] generate(int k) { 125 return new BigFraction[] { 126 BigFraction.ZERO, 127 BigFraction.TWO, 128 new BigFraction(2 * k)}; 129 } 130 }); 131 } 132 133 /** 134 * Create a Laguerre polynomial. 135 * <p><a href="http://mathworld.wolfram.com/LaguerrePolynomial.html">Laguerre 136 * polynomials</a> are orthogonal polynomials. 137 * They can be defined by the following recurrence relations: 138 * <pre> 139 * L<sub>0</sub>(X) = 1 140 * L<sub>1</sub>(X) = 1 - X 141 * (k+1) L<sub>k+1</sub>(X) = (2k + 1 - X) L<sub>k</sub>(X) - k L<sub>k-1</sub>(X) 142 * </pre></p> 143 * @param degree degree of the polynomial 144 * @return Laguerre polynomial of specified degree 145 */ 146 public static PolynomialFunction createLaguerrePolynomial(final int degree) { 147 return buildPolynomial(degree, LAGUERRE_COEFFICIENTS, 148 new RecurrenceCoefficientsGenerator() { 149 /** {@inheritDoc} */ 150 public BigFraction[] generate(int k) { 151 final int kP1 = k + 1; 152 return new BigFraction[] { 153 new BigFraction(2 * k + 1, kP1), 154 new BigFraction(-1, kP1), 155 new BigFraction(k, kP1)}; 156 } 157 }); 158 } 159 160 /** 161 * Create a Legendre polynomial. 162 * <p><a href="http://mathworld.wolfram.com/LegendrePolynomial.html">Legendre 163 * polynomials</a> are orthogonal polynomials. 164 * They can be defined by the following recurrence relations: 165 * <pre> 166 * P<sub>0</sub>(X) = 1 167 * P<sub>1</sub>(X) = X 168 * (k+1) P<sub>k+1</sub>(X) = (2k+1) X P<sub>k</sub>(X) - k P<sub>k-1</sub>(X) 169 * </pre></p> 170 * @param degree degree of the polynomial 171 * @return Legendre polynomial of specified degree 172 */ 173 public static PolynomialFunction createLegendrePolynomial(final int degree) { 174 return buildPolynomial(degree, LEGENDRE_COEFFICIENTS, 175 new RecurrenceCoefficientsGenerator() { 176 /** {@inheritDoc} */ 177 public BigFraction[] generate(int k) { 178 final int kP1 = k + 1; 179 return new BigFraction[] { 180 BigFraction.ZERO, 181 new BigFraction(k + kP1, kP1), 182 new BigFraction(k, kP1)}; 183 } 184 }); 185 } 186 187 /** Get the coefficients array for a given degree. 188 * @param degree degree of the polynomial 189 * @param coefficients list where the computed coefficients are stored 190 * @param generator recurrence coefficients generator 191 * @return coefficients array 192 */ 193 private static PolynomialFunction buildPolynomial(final int degree, 194 final ArrayList<BigFraction> coefficients, 195 final RecurrenceCoefficientsGenerator generator) { 196 197 final int maxDegree = (int) FastMath.floor(FastMath.sqrt(2 * coefficients.size())) - 1; 198 synchronized (PolynomialsUtils.class) { 199 if (degree > maxDegree) { 200 computeUpToDegree(degree, maxDegree, generator, coefficients); 201 } 202 } 203 204 // coefficient for polynomial 0 is l [0] 205 // coefficients for polynomial 1 are l [1] ... l [2] (degrees 0 ... 1) 206 // coefficients for polynomial 2 are l [3] ... l [5] (degrees 0 ... 2) 207 // coefficients for polynomial 3 are l [6] ... l [9] (degrees 0 ... 3) 208 // coefficients for polynomial 4 are l[10] ... l[14] (degrees 0 ... 4) 209 // coefficients for polynomial 5 are l[15] ... l[20] (degrees 0 ... 5) 210 // coefficients for polynomial 6 are l[21] ... l[27] (degrees 0 ... 6) 211 // ... 212 final int start = degree * (degree + 1) / 2; 213 214 final double[] a = new double[degree + 1]; 215 for (int i = 0; i <= degree; ++i) { 216 a[i] = coefficients.get(start + i).doubleValue(); 217 } 218 219 // build the polynomial 220 return new PolynomialFunction(a); 221 222 } 223 224 /** Compute polynomial coefficients up to a given degree. 225 * @param degree maximal degree 226 * @param maxDegree current maximal degree 227 * @param generator recurrence coefficients generator 228 * @param coefficients list where the computed coefficients should be appended 229 */ 230 private static void computeUpToDegree(final int degree, final int maxDegree, 231 final RecurrenceCoefficientsGenerator generator, 232 final ArrayList<BigFraction> coefficients) { 233 234 int startK = (maxDegree - 1) * maxDegree / 2; 235 for (int k = maxDegree; k < degree; ++k) { 236 237 // start indices of two previous polynomials Pk(X) and Pk-1(X) 238 int startKm1 = startK; 239 startK += k; 240 241 // Pk+1(X) = (a[0] + a[1] X) Pk(X) - a[2] Pk-1(X) 242 BigFraction[] ai = generator.generate(k); 243 244 BigFraction ck = coefficients.get(startK); 245 BigFraction ckm1 = coefficients.get(startKm1); 246 247 // degree 0 coefficient 248 coefficients.add(ck.multiply(ai[0]).subtract(ckm1.multiply(ai[2]))); 249 250 // degree 1 to degree k-1 coefficients 251 for (int i = 1; i < k; ++i) { 252 final BigFraction ckPrev = ck; 253 ck = coefficients.get(startK + i); 254 ckm1 = coefficients.get(startKm1 + i); 255 coefficients.add(ck.multiply(ai[0]).add(ckPrev.multiply(ai[1])).subtract(ckm1.multiply(ai[2]))); 256 } 257 258 // degree k coefficient 259 final BigFraction ckPrev = ck; 260 ck = coefficients.get(startK + k); 261 coefficients.add(ck.multiply(ai[0]).add(ckPrev.multiply(ai[1]))); 262 263 // degree k+1 coefficient 264 coefficients.add(ck.multiply(ai[1])); 265 266 } 267 268 } 269 270 /** Interface for recurrence coefficients generation. */ 271 private static interface RecurrenceCoefficientsGenerator { 272 /** 273 * Generate recurrence coefficients. 274 * @param k highest degree of the polynomials used in the recurrence 275 * @return an array of three coefficients such that 276 * P<sub>k+1</sub>(X) = (a[0] + a[1] X) P<sub>k</sub>(X) - a[2] P<sub>k-1</sub>(X) 277 */ 278 BigFraction[] generate(int k); 279 } 280 281} 282