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.math3.analysis.polynomials;
018
019 import java.util.Arrays;
020
021 import org.apache.commons.math3.util.MathArrays;
022 import org.apache.commons.math3.analysis.DifferentiableUnivariateFunction;
023 import org.apache.commons.math3.analysis.UnivariateFunction;
024 import org.apache.commons.math3.analysis.differentiation.DerivativeStructure;
025 import org.apache.commons.math3.analysis.differentiation.UnivariateDifferentiableFunction;
026 import org.apache.commons.math3.exception.NonMonotonicSequenceException;
027 import org.apache.commons.math3.exception.OutOfRangeException;
028 import org.apache.commons.math3.exception.NumberIsTooSmallException;
029 import org.apache.commons.math3.exception.DimensionMismatchException;
030 import org.apache.commons.math3.exception.NullArgumentException;
031 import org.apache.commons.math3.exception.util.LocalizedFormats;
032
033 /**
034 * Represents a polynomial spline function.
035 * <p>
036 * A <strong>polynomial spline function</strong> consists of a set of
037 * <i>interpolating polynomials</i> and an ascending array of domain
038 * <i>knot points</i>, determining the intervals over which the spline function
039 * is defined by the constituent polynomials. The polynomials are assumed to
040 * have been computed to match the values of another function at the knot
041 * points. The value consistency constraints are not currently enforced by
042 * <code>PolynomialSplineFunction</code> itself, but are assumed to hold among
043 * the polynomials and knot points passed to the constructor.</p>
044 * <p>
045 * N.B.: The polynomials in the <code>polynomials</code> property must be
046 * centered on the knot points to compute the spline function values.
047 * See below.</p>
048 * <p>
049 * The domain of the polynomial spline function is
050 * <code>[smallest knot, largest knot]</code>. Attempts to evaluate the
051 * function at values outside of this range generate IllegalArgumentExceptions.
052 * </p>
053 * <p>
054 * The value of the polynomial spline function for an argument <code>x</code>
055 * is computed as follows:
056 * <ol>
057 * <li>The knot array is searched to find the segment to which <code>x</code>
058 * belongs. If <code>x</code> is less than the smallest knot point or greater
059 * than the largest one, an <code>IllegalArgumentException</code>
060 * is thrown.</li>
061 * <li> Let <code>j</code> be the index of the largest knot point that is less
062 * than or equal to <code>x</code>. The value returned is <br>
063 * <code>polynomials[j](x - knot[j])</code></li></ol></p>
064 *
065 * @version $Id: PolynomialSplineFunction.java 1455194 2013-03-11 15:45:54Z luc $
066 */
067 public class PolynomialSplineFunction implements UnivariateDifferentiableFunction, DifferentiableUnivariateFunction {
068 /**
069 * Spline segment interval delimiters (knots).
070 * Size is n + 1 for n segments.
071 */
072 private final double knots[];
073 /**
074 * The polynomial functions that make up the spline. The first element
075 * determines the value of the spline over the first subinterval, the
076 * second over the second, etc. Spline function values are determined by
077 * evaluating these functions at {@code (x - knot[i])} where i is the
078 * knot segment to which x belongs.
079 */
080 private final PolynomialFunction polynomials[];
081 /**
082 * Number of spline segments. It is equal to the number of polynomials and
083 * to the number of partition points - 1.
084 */
085 private final int n;
086
087
088 /**
089 * Construct a polynomial spline function with the given segment delimiters
090 * and interpolating polynomials.
091 * The constructor copies both arrays and assigns the copies to the knots
092 * and polynomials properties, respectively.
093 *
094 * @param knots Spline segment interval delimiters.
095 * @param polynomials Polynomial functions that make up the spline.
096 * @throws NullArgumentException if either of the input arrays is {@code null}.
097 * @throws NumberIsTooSmallException if knots has length less than 2.
098 * @throws DimensionMismatchException if {@code polynomials.length != knots.length - 1}.
099 * @throws NonMonotonicSequenceException if the {@code knots} array is not strictly increasing.
100 *
101 */
102 public PolynomialSplineFunction(double knots[], PolynomialFunction polynomials[])
103 throws NullArgumentException, NumberIsTooSmallException,
104 DimensionMismatchException, NonMonotonicSequenceException{
105 if (knots == null ||
106 polynomials == null) {
107 throw new NullArgumentException();
108 }
109 if (knots.length < 2) {
110 throw new NumberIsTooSmallException(LocalizedFormats.NOT_ENOUGH_POINTS_IN_SPLINE_PARTITION,
111 2, knots.length, false);
112 }
113 if (knots.length - 1 != polynomials.length) {
114 throw new DimensionMismatchException(polynomials.length, knots.length);
115 }
116 MathArrays.checkOrder(knots);
117
118 this.n = knots.length -1;
119 this.knots = new double[n + 1];
120 System.arraycopy(knots, 0, this.knots, 0, n + 1);
121 this.polynomials = new PolynomialFunction[n];
122 System.arraycopy(polynomials, 0, this.polynomials, 0, n);
123 }
124
125 /**
126 * Compute the value for the function.
127 * See {@link PolynomialSplineFunction} for details on the algorithm for
128 * computing the value of the function.
129 *
130 * @param v Point for which the function value should be computed.
131 * @return the value.
132 * @throws OutOfRangeException if {@code v} is outside of the domain of the
133 * spline function (smaller than the smallest knot point or larger than the
134 * largest knot point).
135 */
136 public double value(double v) {
137 if (v < knots[0] || v > knots[n]) {
138 throw new OutOfRangeException(v, knots[0], knots[n]);
139 }
140 int i = Arrays.binarySearch(knots, v);
141 if (i < 0) {
142 i = -i - 2;
143 }
144 // This will handle the case where v is the last knot value
145 // There are only n-1 polynomials, so if v is the last knot
146 // then we will use the last polynomial to calculate the value.
147 if ( i >= polynomials.length ) {
148 i--;
149 }
150 return polynomials[i].value(v - knots[i]);
151 }
152
153 /**
154 * Get the derivative of the polynomial spline function.
155 *
156 * @return the derivative function.
157 */
158 public UnivariateFunction derivative() {
159 return polynomialSplineDerivative();
160 }
161
162 /**
163 * Get the derivative of the polynomial spline function.
164 *
165 * @return the derivative function.
166 */
167 public PolynomialSplineFunction polynomialSplineDerivative() {
168 PolynomialFunction derivativePolynomials[] = new PolynomialFunction[n];
169 for (int i = 0; i < n; i++) {
170 derivativePolynomials[i] = polynomials[i].polynomialDerivative();
171 }
172 return new PolynomialSplineFunction(knots, derivativePolynomials);
173 }
174
175
176 /** {@inheritDoc}
177 * @since 3.1
178 */
179 public DerivativeStructure value(final DerivativeStructure t) {
180 final double t0 = t.getValue();
181 if (t0 < knots[0] || t0 > knots[n]) {
182 throw new OutOfRangeException(t0, knots[0], knots[n]);
183 }
184 int i = Arrays.binarySearch(knots, t0);
185 if (i < 0) {
186 i = -i - 2;
187 }
188 // This will handle the case where t is the last knot value
189 // There are only n-1 polynomials, so if t is the last knot
190 // then we will use the last polynomial to calculate the value.
191 if ( i >= polynomials.length ) {
192 i--;
193 }
194 return polynomials[i].value(t.subtract(knots[i]));
195 }
196
197 /**
198 * Get the number of spline segments.
199 * It is also the number of polynomials and the number of knot points - 1.
200 *
201 * @return the number of spline segments.
202 */
203 public int getN() {
204 return n;
205 }
206
207 /**
208 * Get a copy of the interpolating polynomials array.
209 * It returns a fresh copy of the array. Changes made to the copy will
210 * not affect the polynomials property.
211 *
212 * @return the interpolating polynomials.
213 */
214 public PolynomialFunction[] getPolynomials() {
215 PolynomialFunction p[] = new PolynomialFunction[n];
216 System.arraycopy(polynomials, 0, p, 0, n);
217 return p;
218 }
219
220 /**
221 * Get an array copy of the knot points.
222 * It returns a fresh copy of the array. Changes made to the copy
223 * will not affect the knots property.
224 *
225 * @return the knot points.
226 */
227 public double[] getKnots() {
228 double out[] = new double[n + 1];
229 System.arraycopy(knots, 0, out, 0, n + 1);
230 return out;
231 }
232 }