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.MathRuntimeException;
023    import org.apache.commons.math.analysis.DifferentiableUnivariateRealFunction;
024    import org.apache.commons.math.analysis.UnivariateRealFunction;
025    
026    /**
027     * Immutable representation of a real polynomial function with real coefficients.
028     * <p>
029     * <a href="http://mathworld.wolfram.com/HornersMethod.html">Horner's Method</a>
030     *  is used to evaluate the function.</p>
031     *
032     * @version $Revision: 922714 $ $Date: 2010-03-13 20:35:14 -0500 (Sat, 13 Mar 2010) $
033     */
034    public class PolynomialFunction implements DifferentiableUnivariateRealFunction, Serializable {
035    
036        /** Message for empty coefficients array. */
037        private static final String EMPTY_ARRAY_MESSAGE =
038            "empty polynomials coefficients array";
039    
040        /**
041         * Serialization identifier
042         */
043        private static final long serialVersionUID = -7726511984200295583L;
044    
045        /**
046         * The coefficients of the polynomial, ordered by degree -- i.e.,
047         * coefficients[0] is the constant term and coefficients[n] is the
048         * coefficient of x^n where n is the degree of the polynomial.
049         */
050        private final double coefficients[];
051    
052        /**
053         * Construct a polynomial with the given coefficients.  The first element
054         * of the coefficients array is the constant term.  Higher degree
055         * coefficients follow in sequence.  The degree of the resulting polynomial
056         * is the index of the last non-null element of the array, or 0 if all elements
057         * are null.
058         * <p>
059         * The constructor makes a copy of the input array and assigns the copy to
060         * the coefficients property.</p>
061         *
062         * @param c polynomial coefficients
063         * @throws NullPointerException if c is null
064         * @throws IllegalArgumentException if c is empty
065         */
066        public PolynomialFunction(double c[]) {
067            super();
068            if (c.length < 1) {
069                throw MathRuntimeException.createIllegalArgumentException(EMPTY_ARRAY_MESSAGE);
070            }
071            int l = c.length;
072            while ((l > 1) && (c[l - 1] == 0)) {
073                --l;
074            }
075            this.coefficients = new double[l];
076            System.arraycopy(c, 0, this.coefficients, 0, l);
077        }
078    
079        /**
080         * Compute the value of the function for the given argument.
081         * <p>
082         *  The value returned is <br>
083         *   <code>coefficients[n] * x^n + ... + coefficients[1] * x  + coefficients[0]</code>
084         * </p>
085         *
086         * @param x the argument for which the function value should be computed
087         * @return the value of the polynomial at the given point
088         * @see UnivariateRealFunction#value(double)
089         */
090        public double value(double x) {
091           return evaluate(coefficients, x);
092        }
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  the coefficients of the polynomial to evaluate
121         * @param argument  the input value
122         * @return  the value of the polynomial
123         * @throws IllegalArgumentException if coefficients is empty
124         * @throws NullPointerException if coefficients is null
125         */
126        protected static double evaluate(double[] coefficients, double argument) {
127            int n = coefficients.length;
128            if (n < 1) {
129                throw MathRuntimeException.createIllegalArgumentException(EMPTY_ARRAY_MESSAGE);
130            }
131            double result = coefficients[n - 1];
132            for (int j = n -2; j >=0; j--) {
133                result = argument * result + coefficients[j];
134            }
135            return result;
136        }
137    
138        /**
139         * Add a polynomial to the instance.
140         * @param p polynomial to add
141         * @return a new polynomial which is the sum of the instance and p
142         */
143        public PolynomialFunction add(final PolynomialFunction p) {
144    
145            // identify the lowest degree polynomial
146            final int lowLength  = Math.min(coefficients.length, p.coefficients.length);
147            final int highLength = Math.max(coefficients.length, p.coefficients.length);
148    
149            // build the coefficients array
150            double[] newCoefficients = new double[highLength];
151            for (int i = 0; i < lowLength; ++i) {
152                newCoefficients[i] = coefficients[i] + p.coefficients[i];
153            }
154            System.arraycopy((coefficients.length < p.coefficients.length) ?
155                             p.coefficients : coefficients,
156                             lowLength,
157                             newCoefficients, lowLength,
158                             highLength - lowLength);
159    
160            return new PolynomialFunction(newCoefficients);
161    
162        }
163    
164        /**
165         * Subtract a polynomial from the instance.
166         * @param p polynomial to subtract
167         * @return a new polynomial which is the difference the instance minus p
168         */
169        public PolynomialFunction subtract(final PolynomialFunction p) {
170    
171            // identify the lowest degree polynomial
172            int lowLength  = Math.min(coefficients.length, p.coefficients.length);
173            int highLength = Math.max(coefficients.length, p.coefficients.length);
174    
175            // build the coefficients array
176            double[] newCoefficients = new double[highLength];
177            for (int i = 0; i < lowLength; ++i) {
178                newCoefficients[i] = coefficients[i] - p.coefficients[i];
179            }
180            if (coefficients.length < p.coefficients.length) {
181                for (int i = lowLength; i < highLength; ++i) {
182                    newCoefficients[i] = -p.coefficients[i];
183                }
184            } else {
185                System.arraycopy(coefficients, lowLength, newCoefficients, lowLength,
186                                 highLength - lowLength);
187            }
188    
189            return new PolynomialFunction(newCoefficients);
190    
191        }
192    
193        /**
194         * Negate the instance.
195         * @return a new polynomial
196         */
197        public PolynomialFunction negate() {
198            double[] newCoefficients = new double[coefficients.length];
199            for (int i = 0; i < coefficients.length; ++i) {
200                newCoefficients[i] = -coefficients[i];
201            }
202            return new PolynomialFunction(newCoefficients);
203        }
204    
205        /**
206         * Multiply the instance by a polynomial.
207         * @param p polynomial to multiply by
208         * @return a new polynomial
209         */
210        public PolynomialFunction multiply(final PolynomialFunction p) {
211    
212            double[] newCoefficients = new double[coefficients.length + p.coefficients.length - 1];
213    
214            for (int i = 0; i < newCoefficients.length; ++i) {
215                newCoefficients[i] = 0.0;
216                for (int j = Math.max(0, i + 1 - p.coefficients.length);
217                     j < Math.min(coefficients.length, i + 1);
218                     ++j) {
219                    newCoefficients[i] += coefficients[j] * p.coefficients[i-j];
220                }
221            }
222    
223            return new PolynomialFunction(newCoefficients);
224    
225        }
226    
227        /**
228         * Returns the coefficients of the derivative of the polynomial with the given coefficients.
229         *
230         * @param coefficients  the coefficients of the polynomial to differentiate
231         * @return the coefficients of the derivative or null if coefficients has length 1.
232         * @throws IllegalArgumentException if coefficients is empty
233         * @throws NullPointerException if coefficients is null
234         */
235        protected static double[] differentiate(double[] coefficients) {
236            int n = coefficients.length;
237            if (n < 1) {
238                throw MathRuntimeException.createIllegalArgumentException(EMPTY_ARRAY_MESSAGE);
239            }
240            if (n == 1) {
241                return new double[]{0};
242            }
243            double[] result = new double[n - 1];
244            for (int i = n - 1; i  > 0; i--) {
245                result[i - 1] = i * coefficients[i];
246            }
247            return result;
248        }
249    
250        /**
251         * Returns the derivative as a PolynomialRealFunction
252         *
253         * @return  the derivative polynomial
254         */
255        public PolynomialFunction polynomialDerivative() {
256            return new PolynomialFunction(differentiate(coefficients));
257        }
258    
259        /**
260         * Returns the derivative as a UnivariateRealFunction
261         *
262         * @return  the derivative function
263         */
264        public UnivariateRealFunction derivative() {
265            return polynomialDerivative();
266        }
267    
268        /** Returns a string representation of the polynomial.
269    
270         * <p>The representation is user oriented. Terms are displayed lowest
271         * degrees first. The multiplications signs, coefficients equals to
272         * one and null terms are not displayed (except if the polynomial is 0,
273         * in which case the 0 constant term is displayed). Addition of terms
274         * with negative coefficients are replaced by subtraction of terms
275         * with positive coefficients except for the first displayed term
276         * (i.e. we display <code>-3</code> for a constant negative polynomial,
277         * but <code>1 - 3 x + x^2</code> if the negative coefficient is not
278         * the first one displayed).</p>
279    
280         * @return a string representation of the polynomial
281    
282         */
283        @Override
284         public String toString() {
285    
286           StringBuffer s = new StringBuffer();
287           if (coefficients[0] == 0.0) {
288             if (coefficients.length == 1) {
289               return "0";
290             }
291           } else {
292             s.append(Double.toString(coefficients[0]));
293           }
294    
295           for (int i = 1; i < coefficients.length; ++i) {
296    
297             if (coefficients[i] != 0) {
298    
299               if (s.length() > 0) {
300                 if (coefficients[i] < 0) {
301                   s.append(" - ");
302                 } else {
303                   s.append(" + ");
304                 }
305               } else {
306                 if (coefficients[i] < 0) {
307                   s.append("-");
308                 }
309               }
310    
311               double absAi = Math.abs(coefficients[i]);
312               if ((absAi - 1) != 0) {
313                 s.append(Double.toString(absAi));
314                 s.append(' ');
315               }
316    
317               s.append("x");
318               if (i > 1) {
319                 s.append('^');
320                 s.append(Integer.toString(i));
321               }
322             }
323    
324           }
325    
326           return s.toString();
327    
328         }
329    
330        /** {@inheritDoc} */
331        @Override
332        public int hashCode() {
333            final int prime = 31;
334            int result = 1;
335            result = prime * result + Arrays.hashCode(coefficients);
336            return result;
337        }
338    
339        /** {@inheritDoc} */
340        @Override
341        public boolean equals(Object obj) {
342            if (this == obj)
343                return true;
344            if (!(obj instanceof PolynomialFunction))
345                return false;
346            PolynomialFunction other = (PolynomialFunction) obj;
347            if (!Arrays.equals(coefficients, other.coefficients))
348                return false;
349            return true;
350        }
351    
352    }