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.interpolation;
18  
19  import static org.junit.Assert.assertEquals;
20  import static org.junit.Assert.assertTrue;
21  
22  import org.apache.commons.math4.legacy.analysis.UnivariateFunction;
23  import org.apache.commons.statistics.distribution.ContinuousDistribution;
24  import org.apache.commons.statistics.distribution.UniformContinuousDistribution;
25  import org.apache.commons.math4.legacy.exception.DimensionMismatchException;
26  import org.apache.commons.math4.legacy.exception.NonMonotonicSequenceException;
27  import org.apache.commons.math4.legacy.exception.NullArgumentException;
28  import org.apache.commons.math4.legacy.exception.NumberIsTooSmallException;
29  import org.apache.commons.math4.legacy.analysis.polynomials.PolynomialSplineFunction;
30  import org.apache.commons.rng.UniformRandomProvider;
31  import org.apache.commons.rng.simple.RandomSource;
32  import org.apache.commons.math4.core.jdkmath.JdkMath;
33  import org.apache.commons.numbers.core.Precision;
34  import org.junit.Assert;
35  import org.junit.Test;
36  import org.junit.Ignore;
37  
38  public class AkimaSplineInterpolatorTest {
39      @Test
40      public void testIllegalArguments() {
41          // Data set arrays of different size.
42          UnivariateInterpolator i = new AkimaSplineInterpolator();
43  
44          try {
45              double yval[] = {0.0, 1.0, 2.0, 3.0, 4.0};
46              i.interpolate(null, yval);
47              Assert.fail("Failed to detect x null pointer");
48          } catch (NullArgumentException iae) {
49              // Expected.
50          }
51  
52          try {
53              double xval[] = {0.0, 1.0, 2.0, 3.0, 4.0};
54              i.interpolate(xval, null);
55              Assert.fail("Failed to detect y null pointer");
56          } catch (NullArgumentException iae) {
57              // Expected.
58          }
59  
60          try {
61              double xval[] = {0.0, 1.0, 2.0, 3.0};
62              double yval[] = {0.0, 1.0, 2.0, 3.0};
63              i.interpolate(xval, yval);
64              Assert.fail("Failed to detect insufficient data");
65          } catch (NumberIsTooSmallException iae) {
66              // Expected.
67          }
68  
69          try {
70              double xval[] = {0.0, 1.0, 2.0, 3.0, 4.0};
71              double yval[] = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0};
72              i.interpolate(xval, yval);
73              Assert.fail("Failed to detect data set array with different sizes.");
74          } catch (DimensionMismatchException iae) {
75              // Expected.
76          }
77  
78          // X values not sorted.
79          try {
80              double xval[] = {0.0, 1.0, 0.5, 7.0, 3.5, 2.2, 8.0};
81              double yval[] = {0.0, 1.0, 2.0, 3.0, 4.0, 5.0, 6.0};
82              i.interpolate(xval, yval);
83              Assert.fail("Failed to detect unsorted arguments.");
84          } catch (NonMonotonicSequenceException iae) {
85              // Expected.
86          }
87      }
88  
89      /*
90       * Interpolate a straight line. <p> y = 2 x - 5 <p> Tolerances determined by performing same calculation using
91       * Math.NET over ten runs of 100 random number draws for the same function over the same span with the same number
92       * of elements
93       */
94      @Test
95      public void testInterpolateLine() {
96          final int numberOfElements = 10;
97          final double minimumX = -10;
98          final double maximumX = 10;
99          final int numberOfSamples = 100;
100         final double interpolationTolerance = 1e-15;
101         final double maxTolerance = 1e-15;
102 
103         UnivariateFunction f = new UnivariateFunction() {
104             @Override
105             public double value(double x) {
106                 return 2 * x - 5;
107             }
108         };
109 
110         testInterpolation( minimumX, maximumX, numberOfElements, numberOfSamples, f, interpolationTolerance,
111                            maxTolerance );
112     }
113 
114     /*
115      * Interpolate a straight line. <p> y = 3 x<sup>2</sup> - 5 x + 7 <p> Tolerances determined by performing same
116      * calculation using Math.NET over ten runs of 100 random number draws for the same function over the same span with
117      * the same number of elements
118      */
119 
120     @Test
121     public void testInterpolateParabola() {
122         final int numberOfElements = 10;
123         final double minimumX = -10;
124         final double maximumX = 10;
125         final int numberOfSamples = 100;
126         final double interpolationTolerance = 7e-15;
127         final double maxTolerance = 6e-14;
128 
129         UnivariateFunction f = new UnivariateFunction() {
130             @Override
131             public double value(double x) {
132                 return (3 * x * x) - (5 * x) + 7;
133             }
134         };
135 
136         testInterpolation( minimumX, maximumX, numberOfElements, numberOfSamples, f, interpolationTolerance,
137                            maxTolerance );
138     }
139 
140     /*
141      * Interpolate a straight line. <p> y = 3 x<sup>3</sup> - 0.5 x<sup>2</sup> + x - 1 <p> Tolerances determined by
142      * performing same calculation using Math.NET over ten runs of 100 random number draws for the same function over
143      * the same span with the same number of elements
144      */
145     @Test
146     public void testInterpolateCubic() {
147         final int numberOfElements = 10;
148         final double minimumX = -3;
149         final double maximumX = 3;
150         final int numberOfSamples = 100;
151         final double interpolationTolerance = 0.37;
152         final double maxTolerance = 3.8;
153 
154         UnivariateFunction f = new UnivariateFunction() {
155             @Override
156             public double value(double x) {
157                 return (3 * x * x * x) - (0.5 * x * x) + (1 * x) - 1;
158             }
159         };
160 
161         testInterpolation( minimumX, maximumX, numberOfElements, numberOfSamples, f, interpolationTolerance,
162                            maxTolerance );
163     }
164 
165     // Test currently fails but it is not clear whether
166     //   https://issues.apache.org/jira/browse/MATH-1635
167     // actually describes a bug, or a limitation of the algorithm.
168     @Ignore
169     @Test
170     public void testMath1635() {
171         final double[] x = {
172             5994, 6005, 6555, 6588, 6663,
173             6760, 6770, 6792, 6856, 6964,
174             7028, 7233, 7426, 7469, 7619,
175             7910, 8038, 8178, 8414, 8747,
176             8983, 9316, 9864, 9875
177         };
178 
179         final double[] y = {
180             3.0, 2.0, 2.0, 2.0, 2.0,
181             2.0, 2.0, 2.0, 2.0, 2.0,
182             2.0, 2.0, 2.0, 2.0, 2.0,
183             2.0, 2.0, 2.0, 2.0, 2.0,
184             2.0, 2.0, 2.0, 3.0
185         };
186 
187         final AkimaSplineInterpolator interpolator = new AkimaSplineInterpolator(true);
188         final PolynomialSplineFunction interpolate = interpolator.interpolate(x, y);
189         final double value = interpolate.value(9584);
190         final double expected = 2;
191         Assert.assertEquals(expected, value, 1e-4);
192     }
193 
194     @Test
195     public void testOriginalVsModified() {
196         final UnivariateFunction f = new UnivariateFunction() {
197             @Override
198             public double value(double x) {
199                 return x < -1 ? -1 :
200                     x < 1 ? x : 1;
201             }
202         };
203 
204         final double[] xS = new double[] {-1, 0, 1, 2, 3 };
205         final double[] yS = new double[xS.length];
206 
207         for (int i = 0; i < xS.length; i++) {
208             yS[i] = f.value(xS[i]);
209         }
210 
211         final UnivariateFunction iOriginal = new AkimaSplineInterpolator(false).interpolate(xS, yS);
212         final UnivariateFunction iModified = new AkimaSplineInterpolator(true).interpolate(xS, yS);
213 
214         final int n = 100;
215         final double delta = 1d / n;
216         for (int i = 1; i < n - 1; i++) {
217             final double x = 2 - i * delta;
218 
219             final double value = f.value(x);
220             final double diffOriginal = Math.abs(iOriginal.value(x) - value);
221             final double diffModified = Math.abs(iModified.value(x) - value);
222 
223             // In interval (1, 2), the modified algorithm eliminates interpolation artefacts.
224             Assert.assertTrue(diffOriginal > 0);
225             Assert.assertEquals(0d, diffModified, 0d);
226         }
227     }
228 
229     private void testInterpolation( double minimumX, double maximumX, int numberOfElements, int numberOfSamples,
230                                     UnivariateFunction f, double tolerance, double maxTolerance ) {
231         double expected;
232         double actual;
233         double currentX;
234         final double delta = ( maximumX - minimumX ) / ( (double) numberOfElements );
235         double xValues[] = new double[numberOfElements];
236         double yValues[] = new double[numberOfElements];
237 
238         for (int i = 0; i < numberOfElements; i++) {
239             xValues[i] = minimumX + delta * (double) i;
240             yValues[i] = f.value(xValues[i]);
241         }
242 
243         UnivariateFunction interpolation = new AkimaSplineInterpolator().interpolate( xValues, yValues );
244 
245         for (int i = 0; i < numberOfElements; i++) {
246             currentX = xValues[i];
247             expected = f.value(currentX);
248             actual = interpolation.value( currentX );
249             assertTrue( Precision.equals( expected, actual ) );
250         }
251 
252         final UniformRandomProvider rng = RandomSource.WELL_19937_C.create(1234567L); // "tol" depends on the seed.
253         final ContinuousDistribution.Sampler distX =
254             UniformContinuousDistribution.of(xValues[0], xValues[xValues.length - 1]).createSampler(rng);
255 
256         double sumError = 0;
257         for (int i = 0; i < numberOfSamples; i++) {
258             currentX = distX.sample();
259             expected = f.value(currentX);
260             actual = interpolation.value( currentX );
261             sumError += JdkMath.abs( actual - expected );
262             assertEquals( expected, actual, maxTolerance );
263         }
264 
265         assertEquals( 0.0, sumError / (double) numberOfSamples, tolerance );
266     }
267 }