From 0a5cd11327d50e5906fb4dc08bce5baea6b2d247 Mon Sep 17 00:00:00 2001 From: Thomas Neidhart Date: Wed, 25 Feb 2015 23:02:30 +0100 Subject: [PATCH] Remove deprecated interpolation and fitter classes. --- .../BicubicSplineInterpolatingFunction.java | 638 ----------------- .../BicubicSplineInterpolator.java | 176 ----- ...ngPolynomialBicubicSplineInterpolator.java | 171 ----- .../TricubicSplineInterpolatingFunction.java | 482 ------------- .../TricubicSplineInterpolator.java | 201 ------ .../math4/analysis/solvers/NewtonSolver.java | 92 --- .../commons/math4/fitting/CurveFitter.java | 233 ------ .../commons/math4/fitting/GaussianFitter.java | 365 ---------- .../commons/math4/fitting/HarmonicFitter.java | 384 ---------- .../math4/fitting/PolynomialFitter.java | 72 -- ...icubicSplineInterpolatingFunctionTest.java | 670 ------------------ .../BicubicSplineInterpolatorTest.java | 186 ----- ...lynomialBicubicSplineInterpolatorTest.java | 181 ----- ...icubicSplineInterpolatingFunctionTest.java | 545 -------------- .../TricubicSplineInterpolatorTest.java | 214 ------ .../analysis/solvers/NewtonSolverTest.java | 111 --- .../math4/fitting/CurveFitterTest.java | 143 ---- .../math4/fitting/GaussianFitterTest.java | 364 ---------- .../math4/fitting/HarmonicFitterTest.java | 187 ----- .../math4/fitting/PolynomialFitterTest.java | 288 -------- 20 files changed, 5703 deletions(-) delete mode 100644 src/main/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatingFunction.java delete mode 100644 src/main/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolator.java delete mode 100644 src/main/java/org/apache/commons/math4/analysis/interpolation/SmoothingPolynomialBicubicSplineInterpolator.java delete mode 100644 src/main/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatingFunction.java delete mode 100644 src/main/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolator.java delete mode 100644 src/main/java/org/apache/commons/math4/analysis/solvers/NewtonSolver.java delete mode 100644 src/main/java/org/apache/commons/math4/fitting/CurveFitter.java delete mode 100644 src/main/java/org/apache/commons/math4/fitting/GaussianFitter.java delete mode 100644 src/main/java/org/apache/commons/math4/fitting/HarmonicFitter.java delete mode 100644 src/main/java/org/apache/commons/math4/fitting/PolynomialFitter.java delete mode 100644 src/test/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatingFunctionTest.java delete mode 100644 src/test/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatorTest.java delete mode 100644 src/test/java/org/apache/commons/math4/analysis/interpolation/SmoothingPolynomialBicubicSplineInterpolatorTest.java delete mode 100644 src/test/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatingFunctionTest.java delete mode 100644 src/test/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatorTest.java delete mode 100644 src/test/java/org/apache/commons/math4/analysis/solvers/NewtonSolverTest.java delete mode 100644 src/test/java/org/apache/commons/math4/fitting/CurveFitterTest.java delete mode 100644 src/test/java/org/apache/commons/math4/fitting/GaussianFitterTest.java delete mode 100644 src/test/java/org/apache/commons/math4/fitting/HarmonicFitterTest.java delete mode 100644 src/test/java/org/apache/commons/math4/fitting/PolynomialFitterTest.java diff --git a/src/main/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatingFunction.java b/src/main/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatingFunction.java deleted file mode 100644 index e9f7e190e..000000000 --- a/src/main/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatingFunction.java +++ /dev/null @@ -1,638 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import java.util.Arrays; - -import org.apache.commons.math4.analysis.BivariateFunction; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.NoDataException; -import org.apache.commons.math4.exception.NonMonotonicSequenceException; -import org.apache.commons.math4.exception.OutOfRangeException; -import org.apache.commons.math4.util.MathArrays; - -/** - * Function that implements the - * - * bicubic spline interpolation. Due to numerical accuracy issues this should not - * be used. - * - * @since 2.1 - * @deprecated as of 3.4 replaced by - * {@link org.apache.commons.math4.analysis.interpolation.PiecewiseBicubicSplineInterpolatingFunction} - */ -@Deprecated -public class BicubicSplineInterpolatingFunction - implements BivariateFunction { - /** Number of coefficients. */ - private static final int NUM_COEFF = 16; - /** - * Matrix to compute the spline coefficients from the function values - * and function derivatives values - */ - private static final double[][] AINV = { - { 1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0 }, - { -3,3,0,0,-2,-1,0,0,0,0,0,0,0,0,0,0 }, - { 2,-2,0,0,1,1,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0 }, - { 0,0,0,0,0,0,0,0,-3,3,0,0,-2,-1,0,0 }, - { 0,0,0,0,0,0,0,0,2,-2,0,0,1,1,0,0 }, - { -3,0,3,0,0,0,0,0,-2,0,-1,0,0,0,0,0 }, - { 0,0,0,0,-3,0,3,0,0,0,0,0,-2,0,-1,0 }, - { 9,-9,-9,9,6,3,-6,-3,6,-6,3,-3,4,2,2,1 }, - { -6,6,6,-6,-3,-3,3,3,-4,4,-2,2,-2,-2,-1,-1 }, - { 2,0,-2,0,0,0,0,0,1,0,1,0,0,0,0,0 }, - { 0,0,0,0,2,0,-2,0,0,0,0,0,1,0,1,0 }, - { -6,6,6,-6,-4,-2,4,2,-3,3,-3,3,-2,-1,-2,-1 }, - { 4,-4,-4,4,2,2,-2,-2,2,-2,2,-2,1,1,1,1 } - }; - - /** Samples x-coordinates */ - private final double[] xval; - /** Samples y-coordinates */ - private final double[] yval; - /** Set of cubic splines patching the whole data grid */ - private final BicubicSplineFunction[][] splines; - /** - * Partial derivatives. - * The value of the first index determines the kind of derivatives: - * 0 = first partial derivatives wrt x - * 1 = first partial derivatives wrt y - * 2 = second partial derivatives wrt x - * 3 = second partial derivatives wrt y - * 4 = cross partial derivatives - */ - private final BivariateFunction[][][] partialDerivatives; - - /** - * @param x Sample values of the x-coordinate, in increasing order. - * @param y Sample values of the y-coordinate, in increasing order. - * @param f Values of the function on every grid point. - * @param dFdX Values of the partial derivative of function with respect - * to x on every grid point. - * @param dFdY Values of the partial derivative of function with respect - * to y on every grid point. - * @param d2FdXdY Values of the cross partial derivative of function on - * every grid point. - * @throws DimensionMismatchException if the various arrays do not contain - * the expected number of elements. - * @throws NonMonotonicSequenceException if {@code x} or {@code y} are - * not strictly increasing. - * @throws NoDataException if any of the arrays has zero length. - */ - public BicubicSplineInterpolatingFunction(double[] x, - double[] y, - double[][] f, - double[][] dFdX, - double[][] dFdY, - double[][] d2FdXdY) - throws DimensionMismatchException, - NoDataException, - NonMonotonicSequenceException { - this(x, y, f, dFdX, dFdY, d2FdXdY, false); - } - - /** - * @param x Sample values of the x-coordinate, in increasing order. - * @param y Sample values of the y-coordinate, in increasing order. - * @param f Values of the function on every grid point. - * @param dFdX Values of the partial derivative of function with respect - * to x on every grid point. - * @param dFdY Values of the partial derivative of function with respect - * to y on every grid point. - * @param d2FdXdY Values of the cross partial derivative of function on - * every grid point. - * @param initializeDerivatives Whether to initialize the internal data - * needed for calling any of the methods that compute the partial derivatives - * this function. - * @throws DimensionMismatchException if the various arrays do not contain - * the expected number of elements. - * @throws NonMonotonicSequenceException if {@code x} or {@code y} are - * not strictly increasing. - * @throws NoDataException if any of the arrays has zero length. - * - * @see #partialDerivativeX(double,double) - * @see #partialDerivativeY(double,double) - * @see #partialDerivativeXX(double,double) - * @see #partialDerivativeYY(double,double) - * @see #partialDerivativeXY(double,double) - */ - public BicubicSplineInterpolatingFunction(double[] x, - double[] y, - double[][] f, - double[][] dFdX, - double[][] dFdY, - double[][] d2FdXdY, - boolean initializeDerivatives) - throws DimensionMismatchException, - NoDataException, - NonMonotonicSequenceException { - final int xLen = x.length; - final int yLen = y.length; - - if (xLen == 0 || yLen == 0 || f.length == 0 || f[0].length == 0) { - throw new NoDataException(); - } - if (xLen != f.length) { - throw new DimensionMismatchException(xLen, f.length); - } - if (xLen != dFdX.length) { - throw new DimensionMismatchException(xLen, dFdX.length); - } - if (xLen != dFdY.length) { - throw new DimensionMismatchException(xLen, dFdY.length); - } - if (xLen != d2FdXdY.length) { - throw new DimensionMismatchException(xLen, d2FdXdY.length); - } - - MathArrays.checkOrder(x); - MathArrays.checkOrder(y); - - xval = x.clone(); - yval = y.clone(); - - final int lastI = xLen - 1; - final int lastJ = yLen - 1; - splines = new BicubicSplineFunction[lastI][lastJ]; - - for (int i = 0; i < lastI; i++) { - if (f[i].length != yLen) { - throw new DimensionMismatchException(f[i].length, yLen); - } - if (dFdX[i].length != yLen) { - throw new DimensionMismatchException(dFdX[i].length, yLen); - } - if (dFdY[i].length != yLen) { - throw new DimensionMismatchException(dFdY[i].length, yLen); - } - if (d2FdXdY[i].length != yLen) { - throw new DimensionMismatchException(d2FdXdY[i].length, yLen); - } - final int ip1 = i + 1; - for (int j = 0; j < lastJ; j++) { - final int jp1 = j + 1; - final double[] beta = new double[] { - f[i][j], f[ip1][j], f[i][jp1], f[ip1][jp1], - dFdX[i][j], dFdX[ip1][j], dFdX[i][jp1], dFdX[ip1][jp1], - dFdY[i][j], dFdY[ip1][j], dFdY[i][jp1], dFdY[ip1][jp1], - d2FdXdY[i][j], d2FdXdY[ip1][j], d2FdXdY[i][jp1], d2FdXdY[ip1][jp1] - }; - - splines[i][j] = new BicubicSplineFunction(computeSplineCoefficients(beta), - initializeDerivatives); - } - } - - if (initializeDerivatives) { - // Compute all partial derivatives. - partialDerivatives = new BivariateFunction[5][lastI][lastJ]; - - for (int i = 0; i < lastI; i++) { - for (int j = 0; j < lastJ; j++) { - final BicubicSplineFunction bcs = splines[i][j]; - partialDerivatives[0][i][j] = bcs.partialDerivativeX(); - partialDerivatives[1][i][j] = bcs.partialDerivativeY(); - partialDerivatives[2][i][j] = bcs.partialDerivativeXX(); - partialDerivatives[3][i][j] = bcs.partialDerivativeYY(); - partialDerivatives[4][i][j] = bcs.partialDerivativeXY(); - } - } - } else { - // Partial derivative methods cannot be used. - partialDerivatives = null; - } - } - - /** - * {@inheritDoc} - */ - public double value(double x, double y) - throws OutOfRangeException { - final int i = searchIndex(x, xval); - final int j = searchIndex(y, yval); - - final double xN = (x - xval[i]) / (xval[i + 1] - xval[i]); - final double yN = (y - yval[j]) / (yval[j + 1] - yval[j]); - - return splines[i][j].value(xN, yN); - } - - /** - * Indicates whether a point is within the interpolation range. - * - * @param x First coordinate. - * @param y Second coordinate. - * @return {@code true} if (x, y) is a valid point. - * @since 3.3 - */ - public boolean isValidPoint(double x, double y) { - if (x < xval[0] || - x > xval[xval.length - 1] || - y < yval[0] || - y > yval[yval.length - 1]) { - return false; - } else { - return true; - } - } - - /** - * @param x x-coordinate. - * @param y y-coordinate. - * @return the value at point (x, y) of the first partial derivative with - * respect to x. - * @throws OutOfRangeException if {@code x} (resp. {@code y}) is outside - * the range defined by the boundary values of {@code xval} (resp. - * {@code yval}). - * @throws NullPointerException if the internal data were not initialized - * (cf. {@link #BicubicSplineInterpolatingFunction(double[],double[],double[][], - * double[][],double[][],double[][],boolean) constructor}). - */ - public double partialDerivativeX(double x, double y) - throws OutOfRangeException { - return partialDerivative(0, x, y); - } - /** - * @param x x-coordinate. - * @param y y-coordinate. - * @return the value at point (x, y) of the first partial derivative with - * respect to y. - * @throws OutOfRangeException if {@code x} (resp. {@code y}) is outside - * the range defined by the boundary values of {@code xval} (resp. - * {@code yval}). - * @throws NullPointerException if the internal data were not initialized - * (cf. {@link #BicubicSplineInterpolatingFunction(double[],double[],double[][], - * double[][],double[][],double[][],boolean) constructor}). - */ - public double partialDerivativeY(double x, double y) - throws OutOfRangeException { - return partialDerivative(1, x, y); - } - /** - * @param x x-coordinate. - * @param y y-coordinate. - * @return the value at point (x, y) of the second partial derivative with - * respect to x. - * @throws OutOfRangeException if {@code x} (resp. {@code y}) is outside - * the range defined by the boundary values of {@code xval} (resp. - * {@code yval}). - * @throws NullPointerException if the internal data were not initialized - * (cf. {@link #BicubicSplineInterpolatingFunction(double[],double[],double[][], - * double[][],double[][],double[][],boolean) constructor}). - */ - public double partialDerivativeXX(double x, double y) - throws OutOfRangeException { - return partialDerivative(2, x, y); - } - /** - * @param x x-coordinate. - * @param y y-coordinate. - * @return the value at point (x, y) of the second partial derivative with - * respect to y. - * @throws OutOfRangeException if {@code x} (resp. {@code y}) is outside - * the range defined by the boundary values of {@code xval} (resp. - * {@code yval}). - * @throws NullPointerException if the internal data were not initialized - * (cf. {@link #BicubicSplineInterpolatingFunction(double[],double[],double[][], - * double[][],double[][],double[][],boolean) constructor}). - */ - public double partialDerivativeYY(double x, double y) - throws OutOfRangeException { - return partialDerivative(3, x, y); - } - /** - * @param x x-coordinate. - * @param y y-coordinate. - * @return the value at point (x, y) of the second partial cross-derivative. - * @throws OutOfRangeException if {@code x} (resp. {@code y}) is outside - * the range defined by the boundary values of {@code xval} (resp. - * {@code yval}). - * @throws NullPointerException if the internal data were not initialized - * (cf. {@link #BicubicSplineInterpolatingFunction(double[],double[],double[][], - * double[][],double[][],double[][],boolean) constructor}). - */ - public double partialDerivativeXY(double x, double y) - throws OutOfRangeException { - return partialDerivative(4, x, y); - } - - /** - * @param which First index in {@link #partialDerivatives}. - * @param x x-coordinate. - * @param y y-coordinate. - * @return the value at point (x, y) of the selected partial derivative. - * @throws OutOfRangeException if {@code x} (resp. {@code y}) is outside - * the range defined by the boundary values of {@code xval} (resp. - * {@code yval}). - * @throws NullPointerException if the internal data were not initialized - * (cf. {@link #BicubicSplineInterpolatingFunction(double[],double[],double[][], - * double[][],double[][],double[][],boolean) constructor}). - */ - private double partialDerivative(int which, double x, double y) - throws OutOfRangeException { - final int i = searchIndex(x, xval); - final int j = searchIndex(y, yval); - - final double xN = (x - xval[i]) / (xval[i + 1] - xval[i]); - final double yN = (y - yval[j]) / (yval[j + 1] - yval[j]); - - return partialDerivatives[which][i][j].value(xN, yN); - } - - /** - * @param c Coordinate. - * @param val Coordinate samples. - * @return the index in {@code val} corresponding to the interval - * containing {@code c}. - * @throws OutOfRangeException if {@code c} is out of the - * range defined by the boundary values of {@code val}. - */ - private int searchIndex(double c, double[] val) { - final int r = Arrays.binarySearch(val, c); - - if (r == -1 || - r == -val.length - 1) { - throw new OutOfRangeException(c, val[0], val[val.length - 1]); - } - - if (r < 0) { - // "c" in within an interpolation sub-interval: Return the - // index of the sample at the lower end of the sub-interval. - return -r - 2; - } - final int last = val.length - 1; - if (r == last) { - // "c" is the last sample of the range: Return the index - // of the sample at the lower end of the last sub-interval. - return last - 1; - } - - // "c" is another sample point. - return r; - } - - /** - * Compute the spline coefficients from the list of function values and - * function partial derivatives values at the four corners of a grid - * element. They must be specified in the following order: - * - * where the subscripts indicate the partial derivative with respect to - * the corresponding variable(s). - * - * @param beta List of function values and function partial derivatives - * values. - * @return the spline coefficients. - */ - private double[] computeSplineCoefficients(double[] beta) { - final double[] a = new double[NUM_COEFF]; - - for (int i = 0; i < NUM_COEFF; i++) { - double result = 0; - final double[] row = AINV[i]; - for (int j = 0; j < NUM_COEFF; j++) { - result += row[j] * beta[j]; - } - a[i] = result; - } - - return a; - } -} - -/** - * 2D-spline function. - * - */ -class BicubicSplineFunction implements BivariateFunction { - /** Number of points. */ - private static final short N = 4; - /** Coefficients */ - private final double[][] a; - /** First partial derivative along x. */ - private final BivariateFunction partialDerivativeX; - /** First partial derivative along y. */ - private final BivariateFunction partialDerivativeY; - /** Second partial derivative along x. */ - private final BivariateFunction partialDerivativeXX; - /** Second partial derivative along y. */ - private final BivariateFunction partialDerivativeYY; - /** Second crossed partial derivative. */ - private final BivariateFunction partialDerivativeXY; - - /** - * Simple constructor. - * - * @param coeff Spline coefficients. - */ - public BicubicSplineFunction(double[] coeff) { - this(coeff, false); - } - - /** - * Simple constructor. - * - * @param coeff Spline coefficients. - * @param initializeDerivatives Whether to initialize the internal data - * needed for calling any of the methods that compute the partial derivatives - * this function. - */ - public BicubicSplineFunction(double[] coeff, - boolean initializeDerivatives) { - a = new double[N][N]; - for (int i = 0; i < N; i++) { - for (int j = 0; j < N; j++) { - a[i][j] = coeff[i * N + j]; - } - } - - if (initializeDerivatives) { - // Compute all partial derivatives functions. - final double[][] aX = new double[N][N]; - final double[][] aY = new double[N][N]; - final double[][] aXX = new double[N][N]; - final double[][] aYY = new double[N][N]; - final double[][] aXY = new double[N][N]; - - for (int i = 0; i < N; i++) { - for (int j = 0; j < N; j++) { - final double c = a[i][j]; - aX[i][j] = i * c; - aY[i][j] = j * c; - aXX[i][j] = (i - 1) * aX[i][j]; - aYY[i][j] = (j - 1) * aY[i][j]; - aXY[i][j] = j * aX[i][j]; - } - } - - partialDerivativeX = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double[] pX = {0, 1, x, x2}; - - final double y2 = y * y; - final double y3 = y2 * y; - final double[] pY = {1, y, y2, y3}; - - return apply(pX, pY, aX); - } - }; - partialDerivativeY = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double x3 = x2 * x; - final double[] pX = {1, x, x2, x3}; - - final double y2 = y * y; - final double[] pY = {0, 1, y, y2}; - - return apply(pX, pY, aY); - } - }; - partialDerivativeXX = new BivariateFunction() { - public double value(double x, double y) { - final double[] pX = {0, 0, 1, x}; - - final double y2 = y * y; - final double y3 = y2 * y; - final double[] pY = {1, y, y2, y3}; - - return apply(pX, pY, aXX); - } - }; - partialDerivativeYY = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double x3 = x2 * x; - final double[] pX = {1, x, x2, x3}; - - final double[] pY = {0, 0, 1, y}; - - return apply(pX, pY, aYY); - } - }; - partialDerivativeXY = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double[] pX = {0, 1, x, x2}; - - final double y2 = y * y; - final double[] pY = {0, 1, y, y2}; - - return apply(pX, pY, aXY); - } - }; - } else { - partialDerivativeX = null; - partialDerivativeY = null; - partialDerivativeXX = null; - partialDerivativeYY = null; - partialDerivativeXY = null; - } - } - - /** - * {@inheritDoc} - */ - public double value(double x, double y) { - if (x < 0 || x > 1) { - throw new OutOfRangeException(x, 0, 1); - } - if (y < 0 || y > 1) { - throw new OutOfRangeException(y, 0, 1); - } - - final double x2 = x * x; - final double x3 = x2 * x; - final double[] pX = {1, x, x2, x3}; - - final double y2 = y * y; - final double y3 = y2 * y; - final double[] pY = {1, y, y2, y3}; - - return apply(pX, pY, a); - } - - /** - * Compute the value of the bicubic polynomial. - * - * @param pX Powers of the x-coordinate. - * @param pY Powers of the y-coordinate. - * @param coeff Spline coefficients. - * @return the interpolated value. - */ - private double apply(double[] pX, double[] pY, double[][] coeff) { - double result = 0; - for (int i = 0; i < N; i++) { - for (int j = 0; j < N; j++) { - result += coeff[i][j] * pX[i] * pY[j]; - } - } - - return result; - } - - /** - * @return the partial derivative wrt {@code x}. - */ - public BivariateFunction partialDerivativeX() { - return partialDerivativeX; - } - /** - * @return the partial derivative wrt {@code y}. - */ - public BivariateFunction partialDerivativeY() { - return partialDerivativeY; - } - /** - * @return the second partial derivative wrt {@code x}. - */ - public BivariateFunction partialDerivativeXX() { - return partialDerivativeXX; - } - /** - * @return the second partial derivative wrt {@code y}. - */ - public BivariateFunction partialDerivativeYY() { - return partialDerivativeYY; - } - /** - * @return the second partial cross-derivative. - */ - public BivariateFunction partialDerivativeXY() { - return partialDerivativeXY; - } -} diff --git a/src/main/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolator.java b/src/main/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolator.java deleted file mode 100644 index 53e726f18..000000000 --- a/src/main/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolator.java +++ /dev/null @@ -1,176 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.analysis.UnivariateFunction; -import org.apache.commons.math4.analysis.polynomials.PolynomialSplineFunction; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.NoDataException; -import org.apache.commons.math4.exception.NonMonotonicSequenceException; -import org.apache.commons.math4.exception.NumberIsTooSmallException; -import org.apache.commons.math4.util.MathArrays; - -/** - * Generates a bicubic interpolating function. Due to numerical accuracy issues this should not - * be used. - * - * @since 2.2 - * @deprecated as of 3.4 replaced by {@link org.apache.commons.math4.analysis.interpolation.PiecewiseBicubicSplineInterpolator} - */ -@Deprecated -public class BicubicSplineInterpolator - implements BivariateGridInterpolator { - /** Whether to initialize internal data used to compute the analytical - derivatives of the splines. */ - private final boolean initializeDerivatives; - - /** - * Default constructor. - * The argument {@link #BicubicSplineInterpolator(boolean) initializeDerivatives} - * is set to {@code false}. - */ - public BicubicSplineInterpolator() { - this(false); - } - - /** - * Creates an interpolator. - * - * @param initializeDerivatives Whether to initialize the internal data - * needed for calling any of the methods that compute the partial derivatives - * of the {@link BicubicSplineInterpolatingFunction function} returned from - * the call to {@link #interpolate(double[],double[],double[][]) interpolate}. - */ - public BicubicSplineInterpolator(boolean initializeDerivatives) { - this.initializeDerivatives = initializeDerivatives; - } - - /** - * {@inheritDoc} - */ - public BicubicSplineInterpolatingFunction interpolate(final double[] xval, - final double[] yval, - final double[][] fval) - throws NoDataException, DimensionMismatchException, - NonMonotonicSequenceException, NumberIsTooSmallException { - if (xval.length == 0 || yval.length == 0 || fval.length == 0) { - throw new NoDataException(); - } - if (xval.length != fval.length) { - throw new DimensionMismatchException(xval.length, fval.length); - } - - MathArrays.checkOrder(xval); - MathArrays.checkOrder(yval); - - final int xLen = xval.length; - final int yLen = yval.length; - - // Samples (first index is y-coordinate, i.e. subarray variable is x) - // 0 <= i < xval.length - // 0 <= j < yval.length - // fX[j][i] = f(xval[i], yval[j]) - final double[][] fX = new double[yLen][xLen]; - for (int i = 0; i < xLen; i++) { - if (fval[i].length != yLen) { - throw new DimensionMismatchException(fval[i].length, yLen); - } - - for (int j = 0; j < yLen; j++) { - fX[j][i] = fval[i][j]; - } - } - - final SplineInterpolator spInterpolator = new SplineInterpolator(); - - // For each line y[j] (0 <= j < yLen), construct a 1D spline with - // respect to variable x - final PolynomialSplineFunction[] ySplineX = new PolynomialSplineFunction[yLen]; - for (int j = 0; j < yLen; j++) { - ySplineX[j] = spInterpolator.interpolate(xval, fX[j]); - } - - // For each line x[i] (0 <= i < xLen), construct a 1D spline with - // respect to variable y generated by array fY_1[i] - final PolynomialSplineFunction[] xSplineY = new PolynomialSplineFunction[xLen]; - for (int i = 0; i < xLen; i++) { - xSplineY[i] = spInterpolator.interpolate(yval, fval[i]); - } - - // Partial derivatives with respect to x at the grid knots - final double[][] dFdX = new double[xLen][yLen]; - for (int j = 0; j < yLen; j++) { - final UnivariateFunction f = ySplineX[j].derivative(); - for (int i = 0; i < xLen; i++) { - dFdX[i][j] = f.value(xval[i]); - } - } - - // Partial derivatives with respect to y at the grid knots - final double[][] dFdY = new double[xLen][yLen]; - for (int i = 0; i < xLen; i++) { - final UnivariateFunction f = xSplineY[i].derivative(); - for (int j = 0; j < yLen; j++) { - dFdY[i][j] = f.value(yval[j]); - } - } - - // Cross partial derivatives - final double[][] d2FdXdY = new double[xLen][yLen]; - for (int i = 0; i < xLen ; i++) { - final int nI = nextIndex(i, xLen); - final int pI = previousIndex(i); - for (int j = 0; j < yLen; j++) { - final int nJ = nextIndex(j, yLen); - final int pJ = previousIndex(j); - d2FdXdY[i][j] = (fval[nI][nJ] - fval[nI][pJ] - - fval[pI][nJ] + fval[pI][pJ]) / - ((xval[nI] - xval[pI]) * (yval[nJ] - yval[pJ])); - } - } - - // Create the interpolating splines - return new BicubicSplineInterpolatingFunction(xval, yval, fval, - dFdX, dFdY, d2FdXdY, - initializeDerivatives); - } - - /** - * Computes the next index of an array, clipping if necessary. - * It is assumed (but not checked) that {@code i >= 0}. - * - * @param i Index. - * @param max Upper limit of the array. - * @return the next index. - */ - private int nextIndex(int i, int max) { - final int index = i + 1; - return index < max ? index : index - 1; - } - /** - * Computes the previous index of an array, clipping if necessary. - * It is assumed (but not checked) that {@code i} is smaller than the size - * of the array. - * - * @param i Index. - * @return the previous index. - */ - private int previousIndex(int i) { - final int index = i - 1; - return index >= 0 ? index : 0; - } -} diff --git a/src/main/java/org/apache/commons/math4/analysis/interpolation/SmoothingPolynomialBicubicSplineInterpolator.java b/src/main/java/org/apache/commons/math4/analysis/interpolation/SmoothingPolynomialBicubicSplineInterpolator.java deleted file mode 100644 index 243da0c25..000000000 --- a/src/main/java/org/apache/commons/math4/analysis/interpolation/SmoothingPolynomialBicubicSplineInterpolator.java +++ /dev/null @@ -1,171 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.analysis.polynomials.PolynomialFunction; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.NoDataException; -import org.apache.commons.math4.exception.NonMonotonicSequenceException; -import org.apache.commons.math4.exception.NotPositiveException; -import org.apache.commons.math4.exception.NullArgumentException; -import org.apache.commons.math4.fitting.PolynomialFitter; -import org.apache.commons.math4.optim.SimpleVectorValueChecker; -import org.apache.commons.math4.optim.nonlinear.vector.jacobian.GaussNewtonOptimizer; -import org.apache.commons.math4.util.MathArrays; -import org.apache.commons.math4.util.Precision; - -/** - * Generates a bicubic interpolation function. - * Prior to generating the interpolating function, the input is smoothed using - * polynomial fitting. - * - * @since 2.2 - * @deprecated To be removed in 4.0 (see MATH-1166). - */ -@Deprecated -public class SmoothingPolynomialBicubicSplineInterpolator - extends BicubicSplineInterpolator { - /** Fitter for x. */ - private final PolynomialFitter xFitter; - /** Degree of the fitting polynomial. */ - private final int xDegree; - /** Fitter for y. */ - private final PolynomialFitter yFitter; - /** Degree of the fitting polynomial. */ - private final int yDegree; - - /** - * Default constructor. The degree of the fitting polynomials is set to 3. - */ - public SmoothingPolynomialBicubicSplineInterpolator() { - this(3); - } - - /** - * @param degree Degree of the polynomial fitting functions. - * @exception NotPositiveException if degree is not positive - */ - public SmoothingPolynomialBicubicSplineInterpolator(int degree) - throws NotPositiveException { - this(degree, degree); - } - - /** - * @param xDegree Degree of the polynomial fitting functions along the - * x-dimension. - * @param yDegree Degree of the polynomial fitting functions along the - * y-dimension. - * @exception NotPositiveException if degrees are not positive - */ - public SmoothingPolynomialBicubicSplineInterpolator(int xDegree, int yDegree) - throws NotPositiveException { - if (xDegree < 0) { - throw new NotPositiveException(xDegree); - } - if (yDegree < 0) { - throw new NotPositiveException(yDegree); - } - this.xDegree = xDegree; - this.yDegree = yDegree; - - final double safeFactor = 1e2; - final SimpleVectorValueChecker checker - = new SimpleVectorValueChecker(safeFactor * Precision.EPSILON, - safeFactor * Precision.SAFE_MIN); - xFitter = new PolynomialFitter(new GaussNewtonOptimizer(false, checker)); - yFitter = new PolynomialFitter(new GaussNewtonOptimizer(false, checker)); - } - - /** - * {@inheritDoc} - */ - @Override - public BicubicSplineInterpolatingFunction interpolate(final double[] xval, - final double[] yval, - final double[][] fval) - throws NoDataException, NullArgumentException, - DimensionMismatchException, NonMonotonicSequenceException { - if (xval.length == 0 || yval.length == 0 || fval.length == 0) { - throw new NoDataException(); - } - if (xval.length != fval.length) { - throw new DimensionMismatchException(xval.length, fval.length); - } - - final int xLen = xval.length; - final int yLen = yval.length; - - for (int i = 0; i < xLen; i++) { - if (fval[i].length != yLen) { - throw new DimensionMismatchException(fval[i].length, yLen); - } - } - - MathArrays.checkOrder(xval); - MathArrays.checkOrder(yval); - - // For each line y[j] (0 <= j < yLen), construct a polynomial, with - // respect to variable x, fitting array fval[][j] - final PolynomialFunction[] yPolyX = new PolynomialFunction[yLen]; - for (int j = 0; j < yLen; j++) { - xFitter.clearObservations(); - for (int i = 0; i < xLen; i++) { - xFitter.addObservedPoint(1, xval[i], fval[i][j]); - } - - // Initial guess for the fit is zero for each coefficients (of which - // there are "xDegree" + 1). - yPolyX[j] = new PolynomialFunction(xFitter.fit(new double[xDegree + 1])); - } - - // For every knot (xval[i], yval[j]) of the grid, calculate corrected - // values fval_1 - final double[][] fval_1 = new double[xLen][yLen]; - for (int j = 0; j < yLen; j++) { - final PolynomialFunction f = yPolyX[j]; - for (int i = 0; i < xLen; i++) { - fval_1[i][j] = f.value(xval[i]); - } - } - - // For each line x[i] (0 <= i < xLen), construct a polynomial, with - // respect to variable y, fitting array fval_1[i][] - final PolynomialFunction[] xPolyY = new PolynomialFunction[xLen]; - for (int i = 0; i < xLen; i++) { - yFitter.clearObservations(); - for (int j = 0; j < yLen; j++) { - yFitter.addObservedPoint(1, yval[j], fval_1[i][j]); - } - - // Initial guess for the fit is zero for each coefficients (of which - // there are "yDegree" + 1). - xPolyY[i] = new PolynomialFunction(yFitter.fit(new double[yDegree + 1])); - } - - // For every knot (xval[i], yval[j]) of the grid, calculate corrected - // values fval_2 - final double[][] fval_2 = new double[xLen][yLen]; - for (int i = 0; i < xLen; i++) { - final PolynomialFunction f = xPolyY[i]; - for (int j = 0; j < yLen; j++) { - fval_2[i][j] = f.value(yval[j]); - } - } - - return super.interpolate(xval, yval, fval_2); - } -} diff --git a/src/main/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatingFunction.java b/src/main/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatingFunction.java deleted file mode 100644 index fa5f76ce1..000000000 --- a/src/main/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatingFunction.java +++ /dev/null @@ -1,482 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.analysis.TrivariateFunction; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.NoDataException; -import org.apache.commons.math4.exception.NonMonotonicSequenceException; -import org.apache.commons.math4.exception.OutOfRangeException; -import org.apache.commons.math4.util.MathArrays; - -/** - * Function that implements the - * - * tricubic spline interpolation, as proposed in - * - * Tricubic interpolation in three dimensions
- * F. Lekien and J. Marsden
- * Int. J. Numer. Meth. Engng 2005; 63:455-471 - * - * - * @since 2.2 - * @deprecated To be removed in 4.0 (see MATH-1166). - */ -@Deprecated -public class TricubicSplineInterpolatingFunction - implements TrivariateFunction { - /** - * Matrix to compute the spline coefficients from the function values - * and function derivatives values - */ - private static final double[][] AINV = { - { 1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { -3,3,0,0,0,0,0,0,-2,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 2,-2,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-3,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { -3,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,-3,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 9,-9,-9,9,0,0,0,0,6,3,-6,-3,0,0,0,0,6,-6,3,-3,0,0,0,0,0,0,0,0,0,0,0,0,4,2,2,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { -6,6,6,-6,0,0,0,0,-3,-3,3,3,0,0,0,0,-4,4,-2,2,0,0,0,0,0,0,0,0,0,0,0,0,-2,-2,-1,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 2,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,2,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { -6,6,6,-6,0,0,0,0,-4,-2,4,2,0,0,0,0,-3,3,-3,3,0,0,0,0,0,0,0,0,0,0,0,0,-2,-1,-2,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 4,-4,-4,4,0,0,0,0,2,2,-2,-2,0,0,0,0,2,-2,2,-2,0,0,0,0,0,0,0,0,0,0,0,0,1,1,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-3,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-3,3,0,0,0,0,0,0,-2,-1,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,-2,0,0,0,0,0,0,1,1,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-3,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,-1,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-3,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,-1,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,9,-9,-9,9,0,0,0,0,0,0,0,0,0,0,0,0,6,3,-6,-3,0,0,0,0,6,-6,3,-3,0,0,0,0,4,2,2,1,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-6,6,6,-6,0,0,0,0,0,0,0,0,0,0,0,0,-3,-3,3,3,0,0,0,0,-4,4,-2,2,0,0,0,0,-2,-2,-1,-1,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,1,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-6,6,6,-6,0,0,0,0,0,0,0,0,0,0,0,0,-4,-2,4,2,0,0,0,0,-3,3,-3,3,0,0,0,0,-2,-1,-2,-1,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,-4,-4,4,0,0,0,0,0,0,0,0,0,0,0,0,2,2,-2,-2,0,0,0,0,2,-2,2,-2,0,0,0,0,1,1,1,1,0,0,0,0 }, - {-3,0,0,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,0,0,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,-3,0,0,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,0,0,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 9,-9,0,0,-9,9,0,0,6,3,0,0,-6,-3,0,0,0,0,0,0,0,0,0,0,6,-6,0,0,3,-3,0,0,0,0,0,0,0,0,0,0,4,2,0,0,2,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { -6,6,0,0,6,-6,0,0,-3,-3,0,0,3,3,0,0,0,0,0,0,0,0,0,0,-4,4,0,0,-2,2,0,0,0,0,0,0,0,0,0,0,-2,-2,0,0,-1,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-3,0,0,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,0,0,-1,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-3,0,0,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,0,0,-1,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,9,-9,0,0,-9,9,0,0,0,0,0,0,0,0,0,0,6,3,0,0,-6,-3,0,0,0,0,0,0,0,0,0,0,6,-6,0,0,3,-3,0,0,4,2,0,0,2,1,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-6,6,0,0,6,-6,0,0,0,0,0,0,0,0,0,0,-3,-3,0,0,3,3,0,0,0,0,0,0,0,0,0,0,-4,4,0,0,-2,2,0,0,-2,-2,0,0,-1,-1,0,0 }, - { 9,0,-9,0,-9,0,9,0,0,0,0,0,0,0,0,0,6,0,3,0,-6,0,-3,0,6,0,-6,0,3,0,-3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,0,2,0,2,0,1,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,9,0,-9,0,-9,0,9,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,6,0,3,0,-6,0,-3,0,6,0,-6,0,3,0,-3,0,0,0,0,0,0,0,0,0,4,0,2,0,2,0,1,0 }, - { -27,27,27,-27,27,-27,-27,27,-18,-9,18,9,18,9,-18,-9,-18,18,-9,9,18,-18,9,-9,-18,18,18,-18,-9,9,9,-9,-12,-6,-6,-3,12,6,6,3,-12,-6,12,6,-6,-3,6,3,-12,12,-6,6,-6,6,-3,3,-8,-4,-4,-2,-4,-2,-2,-1 }, - { 18,-18,-18,18,-18,18,18,-18,9,9,-9,-9,-9,-9,9,9,12,-12,6,-6,-12,12,-6,6,12,-12,-12,12,6,-6,-6,6,6,6,3,3,-6,-6,-3,-3,6,6,-6,-6,3,3,-3,-3,8,-8,4,-4,4,-4,2,-2,4,4,2,2,2,2,1,1 }, - { -6,0,6,0,6,0,-6,0,0,0,0,0,0,0,0,0,-3,0,-3,0,3,0,3,0,-4,0,4,0,-2,0,2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,-2,0,-1,0,-1,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,-6,0,6,0,6,0,-6,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-3,0,-3,0,3,0,3,0,-4,0,4,0,-2,0,2,0,0,0,0,0,0,0,0,0,-2,0,-2,0,-1,0,-1,0 }, - { 18,-18,-18,18,-18,18,18,-18,12,6,-12,-6,-12,-6,12,6,9,-9,9,-9,-9,9,-9,9,12,-12,-12,12,6,-6,-6,6,6,3,6,3,-6,-3,-6,-3,8,4,-8,-4,4,2,-4,-2,6,-6,6,-6,3,-3,3,-3,4,2,4,2,2,1,2,1 }, - { -12,12,12,-12,12,-12,-12,12,-6,-6,6,6,6,6,-6,-6,-6,6,-6,6,6,-6,6,-6,-8,8,8,-8,-4,4,4,-4,-3,-3,-3,-3,3,3,3,3,-4,-4,4,4,-2,-2,2,2,-4,4,-4,4,-2,2,-2,2,-2,-2,-2,-2,-1,-1,-1,-1 }, - { 2,0,0,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,2,0,0,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { -6,6,0,0,6,-6,0,0,-4,-2,0,0,4,2,0,0,0,0,0,0,0,0,0,0,-3,3,0,0,-3,3,0,0,0,0,0,0,0,0,0,0,-2,-1,0,0,-2,-1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 4,-4,0,0,-4,4,0,0,2,2,0,0,-2,-2,0,0,0,0,0,0,0,0,0,0,2,-2,0,0,2,-2,0,0,0,0,0,0,0,0,0,0,1,1,0,0,1,1,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,0,0,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,0,0,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,0,0,1,0,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-6,6,0,0,6,-6,0,0,0,0,0,0,0,0,0,0,-4,-2,0,0,4,2,0,0,0,0,0,0,0,0,0,0,-3,3,0,0,-3,3,0,0,-2,-1,0,0,-2,-1,0,0 }, - { 0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,4,-4,0,0,-4,4,0,0,0,0,0,0,0,0,0,0,2,2,0,0,-2,-2,0,0,0,0,0,0,0,0,0,0,2,-2,0,0,2,-2,0,0,1,1,0,0,1,1,0,0 }, - { -6,0,6,0,6,0,-6,0,0,0,0,0,0,0,0,0,-4,0,-2,0,4,0,2,0,-3,0,3,0,-3,0,3,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-2,0,-1,0,-2,0,-1,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,-6,0,6,0,6,0,-6,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,-4,0,-2,0,4,0,2,0,-3,0,3,0,-3,0,3,0,0,0,0,0,0,0,0,0,-2,0,-1,0,-2,0,-1,0 }, - { 18,-18,-18,18,-18,18,18,-18,12,6,-12,-6,-12,-6,12,6,12,-12,6,-6,-12,12,-6,6,9,-9,-9,9,9,-9,-9,9,8,4,4,2,-8,-4,-4,-2,6,3,-6,-3,6,3,-6,-3,6,-6,3,-3,6,-6,3,-3,4,2,2,1,4,2,2,1 }, - { -12,12,12,-12,12,-12,-12,12,-6,-6,6,6,6,6,-6,-6,-8,8,-4,4,8,-8,4,-4,-6,6,6,-6,-6,6,6,-6,-4,-4,-2,-2,4,4,2,2,-3,-3,3,3,-3,-3,3,3,-4,4,-2,2,-4,4,-2,2,-2,-2,-1,-1,-2,-2,-1,-1 }, - { 4,0,-4,0,-4,0,4,0,0,0,0,0,0,0,0,0,2,0,2,0,-2,0,-2,0,2,0,-2,0,2,0,-2,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,1,0,1,0,1,0,1,0,0,0,0,0,0,0,0,0 }, - { 0,0,0,0,0,0,0,0,4,0,-4,0,-4,0,4,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,0,2,0,2,0,-2,0,-2,0,2,0,-2,0,2,0,-2,0,0,0,0,0,0,0,0,0,1,0,1,0,1,0,1,0 }, - { -12,12,12,-12,12,-12,-12,12,-8,-4,8,4,8,4,-8,-4,-6,6,-6,6,6,-6,6,-6,-6,6,6,-6,-6,6,6,-6,-4,-2,-4,-2,4,2,4,2,-4,-2,4,2,-4,-2,4,2,-3,3,-3,3,-3,3,-3,3,-2,-1,-2,-1,-2,-1,-2,-1 }, - { 8,-8,-8,8,-8,8,8,-8,4,4,-4,-4,-4,-4,4,4,4,-4,4,-4,-4,4,-4,4,4,-4,-4,4,4,-4,-4,4,2,2,2,2,-2,-2,-2,-2,2,2,-2,-2,2,2,-2,-2,2,-2,2,-2,2,-2,2,-2,1,1,1,1,1,1,1,1 } - }; - - /** Samples x-coordinates */ - private final double[] xval; - /** Samples y-coordinates */ - private final double[] yval; - /** Samples z-coordinates */ - private final double[] zval; - /** Set of cubic splines pacthing the whole data grid */ - private final TricubicSplineFunction[][][] splines; - - /** - * @param x Sample values of the x-coordinate, in increasing order. - * @param y Sample values of the y-coordinate, in increasing order. - * @param z Sample values of the y-coordinate, in increasing order. - * @param f Values of the function on every grid point. - * @param dFdX Values of the partial derivative of function with respect to x on every grid point. - * @param dFdY Values of the partial derivative of function with respect to y on every grid point. - * @param dFdZ Values of the partial derivative of function with respect to z on every grid point. - * @param d2FdXdY Values of the cross partial derivative of function on every grid point. - * @param d2FdXdZ Values of the cross partial derivative of function on every grid point. - * @param d2FdYdZ Values of the cross partial derivative of function on every grid point. - * @param d3FdXdYdZ Values of the cross partial derivative of function on every grid point. - * @throws NoDataException if any of the arrays has zero length. - * @throws DimensionMismatchException if the various arrays do not contain the expected number of elements. - * @throws NonMonotonicSequenceException if {@code x}, {@code y} or {@code z} are not strictly increasing. - */ - public TricubicSplineInterpolatingFunction(double[] x, - double[] y, - double[] z, - double[][][] f, - double[][][] dFdX, - double[][][] dFdY, - double[][][] dFdZ, - double[][][] d2FdXdY, - double[][][] d2FdXdZ, - double[][][] d2FdYdZ, - double[][][] d3FdXdYdZ) - throws NoDataException, - DimensionMismatchException, - NonMonotonicSequenceException { - final int xLen = x.length; - final int yLen = y.length; - final int zLen = z.length; - - if (xLen == 0 || yLen == 0 || z.length == 0 || f.length == 0 || f[0].length == 0) { - throw new NoDataException(); - } - if (xLen != f.length) { - throw new DimensionMismatchException(xLen, f.length); - } - if (xLen != dFdX.length) { - throw new DimensionMismatchException(xLen, dFdX.length); - } - if (xLen != dFdY.length) { - throw new DimensionMismatchException(xLen, dFdY.length); - } - if (xLen != dFdZ.length) { - throw new DimensionMismatchException(xLen, dFdZ.length); - } - if (xLen != d2FdXdY.length) { - throw new DimensionMismatchException(xLen, d2FdXdY.length); - } - if (xLen != d2FdXdZ.length) { - throw new DimensionMismatchException(xLen, d2FdXdZ.length); - } - if (xLen != d2FdYdZ.length) { - throw new DimensionMismatchException(xLen, d2FdYdZ.length); - } - if (xLen != d3FdXdYdZ.length) { - throw new DimensionMismatchException(xLen, d3FdXdYdZ.length); - } - - MathArrays.checkOrder(x); - MathArrays.checkOrder(y); - MathArrays.checkOrder(z); - - xval = x.clone(); - yval = y.clone(); - zval = z.clone(); - - final int lastI = xLen - 1; - final int lastJ = yLen - 1; - final int lastK = zLen - 1; - splines = new TricubicSplineFunction[lastI][lastJ][lastK]; - - for (int i = 0; i < lastI; i++) { - if (f[i].length != yLen) { - throw new DimensionMismatchException(f[i].length, yLen); - } - if (dFdX[i].length != yLen) { - throw new DimensionMismatchException(dFdX[i].length, yLen); - } - if (dFdY[i].length != yLen) { - throw new DimensionMismatchException(dFdY[i].length, yLen); - } - if (dFdZ[i].length != yLen) { - throw new DimensionMismatchException(dFdZ[i].length, yLen); - } - if (d2FdXdY[i].length != yLen) { - throw new DimensionMismatchException(d2FdXdY[i].length, yLen); - } - if (d2FdXdZ[i].length != yLen) { - throw new DimensionMismatchException(d2FdXdZ[i].length, yLen); - } - if (d2FdYdZ[i].length != yLen) { - throw new DimensionMismatchException(d2FdYdZ[i].length, yLen); - } - if (d3FdXdYdZ[i].length != yLen) { - throw new DimensionMismatchException(d3FdXdYdZ[i].length, yLen); - } - - final int ip1 = i + 1; - for (int j = 0; j < lastJ; j++) { - if (f[i][j].length != zLen) { - throw new DimensionMismatchException(f[i][j].length, zLen); - } - if (dFdX[i][j].length != zLen) { - throw new DimensionMismatchException(dFdX[i][j].length, zLen); - } - if (dFdY[i][j].length != zLen) { - throw new DimensionMismatchException(dFdY[i][j].length, zLen); - } - if (dFdZ[i][j].length != zLen) { - throw new DimensionMismatchException(dFdZ[i][j].length, zLen); - } - if (d2FdXdY[i][j].length != zLen) { - throw new DimensionMismatchException(d2FdXdY[i][j].length, zLen); - } - if (d2FdXdZ[i][j].length != zLen) { - throw new DimensionMismatchException(d2FdXdZ[i][j].length, zLen); - } - if (d2FdYdZ[i][j].length != zLen) { - throw new DimensionMismatchException(d2FdYdZ[i][j].length, zLen); - } - if (d3FdXdYdZ[i][j].length != zLen) { - throw new DimensionMismatchException(d3FdXdYdZ[i][j].length, zLen); - } - - final int jp1 = j + 1; - for (int k = 0; k < lastK; k++) { - final int kp1 = k + 1; - - final double[] beta = new double[] { - f[i][j][k], f[ip1][j][k], - f[i][jp1][k], f[ip1][jp1][k], - f[i][j][kp1], f[ip1][j][kp1], - f[i][jp1][kp1], f[ip1][jp1][kp1], - - dFdX[i][j][k], dFdX[ip1][j][k], - dFdX[i][jp1][k], dFdX[ip1][jp1][k], - dFdX[i][j][kp1], dFdX[ip1][j][kp1], - dFdX[i][jp1][kp1], dFdX[ip1][jp1][kp1], - - dFdY[i][j][k], dFdY[ip1][j][k], - dFdY[i][jp1][k], dFdY[ip1][jp1][k], - dFdY[i][j][kp1], dFdY[ip1][j][kp1], - dFdY[i][jp1][kp1], dFdY[ip1][jp1][kp1], - - dFdZ[i][j][k], dFdZ[ip1][j][k], - dFdZ[i][jp1][k], dFdZ[ip1][jp1][k], - dFdZ[i][j][kp1], dFdZ[ip1][j][kp1], - dFdZ[i][jp1][kp1], dFdZ[ip1][jp1][kp1], - - d2FdXdY[i][j][k], d2FdXdY[ip1][j][k], - d2FdXdY[i][jp1][k], d2FdXdY[ip1][jp1][k], - d2FdXdY[i][j][kp1], d2FdXdY[ip1][j][kp1], - d2FdXdY[i][jp1][kp1], d2FdXdY[ip1][jp1][kp1], - - d2FdXdZ[i][j][k], d2FdXdZ[ip1][j][k], - d2FdXdZ[i][jp1][k], d2FdXdZ[ip1][jp1][k], - d2FdXdZ[i][j][kp1], d2FdXdZ[ip1][j][kp1], - d2FdXdZ[i][jp1][kp1], d2FdXdZ[ip1][jp1][kp1], - - d2FdYdZ[i][j][k], d2FdYdZ[ip1][j][k], - d2FdYdZ[i][jp1][k], d2FdYdZ[ip1][jp1][k], - d2FdYdZ[i][j][kp1], d2FdYdZ[ip1][j][kp1], - d2FdYdZ[i][jp1][kp1], d2FdYdZ[ip1][jp1][kp1], - - d3FdXdYdZ[i][j][k], d3FdXdYdZ[ip1][j][k], - d3FdXdYdZ[i][jp1][k], d3FdXdYdZ[ip1][jp1][k], - d3FdXdYdZ[i][j][kp1], d3FdXdYdZ[ip1][j][kp1], - d3FdXdYdZ[i][jp1][kp1], d3FdXdYdZ[ip1][jp1][kp1], - }; - - splines[i][j][k] = new TricubicSplineFunction(computeSplineCoefficients(beta)); - } - } - } - } - - /** - * {@inheritDoc} - * - * @throws OutOfRangeException if any of the variables is outside its interpolation range. - */ - public double value(double x, double y, double z) - throws OutOfRangeException { - final int i = searchIndex(x, xval); - if (i == -1) { - throw new OutOfRangeException(x, xval[0], xval[xval.length - 1]); - } - final int j = searchIndex(y, yval); - if (j == -1) { - throw new OutOfRangeException(y, yval[0], yval[yval.length - 1]); - } - final int k = searchIndex(z, zval); - if (k == -1) { - throw new OutOfRangeException(z, zval[0], zval[zval.length - 1]); - } - - final double xN = (x - xval[i]) / (xval[i + 1] - xval[i]); - final double yN = (y - yval[j]) / (yval[j + 1] - yval[j]); - final double zN = (z - zval[k]) / (zval[k + 1] - zval[k]); - - return splines[i][j][k].value(xN, yN, zN); - } - - /** - * @param c Coordinate. - * @param val Coordinate samples. - * @return the index in {@code val} corresponding to the interval containing {@code c}, or {@code -1} - * if {@code c} is out of the range defined by the end values of {@code val}. - */ - private int searchIndex(double c, double[] val) { - if (c < val[0]) { - return -1; - } - - final int max = val.length; - for (int i = 1; i < max; i++) { - if (c <= val[i]) { - return i - 1; - } - } - - return -1; - } - - /** - * Compute the spline coefficients from the list of function values and - * function partial derivatives values at the four corners of a grid - * element. They must be specified in the following order: - * - * where the subscripts indicate the partial derivative with respect to - * the corresponding variable(s). - * - * @param beta List of function values and function partial derivatives values. - * @return the spline coefficients. - */ - private double[] computeSplineCoefficients(double[] beta) { - final int sz = 64; - final double[] a = new double[sz]; - - for (int i = 0; i < sz; i++) { - double result = 0; - final double[] row = AINV[i]; - for (int j = 0; j < sz; j++) { - result += row[j] * beta[j]; - } - a[i] = result; - } - - return a; - } -} - -/** - * 3D-spline function. - * - */ -class TricubicSplineFunction - implements TrivariateFunction { - /** Number of points. */ - private static final short N = 4; - /** Coefficients */ - private final double[][][] a = new double[N][N][N]; - - /** - * @param aV List of spline coefficients. - */ - public TricubicSplineFunction(double[] aV) { - for (int i = 0; i < N; i++) { - for (int j = 0; j < N; j++) { - for (int k = 0; k < N; k++) { - a[i][j][k] = aV[i + N * (j + N * k)]; - } - } - } - } - - /** - * @param x x-coordinate of the interpolation point. - * @param y y-coordinate of the interpolation point. - * @param z z-coordinate of the interpolation point. - * @return the interpolated value. - * @throws OutOfRangeException if {@code x}, {@code y} or - * {@code z} are not in the interval {@code [0, 1]}. - */ - public double value(double x, double y, double z) - throws OutOfRangeException { - if (x < 0 || x > 1) { - throw new OutOfRangeException(x, 0, 1); - } - if (y < 0 || y > 1) { - throw new OutOfRangeException(y, 0, 1); - } - if (z < 0 || z > 1) { - throw new OutOfRangeException(z, 0, 1); - } - - final double x2 = x * x; - final double x3 = x2 * x; - final double[] pX = { 1, x, x2, x3 }; - - final double y2 = y * y; - final double y3 = y2 * y; - final double[] pY = { 1, y, y2, y3 }; - - final double z2 = z * z; - final double z3 = z2 * z; - final double[] pZ = { 1, z, z2, z3 }; - - double result = 0; - for (int i = 0; i < N; i++) { - for (int j = 0; j < N; j++) { - for (int k = 0; k < N; k++) { - result += a[i][j][k] * pX[i] * pY[j] * pZ[k]; - } - } - } - - return result; - } -} diff --git a/src/main/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolator.java b/src/main/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolator.java deleted file mode 100644 index c068f74a3..000000000 --- a/src/main/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolator.java +++ /dev/null @@ -1,201 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.NoDataException; -import org.apache.commons.math4.exception.NonMonotonicSequenceException; -import org.apache.commons.math4.exception.NumberIsTooSmallException; -import org.apache.commons.math4.util.MathArrays; - -/** - * Generates a tricubic interpolating function. - * - * @since 2.2 - * @deprecated To be removed in 4.0 (see MATH-1166). - */ -@Deprecated -public class TricubicSplineInterpolator - implements TrivariateGridInterpolator { - /** - * {@inheritDoc} - */ - public TricubicSplineInterpolatingFunction interpolate(final double[] xval, - final double[] yval, - final double[] zval, - final double[][][] fval) - throws NoDataException, NumberIsTooSmallException, - DimensionMismatchException, NonMonotonicSequenceException { - if (xval.length == 0 || yval.length == 0 || zval.length == 0 || fval.length == 0) { - throw new NoDataException(); - } - if (xval.length != fval.length) { - throw new DimensionMismatchException(xval.length, fval.length); - } - - MathArrays.checkOrder(xval); - MathArrays.checkOrder(yval); - MathArrays.checkOrder(zval); - - final int xLen = xval.length; - final int yLen = yval.length; - final int zLen = zval.length; - - // Samples, re-ordered as (z, x, y) and (y, z, x) tuplets - // fvalXY[k][i][j] = f(xval[i], yval[j], zval[k]) - // fvalZX[j][k][i] = f(xval[i], yval[j], zval[k]) - final double[][][] fvalXY = new double[zLen][xLen][yLen]; - final double[][][] fvalZX = new double[yLen][zLen][xLen]; - for (int i = 0; i < xLen; i++) { - if (fval[i].length != yLen) { - throw new DimensionMismatchException(fval[i].length, yLen); - } - - for (int j = 0; j < yLen; j++) { - if (fval[i][j].length != zLen) { - throw new DimensionMismatchException(fval[i][j].length, zLen); - } - - for (int k = 0; k < zLen; k++) { - final double v = fval[i][j][k]; - fvalXY[k][i][j] = v; - fvalZX[j][k][i] = v; - } - } - } - - final BicubicSplineInterpolator bsi = new BicubicSplineInterpolator(true); - - // For each line x[i] (0 <= i < xLen), construct a 2D spline in y and z - final BicubicSplineInterpolatingFunction[] xSplineYZ - = new BicubicSplineInterpolatingFunction[xLen]; - for (int i = 0; i < xLen; i++) { - xSplineYZ[i] = bsi.interpolate(yval, zval, fval[i]); - } - - // For each line y[j] (0 <= j < yLen), construct a 2D spline in z and x - final BicubicSplineInterpolatingFunction[] ySplineZX - = new BicubicSplineInterpolatingFunction[yLen]; - for (int j = 0; j < yLen; j++) { - ySplineZX[j] = bsi.interpolate(zval, xval, fvalZX[j]); - } - - // For each line z[k] (0 <= k < zLen), construct a 2D spline in x and y - final BicubicSplineInterpolatingFunction[] zSplineXY - = new BicubicSplineInterpolatingFunction[zLen]; - for (int k = 0; k < zLen; k++) { - zSplineXY[k] = bsi.interpolate(xval, yval, fvalXY[k]); - } - - // Partial derivatives wrt x and wrt y - final double[][][] dFdX = new double[xLen][yLen][zLen]; - final double[][][] dFdY = new double[xLen][yLen][zLen]; - final double[][][] d2FdXdY = new double[xLen][yLen][zLen]; - for (int k = 0; k < zLen; k++) { - final BicubicSplineInterpolatingFunction f = zSplineXY[k]; - for (int i = 0; i < xLen; i++) { - final double x = xval[i]; - for (int j = 0; j < yLen; j++) { - final double y = yval[j]; - dFdX[i][j][k] = f.partialDerivativeX(x, y); - dFdY[i][j][k] = f.partialDerivativeY(x, y); - d2FdXdY[i][j][k] = f.partialDerivativeXY(x, y); - } - } - } - - // Partial derivatives wrt y and wrt z - final double[][][] dFdZ = new double[xLen][yLen][zLen]; - final double[][][] d2FdYdZ = new double[xLen][yLen][zLen]; - for (int i = 0; i < xLen; i++) { - final BicubicSplineInterpolatingFunction f = xSplineYZ[i]; - for (int j = 0; j < yLen; j++) { - final double y = yval[j]; - for (int k = 0; k < zLen; k++) { - final double z = zval[k]; - dFdZ[i][j][k] = f.partialDerivativeY(y, z); - d2FdYdZ[i][j][k] = f.partialDerivativeXY(y, z); - } - } - } - - // Partial derivatives wrt x and wrt z - final double[][][] d2FdZdX = new double[xLen][yLen][zLen]; - for (int j = 0; j < yLen; j++) { - final BicubicSplineInterpolatingFunction f = ySplineZX[j]; - for (int k = 0; k < zLen; k++) { - final double z = zval[k]; - for (int i = 0; i < xLen; i++) { - final double x = xval[i]; - d2FdZdX[i][j][k] = f.partialDerivativeXY(z, x); - } - } - } - - // Third partial cross-derivatives - final double[][][] d3FdXdYdZ = new double[xLen][yLen][zLen]; - for (int i = 0; i < xLen ; i++) { - final int nI = nextIndex(i, xLen); - final int pI = previousIndex(i); - for (int j = 0; j < yLen; j++) { - final int nJ = nextIndex(j, yLen); - final int pJ = previousIndex(j); - for (int k = 0; k < zLen; k++) { - final int nK = nextIndex(k, zLen); - final int pK = previousIndex(k); - - // XXX Not sure about this formula - d3FdXdYdZ[i][j][k] = (fval[nI][nJ][nK] - fval[nI][pJ][nK] - - fval[pI][nJ][nK] + fval[pI][pJ][nK] - - fval[nI][nJ][pK] + fval[nI][pJ][pK] + - fval[pI][nJ][pK] - fval[pI][pJ][pK]) / - ((xval[nI] - xval[pI]) * (yval[nJ] - yval[pJ]) * (zval[nK] - zval[pK])) ; - } - } - } - - // Create the interpolating splines - return new TricubicSplineInterpolatingFunction(xval, yval, zval, fval, - dFdX, dFdY, dFdZ, - d2FdXdY, d2FdZdX, d2FdYdZ, - d3FdXdYdZ); - } - - /** - * Compute the next index of an array, clipping if necessary. - * It is assumed (but not checked) that {@code i} is larger than or equal to 0}. - * - * @param i Index - * @param max Upper limit of the array - * @return the next index - */ - private int nextIndex(int i, int max) { - final int index = i + 1; - return index < max ? index : index - 1; - } - /** - * Compute the previous index of an array, clipping if necessary. - * It is assumed (but not checked) that {@code i} is smaller than the size of the array. - * - * @param i Index - * @return the previous index - */ - private int previousIndex(int i) { - final int index = i - 1; - return index >= 0 ? index : 0; - } -} diff --git a/src/main/java/org/apache/commons/math4/analysis/solvers/NewtonSolver.java b/src/main/java/org/apache/commons/math4/analysis/solvers/NewtonSolver.java deleted file mode 100644 index f3770308e..000000000 --- a/src/main/java/org/apache/commons/math4/analysis/solvers/NewtonSolver.java +++ /dev/null @@ -1,92 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ - -package org.apache.commons.math4.analysis.solvers; - -import org.apache.commons.math4.analysis.DifferentiableUnivariateFunction; -import org.apache.commons.math4.exception.TooManyEvaluationsException; -import org.apache.commons.math4.util.FastMath; - -/** - * Implements - * Newton's Method for finding zeros of real univariate functions. - *

- * The function should be continuous but not necessarily smooth.

- * - * @deprecated as of 3.1, replaced by {@link NewtonRaphsonSolver} - */ -@Deprecated -public class NewtonSolver extends AbstractDifferentiableUnivariateSolver { - /** Default absolute accuracy. */ - private static final double DEFAULT_ABSOLUTE_ACCURACY = 1e-6; - - /** - * Construct a solver. - */ - public NewtonSolver() { - this(DEFAULT_ABSOLUTE_ACCURACY); - } - /** - * Construct a solver. - * - * @param absoluteAccuracy Absolute accuracy. - */ - public NewtonSolver(double absoluteAccuracy) { - super(absoluteAccuracy); - } - - /** - * Find a zero near the midpoint of {@code min} and {@code max}. - * - * @param f Function to solve. - * @param min Lower bound for the interval. - * @param max Upper bound for the interval. - * @param maxEval Maximum number of evaluations. - * @return the value where the function is zero. - * @throws org.apache.commons.math4.exception.TooManyEvaluationsException - * if the maximum evaluation count is exceeded. - * @throws org.apache.commons.math4.exception.NumberIsTooLargeException - * if {@code min >= max}. - */ - @Override - public double solve(int maxEval, final DifferentiableUnivariateFunction f, - final double min, final double max) - throws TooManyEvaluationsException { - return super.solve(maxEval, f, UnivariateSolverUtils.midpoint(min, max)); - } - - /** - * {@inheritDoc} - */ - @Override - protected double doSolve() - throws TooManyEvaluationsException { - final double startValue = getStartValue(); - final double absoluteAccuracy = getAbsoluteAccuracy(); - - double x0 = startValue; - double x1; - while (true) { - x1 = x0 - (computeObjectiveValue(x0) / computeDerivativeObjectiveValue(x0)); - if (FastMath.abs(x1 - x0) <= absoluteAccuracy) { - return x1; - } - - x0 = x1; - } - } -} diff --git a/src/main/java/org/apache/commons/math4/fitting/CurveFitter.java b/src/main/java/org/apache/commons/math4/fitting/CurveFitter.java deleted file mode 100644 index 8cce42659..000000000 --- a/src/main/java/org/apache/commons/math4/fitting/CurveFitter.java +++ /dev/null @@ -1,233 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.fitting; - -import java.util.ArrayList; -import java.util.List; - -import org.apache.commons.math4.analysis.MultivariateMatrixFunction; -import org.apache.commons.math4.analysis.MultivariateVectorFunction; -import org.apache.commons.math4.analysis.ParametricUnivariateFunction; -import org.apache.commons.math4.optim.InitialGuess; -import org.apache.commons.math4.optim.MaxEval; -import org.apache.commons.math4.optim.PointVectorValuePair; -import org.apache.commons.math4.optim.nonlinear.vector.ModelFunction; -import org.apache.commons.math4.optim.nonlinear.vector.ModelFunctionJacobian; -import org.apache.commons.math4.optim.nonlinear.vector.MultivariateVectorOptimizer; -import org.apache.commons.math4.optim.nonlinear.vector.Target; -import org.apache.commons.math4.optim.nonlinear.vector.Weight; - -/** - * Fitter for parametric univariate real functions y = f(x). - *
- * When a univariate real function y = f(x) does depend on some - * unknown parameters p0, p1 ... pn-1, - * this class can be used to find these parameters. It does this - * by fitting the curve so it remains very close to a set of - * observed points (x0, y0), (x1, - * y1) ... (xk-1, yk-1). This fitting - * is done by finding the parameters values that minimizes the objective - * function ∑(yi-f(xi))2. This is - * really a least squares problem. - * - * @param Function to use for the fit. - * - * @since 2.0 - * @deprecated As of 3.3. Please use {@link AbstractCurveFitter} and - * {@link WeightedObservedPoints} instead. - */ -@Deprecated -public class CurveFitter { - /** Optimizer to use for the fitting. */ - private final MultivariateVectorOptimizer optimizer; - /** Observed points. */ - private final List observations; - - /** - * Simple constructor. - * - * @param optimizer Optimizer to use for the fitting. - * @since 3.1 - */ - public CurveFitter(final MultivariateVectorOptimizer optimizer) { - this.optimizer = optimizer; - observations = new ArrayList(); - } - - /** Add an observed (x,y) point to the sample with unit weight. - *

Calling this method is equivalent to call - * {@code addObservedPoint(1.0, x, y)}.

- * @param x abscissa of the point - * @param y observed value of the point at x, after fitting we should - * have f(x) as close as possible to this value - * @see #addObservedPoint(double, double, double) - * @see #addObservedPoint(WeightedObservedPoint) - * @see #getObservations() - */ - public void addObservedPoint(double x, double y) { - addObservedPoint(1.0, x, y); - } - - /** Add an observed weighted (x,y) point to the sample. - * @param weight weight of the observed point in the fit - * @param x abscissa of the point - * @param y observed value of the point at x, after fitting we should - * have f(x) as close as possible to this value - * @see #addObservedPoint(double, double) - * @see #addObservedPoint(WeightedObservedPoint) - * @see #getObservations() - */ - public void addObservedPoint(double weight, double x, double y) { - observations.add(new WeightedObservedPoint(weight, x, y)); - } - - /** Add an observed weighted (x,y) point to the sample. - * @param observed observed point to add - * @see #addObservedPoint(double, double) - * @see #addObservedPoint(double, double, double) - * @see #getObservations() - */ - public void addObservedPoint(WeightedObservedPoint observed) { - observations.add(observed); - } - - /** Get the observed points. - * @return observed points - * @see #addObservedPoint(double, double) - * @see #addObservedPoint(double, double, double) - * @see #addObservedPoint(WeightedObservedPoint) - */ - public WeightedObservedPoint[] getObservations() { - return observations.toArray(new WeightedObservedPoint[observations.size()]); - } - - /** - * Remove all observations. - */ - public void clearObservations() { - observations.clear(); - } - - /** - * Fit a curve. - * This method compute the coefficients of the curve that best - * fit the sample of observed points previously given through calls - * to the {@link #addObservedPoint(WeightedObservedPoint) - * addObservedPoint} method. - * - * @param f parametric function to fit. - * @param initialGuess first guess of the function parameters. - * @return the fitted parameters. - * @throws org.apache.commons.math4.exception.DimensionMismatchException - * if the start point dimension is wrong. - */ - public double[] fit(T f, final double[] initialGuess) { - return fit(Integer.MAX_VALUE, f, initialGuess); - } - - /** - * Fit a curve. - * This method compute the coefficients of the curve that best - * fit the sample of observed points previously given through calls - * to the {@link #addObservedPoint(WeightedObservedPoint) - * addObservedPoint} method. - * - * @param f parametric function to fit. - * @param initialGuess first guess of the function parameters. - * @param maxEval Maximum number of function evaluations. - * @return the fitted parameters. - * @throws org.apache.commons.math4.exception.TooManyEvaluationsException - * if the number of allowed evaluations is exceeded. - * @throws org.apache.commons.math4.exception.DimensionMismatchException - * if the start point dimension is wrong. - * @since 3.0 - */ - public double[] fit(int maxEval, T f, - final double[] initialGuess) { - // Prepare least squares problem. - double[] target = new double[observations.size()]; - double[] weights = new double[observations.size()]; - int i = 0; - for (WeightedObservedPoint point : observations) { - target[i] = point.getY(); - weights[i] = point.getWeight(); - ++i; - } - - // Input to the optimizer: the model and its Jacobian. - final TheoreticalValuesFunction model = new TheoreticalValuesFunction(f); - - // Perform the fit. - final PointVectorValuePair optimum - = optimizer.optimize(new MaxEval(maxEval), - model.getModelFunction(), - model.getModelFunctionJacobian(), - new Target(target), - new Weight(weights), - new InitialGuess(initialGuess)); - // Extract the coefficients. - return optimum.getPointRef(); - } - - /** Vectorial function computing function theoretical values. */ - private class TheoreticalValuesFunction { - /** Function to fit. */ - private final ParametricUnivariateFunction f; - - /** - * @param f function to fit. - */ - public TheoreticalValuesFunction(final ParametricUnivariateFunction f) { - this.f = f; - } - - /** - * @return the model function values. - */ - public ModelFunction getModelFunction() { - return new ModelFunction(new MultivariateVectorFunction() { - /** {@inheritDoc} */ - public double[] value(double[] point) { - // compute the residuals - final double[] values = new double[observations.size()]; - int i = 0; - for (WeightedObservedPoint observed : observations) { - values[i++] = f.value(observed.getX(), point); - } - - return values; - } - }); - } - - /** - * @return the model function Jacobian. - */ - public ModelFunctionJacobian getModelFunctionJacobian() { - return new ModelFunctionJacobian(new MultivariateMatrixFunction() { - public double[][] value(double[] point) { - final double[][] jacobian = new double[observations.size()][]; - int i = 0; - for (WeightedObservedPoint observed : observations) { - jacobian[i++] = f.gradient(observed.getX(), point); - } - return jacobian; - } - }); - } - } -} diff --git a/src/main/java/org/apache/commons/math4/fitting/GaussianFitter.java b/src/main/java/org/apache/commons/math4/fitting/GaussianFitter.java deleted file mode 100644 index 285a4678f..000000000 --- a/src/main/java/org/apache/commons/math4/fitting/GaussianFitter.java +++ /dev/null @@ -1,365 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.fitting; - -import java.util.Arrays; -import java.util.Comparator; - -import org.apache.commons.math4.analysis.function.Gaussian; -import org.apache.commons.math4.exception.NotStrictlyPositiveException; -import org.apache.commons.math4.exception.NullArgumentException; -import org.apache.commons.math4.exception.NumberIsTooSmallException; -import org.apache.commons.math4.exception.OutOfRangeException; -import org.apache.commons.math4.exception.ZeroException; -import org.apache.commons.math4.exception.util.LocalizedFormats; -import org.apache.commons.math4.optim.nonlinear.vector.MultivariateVectorOptimizer; -import org.apache.commons.math4.util.FastMath; - -/** - * Fits points to a {@link - * org.apache.commons.math4.analysis.function.Gaussian.Parametric Gaussian} function. - *

- * Usage example: - * 

- *   GaussianFitter fitter = new GaussianFitter(
- *     new LevenbergMarquardtOptimizer());
- *   fitter.addObservedPoint(4.0254623,  531026.0);
- *   fitter.addObservedPoint(4.03128248, 984167.0);
- *   fitter.addObservedPoint(4.03839603, 1887233.0);
- *   fitter.addObservedPoint(4.04421621, 2687152.0);
- *   fitter.addObservedPoint(4.05132976, 3461228.0);
- *   fitter.addObservedPoint(4.05326982, 3580526.0);
- *   fitter.addObservedPoint(4.05779662, 3439750.0);
- *   fitter.addObservedPoint(4.0636168,  2877648.0);
- *   fitter.addObservedPoint(4.06943698, 2175960.0);
- *   fitter.addObservedPoint(4.07525716, 1447024.0);
- *   fitter.addObservedPoint(4.08237071, 717104.0);
- *   fitter.addObservedPoint(4.08366408, 620014.0);
- *   double[] parameters = fitter.fit();
- *
- * - * @since 2.2 - * @deprecated As of 3.3. Please use {@link GaussianCurveFitter} and - * {@link WeightedObservedPoints} instead. - */ -@Deprecated -public class GaussianFitter extends CurveFitter { - /** - * Constructs an instance using the specified optimizer. - * - * @param optimizer Optimizer to use for the fitting. - */ - public GaussianFitter(MultivariateVectorOptimizer optimizer) { - super(optimizer); - } - - /** - * Fits a Gaussian function to the observed points. - * - * @param initialGuess First guess values in the following order: - *
    - *
  • Norm
    • - *
  • Mean
    • - *
  • Sigma
    • - *
- * @return the parameters of the Gaussian function that best fits the - * observed points (in the same order as above). - * @since 3.0 - */ - public double[] fit(double[] initialGuess) { - final Gaussian.Parametric f = new Gaussian.Parametric() { - @Override - public double value(double x, double ... p) { - double v = Double.POSITIVE_INFINITY; - try { - v = super.value(x, p); - } catch (NotStrictlyPositiveException e) { // NOPMD - // Do nothing. - } - return v; - } - - @Override - public double[] gradient(double x, double ... p) { - double[] v = { Double.POSITIVE_INFINITY, - Double.POSITIVE_INFINITY, - Double.POSITIVE_INFINITY }; - try { - v = super.gradient(x, p); - } catch (NotStrictlyPositiveException e) { // NOPMD - // Do nothing. - } - return v; - } - }; - - return fit(f, initialGuess); - } - - /** - * Fits a Gaussian function to the observed points. - * - * @return the parameters of the Gaussian function that best fits the - * observed points (in the same order as above). - */ - public double[] fit() { - final double[] guess = (new ParameterGuesser(getObservations())).guess(); - return fit(guess); - } - - /** - * Guesses the parameters {@code norm}, {@code mean}, and {@code sigma} - * of a {@link org.apache.commons.math4.analysis.function.Gaussian.Parametric} - * based on the specified observed points. - */ - public static class ParameterGuesser { - /** Normalization factor. */ - private final double norm; - /** Mean. */ - private final double mean; - /** Standard deviation. */ - private final double sigma; - - /** - * Constructs instance with the specified observed points. - * - * @param observations Observed points from which to guess the - * parameters of the Gaussian. - * @throws NullArgumentException if {@code observations} is - * {@code null}. - * @throws NumberIsTooSmallException if there are less than 3 - * observations. - */ - public ParameterGuesser(WeightedObservedPoint[] observations) { - if (observations == null) { - throw new NullArgumentException(LocalizedFormats.INPUT_ARRAY); - } - if (observations.length < 3) { - throw new NumberIsTooSmallException(observations.length, 3, true); - } - - final WeightedObservedPoint[] sorted = sortObservations(observations); - final double[] params = basicGuess(sorted); - - norm = params[0]; - mean = params[1]; - sigma = params[2]; - } - - /** - * Gets an estimation of the parameters. - * - * @return the guessed parameters, in the following order: - *
    - *
  • Normalization factor
    • - *
  • Mean
    • - *
  • Standard deviation
    • - *
- */ - public double[] guess() { - return new double[] { norm, mean, sigma }; - } - - /** - * Sort the observations. - * - * @param unsorted Input observations. - * @return the input observations, sorted. - */ - private WeightedObservedPoint[] sortObservations(WeightedObservedPoint[] unsorted) { - final WeightedObservedPoint[] observations = unsorted.clone(); - final Comparator cmp - = new Comparator() { - public int compare(WeightedObservedPoint p1, - WeightedObservedPoint p2) { - if (p1 == null && p2 == null) { - return 0; - } - if (p1 == null) { - return -1; - } - if (p2 == null) { - return 1; - } - if (p1.getX() < p2.getX()) { - return -1; - } - if (p1.getX() > p2.getX()) { - return 1; - } - if (p1.getY() < p2.getY()) { - return -1; - } - if (p1.getY() > p2.getY()) { - return 1; - } - if (p1.getWeight() < p2.getWeight()) { - return -1; - } - if (p1.getWeight() > p2.getWeight()) { - return 1; - } - return 0; - } - }; - - Arrays.sort(observations, cmp); - return observations; - } - - /** - * Guesses the parameters based on the specified observed points. - * - * @param points Observed points, sorted. - * @return the guessed parameters (normalization factor, mean and - * sigma). - */ - private double[] basicGuess(WeightedObservedPoint[] points) { - final int maxYIdx = findMaxY(points); - final double n = points[maxYIdx].getY(); - final double m = points[maxYIdx].getX(); - - double fwhmApprox; - try { - final double halfY = n + ((m - n) / 2); - final double fwhmX1 = interpolateXAtY(points, maxYIdx, -1, halfY); - final double fwhmX2 = interpolateXAtY(points, maxYIdx, 1, halfY); - fwhmApprox = fwhmX2 - fwhmX1; - } catch (OutOfRangeException e) { - // TODO: Exceptions should not be used for flow control. - fwhmApprox = points[points.length - 1].getX() - points[0].getX(); - } - final double s = fwhmApprox / (2 * FastMath.sqrt(2 * FastMath.log(2))); - - return new double[] { n, m, s }; - } - - /** - * Finds index of point in specified points with the largest Y. - * - * @param points Points to search. - * @return the index in specified points array. - */ - private int findMaxY(WeightedObservedPoint[] points) { - int maxYIdx = 0; - for (int i = 1; i < points.length; i++) { - if (points[i].getY() > points[maxYIdx].getY()) { - maxYIdx = i; - } - } - return maxYIdx; - } - - /** - * Interpolates using the specified points to determine X at the - * specified Y. - * - * @param points Points to use for interpolation. - * @param startIdx Index within points from which to start the search for - * interpolation bounds points. - * @param idxStep Index step for searching interpolation bounds points. - * @param y Y value for which X should be determined. - * @return the value of X for the specified Y. - * @throws ZeroException if {@code idxStep} is 0. - * @throws OutOfRangeException if specified {@code y} is not within the - * range of the specified {@code points}. - */ - private double interpolateXAtY(WeightedObservedPoint[] points, - int startIdx, - int idxStep, - double y) - throws OutOfRangeException { - if (idxStep == 0) { - throw new ZeroException(); - } - final WeightedObservedPoint[] twoPoints - = getInterpolationPointsForY(points, startIdx, idxStep, y); - final WeightedObservedPoint p1 = twoPoints[0]; - final WeightedObservedPoint p2 = twoPoints[1]; - if (p1.getY() == y) { - return p1.getX(); - } - if (p2.getY() == y) { - return p2.getX(); - } - return p1.getX() + (((y - p1.getY()) * (p2.getX() - p1.getX())) / - (p2.getY() - p1.getY())); - } - - /** - * Gets the two bounding interpolation points from the specified points - * suitable for determining X at the specified Y. - * - * @param points Points to use for interpolation. - * @param startIdx Index within points from which to start search for - * interpolation bounds points. - * @param idxStep Index step for search for interpolation bounds points. - * @param y Y value for which X should be determined. - * @return the array containing two points suitable for determining X at - * the specified Y. - * @throws ZeroException if {@code idxStep} is 0. - * @throws OutOfRangeException if specified {@code y} is not within the - * range of the specified {@code points}. - */ - private WeightedObservedPoint[] getInterpolationPointsForY(WeightedObservedPoint[] points, - int startIdx, - int idxStep, - double y) - throws OutOfRangeException { - if (idxStep == 0) { - throw new ZeroException(); - } - for (int i = startIdx; - idxStep < 0 ? i + idxStep >= 0 : i + idxStep < points.length; - i += idxStep) { - final WeightedObservedPoint p1 = points[i]; - final WeightedObservedPoint p2 = points[i + idxStep]; - if (isBetween(y, p1.getY(), p2.getY())) { - if (idxStep < 0) { - return new WeightedObservedPoint[] { p2, p1 }; - } else { - return new WeightedObservedPoint[] { p1, p2 }; - } - } - } - - // Boundaries are replaced by dummy values because the raised - // exception is caught and the message never displayed. - // TODO: Exceptions should not be used for flow control. - throw new OutOfRangeException(y, - Double.NEGATIVE_INFINITY, - Double.POSITIVE_INFINITY); - } - - /** - * Determines whether a value is between two other values. - * - * @param value Value to test whether it is between {@code boundary1} - * and {@code boundary2}. - * @param boundary1 One end of the range. - * @param boundary2 Other end of the range. - * @return {@code true} if {@code value} is between {@code boundary1} and - * {@code boundary2} (inclusive), {@code false} otherwise. - */ - private boolean isBetween(double value, - double boundary1, - double boundary2) { - return (value >= boundary1 && value <= boundary2) || - (value >= boundary2 && value <= boundary1); - } - } -} diff --git a/src/main/java/org/apache/commons/math4/fitting/HarmonicFitter.java b/src/main/java/org/apache/commons/math4/fitting/HarmonicFitter.java deleted file mode 100644 index e74d0ac30..000000000 --- a/src/main/java/org/apache/commons/math4/fitting/HarmonicFitter.java +++ /dev/null @@ -1,384 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.fitting; - -import org.apache.commons.math4.analysis.function.HarmonicOscillator; -import org.apache.commons.math4.exception.MathIllegalStateException; -import org.apache.commons.math4.exception.NumberIsTooSmallException; -import org.apache.commons.math4.exception.ZeroException; -import org.apache.commons.math4.exception.util.LocalizedFormats; -import org.apache.commons.math4.optim.nonlinear.vector.MultivariateVectorOptimizer; -import org.apache.commons.math4.util.FastMath; - -/** - * Class that implements a curve fitting specialized for sinusoids. - * - * Harmonic fitting is a very simple case of curve fitting. The - * estimated coefficients are the amplitude a, the pulsation ω and - * the phase φ: f (t) = a cos (ω t + φ). They are - * searched by a least square estimator initialized with a rough guess - * based on integrals. - * - * @since 2.0 - * @deprecated As of 3.3. Please use {@link HarmonicCurveFitter} and - * {@link WeightedObservedPoints} instead. - */ -@Deprecated -public class HarmonicFitter extends CurveFitter { - /** - * Simple constructor. - * @param optimizer Optimizer to use for the fitting. - */ - public HarmonicFitter(final MultivariateVectorOptimizer optimizer) { - super(optimizer); - } - - /** - * Fit an harmonic function to the observed points. - * - * @param initialGuess First guess values in the following order: - *
    - *
  • Amplitude
    • - *
  • Angular frequency
    • - *
  • Phase
    • - *
- * @return the parameters of the harmonic function that best fits the - * observed points (in the same order as above). - */ - public double[] fit(double[] initialGuess) { - return fit(new HarmonicOscillator.Parametric(), initialGuess); - } - - /** - * Fit an harmonic function to the observed points. - * An initial guess will be automatically computed. - * - * @return the parameters of the harmonic function that best fits the - * observed points (see the other {@link #fit(double[]) fit} method. - * @throws NumberIsTooSmallException if the sample is too short for the - * the first guess to be computed. - * @throws ZeroException if the first guess cannot be computed because - * the abscissa range is zero. - */ - public double[] fit() { - return fit((new ParameterGuesser(getObservations())).guess()); - } - - /** - * This class guesses harmonic coefficients from a sample. - *

The algorithm used to guess the coefficients is as follows:

- * - *

We know f (t) at some sampling points ti and want to find a, - * ω and φ such that f (t) = a cos (ω t + φ). - *

- * - *

From the analytical expression, we can compute two primitives : - * 

-     *     If2  (t) = ∫ f2  = a2 × [t + S (t)] / 2
-     *     If'2 (t) = ∫ f'2 = a2 ω2 × [t - S (t)] / 2
-     *     where S (t) = sin (2 (ω t + φ)) / (2 ω)
-     *
- *

- * - *

We can remove S between these expressions : - * 

-     *     If'2 (t) = a2 ω2 t - ω2 If2 (t)
-     *
- *

- * - *

The preceding expression shows that If'2 (t) is a linear - * combination of both t and If2 (t): If'2 (t) = A × t + B × If2 (t) - *

- * - *

From the primitive, we can deduce the same form for definite - * integrals between t1 and ti for each ti : - * 

-     *   If2 (ti) - If2 (t1) = A × (ti - t1) + B × (If2 (ti) - If2 (t1))
-     *
- *

- * - *

We can find the coefficients A and B that best fit the sample - * to this linear expression by computing the definite integrals for - * each sample points. - *

- * - *

For a bilinear expression z (xi, yi) = A × xi + B × yi, the - * coefficients A and B that minimize a least square criterion - * ∑ (zi - z (xi, yi))2 are given by these expressions:

- * 
-     *
-     *         ∑yiyi ∑xizi - ∑xiyi ∑yizi
-     *     A = ------------------------
-     *         ∑xixi ∑yiyi - ∑xiyi ∑xiyi
-     *
-     *         ∑xixi ∑yizi - ∑xiyi ∑xizi
-     *     B = ------------------------
-     *         ∑xixi ∑yiyi - ∑xiyi ∑xiyi
-     *
- *

- * - * - *

In fact, we can assume both a and ω are positive and - * compute them directly, knowing that A = a2 ω2 and that - * B = - ω2. The complete algorithm is therefore:

- * 
-     *
-     * for each ti from t1 to tn-1, compute:
-     *   f  (ti)
-     *   f' (ti) = (f (ti+1) - f(ti-1)) / (ti+1 - ti-1)
-     *   xi = ti - t1
-     *   yi = ∫ f2 from t1 to ti
-     *   zi = ∫ f'2 from t1 to ti
-     *   update the sums ∑xixi, ∑yiyi, ∑xiyi, ∑xizi and ∑yizi
-     * end for
-     *
-     *            |--------------------------
-     *         \  | ∑yiyi ∑xizi - ∑xiyi ∑yizi
-     * a     =  \ | ------------------------
-     *           \| ∑xiyi ∑xizi - ∑xixi ∑yizi
-     *
-     *
-     *            |--------------------------
-     *         \  | ∑xiyi ∑xizi - ∑xixi ∑yizi
-     * ω     =  \ | ------------------------
-     *           \| ∑xixi ∑yiyi - ∑xiyi ∑xiyi
-     *
-     *
- *

- * - *

Once we know ω, we can compute: - * 

-     *    fc = ω f (t) cos (ω t) - f' (t) sin (ω t)
-     *    fs = ω f (t) sin (ω t) + f' (t) cos (ω t)
-     *
- *

- * - *

It appears that fc = a ω cos (φ) and - * fs = -a ω sin (φ), so we can use these - * expressions to compute φ. The best estimate over the sample is - * given by averaging these expressions. - *

- * - *

Since integrals and means are involved in the preceding - * estimations, these operations run in O(n) time, where n is the - * number of measurements.

- */ - public static class ParameterGuesser { - /** Amplitude. */ - private final double a; - /** Angular frequency. */ - private final double omega; - /** Phase. */ - private final double phi; - - /** - * Simple constructor. - * - * @param observations Sampled observations. - * @throws NumberIsTooSmallException if the sample is too short. - * @throws ZeroException if the abscissa range is zero. - * @throws MathIllegalStateException when the guessing procedure cannot - * produce sensible results. - */ - public ParameterGuesser(WeightedObservedPoint[] observations) { - if (observations.length < 4) { - throw new NumberIsTooSmallException(LocalizedFormats.INSUFFICIENT_OBSERVED_POINTS_IN_SAMPLE, - observations.length, 4, true); - } - - final WeightedObservedPoint[] sorted = sortObservations(observations); - - final double aOmega[] = guessAOmega(sorted); - a = aOmega[0]; - omega = aOmega[1]; - - phi = guessPhi(sorted); - } - - /** - * Gets an estimation of the parameters. - * - * @return the guessed parameters, in the following order: - *
    - *
  • Amplitude
    • - *
  • Angular frequency
    • - *
  • Phase
    • - *
- */ - public double[] guess() { - return new double[] { a, omega, phi }; - } - - /** - * Sort the observations with respect to the abscissa. - * - * @param unsorted Input observations. - * @return the input observations, sorted. - */ - private WeightedObservedPoint[] sortObservations(WeightedObservedPoint[] unsorted) { - final WeightedObservedPoint[] observations = unsorted.clone(); - - // Since the samples are almost always already sorted, this - // method is implemented as an insertion sort that reorders the - // elements in place. Insertion sort is very efficient in this case. - WeightedObservedPoint curr = observations[0]; - for (int j = 1; j < observations.length; ++j) { - WeightedObservedPoint prec = curr; - curr = observations[j]; - if (curr.getX() < prec.getX()) { - // the current element should be inserted closer to the beginning - int i = j - 1; - WeightedObservedPoint mI = observations[i]; - while ((i >= 0) && (curr.getX() < mI.getX())) { - observations[i + 1] = mI; - if (i-- != 0) { - mI = observations[i]; - } - } - observations[i + 1] = curr; - curr = observations[j]; - } - } - - return observations; - } - - /** - * Estimate a first guess of the amplitude and angular frequency. - * This method assumes that the {@link #sortObservations(WeightedObservedPoint[])} method - * has been called previously. - * - * @param observations Observations, sorted w.r.t. abscissa. - * @throws ZeroException if the abscissa range is zero. - * @throws MathIllegalStateException when the guessing procedure cannot - * produce sensible results. - * @return the guessed amplitude (at index 0) and circular frequency - * (at index 1). - */ - private double[] guessAOmega(WeightedObservedPoint[] observations) { - final double[] aOmega = new double[2]; - - // initialize the sums for the linear model between the two integrals - double sx2 = 0; - double sy2 = 0; - double sxy = 0; - double sxz = 0; - double syz = 0; - - double currentX = observations[0].getX(); - double currentY = observations[0].getY(); - double f2Integral = 0; - double fPrime2Integral = 0; - final double startX = currentX; - for (int i = 1; i < observations.length; ++i) { - // one step forward - final double previousX = currentX; - final double previousY = currentY; - currentX = observations[i].getX(); - currentY = observations[i].getY(); - - // update the integrals of f2 and f'2 - // considering a linear model for f (and therefore constant f') - final double dx = currentX - previousX; - final double dy = currentY - previousY; - final double f2StepIntegral = - dx * (previousY * previousY + previousY * currentY + currentY * currentY) / 3; - final double fPrime2StepIntegral = dy * dy / dx; - - final double x = currentX - startX; - f2Integral += f2StepIntegral; - fPrime2Integral += fPrime2StepIntegral; - - sx2 += x * x; - sy2 += f2Integral * f2Integral; - sxy += x * f2Integral; - sxz += x * fPrime2Integral; - syz += f2Integral * fPrime2Integral; - } - - // compute the amplitude and pulsation coefficients - double c1 = sy2 * sxz - sxy * syz; - double c2 = sxy * sxz - sx2 * syz; - double c3 = sx2 * sy2 - sxy * sxy; - if ((c1 / c2 < 0) || (c2 / c3 < 0)) { - final int last = observations.length - 1; - // Range of the observations, assuming that the - // observations are sorted. - final double xRange = observations[last].getX() - observations[0].getX(); - if (xRange == 0) { - throw new ZeroException(); - } - aOmega[1] = 2 * Math.PI / xRange; - - double yMin = Double.POSITIVE_INFINITY; - double yMax = Double.NEGATIVE_INFINITY; - for (int i = 1; i < observations.length; ++i) { - final double y = observations[i].getY(); - if (y < yMin) { - yMin = y; - } - if (y > yMax) { - yMax = y; - } - } - aOmega[0] = 0.5 * (yMax - yMin); - } else { - if (c2 == 0) { - // In some ill-conditioned cases (cf. MATH-844), the guesser - // procedure cannot produce sensible results. - throw new MathIllegalStateException(LocalizedFormats.ZERO_DENOMINATOR); - } - - aOmega[0] = FastMath.sqrt(c1 / c2); - aOmega[1] = FastMath.sqrt(c2 / c3); - } - - return aOmega; - } - - /** - * Estimate a first guess of the phase. - * - * @param observations Observations, sorted w.r.t. abscissa. - * @return the guessed phase. - */ - private double guessPhi(WeightedObservedPoint[] observations) { - // initialize the means - double fcMean = 0; - double fsMean = 0; - - double currentX = observations[0].getX(); - double currentY = observations[0].getY(); - for (int i = 1; i < observations.length; ++i) { - // one step forward - final double previousX = currentX; - final double previousY = currentY; - currentX = observations[i].getX(); - currentY = observations[i].getY(); - final double currentYPrime = (currentY - previousY) / (currentX - previousX); - - double omegaX = omega * currentX; - double cosine = FastMath.cos(omegaX); - double sine = FastMath.sin(omegaX); - fcMean += omega * currentY * cosine - currentYPrime * sine; - fsMean += omega * currentY * sine + currentYPrime * cosine; - } - - return FastMath.atan2(-fsMean, fcMean); - } - } -} diff --git a/src/main/java/org/apache/commons/math4/fitting/PolynomialFitter.java b/src/main/java/org/apache/commons/math4/fitting/PolynomialFitter.java deleted file mode 100644 index 38ebe9108..000000000 --- a/src/main/java/org/apache/commons/math4/fitting/PolynomialFitter.java +++ /dev/null @@ -1,72 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.fitting; - -import org.apache.commons.math4.analysis.polynomials.PolynomialFunction; -import org.apache.commons.math4.optim.nonlinear.vector.MultivariateVectorOptimizer; - -/** - * Polynomial fitting is a very simple case of {@link CurveFitter curve fitting}. - * The estimated coefficients are the polynomial coefficients (see the - * {@link #fit(double[]) fit} method). - * - * @since 2.0 - * @deprecated As of 3.3. Please use {@link PolynomialCurveFitter} and - * {@link WeightedObservedPoints} instead. - */ -@Deprecated -public class PolynomialFitter extends CurveFitter { - /** - * Simple constructor. - * - * @param optimizer Optimizer to use for the fitting. - */ - public PolynomialFitter(MultivariateVectorOptimizer optimizer) { - super(optimizer); - } - - /** - * Get the coefficients of the polynomial fitting the weighted data points. - * The degree of the fitting polynomial is {@code guess.length - 1}. - * - * @param guess First guess for the coefficients. They must be sorted in - * increasing order of the polynomial's degree. - * @param maxEval Maximum number of evaluations of the polynomial. - * @return the coefficients of the polynomial that best fits the observed points. - * @throws org.apache.commons.math4.exception.TooManyEvaluationsException if - * the number of evaluations exceeds {@code maxEval}. - * @throws org.apache.commons.math4.exception.ConvergenceException - * if the algorithm failed to converge. - */ - public double[] fit(int maxEval, double[] guess) { - return fit(maxEval, new PolynomialFunction.Parametric(), guess); - } - - /** - * Get the coefficients of the polynomial fitting the weighted data points. - * The degree of the fitting polynomial is {@code guess.length - 1}. - * - * @param guess First guess for the coefficients. They must be sorted in - * increasing order of the polynomial's degree. - * @return the coefficients of the polynomial that best fits the observed points. - * @throws org.apache.commons.math4.exception.ConvergenceException - * if the algorithm failed to converge. - */ - public double[] fit(double[] guess) { - return fit(new PolynomialFunction.Parametric(), guess); - } -} diff --git a/src/test/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatingFunctionTest.java b/src/test/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatingFunctionTest.java deleted file mode 100644 index 92bf82b50..000000000 --- a/src/test/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatingFunctionTest.java +++ /dev/null @@ -1,670 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.analysis.BivariateFunction; -import org.apache.commons.math4.analysis.interpolation.BicubicSplineFunction; -import org.apache.commons.math4.analysis.interpolation.BicubicSplineInterpolatingFunction; -import org.apache.commons.math4.distribution.UniformRealDistribution; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.MathIllegalArgumentException; -import org.apache.commons.math4.exception.OutOfRangeException; -import org.apache.commons.math4.random.RandomGenerator; -import org.apache.commons.math4.random.Well19937c; -import org.junit.Assert; -import org.junit.Test; -import org.junit.Ignore; - -/** - * Test case for the bicubic function. - * - */ -public final class BicubicSplineInterpolatingFunctionTest { - /** - * Test preconditions. - */ - @Test - public void testPreconditions() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2.5}; - double[][] zval = new double[xval.length][yval.length]; - - @SuppressWarnings("unused") - BivariateFunction bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, - zval, zval, zval); - - double[] wxval = new double[] {3, 2, 5, 6.5}; - try { - bcf = new BicubicSplineInterpolatingFunction(wxval, yval, zval, zval, zval, zval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - double[] wyval = new double[] {-4, -1, -1, 2.5}; - try { - bcf = new BicubicSplineInterpolatingFunction(xval, wyval, zval, zval, zval, zval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - double[][] wzval = new double[xval.length][yval.length - 1]; - try { - bcf = new BicubicSplineInterpolatingFunction(xval, yval, wzval, zval, zval, zval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, wzval, zval, zval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, zval, wzval, zval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, zval, zval, wzval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - - wzval = new double[xval.length - 1][yval.length]; - try { - bcf = new BicubicSplineInterpolatingFunction(xval, yval, wzval, zval, zval, zval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, wzval, zval, zval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, zval, wzval, zval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, zval, zval, wzval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - } - - /** - * Test for a plane. - *

- * z = 2 x - 3 y + 5 - */ - @Ignore@Test - public void testPlane() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2, 2.5}; - // Function values - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - return 2 * x - 3 * y + 5; - } - }; - double[][] zval = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - zval[i][j] = f.value(xval[i], yval[j]); - } - } - // Partial derivatives with respect to x - double[][] dZdX = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdX[i][j] = 2; - } - } - // Partial derivatives with respect to y - double[][] dZdY = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdY[i][j] = -3; - } - } - // Partial cross-derivatives - double[][] dZdXdY = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdXdY[i][j] = 0; - } - } - - BivariateFunction bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, - dZdX, dZdY, dZdXdY); - double x, y; - double expected, result; - - x = 4; - y = -3; - expected = f.value(x, y); - result = bcf.value(x, y); - Assert.assertEquals("On sample point", - expected, result, 1e-15); - - x = 4.5; - y = -1.5; - expected = f.value(x, y); - result = bcf.value(x, y); - Assert.assertEquals("Half-way between sample points (middle of the patch)", - expected, result, 0.3); - - x = 3.5; - y = -3.5; - expected = f.value(x, y); - result = bcf.value(x, y); - Assert.assertEquals("Half-way between sample points (border of the patch)", - expected, result, 0.3); - } - - /** - * Test for a paraboloid. - *

- * z = 2 x2 - 3 y2 + 4 x y - 5 - */ - @Ignore@Test - public void testParaboloid() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2, 2.5}; - // Function values - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - return 2 * x * x - 3 * y * y + 4 * x * y - 5; - } - }; - double[][] zval = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - zval[i][j] = f.value(xval[i], yval[j]); - } - } - // Partial derivatives with respect to x - double[][] dZdX = new double[xval.length][yval.length]; - BivariateFunction dfdX = new BivariateFunction() { - public double value(double x, double y) { - return 4 * (x + y); - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdX[i][j] = dfdX.value(xval[i], yval[j]); - } - } - // Partial derivatives with respect to y - double[][] dZdY = new double[xval.length][yval.length]; - BivariateFunction dfdY = new BivariateFunction() { - public double value(double x, double y) { - return 4 * x - 6 * y; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdY[i][j] = dfdY.value(xval[i], yval[j]); - } - } - // Partial cross-derivatives - double[][] dZdXdY = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdXdY[i][j] = 4; - } - } - - BivariateFunction bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, - dZdX, dZdY, dZdXdY); - double x, y; - double expected, result; - - x = 4; - y = -3; - expected = f.value(x, y); - result = bcf.value(x, y); - Assert.assertEquals("On sample point", - expected, result, 1e-15); - - x = 4.5; - y = -1.5; - expected = f.value(x, y); - result = bcf.value(x, y); - Assert.assertEquals("Half-way between sample points (middle of the patch)", - expected, result, 2); - - x = 3.5; - y = -3.5; - expected = f.value(x, y); - result = bcf.value(x, y); - Assert.assertEquals("Half-way between sample points (border of the patch)", - expected, result, 2); - } - - /** - * Test for partial derivatives of {@link BicubicSplineFunction}. - *

- * f(x, y) = ΣiΣj (i+1) (j+2) xi yj - */ - @Ignore@Test - public void testSplinePartialDerivatives() { - final int N = 4; - final double[] coeff = new double[16]; - - for (int i = 0; i < N; i++) { - for (int j = 0; j < N; j++) { - coeff[i + N * j] = (i + 1) * (j + 2); - } - } - - final BicubicSplineFunction f = new BicubicSplineFunction(coeff); - BivariateFunction derivative; - final double x = 0.435; - final double y = 0.776; - final double tol = 1e-13; - - derivative = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double y2 = y * y; - final double y3 = y2 * y; - final double yFactor = 2 + 3 * y + 4 * y2 + 5 * y3; - return yFactor * (2 + 6 * x + 12 * x2); - } - }; - Assert.assertEquals("dFdX", derivative.value(x, y), - f.partialDerivativeX().value(x, y), tol); - - derivative = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double x3 = x2 * x; - final double y2 = y * y; - final double xFactor = 1 + 2 * x + 3 * x2 + 4 * x3; - return xFactor * (3 + 8 * y + 15 * y2); - } - }; - Assert.assertEquals("dFdY", derivative.value(x, y), - f.partialDerivativeY().value(x, y), tol); - - derivative = new BivariateFunction() { - public double value(double x, double y) { - final double y2 = y * y; - final double y3 = y2 * y; - final double yFactor = 2 + 3 * y + 4 * y2 + 5 * y3; - return yFactor * (6 + 24 * x); - } - }; - Assert.assertEquals("d2FdX2", derivative.value(x, y), - f.partialDerivativeXX().value(x, y), tol); - - derivative = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double x3 = x2 * x; - final double xFactor = 1 + 2 * x + 3 * x2 + 4 * x3; - return xFactor * (8 + 30 * y); - } - }; - Assert.assertEquals("d2FdY2", derivative.value(x, y), - f.partialDerivativeYY().value(x, y), tol); - - derivative = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double y2 = y * y; - final double yFactor = 3 + 8 * y + 15 * y2; - return yFactor * (2 + 6 * x + 12 * x2); - } - }; - Assert.assertEquals("d2FdXdY", derivative.value(x, y), - f.partialDerivativeXY().value(x, y), tol); - } - - /** - * Test that the partial derivatives computed from a - * {@link BicubicSplineInterpolatingFunction} match the input data. - *

- * f(x, y) = 5 - * - 3 x + 2 y - * - x y + 2 x2 - 3 y2 - * + 4 x2 y - x y2 - 3 x3 + y3 - */ - @Ignore@Test - public void testMatchingPartialDerivatives() { - final int sz = 21; - double[] val = new double[sz]; - // Coordinate values - final double delta = 1d / (sz - 1); - for (int i = 0; i < sz; i++) { - val[i] = i * delta; - } - // Function values - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double x3 = x2 * x; - final double y2 = y * y; - final double y3 = y2 * y; - - return 5 - - 3 * x + 2 * y - - x * y + 2 * x2 - 3 * y2 - + 4 * x2 * y - x * y2 - 3 * x3 + y3; - } - }; - double[][] fval = new double[sz][sz]; - for (int i = 0; i < sz; i++) { - for (int j = 0; j < sz; j++) { - fval[i][j] = f.value(val[i], val[j]); - } - } - // Partial derivatives with respect to x - double[][] dFdX = new double[sz][sz]; - BivariateFunction dfdX = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double y2 = y * y; - return - 3 - y + 4 * x + 8 * x * y - y2 - 9 * x2; - } - }; - for (int i = 0; i < sz; i++) { - for (int j = 0; j < sz; j++) { - dFdX[i][j] = dfdX.value(val[i], val[j]); - } - } - // Partial derivatives with respect to y - double[][] dFdY = new double[sz][sz]; - BivariateFunction dfdY = new BivariateFunction() { - public double value(double x, double y) { - final double x2 = x * x; - final double y2 = y * y; - return 2 - x - 6 * y + 4 * x2 - 2 * x * y + 3 * y2; - } - }; - for (int i = 0; i < sz; i++) { - for (int j = 0; j < sz; j++) { - dFdY[i][j] = dfdY.value(val[i], val[j]); - } - } - // Partial cross-derivatives - double[][] d2FdXdY = new double[sz][sz]; - BivariateFunction d2fdXdY = new BivariateFunction() { - public double value(double x, double y) { - return -1 + 8 * x - 2 * y; - } - }; - for (int i = 0; i < sz; i++) { - for (int j = 0; j < sz; j++) { - d2FdXdY[i][j] = d2fdXdY.value(val[i], val[j]); - } - } - - BicubicSplineInterpolatingFunction bcf - = new BicubicSplineInterpolatingFunction(val, val, fval, dFdX, dFdY, d2FdXdY); - - double x, y; - double expected, result; - - final double tol = 1e-12; - for (int i = 0; i < sz; i++) { - x = val[i]; - for (int j = 0; j < sz; j++) { - y = val[j]; - - expected = dfdX.value(x, y); - result = bcf.partialDerivativeX(x, y); - Assert.assertEquals(x + " " + y + " dFdX", expected, result, tol); - - expected = dfdY.value(x, y); - result = bcf.partialDerivativeY(x, y); - Assert.assertEquals(x + " " + y + " dFdY", expected, result, tol); - - expected = d2fdXdY.value(x, y); - result = bcf.partialDerivativeXY(x, y); - Assert.assertEquals(x + " " + y + " d2FdXdY", expected, result, tol); - } - } - } - - /** - * Interpolating a plane. - *

- * z = 2 x - 3 y + 5 - */ - @Test - public void testInterpolation1() { - final int sz = 21; - double[] xval = new double[sz]; - double[] yval = new double[sz]; - // Coordinate values - final double delta = 1d / (sz - 1); - for (int i = 0; i < sz; i++) { - xval[i] = -1 + 15 * i * delta; - yval[i] = -20 + 30 * i * delta; - } - - // Function values - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - return 2 * x - 3 * y + 5; - } - }; - double[][] zval = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - zval[i][j] = f.value(xval[i], yval[j]); - } - } - // Partial derivatives with respect to x - double[][] dZdX = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdX[i][j] = 2; - } - } - // Partial derivatives with respect to y - double[][] dZdY = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdY[i][j] = -3; - } - } - // Partial cross-derivatives - double[][] dZdXdY = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdXdY[i][j] = 0; - } - } - - final BivariateFunction bcf - = new BicubicSplineInterpolatingFunction(xval, yval, zval, - dZdX, dZdY, dZdXdY); - double x, y; - - final RandomGenerator rng = new Well19937c(1234567L); // "tol" depends on the seed. - final UniformRealDistribution distX - = new UniformRealDistribution(rng, xval[0], xval[xval.length - 1]); - final UniformRealDistribution distY - = new UniformRealDistribution(rng, yval[0], yval[yval.length - 1]); - - final int numSamples = 50; - final double tol = 6; - for (int i = 0; i < numSamples; i++) { - x = distX.sample(); - for (int j = 0; j < numSamples; j++) { - y = distY.sample(); -// System.out.println(x + " " + y + " " + f.value(x, y) + " " + bcf.value(x, y)); - Assert.assertEquals(f.value(x, y), bcf.value(x, y), tol); - } -// System.out.println(); - } - } - - /** - * Interpolating a paraboloid. - *

- * z = 2 x2 - 3 y2 + 4 x y - 5 - */ - @Test - public void testInterpolation2() { - final int sz = 21; - double[] xval = new double[sz]; - double[] yval = new double[sz]; - // Coordinate values - final double delta = 1d / (sz - 1); - for (int i = 0; i < sz; i++) { - xval[i] = -1 + 15 * i * delta; - yval[i] = -20 + 30 * i * delta; - } - - // Function values - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - return 2 * x * x - 3 * y * y + 4 * x * y - 5; - } - }; - double[][] zval = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - zval[i][j] = f.value(xval[i], yval[j]); - } - } - // Partial derivatives with respect to x - double[][] dZdX = new double[xval.length][yval.length]; - BivariateFunction dfdX = new BivariateFunction() { - public double value(double x, double y) { - return 4 * (x + y); - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdX[i][j] = dfdX.value(xval[i], yval[j]); - } - } - // Partial derivatives with respect to y - double[][] dZdY = new double[xval.length][yval.length]; - BivariateFunction dfdY = new BivariateFunction() { - public double value(double x, double y) { - return 4 * x - 6 * y; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdY[i][j] = dfdY.value(xval[i], yval[j]); - } - } - // Partial cross-derivatives - double[][] dZdXdY = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - dZdXdY[i][j] = 4; - } - } - - BivariateFunction bcf = new BicubicSplineInterpolatingFunction(xval, yval, zval, - dZdX, dZdY, dZdXdY); - double x, y; - - final RandomGenerator rng = new Well19937c(1234567L); // "tol" depends on the seed. - final UniformRealDistribution distX - = new UniformRealDistribution(rng, xval[0], xval[xval.length - 1]); - final UniformRealDistribution distY - = new UniformRealDistribution(rng, yval[0], yval[yval.length - 1]); - - final double tol = 224; - for (int i = 0; i < sz; i++) { - x = distX.sample(); - for (int j = 0; j < sz; j++) { - y = distY.sample(); -// System.out.println(x + " " + y + " " + f.value(x, y) + " " + bcf.value(x, y)); - Assert.assertEquals(f.value(x, y), bcf.value(x, y), tol); - } -// System.out.println(); - } - } - - @Test - public void testIsValidPoint() { - final double xMin = -12; - final double xMax = 34; - final double yMin = 5; - final double yMax = 67; - final double[] xval = new double[] { xMin, xMax }; - final double[] yval = new double[] { yMin, yMax }; - final double[][] f = new double[][] { { 1, 2 }, - { 3, 4 } }; - final double[][] dFdX = f; - final double[][] dFdY = f; - final double[][] dFdXdY = f; - - final BicubicSplineInterpolatingFunction bcf - = new BicubicSplineInterpolatingFunction(xval, yval, f, - dFdX, dFdY, dFdXdY); - - double x, y; - - x = xMin; - y = yMin; - Assert.assertTrue(bcf.isValidPoint(x, y)); - // Ensure that no exception is thrown. - bcf.value(x, y); - - x = xMax; - y = yMax; - Assert.assertTrue(bcf.isValidPoint(x, y)); - // Ensure that no exception is thrown. - bcf.value(x, y); - - final double xRange = xMax - xMin; - final double yRange = yMax - yMin; - x = xMin + xRange / 3.4; - y = yMin + yRange / 1.2; - Assert.assertTrue(bcf.isValidPoint(x, y)); - // Ensure that no exception is thrown. - bcf.value(x, y); - - final double small = 1e-8; - x = xMin - small; - y = yMax; - Assert.assertFalse(bcf.isValidPoint(x, y)); - // Ensure that an exception would have been thrown. - try { - bcf.value(x, y); - Assert.fail("OutOfRangeException expected"); - } catch (OutOfRangeException expected) {} - - x = xMin; - y = yMax + small; - Assert.assertFalse(bcf.isValidPoint(x, y)); - // Ensure that an exception would have been thrown. - try { - bcf.value(x, y); - Assert.fail("OutOfRangeException expected"); - } catch (OutOfRangeException expected) {} - } -} diff --git a/src/test/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatorTest.java b/src/test/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatorTest.java deleted file mode 100644 index 91f1f664b..000000000 --- a/src/test/java/org/apache/commons/math4/analysis/interpolation/BicubicSplineInterpolatorTest.java +++ /dev/null @@ -1,186 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.analysis.BivariateFunction; -import org.apache.commons.math4.analysis.interpolation.BicubicSplineInterpolator; -import org.apache.commons.math4.analysis.interpolation.BivariateGridInterpolator; -import org.apache.commons.math4.distribution.UniformRealDistribution; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.MathIllegalArgumentException; -import org.apache.commons.math4.random.RandomGenerator; -import org.apache.commons.math4.random.Well19937c; -import org.junit.Assert; -import org.junit.Test; - -/** - * Test case for the bicubic interpolator. - * - */ -public final class BicubicSplineInterpolatorTest { - /** - * Test preconditions. - */ - @Test - public void testPreconditions() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2.5}; - double[][] zval = new double[xval.length][yval.length]; - - BivariateGridInterpolator interpolator = new BicubicSplineInterpolator(); - - @SuppressWarnings("unused") - BivariateFunction p = interpolator.interpolate(xval, yval, zval); - - double[] wxval = new double[] {3, 2, 5, 6.5}; - try { - p = interpolator.interpolate(wxval, yval, zval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - - double[] wyval = new double[] {-4, -3, -1, -1}; - try { - p = interpolator.interpolate(xval, wyval, zval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - - double[][] wzval = new double[xval.length][yval.length + 1]; - try { - p = interpolator.interpolate(xval, yval, wzval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - wzval = new double[xval.length - 1][yval.length]; - try { - p = interpolator.interpolate(xval, yval, wzval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - } - - /** - * Interpolating a plane. - *

- * z = 2 x - 3 y + 5 - */ - @Test - public void testInterpolation1() { - final int sz = 21; - double[] xval = new double[sz]; - double[] yval = new double[sz]; - // Coordinate values - final double delta = 1d / (sz - 1); - for (int i = 0; i < sz; i++) { - xval[i] = -1 + 15 * i * delta; - yval[i] = -20 + 30 * i * delta; - } - - // Function values - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - return 2 * x - 3 * y + 5; - } - }; - double[][] zval = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - zval[i][j] = f.value(xval[i], yval[j]); - } - } - - BivariateGridInterpolator interpolator = new BicubicSplineInterpolator(); - BivariateFunction p = interpolator.interpolate(xval, yval, zval); - double x, y; - - final RandomGenerator rng = new Well19937c(1234567L); // "tol" depends on the seed. - final UniformRealDistribution distX - = new UniformRealDistribution(rng, xval[0], xval[xval.length - 1]); - final UniformRealDistribution distY - = new UniformRealDistribution(rng, yval[0], yval[yval.length - 1]); - - final int numSamples = 50; - final double tol = 6; - for (int i = 0; i < numSamples; i++) { - x = distX.sample(); - for (int j = 0; j < numSamples; j++) { - y = distY.sample(); -// System.out.println(x + " " + y + " " + f.value(x, y) + " " + p.value(x, y)); - Assert.assertEquals(f.value(x, y), p.value(x, y), tol); - } -// System.out.println(); - } - } - - /** - * Interpolating a paraboloid. - *

- * z = 2 x2 - 3 y2 + 4 x y - 5 - */ - @Test - public void testInterpolation2() { - final int sz = 21; - double[] xval = new double[sz]; - double[] yval = new double[sz]; - // Coordinate values - final double delta = 1d / (sz - 1); - for (int i = 0; i < sz; i++) { - xval[i] = -1 + 15 * i * delta; - yval[i] = -20 + 30 * i * delta; - } - - // Function values - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - return 2 * x * x - 3 * y * y + 4 * x * y - 5; - } - }; - double[][] zval = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - zval[i][j] = f.value(xval[i], yval[j]); - } - } - - BivariateGridInterpolator interpolator = new BicubicSplineInterpolator(); - BivariateFunction p = interpolator.interpolate(xval, yval, zval); - double x, y; - - final RandomGenerator rng = new Well19937c(1234567L); // "tol" depends on the seed. - final UniformRealDistribution distX - = new UniformRealDistribution(rng, xval[0], xval[xval.length - 1]); - final UniformRealDistribution distY - = new UniformRealDistribution(rng, yval[0], yval[yval.length - 1]); - - final int numSamples = 50; - final double tol = 251; - for (int i = 0; i < numSamples; i++) { - x = distX.sample(); - for (int j = 0; j < numSamples; j++) { - y = distY.sample(); -// System.out.println(x + " " + y + " " + f.value(x, y) + " " + p.value(x, y)); - Assert.assertEquals(f.value(x, y), p.value(x, y), tol); - } -// System.out.println(); - } - } -} diff --git a/src/test/java/org/apache/commons/math4/analysis/interpolation/SmoothingPolynomialBicubicSplineInterpolatorTest.java b/src/test/java/org/apache/commons/math4/analysis/interpolation/SmoothingPolynomialBicubicSplineInterpolatorTest.java deleted file mode 100644 index 9cdd88880..000000000 --- a/src/test/java/org/apache/commons/math4/analysis/interpolation/SmoothingPolynomialBicubicSplineInterpolatorTest.java +++ /dev/null @@ -1,181 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.analysis.BivariateFunction; -import org.apache.commons.math4.analysis.interpolation.BivariateGridInterpolator; -import org.apache.commons.math4.analysis.interpolation.SmoothingPolynomialBicubicSplineInterpolator; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.MathIllegalArgumentException; -import org.apache.commons.math4.util.FastMath; -import org.junit.Assert; -import org.junit.Test; - -/** - * Test case for the smoothing bicubic interpolator. - * - */ -public final class SmoothingPolynomialBicubicSplineInterpolatorTest { - /** - * Test preconditions. - */ - @Test - public void testPreconditions() { - double[] xval = new double[] {3, 4, 5, 6.5, 7.5}; - double[] yval = new double[] {-4, -3, -1, 2.5, 3}; - double[][] zval = new double[xval.length][yval.length]; - - BivariateGridInterpolator interpolator = new SmoothingPolynomialBicubicSplineInterpolator(0); - - @SuppressWarnings("unused") - BivariateFunction p = interpolator.interpolate(xval, yval, zval); - - double[] wxval = new double[] {3, 2, 5, 6.5}; - try { - p = interpolator.interpolate(wxval, yval, zval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - - double[] wyval = new double[] {-4, -3, -1, -1}; - try { - p = interpolator.interpolate(xval, wyval, zval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - - double[][] wzval = new double[xval.length][yval.length + 1]; - try { - p = interpolator.interpolate(xval, yval, wzval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - wzval = new double[xval.length - 1][yval.length]; - try { - p = interpolator.interpolate(xval, yval, wzval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - wzval = new double[xval.length][yval.length - 1]; - try { - p = interpolator.interpolate(xval, yval, wzval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - } - - /** - * Test of interpolator for a plane. - *

- * z = 2 x - 3 y + 5 - */ - @Test - public void testPlane() { - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - return 2 * x - 3 * y + 5 - + ((int) (FastMath.abs(5 * x + 3 * y)) % 2 == 0 ? 1 : -1); - } - }; - - BivariateGridInterpolator interpolator = new SmoothingPolynomialBicubicSplineInterpolator(1); - - double[] xval = new double[] {3, 4, 5, 6.5, 7.5}; - double[] yval = new double[] {-4, -3, -1, 2, 2.5, 3.5}; - double[][] zval = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - zval[i][j] = f.value(xval[i], yval[j]); - } - } - - BivariateFunction p = interpolator.interpolate(xval, yval, zval); - double x, y; - double expected, result; - - x = 4; - y = -3; - expected = f.value(x, y); - result = p.value(x, y); - Assert.assertEquals("On sample point", expected, result, 2); - - x = 4.5; - y = -1.5; - expected = f.value(x, y); - result = p.value(x, y); - Assert.assertEquals("half-way between sample points (middle of the patch)", expected, result, 2); - - x = 3.5; - y = -3.5; - expected = f.value(x, y); - result = p.value(x, y); - Assert.assertEquals("half-way between sample points (border of the patch)", expected, result, 2); - } - - /** - * Test of interpolator for a paraboloid. - *

- * z = 2 x2 - 3 y2 + 4 x y - 5 - */ - @Test - public void testParaboloid() { - BivariateFunction f = new BivariateFunction() { - public double value(double x, double y) { - return 2 * x * x - 3 * y * y + 4 * x * y - 5 - + ((int) (FastMath.abs(5 * x + 3 * y)) % 2 == 0 ? 1 : -1); - } - }; - - BivariateGridInterpolator interpolator = new SmoothingPolynomialBicubicSplineInterpolator(4); - - double[] xval = new double[] {3, 4, 5, 6.5, 7.5, 8}; - double[] yval = new double[] {-4, -3, -2, -1, 0.5, 2.5}; - double[][] zval = new double[xval.length][yval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - zval[i][j] = f.value(xval[i], yval[j]); - } - } - - BivariateFunction p = interpolator.interpolate(xval, yval, zval); - double x, y; - double expected, result; - - x = 5; - y = 0.5; - expected = f.value(x, y); - result = p.value(x, y); - Assert.assertEquals("On sample point", expected, result, 2); - - x = 4.5; - y = -1.5; - expected = f.value(x, y); - result = p.value(x, y); - Assert.assertEquals("half-way between sample points (middle of the patch)", expected, result, 2); - - x = 3.5; - y = -3.5; - expected = f.value(x, y); - result = p.value(x, y); - Assert.assertEquals("half-way between sample points (border of the patch)", expected, result, 2); - } -} diff --git a/src/test/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatingFunctionTest.java b/src/test/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatingFunctionTest.java deleted file mode 100644 index 945e9d510..000000000 --- a/src/test/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatingFunctionTest.java +++ /dev/null @@ -1,545 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.analysis.TrivariateFunction; -import org.apache.commons.math4.analysis.interpolation.TricubicSplineInterpolatingFunction; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.MathIllegalArgumentException; -import org.apache.commons.math4.util.FastMath; -import org.junit.Assert; -import org.junit.Test; - -/** - * Test case for the bicubic function. - * - */ -public final class TricubicSplineInterpolatingFunctionTest { - /** - * Test preconditions. - */ - @Test - public void testPreconditions() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2.5}; - double[] zval = new double[] {-12, -8, -5.5, -3, 0, 2.5}; - double[][][] fval = new double[xval.length][yval.length][zval.length]; - - @SuppressWarnings("unused") - TrivariateFunction tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, fval, fval, fval); - - double[] wxval = new double[] {3, 2, 5, 6.5}; - try { - tcf = new TricubicSplineInterpolatingFunction(wxval, yval, zval, - fval, fval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - double[] wyval = new double[] {-4, -1, -1, 2.5}; - try { - tcf = new TricubicSplineInterpolatingFunction(xval, wyval, zval, - fval, fval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - double[] wzval = new double[] {-12, -8, -9, -3, 0, 2.5}; - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, wzval, - fval, fval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - double[][][] wfval = new double[xval.length - 1][yval.length - 1][zval.length]; - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - wfval, fval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, wfval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, wfval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, wfval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - wfval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, wfval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, fval, wfval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, fval, fval, wfval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - wfval = new double[xval.length][yval.length - 1][zval.length]; - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - wfval, fval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, wfval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, wfval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, wfval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - wfval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, wfval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, fval, wfval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, fval, fval, wfval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - wfval = new double[xval.length][yval.length][zval.length - 1]; - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - wfval, fval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, wfval, fval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, wfval, fval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, wfval, - fval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - wfval, fval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, wfval, fval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, fval, wfval, fval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - try { - tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, fval, fval, fval, - fval, fval, fval, wfval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - } - - /** - * Test for a plane. - *

- * f(x, y, z) = 2 x - 3 y - 4 z + 5 - *

- */ - @Test - public void testPlane() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2, 2.5}; - double[] zval = new double[] {-12, -8, -5.5, -3, 0, 2.5}; - - // Function values - TrivariateFunction f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return 2 * x - 3 * y - 4 * z + 5; - } - }; - - double[][][] fval = new double[xval.length][yval.length][zval.length]; - - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - fval[i][j][k] = f.value(xval[i], yval[j], zval[k]); - } - } - } - // Partial derivatives with respect to x - double[][][] dFdX = new double[xval.length][yval.length][zval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - dFdX[i][j][k] = 2; - } - } - } - // Partial derivatives with respect to y - double[][][] dFdY = new double[xval.length][yval.length][zval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - dFdY[i][j][k] = -3; - } - } - } - - // Partial derivatives with respect to z - double[][][] dFdZ = new double[xval.length][yval.length][zval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - dFdZ[i][j][k] = -4; - } - } - } - // Partial cross-derivatives - double[][][] d2FdXdY = new double[xval.length][yval.length][zval.length]; - double[][][] d2FdXdZ = new double[xval.length][yval.length][zval.length]; - double[][][] d2FdYdZ = new double[xval.length][yval.length][zval.length]; - double[][][] d3FdXdYdZ = new double[xval.length][yval.length][zval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - d2FdXdY[i][j][k] = 0; - d2FdXdZ[i][j][k] = 0; - d2FdYdZ[i][j][k] = 0; - d3FdXdYdZ[i][j][k] = 0; - } - } - } - - TrivariateFunction tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, dFdX, dFdY, dFdZ, - d2FdXdY, d2FdXdZ, d2FdYdZ, - d3FdXdYdZ); - double x, y, z; - double expected, result; - - x = 4; - y = -3; - z = 0; - expected = f.value(x, y, z); - result = tcf.value(x, y, z); - Assert.assertEquals("On sample point", - expected, result, 1e-15); - - x = 4.5; - y = -1.5; - z = -4.25; - expected = f.value(x, y, z); - result = tcf.value(x, y, z); - Assert.assertEquals("Half-way between sample points (middle of the patch)", - expected, result, 0.3); - - x = 3.5; - y = -3.5; - z = -10; - expected = f.value(x, y, z); - result = tcf.value(x, y, z); - Assert.assertEquals("Half-way between sample points (border of the patch)", - expected, result, 0.3); - } - - /** - * Sine wave. - *

- * f(x, y, z) = a cos [ω z - ky x - ky y] - *

- * with A = 0.2, ω = 0.5, kx = 2, ky = 1. - */ - @Test - public void testWave() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2, 2.5}; - double[] zval = new double[] {-12, -8, -5.5, -3, 0, 4}; - - final double a = 0.2; - final double omega = 0.5; - final double kx = 2; - final double ky = 1; - - // Function values - TrivariateFunction f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return a * FastMath.cos(omega * z - kx * x - ky * y); - } - }; - - double[][][] fval = new double[xval.length][yval.length][zval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - fval[i][j][k] = f.value(xval[i], yval[j], zval[k]); - } - } - } - - // Partial derivatives with respect to x - double[][][] dFdX = new double[xval.length][yval.length][zval.length]; - TrivariateFunction dFdX_f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return a * FastMath.sin(omega * z - kx * x - ky * y) * kx; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - dFdX[i][j][k] = dFdX_f.value(xval[i], yval[j], zval[k]); - } - } - } - - // Partial derivatives with respect to y - double[][][] dFdY = new double[xval.length][yval.length][zval.length]; - TrivariateFunction dFdY_f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return a * FastMath.sin(omega * z - kx * x - ky * y) * ky; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - dFdY[i][j][k] = dFdY_f.value(xval[i], yval[j], zval[k]); - } - } - } - - // Partial derivatives with respect to z - double[][][] dFdZ = new double[xval.length][yval.length][zval.length]; - TrivariateFunction dFdZ_f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return -a * FastMath.sin(omega * z - kx * x - ky * y) * omega; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - dFdZ[i][j][k] = dFdZ_f.value(xval[i], yval[j], zval[k]); - } - } - } - - // Partial second derivatives w.r.t. (x, y) - double[][][] d2FdXdY = new double[xval.length][yval.length][zval.length]; - TrivariateFunction d2FdXdY_f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return -a * FastMath.cos(omega * z - kx * x - ky * y) * kx * ky; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - d2FdXdY[i][j][k] = d2FdXdY_f.value(xval[i], yval[j], zval[k]); - } - } - } - - // Partial second derivatives w.r.t. (x, z) - double[][][] d2FdXdZ = new double[xval.length][yval.length][zval.length]; - TrivariateFunction d2FdXdZ_f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return a * FastMath.cos(omega * z - kx * x - ky * y) * kx * omega; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - d2FdXdZ[i][j][k] = d2FdXdZ_f.value(xval[i], yval[j], zval[k]); - } - } - } - - // Partial second derivatives w.r.t. (y, z) - double[][][] d2FdYdZ = new double[xval.length][yval.length][zval.length]; - TrivariateFunction d2FdYdZ_f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return a * FastMath.cos(omega * z - kx * x - ky * y) * ky * omega; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - d2FdYdZ[i][j][k] = d2FdYdZ_f.value(xval[i], yval[j], zval[k]); - } - } - } - - // Partial third derivatives - double[][][] d3FdXdYdZ = new double[xval.length][yval.length][zval.length]; - TrivariateFunction d3FdXdYdZ_f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return a * FastMath.sin(omega * z - kx * x - ky * y) * kx * ky * omega; - } - }; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - d3FdXdYdZ[i][j][k] = d3FdXdYdZ_f.value(xval[i], yval[j], zval[k]); - } - } - } - - TrivariateFunction tcf = new TricubicSplineInterpolatingFunction(xval, yval, zval, - fval, dFdX, dFdY, dFdZ, - d2FdXdY, d2FdXdZ, d2FdYdZ, - d3FdXdYdZ); - double x, y, z; - double expected, result; - - x = 4; - y = -3; - z = 0; - expected = f.value(x, y, z); - result = tcf.value(x, y, z); - Assert.assertEquals("On sample point", - expected, result, 1e-14); - - x = 4.5; - y = -1.5; - z = -4.25; - expected = f.value(x, y, z); - result = tcf.value(x, y, z); - Assert.assertEquals("Half-way between sample points (middle of the patch)", - expected, result, 0.1); - - x = 3.5; - y = -3.5; - z = -10; - expected = f.value(x, y, z); - result = tcf.value(x, y, z); - Assert.assertEquals("Half-way between sample points (border of the patch)", - expected, result, 0.1); - } -} diff --git a/src/test/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatorTest.java b/src/test/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatorTest.java deleted file mode 100644 index 83023d3ef..000000000 --- a/src/test/java/org/apache/commons/math4/analysis/interpolation/TricubicSplineInterpolatorTest.java +++ /dev/null @@ -1,214 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.interpolation; - -import org.apache.commons.math4.analysis.TrivariateFunction; -import org.apache.commons.math4.analysis.interpolation.TricubicSplineInterpolator; -import org.apache.commons.math4.analysis.interpolation.TrivariateGridInterpolator; -import org.apache.commons.math4.exception.DimensionMismatchException; -import org.apache.commons.math4.exception.MathIllegalArgumentException; -import org.apache.commons.math4.util.FastMath; -import org.junit.Assert; -import org.junit.Test; -import org.junit.Ignore; - -/** - * Test case for the tricubic interpolator. - * - */ -public final class TricubicSplineInterpolatorTest { - /** - * Test preconditions. - */ - @Test - public void testPreconditions() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2.5}; - double[] zval = new double[] {-12, -8, -5.5, -3, 0, 2.5}; - double[][][] fval = new double[xval.length][yval.length][zval.length]; - - TrivariateGridInterpolator interpolator = new TricubicSplineInterpolator(); - - @SuppressWarnings("unused") - TrivariateFunction p = interpolator.interpolate(xval, yval, zval, fval); - - double[] wxval = new double[] {3, 2, 5, 6.5}; - try { - p = interpolator.interpolate(wxval, yval, zval, fval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - - double[] wyval = new double[] {-4, -3, -1, -1}; - try { - p = interpolator.interpolate(xval, wyval, zval, fval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - - double[] wzval = new double[] {-12, -8, -5.5, -3, -4, 2.5}; - try { - p = interpolator.interpolate(xval, yval, wzval, fval); - Assert.fail("an exception should have been thrown"); - } catch (MathIllegalArgumentException e) { - // Expected - } - - double[][][] wfval = new double[xval.length][yval.length + 1][zval.length]; - try { - p = interpolator.interpolate(xval, yval, zval, wfval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - wfval = new double[xval.length - 1][yval.length][zval.length]; - try { - p = interpolator.interpolate(xval, yval, zval, wfval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - wfval = new double[xval.length][yval.length][zval.length - 1]; - try { - p = interpolator.interpolate(xval, yval, zval, wfval); - Assert.fail("an exception should have been thrown"); - } catch (DimensionMismatchException e) { - // Expected - } - } - - /** - * Test of interpolator for a plane. - *

- * f(x, y, z) = 2 x - 3 y - z + 5 - */ - @Ignore@Test - public void testPlane() { - TrivariateFunction f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return 2 * x - 3 * y - z + 5; - } - }; - - TrivariateGridInterpolator interpolator = new TricubicSplineInterpolator(); - - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2, 2.5}; - double[] zval = new double[] {-12, -8, -5.5, -3, 0, 2.5}; - double[][][] fval = new double[xval.length][yval.length][zval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - fval[i][j][k] = f.value(xval[i], yval[j], zval[k]); - } - } - } - - TrivariateFunction p = interpolator.interpolate(xval, yval, zval, fval); - double x, y, z; - double expected, result; - - x = 4; - y = -3; - z = 0; - expected = f.value(x, y, z); - result = p.value(x, y, z); - Assert.assertEquals("On sample point", expected, result, 1e-15); - - x = 4.5; - y = -1.5; - z = -4.25; - expected = f.value(x, y, z); - result = p.value(x, y, z); - Assert.assertEquals("half-way between sample points (middle of the patch)", expected, result, 0.3); - - x = 3.5; - y = -3.5; - z = -10; - expected = f.value(x, y, z); - result = p.value(x, y, z); - Assert.assertEquals("half-way between sample points (border of the patch)", expected, result, 0.3); - } - - /** - * Test of interpolator for a sine wave. - *

- *

- * f(x, y, z) = a cos [ω z - ky x - ky y] - *

- * with A = 0.2, ω = 0.5, kx = 2, ky = 1. - */ - @Ignore@Test - public void testWave() { - double[] xval = new double[] {3, 4, 5, 6.5}; - double[] yval = new double[] {-4, -3, -1, 2, 2.5}; - double[] zval = new double[] {-12, -8, -5.5, -3, 0, 4}; - - final double a = 0.2; - final double omega = 0.5; - final double kx = 2; - final double ky = 1; - - // Function values - TrivariateFunction f = new TrivariateFunction() { - public double value(double x, double y, double z) { - return a * FastMath.cos(omega * z - kx * x - ky * y); - } - }; - - double[][][] fval = new double[xval.length][yval.length][zval.length]; - for (int i = 0; i < xval.length; i++) { - for (int j = 0; j < yval.length; j++) { - for (int k = 0; k < zval.length; k++) { - fval[i][j][k] = f.value(xval[i], yval[j], zval[k]); - } - } - } - - TrivariateGridInterpolator interpolator = new TricubicSplineInterpolator(); - - TrivariateFunction p = interpolator.interpolate(xval, yval, zval, fval); - double x, y, z; - double expected, result; - - x = 4; - y = -3; - z = 0; - expected = f.value(x, y, z); - result = p.value(x, y, z); - Assert.assertEquals("On sample point", - expected, result, 1e-12); - - x = 4.5; - y = -1.5; - z = -4.25; - expected = f.value(x, y, z); - result = p.value(x, y, z); - Assert.assertEquals("Half-way between sample points (middle of the patch)", - expected, result, 0.1); - - x = 3.5; - y = -3.5; - z = -10; - expected = f.value(x, y, z); - result = p.value(x, y, z); - Assert.assertEquals("Half-way between sample points (border of the patch)", - expected, result, 0.1); - } -} diff --git a/src/test/java/org/apache/commons/math4/analysis/solvers/NewtonSolverTest.java b/src/test/java/org/apache/commons/math4/analysis/solvers/NewtonSolverTest.java deleted file mode 100644 index e9b0896ca..000000000 --- a/src/test/java/org/apache/commons/math4/analysis/solvers/NewtonSolverTest.java +++ /dev/null @@ -1,111 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.analysis.solvers; - -import org.apache.commons.math4.analysis.DifferentiableUnivariateFunction; -import org.apache.commons.math4.analysis.QuinticFunction; -import org.apache.commons.math4.analysis.UnivariateFunction; -import org.apache.commons.math4.analysis.differentiation.DerivativeStructure; -import org.apache.commons.math4.analysis.differentiation.UnivariateDifferentiableFunction; -import org.apache.commons.math4.analysis.function.Sin; -import org.apache.commons.math4.analysis.solvers.NewtonSolver; -import org.apache.commons.math4.util.FastMath; -import org.junit.Assert; -import org.junit.Test; - - -/** - * @deprecated - */ -@Deprecated -public final class NewtonSolverTest { - /** - * - */ - @Test - public void testSinZero() { - DifferentiableUnivariateFunction f = new Sin(); - double result; - - NewtonSolver solver = new NewtonSolver(); - result = solver.solve(100, f, 3, 4); - Assert.assertEquals(result, FastMath.PI, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 1, 4); - Assert.assertEquals(result, FastMath.PI, solver.getAbsoluteAccuracy()); - - Assert.assertTrue(solver.getEvaluations() > 0); - } - - /** - * - */ - @Test - public void testQuinticZero() { - final UnivariateDifferentiableFunction q = new QuinticFunction(); - DifferentiableUnivariateFunction f = new DifferentiableUnivariateFunction() { - - public double value(double x) { - return q.value(x); - } - - public UnivariateFunction derivative() { - return new UnivariateFunction() { - public double value(double x) { - return q.value(new DerivativeStructure(1, 1, 0, x)).getPartialDerivative(1); - } - }; - } - - }; - double result; - - NewtonSolver solver = new NewtonSolver(); - result = solver.solve(100, f, -0.2, 0.2); - Assert.assertEquals(result, 0, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, -0.1, 0.3); - Assert.assertEquals(result, 0, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, -0.3, 0.45); - Assert.assertEquals(result, 0, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 0.3, 0.7); - Assert.assertEquals(result, 0.5, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 0.2, 0.6); - Assert.assertEquals(result, 0.5, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 0.05, 0.95); - Assert.assertEquals(result, 0.5, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 0.85, 1.25); - Assert.assertEquals(result, 1.0, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 0.8, 1.2); - Assert.assertEquals(result, 1.0, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 0.85, 1.75); - Assert.assertEquals(result, 1.0, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 0.55, 1.45); - Assert.assertEquals(result, 1.0, solver.getAbsoluteAccuracy()); - - result = solver.solve(100, f, 0.85, 5); - Assert.assertEquals(result, 1.0, solver.getAbsoluteAccuracy()); - } -} diff --git a/src/test/java/org/apache/commons/math4/fitting/CurveFitterTest.java b/src/test/java/org/apache/commons/math4/fitting/CurveFitterTest.java deleted file mode 100644 index 18915e074..000000000 --- a/src/test/java/org/apache/commons/math4/fitting/CurveFitterTest.java +++ /dev/null @@ -1,143 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.fitting; - -import org.apache.commons.math4.analysis.ParametricUnivariateFunction; -import org.apache.commons.math4.fitting.CurveFitter; -import org.apache.commons.math4.optim.nonlinear.vector.jacobian.LevenbergMarquardtOptimizer; -import org.apache.commons.math4.util.FastMath; -import org.junit.Assert; -import org.junit.Test; - -@Deprecated -public class CurveFitterTest { - @Test - public void testMath303() { - LevenbergMarquardtOptimizer optimizer = new LevenbergMarquardtOptimizer(); - CurveFitter fitter = new CurveFitter(optimizer); - fitter.addObservedPoint(2.805d, 0.6934785852953367d); - fitter.addObservedPoint(2.74333333333333d, 0.6306772025518496d); - fitter.addObservedPoint(1.655d, 0.9474675497289684); - fitter.addObservedPoint(1.725d, 0.9013594835804194d); - - ParametricUnivariateFunction sif = new SimpleInverseFunction(); - - double[] initialguess1 = new double[1]; - initialguess1[0] = 1.0d; - Assert.assertEquals(1, fitter.fit(sif, initialguess1).length); - - double[] initialguess2 = new double[2]; - initialguess2[0] = 1.0d; - initialguess2[1] = .5d; - Assert.assertEquals(2, fitter.fit(sif, initialguess2).length); - } - - @Test - public void testMath304() { - LevenbergMarquardtOptimizer optimizer = new LevenbergMarquardtOptimizer(); - CurveFitter fitter = new CurveFitter(optimizer); - fitter.addObservedPoint(2.805d, 0.6934785852953367d); - fitter.addObservedPoint(2.74333333333333d, 0.6306772025518496d); - fitter.addObservedPoint(1.655d, 0.9474675497289684); - fitter.addObservedPoint(1.725d, 0.9013594835804194d); - - ParametricUnivariateFunction sif = new SimpleInverseFunction(); - - double[] initialguess1 = new double[1]; - initialguess1[0] = 1.0d; - Assert.assertEquals(1.6357215104109237, fitter.fit(sif, initialguess1)[0], 1.0e-14); - - double[] initialguess2 = new double[1]; - initialguess2[0] = 10.0d; - Assert.assertEquals(1.6357215104109237, fitter.fit(sif, initialguess1)[0], 1.0e-14); - } - - @Test - public void testMath372() { - LevenbergMarquardtOptimizer optimizer = new LevenbergMarquardtOptimizer(); - CurveFitter curveFitter = new CurveFitter(optimizer); - - curveFitter.addObservedPoint( 15, 4443); - curveFitter.addObservedPoint( 31, 8493); - curveFitter.addObservedPoint( 62, 17586); - curveFitter.addObservedPoint(125, 30582); - curveFitter.addObservedPoint(250, 45087); - curveFitter.addObservedPoint(500, 50683); - - ParametricUnivariateFunction f = new ParametricUnivariateFunction() { - public double value(double x, double ... parameters) { - double a = parameters[0]; - double b = parameters[1]; - double c = parameters[2]; - double d = parameters[3]; - - return d + ((a - d) / (1 + FastMath.pow(x / c, b))); - } - - public double[] gradient(double x, double ... parameters) { - double a = parameters[0]; - double b = parameters[1]; - double c = parameters[2]; - double d = parameters[3]; - - double[] gradients = new double[4]; - double den = 1 + FastMath.pow(x / c, b); - - // derivative with respect to a - gradients[0] = 1 / den; - - // derivative with respect to b - // in the reported (invalid) issue, there was a sign error here - gradients[1] = -((a - d) * FastMath.pow(x / c, b) * FastMath.log(x / c)) / (den * den); - - // derivative with respect to c - gradients[2] = (b * FastMath.pow(x / c, b - 1) * (x / (c * c)) * (a - d)) / (den * den); - - // derivative with respect to d - gradients[3] = 1 - (1 / den); - - return gradients; - - } - }; - - double[] initialGuess = new double[] { 1500, 0.95, 65, 35000 }; - double[] estimatedParameters = curveFitter.fit(f, initialGuess); - - Assert.assertEquals( 2411.00, estimatedParameters[0], 500.00); - Assert.assertEquals( 1.62, estimatedParameters[1], 0.04); - Assert.assertEquals( 111.22, estimatedParameters[2], 0.30); - Assert.assertEquals(55347.47, estimatedParameters[3], 300.00); - Assert.assertTrue(optimizer.getRMS() < 600.0); - } - - private static class SimpleInverseFunction implements ParametricUnivariateFunction { - - public double value(double x, double ... parameters) { - return parameters[0] / x + (parameters.length < 2 ? 0 : parameters[1]); - } - - public double[] gradient(double x, double ... doubles) { - double[] gradientVector = new double[doubles.length]; - gradientVector[0] = 1 / x; - if (doubles.length >= 2) { - gradientVector[1] = 1; - } - return gradientVector; - } - } -} diff --git a/src/test/java/org/apache/commons/math4/fitting/GaussianFitterTest.java b/src/test/java/org/apache/commons/math4/fitting/GaussianFitterTest.java deleted file mode 100644 index a7ca9b218..000000000 --- a/src/test/java/org/apache/commons/math4/fitting/GaussianFitterTest.java +++ /dev/null @@ -1,364 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.fitting; - -import org.apache.commons.math4.exception.MathIllegalArgumentException; -import org.apache.commons.math4.fitting.GaussianFitter; -import org.apache.commons.math4.optim.nonlinear.vector.jacobian.LevenbergMarquardtOptimizer; -import org.junit.Assert; -import org.junit.Test; - -/** - * Tests {@link GaussianFitter}. - * - * @since 2.2 - */ -@Deprecated -public class GaussianFitterTest { - /** Good data. */ - protected static final double[][] DATASET1 = new double[][] { - {4.0254623, 531026.0}, - {4.02804905, 664002.0}, - {4.02934242, 787079.0}, - {4.03128248, 984167.0}, - {4.03386923, 1294546.0}, - {4.03580929, 1560230.0}, - {4.03839603, 1887233.0}, - {4.0396894, 2113240.0}, - {4.04162946, 2375211.0}, - {4.04421621, 2687152.0}, - {4.04550958, 2862644.0}, - {4.04744964, 3078898.0}, - {4.05003639, 3327238.0}, - {4.05132976, 3461228.0}, - {4.05326982, 3580526.0}, - {4.05585657, 3576946.0}, - {4.05779662, 3439750.0}, - {4.06038337, 3220296.0}, - {4.06167674, 3070073.0}, - {4.0636168, 2877648.0}, - {4.06620355, 2595848.0}, - {4.06749692, 2390157.0}, - {4.06943698, 2175960.0}, - {4.07202373, 1895104.0}, - {4.0733171, 1687576.0}, - {4.07525716, 1447024.0}, - {4.0778439, 1130879.0}, - {4.07978396, 904900.0}, - {4.08237071, 717104.0}, - {4.08366408, 620014.0} - }; - /** Poor data: right of peak not symmetric with left of peak. */ - protected static final double[][] DATASET2 = new double[][] { - {-20.15, 1523.0}, - {-19.65, 1566.0}, - {-19.15, 1592.0}, - {-18.65, 1927.0}, - {-18.15, 3089.0}, - {-17.65, 6068.0}, - {-17.15, 14239.0}, - {-16.65, 34124.0}, - {-16.15, 64097.0}, - {-15.65, 110352.0}, - {-15.15, 164742.0}, - {-14.65, 209499.0}, - {-14.15, 267274.0}, - {-13.65, 283290.0}, - {-13.15, 275363.0}, - {-12.65, 258014.0}, - {-12.15, 225000.0}, - {-11.65, 200000.0}, - {-11.15, 190000.0}, - {-10.65, 185000.0}, - {-10.15, 180000.0}, - { -9.65, 179000.0}, - { -9.15, 178000.0}, - { -8.65, 177000.0}, - { -8.15, 176000.0}, - { -7.65, 175000.0}, - { -7.15, 174000.0}, - { -6.65, 173000.0}, - { -6.15, 172000.0}, - { -5.65, 171000.0}, - { -5.15, 170000.0} - }; - /** Poor data: long tails. */ - protected static final double[][] DATASET3 = new double[][] { - {-90.15, 1513.0}, - {-80.15, 1514.0}, - {-70.15, 1513.0}, - {-60.15, 1514.0}, - {-50.15, 1513.0}, - {-40.15, 1514.0}, - {-30.15, 1513.0}, - {-20.15, 1523.0}, - {-19.65, 1566.0}, - {-19.15, 1592.0}, - {-18.65, 1927.0}, - {-18.15, 3089.0}, - {-17.65, 6068.0}, - {-17.15, 14239.0}, - {-16.65, 34124.0}, - {-16.15, 64097.0}, - {-15.65, 110352.0}, - {-15.15, 164742.0}, - {-14.65, 209499.0}, - {-14.15, 267274.0}, - {-13.65, 283290.0}, - {-13.15, 275363.0}, - {-12.65, 258014.0}, - {-12.15, 214073.0}, - {-11.65, 182244.0}, - {-11.15, 136419.0}, - {-10.65, 97823.0}, - {-10.15, 58930.0}, - { -9.65, 35404.0}, - { -9.15, 16120.0}, - { -8.65, 9823.0}, - { -8.15, 5064.0}, - { -7.65, 2575.0}, - { -7.15, 1642.0}, - { -6.65, 1101.0}, - { -6.15, 812.0}, - { -5.65, 690.0}, - { -5.15, 565.0}, - { 5.15, 564.0}, - { 15.15, 565.0}, - { 25.15, 564.0}, - { 35.15, 565.0}, - { 45.15, 564.0}, - { 55.15, 565.0}, - { 65.15, 564.0}, - { 75.15, 565.0} - }; - /** Poor data: right of peak is missing. */ - protected static final double[][] DATASET4 = new double[][] { - {-20.15, 1523.0}, - {-19.65, 1566.0}, - {-19.15, 1592.0}, - {-18.65, 1927.0}, - {-18.15, 3089.0}, - {-17.65, 6068.0}, - {-17.15, 14239.0}, - {-16.65, 34124.0}, - {-16.15, 64097.0}, - {-15.65, 110352.0}, - {-15.15, 164742.0}, - {-14.65, 209499.0}, - {-14.15, 267274.0}, - {-13.65, 283290.0} - }; - /** Good data, but few points. */ - protected static final double[][] DATASET5 = new double[][] { - {4.0254623, 531026.0}, - {4.03128248, 984167.0}, - {4.03839603, 1887233.0}, - {4.04421621, 2687152.0}, - {4.05132976, 3461228.0}, - {4.05326982, 3580526.0}, - {4.05779662, 3439750.0}, - {4.0636168, 2877648.0}, - {4.06943698, 2175960.0}, - {4.07525716, 1447024.0}, - {4.08237071, 717104.0}, - {4.08366408, 620014.0} - }; - - /** - * Basic. - */ - @Test - public void testFit01() { - GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - addDatasetToGaussianFitter(DATASET1, fitter); - double[] parameters = fitter.fit(); - - Assert.assertEquals(3496978.1837704973, parameters[0], 1e-4); - Assert.assertEquals(4.054933085999146, parameters[1], 1e-4); - Assert.assertEquals(0.015039355620304326, parameters[2], 1e-4); - } - - /** - * Zero points is not enough observed points. - */ - @Test(expected=MathIllegalArgumentException.class) - public void testFit02() { - GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - fitter.fit(); - } - - /** - * Two points is not enough observed points. - */ - @Test(expected=MathIllegalArgumentException.class) - public void testFit03() { - GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - addDatasetToGaussianFitter(new double[][] { - {4.0254623, 531026.0}, - {4.02804905, 664002.0}}, - fitter); - fitter.fit(); - } - - /** - * Poor data: right of peak not symmetric with left of peak. - */ - @Test - public void testFit04() { - GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - addDatasetToGaussianFitter(DATASET2, fitter); - double[] parameters = fitter.fit(); - - Assert.assertEquals(233003.2967252038, parameters[0], 1e-4); - Assert.assertEquals(-10.654887521095983, parameters[1], 1e-4); - Assert.assertEquals(4.335937353196641, parameters[2], 1e-4); - } - - /** - * Poor data: long tails. - */ - @Test - public void testFit05() { - GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - addDatasetToGaussianFitter(DATASET3, fitter); - double[] parameters = fitter.fit(); - - Assert.assertEquals(283863.81929180305, parameters[0], 1e-4); - Assert.assertEquals(-13.29641995105174, parameters[1], 1e-4); - Assert.assertEquals(1.7297330293549908, parameters[2], 1e-4); - } - - /** - * Poor data: right of peak is missing. - */ - @Test - public void testFit06() { - GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - addDatasetToGaussianFitter(DATASET4, fitter); - double[] parameters = fitter.fit(); - - Assert.assertEquals(285250.66754309234, parameters[0], 1e-4); - Assert.assertEquals(-13.528375695228455, parameters[1], 1e-4); - Assert.assertEquals(1.5204344894331614, parameters[2], 1e-4); - } - - /** - * Basic with smaller dataset. - */ - @Test - public void testFit07() { - GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - addDatasetToGaussianFitter(DATASET5, fitter); - double[] parameters = fitter.fit(); - - Assert.assertEquals(3514384.729342235, parameters[0], 1e-4); - Assert.assertEquals(4.054970307455625, parameters[1], 1e-4); - Assert.assertEquals(0.015029412832160017, parameters[2], 1e-4); - } - - @Test - public void testMath519() { - // The optimizer will try negative sigma values but "GaussianFitter" - // will catch the raised exceptions and return NaN values instead. - - final double[] data = { - 1.1143831578403364E-29, - 4.95281403484594E-28, - 1.1171347211930288E-26, - 1.7044813962636277E-25, - 1.9784716574832164E-24, - 1.8630236407866774E-23, - 1.4820532905097742E-22, - 1.0241963854632831E-21, - 6.275077366673128E-21, - 3.461808994532493E-20, - 1.7407124684715706E-19, - 8.056687953553974E-19, - 3.460193945992071E-18, - 1.3883326374011525E-17, - 5.233894983671116E-17, - 1.8630791465263745E-16, - 6.288759227922111E-16, - 2.0204433920597856E-15, - 6.198768938576155E-15, - 1.821419346860626E-14, - 5.139176445538471E-14, - 1.3956427429045787E-13, - 3.655705706448139E-13, - 9.253753324779779E-13, - 2.267636001476696E-12, - 5.3880460095836855E-12, - 1.2431632654852931E-11 - }; - - GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - for (int i = 0; i < data.length; i++) { - fitter.addObservedPoint(i, data[i]); - } - final double[] p = fitter.fit(); - - Assert.assertEquals(53.1572792, p[1], 1e-7); - Assert.assertEquals(5.75214622, p[2], 1e-8); - } - - @Test - public void testMath798() { - final GaussianFitter fitter = new GaussianFitter(new LevenbergMarquardtOptimizer()); - - // When the data points are not commented out below, the fit stalls. - // This is expected however, since the whole dataset hardly looks like - // a Gaussian. - // When commented out, the fit proceeds fine. - - fitter.addObservedPoint(0.23, 395.0); - //fitter.addObservedPoint(0.68, 0.0); - fitter.addObservedPoint(1.14, 376.0); - //fitter.addObservedPoint(1.59, 0.0); - fitter.addObservedPoint(2.05, 163.0); - //fitter.addObservedPoint(2.50, 0.0); - fitter.addObservedPoint(2.95, 49.0); - //fitter.addObservedPoint(3.41, 0.0); - fitter.addObservedPoint(3.86, 16.0); - //fitter.addObservedPoint(4.32, 0.0); - fitter.addObservedPoint(4.77, 1.0); - - final double[] p = fitter.fit(); - - // Values are copied from a previous run of this test. - Assert.assertEquals(420.8397296167364, p[0], 1e-12); - Assert.assertEquals(0.603770729862231, p[1], 1e-15); - Assert.assertEquals(1.0786447936766612, p[2], 1e-14); - } - - /** - * Adds the specified points to specified GaussianFitter - * instance. - * - * @param points data points where first dimension is a point index and - * second dimension is an array of length two representing the point - * with the first value corresponding to X and the second value - * corresponding to Y - * @param fitter fitter to which the points in points should be - * added as observed points - */ - protected static void addDatasetToGaussianFitter(double[][] points, - GaussianFitter fitter) { - for (int i = 0; i < points.length; i++) { - fitter.addObservedPoint(points[i][0], points[i][1]); - } - } -} diff --git a/src/test/java/org/apache/commons/math4/fitting/HarmonicFitterTest.java b/src/test/java/org/apache/commons/math4/fitting/HarmonicFitterTest.java deleted file mode 100644 index 328d060e5..000000000 --- a/src/test/java/org/apache/commons/math4/fitting/HarmonicFitterTest.java +++ /dev/null @@ -1,187 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.fitting; - -import java.util.Random; - -import org.apache.commons.math4.analysis.function.HarmonicOscillator; -import org.apache.commons.math4.exception.MathIllegalStateException; -import org.apache.commons.math4.exception.NumberIsTooSmallException; -import org.apache.commons.math4.fitting.HarmonicFitter; -import org.apache.commons.math4.fitting.WeightedObservedPoint; -import org.apache.commons.math4.optim.nonlinear.vector.jacobian.LevenbergMarquardtOptimizer; -import org.apache.commons.math4.util.FastMath; -import org.apache.commons.math4.util.MathUtils; -import org.junit.Test; -import org.junit.Assert; - -@Deprecated -public class HarmonicFitterTest { - @Test(expected=NumberIsTooSmallException.class) - public void testPreconditions1() { - HarmonicFitter fitter = - new HarmonicFitter(new LevenbergMarquardtOptimizer()); - - fitter.fit(); - } - - @Test - public void testNoError() { - final double a = 0.2; - final double w = 3.4; - final double p = 4.1; - HarmonicOscillator f = new HarmonicOscillator(a, w, p); - - HarmonicFitter fitter = - new HarmonicFitter(new LevenbergMarquardtOptimizer()); - for (double x = 0.0; x < 1.3; x += 0.01) { - fitter.addObservedPoint(1, x, f.value(x)); - } - - final double[] fitted = fitter.fit(); - Assert.assertEquals(a, fitted[0], 1.0e-13); - Assert.assertEquals(w, fitted[1], 1.0e-13); - Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1e-13); - - HarmonicOscillator ff = new HarmonicOscillator(fitted[0], fitted[1], fitted[2]); - - for (double x = -1.0; x < 1.0; x += 0.01) { - Assert.assertTrue(FastMath.abs(f.value(x) - ff.value(x)) < 1e-13); - } - } - - @Test - public void test1PercentError() { - Random randomizer = new Random(64925784252l); - final double a = 0.2; - final double w = 3.4; - final double p = 4.1; - HarmonicOscillator f = new HarmonicOscillator(a, w, p); - - HarmonicFitter fitter = - new HarmonicFitter(new LevenbergMarquardtOptimizer()); - for (double x = 0.0; x < 10.0; x += 0.1) { - fitter.addObservedPoint(1, x, - f.value(x) + 0.01 * randomizer.nextGaussian()); - } - - final double[] fitted = fitter.fit(); - Assert.assertEquals(a, fitted[0], 7.6e-4); - Assert.assertEquals(w, fitted[1], 2.7e-3); - Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.3e-2); - } - - @Test - public void testTinyVariationsData() { - Random randomizer = new Random(64925784252l); - - HarmonicFitter fitter = - new HarmonicFitter(new LevenbergMarquardtOptimizer()); - for (double x = 0.0; x < 10.0; x += 0.1) { - fitter.addObservedPoint(1, x, 1e-7 * randomizer.nextGaussian()); - } - - fitter.fit(); - // This test serves to cover the part of the code of "guessAOmega" - // when the algorithm using integrals fails. - } - - @Test - public void testInitialGuess() { - Random randomizer = new Random(45314242l); - final double a = 0.2; - final double w = 3.4; - final double p = 4.1; - HarmonicOscillator f = new HarmonicOscillator(a, w, p); - - HarmonicFitter fitter = - new HarmonicFitter(new LevenbergMarquardtOptimizer()); - for (double x = 0.0; x < 10.0; x += 0.1) { - fitter.addObservedPoint(1, x, - f.value(x) + 0.01 * randomizer.nextGaussian()); - } - - final double[] fitted = fitter.fit(new double[] { 0.15, 3.6, 4.5 }); - Assert.assertEquals(a, fitted[0], 1.2e-3); - Assert.assertEquals(w, fitted[1], 3.3e-3); - Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.7e-2); - } - - @Test - public void testUnsorted() { - Random randomizer = new Random(64925784252l); - final double a = 0.2; - final double w = 3.4; - final double p = 4.1; - HarmonicOscillator f = new HarmonicOscillator(a, w, p); - - HarmonicFitter fitter = - new HarmonicFitter(new LevenbergMarquardtOptimizer()); - - // build a regularly spaced array of measurements - int size = 100; - double[] xTab = new double[size]; - double[] yTab = new double[size]; - for (int i = 0; i < size; ++i) { - xTab[i] = 0.1 * i; - yTab[i] = f.value(xTab[i]) + 0.01 * randomizer.nextGaussian(); - } - - // shake it - for (int i = 0; i < size; ++i) { - int i1 = randomizer.nextInt(size); - int i2 = randomizer.nextInt(size); - double xTmp = xTab[i1]; - double yTmp = yTab[i1]; - xTab[i1] = xTab[i2]; - yTab[i1] = yTab[i2]; - xTab[i2] = xTmp; - yTab[i2] = yTmp; - } - - // pass it to the fitter - for (int i = 0; i < size; ++i) { - fitter.addObservedPoint(1, xTab[i], yTab[i]); - } - - final double[] fitted = fitter.fit(); - Assert.assertEquals(a, fitted[0], 7.6e-4); - Assert.assertEquals(w, fitted[1], 3.5e-3); - Assert.assertEquals(p, MathUtils.normalizeAngle(fitted[2], p), 1.5e-2); - } - - @Test(expected=MathIllegalStateException.class) - public void testMath844() { - final double[] y = { 0, 1, 2, 3, 2, 1, - 0, -1, -2, -3, -2, -1, - 0, 1, 2, 3, 2, 1, - 0, -1, -2, -3, -2, -1, - 0, 1, 2, 3, 2, 1, 0 }; - final int len = y.length; - final WeightedObservedPoint[] points = new WeightedObservedPoint[len]; - for (int i = 0; i < len; i++) { - points[i] = new WeightedObservedPoint(1, i, y[i]); - } - - // The guesser fails because the function is far from an harmonic - // function: It is a triangular periodic function with amplitude 3 - // and period 12, and all sample points are taken at integer abscissae - // so function values all belong to the integer subset {-3, -2, -1, 0, - // 1, 2, 3}. - new HarmonicFitter.ParameterGuesser(points); - } -} diff --git a/src/test/java/org/apache/commons/math4/fitting/PolynomialFitterTest.java b/src/test/java/org/apache/commons/math4/fitting/PolynomialFitterTest.java deleted file mode 100644 index ff543ab7b..000000000 --- a/src/test/java/org/apache/commons/math4/fitting/PolynomialFitterTest.java +++ /dev/null @@ -1,288 +0,0 @@ -/* - * Licensed to the Apache Software Foundation (ASF) under one or more - * contributor license agreements. See the NOTICE file distributed with - * this work for additional information regarding copyright ownership. - * The ASF licenses this file to You under the Apache License, Version 2.0 - * (the "License"); you may not use this file except in compliance with - * the License. You may obtain a copy of the License at - * - * http://www.apache.org/licenses/LICENSE-2.0 - * - * Unless required by applicable law or agreed to in writing, software - * distributed under the License is distributed on an "AS IS" BASIS, - * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied. - * See the License for the specific language governing permissions and - * limitations under the License. - */ -package org.apache.commons.math4.fitting; - -import java.util.Random; - -import org.apache.commons.math4.TestUtils; -import org.apache.commons.math4.analysis.polynomials.PolynomialFunction; -import org.apache.commons.math4.analysis.polynomials.PolynomialFunction.Parametric; -import org.apache.commons.math4.distribution.RealDistribution; -import org.apache.commons.math4.distribution.UniformRealDistribution; -import org.apache.commons.math4.exception.ConvergenceException; -import org.apache.commons.math4.exception.TooManyEvaluationsException; -import org.apache.commons.math4.fitting.CurveFitter; -import org.apache.commons.math4.fitting.PolynomialFitter; -import org.apache.commons.math4.optim.SimpleVectorValueChecker; -import org.apache.commons.math4.optim.nonlinear.vector.MultivariateVectorOptimizer; -import org.apache.commons.math4.optim.nonlinear.vector.jacobian.GaussNewtonOptimizer; -import org.apache.commons.math4.optim.nonlinear.vector.jacobian.LevenbergMarquardtOptimizer; -import org.apache.commons.math4.util.FastMath; -import org.junit.Test; -import org.junit.Assert; - -/** - * Test for class {@link CurveFitter} where the function to fit is a - * polynomial. - */ -@Deprecated -public class PolynomialFitterTest { - @Test - public void testFit() { - final RealDistribution rng = new UniformRealDistribution(-100, 100); - rng.reseedRandomGenerator(64925784252L); - - final LevenbergMarquardtOptimizer optim = new LevenbergMarquardtOptimizer(); - final PolynomialFitter fitter = new PolynomialFitter(optim); - final double[] coeff = { 12.9, -3.4, 2.1 }; // 12.9 - 3.4 x + 2.1 x^2 - final PolynomialFunction f = new PolynomialFunction(coeff); - - // Collect data from a known polynomial. - for (int i = 0; i < 100; i++) { - final double x = rng.sample(); - fitter.addObservedPoint(x, f.value(x)); - } - - // Start fit from initial guesses that are far from the optimal values. - final double[] best = fitter.fit(new double[] { -1e-20, 3e15, -5e25 }); - - TestUtils.assertEquals("best != coeff", coeff, best, 1e-12); - } - - @Test - public void testNoError() { - Random randomizer = new Random(64925784252l); - for (int degree = 1; degree < 10; ++degree) { - PolynomialFunction p = buildRandomPolynomial(degree, randomizer); - - PolynomialFitter fitter = new PolynomialFitter(new LevenbergMarquardtOptimizer()); - for (int i = 0; i <= degree; ++i) { - fitter.addObservedPoint(1.0, i, p.value(i)); - } - - final double[] init = new double[degree + 1]; - PolynomialFunction fitted = new PolynomialFunction(fitter.fit(init)); - - for (double x = -1.0; x < 1.0; x += 0.01) { - double error = FastMath.abs(p.value(x) - fitted.value(x)) / - (1.0 + FastMath.abs(p.value(x))); - Assert.assertEquals(0.0, error, 1.0e-6); - } - } - } - - @Test - public void testSmallError() { - Random randomizer = new Random(53882150042l); - double maxError = 0; - for (int degree = 0; degree < 10; ++degree) { - PolynomialFunction p = buildRandomPolynomial(degree, randomizer); - - PolynomialFitter fitter = new PolynomialFitter(new LevenbergMarquardtOptimizer()); - for (double x = -1.0; x < 1.0; x += 0.01) { - fitter.addObservedPoint(1.0, x, - p.value(x) + 0.1 * randomizer.nextGaussian()); - } - - final double[] init = new double[degree + 1]; - PolynomialFunction fitted = new PolynomialFunction(fitter.fit(init)); - - for (double x = -1.0; x < 1.0; x += 0.01) { - double error = FastMath.abs(p.value(x) - fitted.value(x)) / - (1.0 + FastMath.abs(p.value(x))); - maxError = FastMath.max(maxError, error); - Assert.assertTrue(FastMath.abs(error) < 0.1); - } - } - Assert.assertTrue(maxError > 0.01); - } - - @Test - public void testMath798() { - final double tol = 1e-14; - final SimpleVectorValueChecker checker = new SimpleVectorValueChecker(tol, tol); - final double[] init = new double[] { 0, 0 }; - final int maxEval = 3; - - final double[] lm = doMath798(new LevenbergMarquardtOptimizer(checker), maxEval, init); - final double[] gn = doMath798(new GaussNewtonOptimizer(checker), maxEval, init); - - for (int i = 0; i <= 1; i++) { - Assert.assertEquals(lm[i], gn[i], tol); - } - } - - /** - * This test shows that the user can set the maximum number of iterations - * to avoid running for too long. - * But in the test case, the real problem is that the tolerance is way too - * stringent. - */ - @Test(expected=TooManyEvaluationsException.class) - public void testMath798WithToleranceTooLow() { - final double tol = 1e-100; - final SimpleVectorValueChecker checker = new SimpleVectorValueChecker(tol, tol); - final double[] init = new double[] { 0, 0 }; - final int maxEval = 10000; // Trying hard to fit. - - @SuppressWarnings("unused") - final double[] gn = doMath798(new GaussNewtonOptimizer(checker), maxEval, init); - } - - /** - * This test shows that the user can set the maximum number of iterations - * to avoid running for too long. - * Even if the real problem is that the tolerance is way too stringent, it - * is possible to get the best solution so far, i.e. a checker will return - * the point when the maximum iteration count has been reached. - */ - @Test - public void testMath798WithToleranceTooLowButNoException() { - final double tol = 1e-100; - final double[] init = new double[] { 0, 0 }; - final int maxEval = 10000; // Trying hard to fit. - final SimpleVectorValueChecker checker = new SimpleVectorValueChecker(tol, tol, maxEval); - - final double[] lm = doMath798(new LevenbergMarquardtOptimizer(checker), maxEval, init); - final double[] gn = doMath798(new GaussNewtonOptimizer(checker), maxEval, init); - - for (int i = 0; i <= 1; i++) { - Assert.assertEquals(lm[i], gn[i], 1e-15); - } - } - - /** - * @param optimizer Optimizer. - * @param maxEval Maximum number of function evaluations. - * @param init First guess. - * @return the solution found by the given optimizer. - */ - private double[] doMath798(MultivariateVectorOptimizer optimizer, - int maxEval, - double[] init) { - final CurveFitter fitter = new CurveFitter(optimizer); - - fitter.addObservedPoint(-0.2, -7.12442E-13); - fitter.addObservedPoint(-0.199, -4.33397E-13); - fitter.addObservedPoint(-0.198, -2.823E-13); - fitter.addObservedPoint(-0.197, -1.40405E-13); - fitter.addObservedPoint(-0.196, -7.80821E-15); - fitter.addObservedPoint(-0.195, 6.20484E-14); - fitter.addObservedPoint(-0.194, 7.24673E-14); - fitter.addObservedPoint(-0.193, 1.47152E-13); - fitter.addObservedPoint(-0.192, 1.9629E-13); - fitter.addObservedPoint(-0.191, 2.12038E-13); - fitter.addObservedPoint(-0.19, 2.46906E-13); - fitter.addObservedPoint(-0.189, 2.77495E-13); - fitter.addObservedPoint(-0.188, 2.51281E-13); - fitter.addObservedPoint(-0.187, 2.64001E-13); - fitter.addObservedPoint(-0.186, 2.8882E-13); - fitter.addObservedPoint(-0.185, 3.13604E-13); - fitter.addObservedPoint(-0.184, 3.14248E-13); - fitter.addObservedPoint(-0.183, 3.1172E-13); - fitter.addObservedPoint(-0.182, 3.12912E-13); - fitter.addObservedPoint(-0.181, 3.06761E-13); - fitter.addObservedPoint(-0.18, 2.8559E-13); - fitter.addObservedPoint(-0.179, 2.86806E-13); - fitter.addObservedPoint(-0.178, 2.985E-13); - fitter.addObservedPoint(-0.177, 2.67148E-13); - fitter.addObservedPoint(-0.176, 2.94173E-13); - fitter.addObservedPoint(-0.175, 3.27528E-13); - fitter.addObservedPoint(-0.174, 3.33858E-13); - fitter.addObservedPoint(-0.173, 2.97511E-13); - fitter.addObservedPoint(-0.172, 2.8615E-13); - fitter.addObservedPoint(-0.171, 2.84624E-13); - - final double[] coeff = fitter.fit(maxEval, - new PolynomialFunction.Parametric(), - init); - return coeff; - } - - @Test - public void testRedundantSolvable() { - // Levenberg-Marquardt should handle redundant information gracefully - checkUnsolvableProblem(new LevenbergMarquardtOptimizer(), true); - } - - @Test - public void testRedundantUnsolvable() { - // Gauss-Newton should not be able to solve redundant information - checkUnsolvableProblem(new GaussNewtonOptimizer(true, new SimpleVectorValueChecker(1e-15, 1e-15)), false); - } - - @Test - public void testLargeSample() { - Random randomizer = new Random(0x5551480dca5b369bl); - double maxError = 0; - for (int degree = 0; degree < 10; ++degree) { - PolynomialFunction p = buildRandomPolynomial(degree, randomizer); - - PolynomialFitter fitter = new PolynomialFitter(new LevenbergMarquardtOptimizer()); - for (int i = 0; i < 40000; ++i) { - double x = -1.0 + i / 20000.0; - fitter.addObservedPoint(1.0, x, - p.value(x) + 0.1 * randomizer.nextGaussian()); - } - - final double[] init = new double[degree + 1]; - PolynomialFunction fitted = new PolynomialFunction(fitter.fit(init)); - - for (double x = -1.0; x < 1.0; x += 0.01) { - double error = FastMath.abs(p.value(x) - fitted.value(x)) / - (1.0 + FastMath.abs(p.value(x))); - maxError = FastMath.max(maxError, error); - Assert.assertTrue(FastMath.abs(error) < 0.01); - } - } - Assert.assertTrue(maxError > 0.001); - } - - private void checkUnsolvableProblem(MultivariateVectorOptimizer optimizer, - boolean solvable) { - Random randomizer = new Random(1248788532l); - for (int degree = 0; degree < 10; ++degree) { - PolynomialFunction p = buildRandomPolynomial(degree, randomizer); - - PolynomialFitter fitter = new PolynomialFitter(optimizer); - - // reusing the same point over and over again does not bring - // information, the problem cannot be solved in this case for - // degrees greater than 1 (but one point is sufficient for - // degree 0) - for (double x = -1.0; x < 1.0; x += 0.01) { - fitter.addObservedPoint(1.0, 0.0, p.value(0.0)); - } - - try { - final double[] init = new double[degree + 1]; - fitter.fit(init); - Assert.assertTrue(solvable || (degree == 0)); - } catch(ConvergenceException e) { - Assert.assertTrue((! solvable) && (degree > 0)); - } - } - } - - private PolynomialFunction buildRandomPolynomial(int degree, Random randomizer) { - final double[] coefficients = new double[degree + 1]; - for (int i = 0; i <= degree; ++i) { - coefficients[i] = randomizer.nextGaussian(); - } - return new PolynomialFunction(coefficients); - } -}