001/* 002 * Licensed to the Apache Software Foundation (ASF) under one or more 003 * contributor license agreements. See the NOTICE file distributed with 004 * this work for additional information regarding copyright ownership. 005 * The ASF licenses this file to You under the Apache License, Version 2.0 006 * (the "License"); you may not use this file except in compliance with 007 * the License. You may obtain a copy of the License at 008 * 009 * http://www.apache.org/licenses/LICENSE-2.0 010 * 011 * Unless required by applicable law or agreed to in writing, software 012 * distributed under the License is distributed on an "AS IS" BASIS, 013 * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. 014 * See the License for the specific language governing permissions and 015 * limitations under the License. 016 */ 017package org.apache.commons.math4.legacy.analysis.polynomials; 018 019import java.util.Arrays; 020 021import org.apache.commons.math4.legacy.analysis.ParametricUnivariateFunction; 022import org.apache.commons.math4.legacy.analysis.differentiation.DerivativeStructure; 023import org.apache.commons.math4.legacy.analysis.differentiation.UnivariateDifferentiableFunction; 024import org.apache.commons.math4.legacy.exception.NoDataException; 025import org.apache.commons.math4.legacy.exception.NullArgumentException; 026import org.apache.commons.math4.legacy.exception.util.LocalizedFormats; 027import org.apache.commons.math4.core.jdkmath.JdkMath; 028 029/** 030 * Immutable representation of a real polynomial function with real coefficients. 031 * <p> 032 * <a href="http://mathworld.wolfram.com/HornersMethod.html">Horner's Method</a> 033 * is used to evaluate the function.</p> 034 * 035 */ 036public class PolynomialFunction implements UnivariateDifferentiableFunction { 037 /** 038 * The coefficients of the polynomial, ordered by degree -- i.e., 039 * coefficients[0] is the constant term and coefficients[n] is the 040 * coefficient of x^n where n is the degree of the polynomial. 041 */ 042 private final double[] coefficients; 043 044 /** 045 * Construct a polynomial with the given coefficients. The first element 046 * of the coefficients array is the constant term. Higher degree 047 * coefficients follow in sequence. The degree of the resulting polynomial 048 * is the index of the last non-null element of the array, or 0 if all elements 049 * are null. 050 * <p> 051 * The constructor makes a copy of the input array and assigns the copy to 052 * the coefficients property.</p> 053 * 054 * @param c Polynomial coefficients. 055 * @throws NullArgumentException if {@code c} is {@code null}. 056 * @throws NoDataException if {@code c} is empty. 057 */ 058 public PolynomialFunction(double[] c) 059 throws NullArgumentException, NoDataException { 060 super(); 061 NullArgumentException.check(c); 062 int n = c.length; 063 if (n == 0) { 064 throw new NoDataException(LocalizedFormats.EMPTY_POLYNOMIALS_COEFFICIENTS_ARRAY); 065 } 066 while (n > 1 && c[n - 1] == 0) { 067 --n; 068 } 069 this.coefficients = new double[n]; 070 System.arraycopy(c, 0, this.coefficients, 0, n); 071 } 072 073 /** 074 * Compute the value of the function for the given argument. 075 * <p> 076 * The value returned is </p><p> 077 * {@code coefficients[n] * x^n + ... + coefficients[1] * x + coefficients[0]} 078 * </p> 079 * 080 * @param x Argument for which the function value should be computed. 081 * @return the value of the polynomial at the given point. 082 * 083 * @see org.apache.commons.math4.legacy.analysis.UnivariateFunction#value(double) 084 */ 085 @Override 086 public double value(double x) { 087 return evaluate(coefficients, x); 088 } 089 090 /** 091 * Returns the degree of the polynomial. 092 * 093 * @return the degree of the polynomial. 094 */ 095 public int degree() { 096 return coefficients.length - 1; 097 } 098 099 /** 100 * Returns a copy of the coefficients array. 101 * <p> 102 * Changes made to the returned copy will not affect the coefficients of 103 * the polynomial.</p> 104 * 105 * @return a fresh copy of the coefficients array. 106 */ 107 public double[] getCoefficients() { 108 return coefficients.clone(); 109 } 110 111 /** 112 * Uses Horner's Method to evaluate the polynomial with the given coefficients at 113 * the argument. 114 * 115 * @param coefficients Coefficients of the polynomial to evaluate. 116 * @param argument Input value. 117 * @return the value of the polynomial. 118 * @throws NoDataException if {@code coefficients} is empty. 119 * @throws NullArgumentException if {@code coefficients} is {@code null}. 120 */ 121 protected static double evaluate(double[] coefficients, double argument) 122 throws NullArgumentException, NoDataException { 123 NullArgumentException.check(coefficients); 124 int n = coefficients.length; 125 if (n == 0) { 126 throw new NoDataException(LocalizedFormats.EMPTY_POLYNOMIALS_COEFFICIENTS_ARRAY); 127 } 128 double result = coefficients[n - 1]; 129 for (int j = n - 2; j >= 0; j--) { 130 result = argument * result + coefficients[j]; 131 } 132 return result; 133 } 134 135 136 /** {@inheritDoc} 137 * @since 3.1 138 * @throws NoDataException if {@code coefficients} is empty. 139 * @throws NullArgumentException if {@code coefficients} is {@code null}. 140 */ 141 @Override 142 public DerivativeStructure value(final DerivativeStructure t) 143 throws NullArgumentException, NoDataException { 144 NullArgumentException.check(coefficients); 145 int n = coefficients.length; 146 if (n == 0) { 147 throw new NoDataException(LocalizedFormats.EMPTY_POLYNOMIALS_COEFFICIENTS_ARRAY); 148 } 149 DerivativeStructure result = 150 new DerivativeStructure(t.getFreeParameters(), t.getOrder(), coefficients[n - 1]); 151 for (int j = n - 2; j >= 0; j--) { 152 result = result.multiply(t).add(coefficients[j]); 153 } 154 return result; 155 } 156 157 /** 158 * Add a polynomial to the instance. 159 * 160 * @param p Polynomial to add. 161 * @return a new polynomial which is the sum of the instance and {@code p}. 162 */ 163 public PolynomialFunction add(final PolynomialFunction p) { 164 // identify the lowest degree polynomial 165 final int lowLength = JdkMath.min(coefficients.length, p.coefficients.length); 166 final int highLength = JdkMath.max(coefficients.length, p.coefficients.length); 167 168 // build the coefficients array 169 double[] newCoefficients = new double[highLength]; 170 for (int i = 0; i < lowLength; ++i) { 171 newCoefficients[i] = coefficients[i] + p.coefficients[i]; 172 } 173 System.arraycopy((coefficients.length < p.coefficients.length) ? 174 p.coefficients : coefficients, 175 lowLength, 176 newCoefficients, lowLength, 177 highLength - lowLength); 178 179 return new PolynomialFunction(newCoefficients); 180 } 181 182 /** 183 * Subtract a polynomial from the instance. 184 * 185 * @param p Polynomial to subtract. 186 * @return a new polynomial which is the instance minus {@code p}. 187 */ 188 public PolynomialFunction subtract(final PolynomialFunction p) { 189 // identify the lowest degree polynomial 190 int lowLength = JdkMath.min(coefficients.length, p.coefficients.length); 191 int highLength = JdkMath.max(coefficients.length, p.coefficients.length); 192 193 // build the coefficients array 194 double[] newCoefficients = new double[highLength]; 195 for (int i = 0; i < lowLength; ++i) { 196 newCoefficients[i] = coefficients[i] - p.coefficients[i]; 197 } 198 if (coefficients.length < p.coefficients.length) { 199 for (int i = lowLength; i < highLength; ++i) { 200 newCoefficients[i] = -p.coefficients[i]; 201 } 202 } else { 203 System.arraycopy(coefficients, lowLength, newCoefficients, lowLength, 204 highLength - lowLength); 205 } 206 207 return new PolynomialFunction(newCoefficients); 208 } 209 210 /** 211 * Negate the instance. 212 * 213 * @return a new polynomial with all coefficients negated 214 */ 215 public PolynomialFunction negate() { 216 double[] newCoefficients = new double[coefficients.length]; 217 for (int i = 0; i < coefficients.length; ++i) { 218 newCoefficients[i] = -coefficients[i]; 219 } 220 return new PolynomialFunction(newCoefficients); 221 } 222 223 /** 224 * Multiply the instance by a polynomial. 225 * 226 * @param p Polynomial to multiply by. 227 * @return a new polynomial equal to this times {@code p} 228 */ 229 public PolynomialFunction multiply(final PolynomialFunction p) { 230 double[] newCoefficients = new double[coefficients.length + p.coefficients.length - 1]; 231 232 for (int i = 0; i < newCoefficients.length; ++i) { 233 newCoefficients[i] = 0.0; 234 for (int j = JdkMath.max(0, i + 1 - p.coefficients.length); 235 j < JdkMath.min(coefficients.length, i + 1); 236 ++j) { 237 newCoefficients[i] += coefficients[j] * p.coefficients[i-j]; 238 } 239 } 240 241 return new PolynomialFunction(newCoefficients); 242 } 243 244 /** 245 * Returns the coefficients of the derivative of the polynomial with the given coefficients. 246 * 247 * @param coefficients Coefficients of the polynomial to differentiate. 248 * @return the coefficients of the derivative or {@code null} if coefficients has length 1. 249 * @throws NoDataException if {@code coefficients} is empty. 250 * @throws NullArgumentException if {@code coefficients} is {@code null}. 251 */ 252 protected static double[] differentiate(double[] coefficients) 253 throws NullArgumentException, NoDataException { 254 NullArgumentException.check(coefficients); 255 int n = coefficients.length; 256 if (n == 0) { 257 throw new NoDataException(LocalizedFormats.EMPTY_POLYNOMIALS_COEFFICIENTS_ARRAY); 258 } 259 if (n == 1) { 260 return new double[]{0}; 261 } 262 double[] result = new double[n - 1]; 263 for (int i = n - 1; i > 0; i--) { 264 result[i - 1] = i * coefficients[i]; 265 } 266 return result; 267 } 268 269 /** 270 * Returns the derivative as a {@link PolynomialFunction}. 271 * 272 * @return the derivative polynomial. 273 */ 274 public PolynomialFunction polynomialDerivative() { 275 return new PolynomialFunction(differentiate(coefficients)); 276 } 277 278 /** 279 * Returns a string representation of the polynomial. 280 * 281 * <p>The representation is user oriented. Terms are displayed lowest 282 * degrees first. The multiplications signs, coefficients equals to 283 * one and null terms are not displayed (except if the polynomial is 0, 284 * in which case the 0 constant term is displayed). Addition of terms 285 * with negative coefficients are replaced by subtraction of terms 286 * with positive coefficients except for the first displayed term 287 * (i.e. we display <code>-3</code> for a constant negative polynomial, 288 * but <code>1 - 3 x + x^2</code> if the negative coefficient is not 289 * the first one displayed).</p> 290 * 291 * @return a string representation of the polynomial. 292 */ 293 @Override 294 public String toString() { 295 StringBuilder s = new StringBuilder(); 296 if (coefficients[0] == 0.0) { 297 if (coefficients.length == 1) { 298 return "0"; 299 } 300 } else { 301 s.append(toString(coefficients[0])); 302 } 303 304 for (int i = 1; i < coefficients.length; ++i) { 305 if (coefficients[i] != 0) { 306 if (s.length() > 0) { 307 if (coefficients[i] < 0) { 308 s.append(" - "); 309 } else { 310 s.append(" + "); 311 } 312 } else { 313 if (coefficients[i] < 0) { 314 s.append("-"); 315 } 316 } 317 318 double absAi = JdkMath.abs(coefficients[i]); 319 if ((absAi - 1) != 0) { 320 s.append(toString(absAi)); 321 s.append(' '); 322 } 323 324 s.append("x"); 325 if (i > 1) { 326 s.append('^'); 327 s.append(Integer.toString(i)); 328 } 329 } 330 } 331 332 return s.toString(); 333 } 334 335 /** 336 * Creates a string representing a coefficient, removing ".0" endings. 337 * 338 * @param coeff Coefficient. 339 * @return a string representation of {@code coeff}. 340 */ 341 private static String toString(double coeff) { 342 final String c = Double.toString(coeff); 343 if (c.endsWith(".0")) { 344 return c.substring(0, c.length() - 2); 345 } else { 346 return c; 347 } 348 } 349 350 /** {@inheritDoc} */ 351 @Override 352 public int hashCode() { 353 final int prime = 31; 354 int result = 1; 355 result = prime * result + Arrays.hashCode(coefficients); 356 return result; 357 } 358 359 /** {@inheritDoc} */ 360 @Override 361 public boolean equals(Object obj) { 362 if (this == obj) { 363 return true; 364 } 365 if (!(obj instanceof PolynomialFunction)) { 366 return false; 367 } 368 PolynomialFunction other = (PolynomialFunction) obj; 369 return Arrays.equals(coefficients, other.coefficients); 370 } 371 372 /** 373 * Dedicated parametric polynomial class. 374 * 375 * @since 3.0 376 */ 377 public static class Parametric implements ParametricUnivariateFunction { 378 /** {@inheritDoc} */ 379 @Override 380 public double[] gradient(double x, double ... parameters) { 381 final double[] gradient = new double[parameters.length]; 382 double xn = 1.0; 383 for (int i = 0; i < parameters.length; ++i) { 384 gradient[i] = xn; 385 xn *= x; 386 } 387 return gradient; 388 } 389 390 /** {@inheritDoc} */ 391 @Override 392 public double value(final double x, final double ... parameters) 393 throws NoDataException { 394 return PolynomialFunction.evaluate(parameters, x); 395 } 396 } 397}