View Javadoc
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   */
17  package org.apache.commons.math4.legacy.analysis.polynomials;
18  
19  import java.util.Arrays;
20  
21  import org.apache.commons.math4.legacy.analysis.ParametricUnivariateFunction;
22  import org.apache.commons.math4.legacy.analysis.differentiation.DerivativeStructure;
23  import org.apache.commons.math4.legacy.analysis.differentiation.UnivariateDifferentiableFunction;
24  import org.apache.commons.math4.legacy.exception.NoDataException;
25  import org.apache.commons.math4.legacy.exception.NullArgumentException;
26  import org.apache.commons.math4.legacy.exception.util.LocalizedFormats;
27  import org.apache.commons.math4.core.jdkmath.JdkMath;
28  
29  /**
30   * Immutable representation of a real polynomial function with real coefficients.
31   * <p>
32   * <a href="http://mathworld.wolfram.com/HornersMethod.html">Horner's Method</a>
33   * is used to evaluate the function.</p>
34   *
35   */
36  public class PolynomialFunction implements UnivariateDifferentiableFunction {
37      /**
38       * The coefficients of the polynomial, ordered by degree -- i.e.,
39       * coefficients[0] is the constant term and coefficients[n] is the
40       * coefficient of x^n where n is the degree of the polynomial.
41       */
42      private final double[] coefficients;
43  
44      /**
45       * Construct a polynomial with the given coefficients.  The first element
46       * of the coefficients array is the constant term.  Higher degree
47       * coefficients follow in sequence.  The degree of the resulting polynomial
48       * is the index of the last non-null element of the array, or 0 if all elements
49       * are null.
50       * <p>
51       * The constructor makes a copy of the input array and assigns the copy to
52       * the coefficients property.</p>
53       *
54       * @param c Polynomial coefficients.
55       * @throws NullArgumentException if {@code c} is {@code null}.
56       * @throws NoDataException if {@code c} is empty.
57       */
58      public PolynomialFunction(double[] c)
59          throws NullArgumentException, NoDataException {
60          super();
61          NullArgumentException.check(c);
62          int n = c.length;
63          if (n == 0) {
64              throw new NoDataException(LocalizedFormats.EMPTY_POLYNOMIALS_COEFFICIENTS_ARRAY);
65          }
66          while (n > 1 && c[n - 1] == 0) {
67              --n;
68          }
69          this.coefficients = new double[n];
70          System.arraycopy(c, 0, this.coefficients, 0, n);
71      }
72  
73      /**
74       * Compute the value of the function for the given argument.
75       * <p>
76       *  The value returned is </p><p>
77       *  {@code coefficients[n] * x^n + ... + coefficients[1] * x  + coefficients[0]}
78       * </p>
79       *
80       * @param x Argument for which the function value should be computed.
81       * @return the value of the polynomial at the given point.
82       *
83       * @see org.apache.commons.math4.legacy.analysis.UnivariateFunction#value(double)
84       */
85      @Override
86      public double value(double x) {
87         return evaluate(coefficients, x);
88      }
89  
90      /**
91       * Returns the degree of the polynomial.
92       *
93       * @return the degree of the polynomial.
94       */
95      public int degree() {
96          return coefficients.length - 1;
97      }
98  
99      /**
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 }