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