private void testSensitivity(double[] xValues, double[] yValues, double[] keys, double delta) { PiecewisePolynomialWithSensitivityFunction1D func = new PiecewisePolynomialWithSensitivityFunction1D(); PiecewisePolynomialResultsWithSensitivity resultSensi = INTERP.interpolateWithSensitivity(xValues, yValues); DoubleArray[] computedArray = func.nodeSensitivity(resultSensi, keys); for (int i = 0; i < keys.length; ++i) { double base = func.evaluate(resultSensi, keys[i]).get(0); DoubleArray computed = func.nodeSensitivity(resultSensi, keys[i]); assertEquals(computed, computedArray[i]); for (int j = 0; j < yValues.length; ++j) { double[] yValuesBump = Arrays.copyOf(yValues, yValues.length); yValuesBump[j] += delta; PiecewisePolynomialResult resultBump = INTERP.interpolate(xValues, yValuesBump); double expected = (func.evaluate(resultBump, keys[i]).get(0) - base) / delta; assertEquals(computed.get(j), expected, delta); } } } }
final PiecewisePolynomialResultsWithSensitivity resultDw = interps[k].interpolateWithSensitivity(xValues, yValuesDw); final double[] valuesUp = FUNCTION.evaluate(resultUp, xKeys).toArray()[0]; final double[] valuesDw = FUNCTION.evaluate(resultDw, xKeys).toArray()[0]; final double[] diffUp = FUNCTION.differentiate(resultUp, xKeys).toArray()[0]; final double[] diffDw = FUNCTION.differentiate(resultDw, xKeys).toArray()[0]; for (int i = 0; i < 10 * nData; ++i) { final double xKeyUp = xKeys[i] * (1. + EPS); double valueFinite = 0.5 * (valuesUp[i] - valuesDw[i]) / EPS / yValues[j + 1]; double senseFinite = 0.5 * (diffUp[i] - diffDw[i]) / EPS / yValues[j + 1]; final double resNodeSensitivity = FUNCTION.nodeSensitivity(result, xKeys[i]).get(j); final double resNodeSensitivityXkeyUp = FUNCTION.nodeSensitivity(result, xKeyUp).get(j); final double resNodeSensitivityXkeyDw = FUNCTION.nodeSensitivity(result, xKeyDw).get(j); final double senseFiniteXkey = 0.5 * (resNodeSensitivityXkeyUp - resNodeSensitivityXkeyDw) / EPS / xKeys[i]; final double resDiffNodeSensitivity = FUNCTION.differentiateNodeSensitivity(result, xKeys[i]).get(j); assertEquals(valueFinite, resNodeSensitivity, Math.max(Math.abs(yValues[j + 1]) * EPS, EPS)); assertEquals(senseFinite, resDiffNodeSensitivity, Math.max(Math.abs(yValues[j + 1]) * EPS, EPS));
for (int j = 0; j < nKeys; ++j) { double key = xValues[0] + interval * j; InterpolatorTestUtil.assertRelative("clampedTest", FUNC.evaluate(resultBase, key).get(0), FUNC.evaluate(result, key) .get(0), EPS); InterpolatorTestUtil.assertArrayRelative("clampedTest", FUNC.nodeSensitivity(resultBase, key).toArray(), FUNC.nodeSensitivity(result, key).toArray(), EPS);
@Override protected DoubleArray doParameterSensitivity(double xValue) { ArgChecker.isTrue(Math.abs(xValue) > SMALL, "magnitude of xValue must not be small"); DoubleArray resSense = FUNCTION.nodeSensitivity(polySens.get(), xValue); return resSense.multipliedBy(DoubleArray.of(resSense.size(), i -> xValues[i] / xValue)); }
@Override protected double doInterpolate(double xValue) { ArgChecker.isTrue(Math.abs(xValue) > SMALL, "magnitude of xValue must not be small"); double resValue = FUNCTION.evaluate(poly, xValue).get(0); return resValue / xValue; }
public void baseInterpolationTest() { int nExamples = Y.length; int n = XX.length; for (int example = 0; example < nExamples; example++) { PiecewisePolynomialResult pp = BASE.interpolate(X, Y[example]); BoundCurveInterpolator bound = PCHIP.bind( DoubleArray.ofUnsafe(X), DoubleArray.ofUnsafe(Y[example]), INTERPOLATOR, INTERPOLATOR); for (int i = 0; i < n; i++) { double computedValue = bound.interpolate(XX[i]); double expectedValue = PPVAL.evaluate(pp, XX[i]).get(0); assertEquals(computedValue, expectedValue, 1e-14); double computedDerivative = bound.firstDerivative(XX[i]); double expectedDerivative = PPVAL.differentiate(pp, XX[i]).get(0); assertEquals(computedDerivative, expectedDerivative, 1e-14); } } }
final PiecewisePolynomialResult result = _method.interpolate(xValues, yValues); ArgChecker.isTrue(result.getOrder() == 4, "Primary interpolant is not cubic"); final double[] initialFirst = _function.differentiate(result, xValues).rowArray(0); firstWithSensitivity[0] = DoubleArray.copyOf(firstDerivativeCalculator(intervals, slopes, initialFirst)); final double[] initialFirstUp = _function.differentiate(_method.interpolate(xValues, yValuesUp), xValues).rowArray(0); final double[] initialFirstDw = _function.differentiate(_method.interpolate(xValues, yValuesDw), xValues).rowArray(0); final double[] firstUp = firstDerivativeCalculator(intervals, slopesUp, initialFirstUp); final double[] firstDw = firstDerivativeCalculator(intervals, slopesDw, initialFirstDw); ArgChecker.isTrue(resultWithSensitivity.getOrder() == 4, "Primary interpolant is not cubic"); final double[] initialFirst = _function.differentiate(resultWithSensitivity, xValues).rowArray(0); final DoubleArray[] initialFirstSense = _function.differentiateNodeSensitivity(resultWithSensitivity, xValues); firstWithSensitivity = firstDerivativeWithSensitivityCalculator(intervals, slopes, slopesSensitivity, slopesAbsSensitivity, initialFirst, initialFirstSense);
final double[] initialFirst = _function.differentiate(result, xValuesSrt).rowArray(0); final double[] first = firstDerivativeCalculator(intervals, slopes, initialFirst); final double[][] coefs = _solver.solve(yValuesSrt, intervals, slopes, first);
final PiecewisePolynomialResultsWithSensitivity resultDw = interps[k].interpolateWithSensitivity(xValues, yValuesDw); final double[] diffUp = FUNCTION.differentiateTwice(resultUp, xKeys).toArray()[0]; final double[] diffDw = FUNCTION.differentiateTwice(resultDw, xKeys).toArray()[0]; for (int i = 0; i < 10 * nData; ++i) { final double xKeyUp = xKeys[i] * (1. + EPS); final double resdiffNodeSensitivityXkeyUp = FUNCTION.differentiateNodeSensitivity(result, xKeyUp).get(j); final double resdiffNodeSensitivityXkeyDw = FUNCTION.differentiateNodeSensitivity(result, xKeyDw).get(j); final double senseFiniteXkey = 0.5 * (resdiffNodeSensitivityXkeyUp - resdiffNodeSensitivityXkeyDw) / EPS / xKeys[i]; final double resDiffTwiceNodeSensitivity = FUNCTION.differentiateTwiceNodeSensitivity(result, xKeys[i]).get(j);
final PiecewisePolynomialResultsWithSensitivity resultDw = interps[k].interpolateWithSensitivity(xValues, yValuesDw); final double[] valuesUp = FUNCTION.evaluate(resultUp, xKeys).rowArray(0); final double[] valuesDw = FUNCTION.evaluate(resultDw, xKeys).rowArray(0); final double[] diffUp = FUNCTION.differentiate(resultUp, xKeys).rowArray(0); final double[] diffDw = FUNCTION.differentiate(resultDw, xKeys).rowArray(0); for (int i = 0; i < 10 * nData; ++i) { final double xKeyUp = xKeys[i] * (1. + EPS); double valueFinite = 0.5 * (valuesUp[i] - valuesDw[i]) / EPS / yValues[j]; double senseFinite = 0.5 * (diffUp[i] - diffDw[i]) / EPS / yValues[j]; final double resNodeSensitivity = FUNCTION.nodeSensitivity(result, xKeys[i]).get(j); final double resNodeSensitivityXkeyUp = FUNCTION.nodeSensitivity(result, xKeyUp).get(j); final double resNodeSensitivityXkeyDw = FUNCTION.nodeSensitivity(result, xKeyDw).get(j); final double senseFiniteXkey = 0.5 * (resNodeSensitivityXkeyUp - resNodeSensitivityXkeyDw) / EPS / xKeys[i]; final double resDiffNodeSensitivity = FUNCTION.differentiateNodeSensitivity(result, xKeys[i]).get(j); assertEquals(valueFinite, resNodeSensitivity, Math.max(Math.abs(yValues[j]) * EPS, EPS)); assertEquals(senseFinite, resDiffNodeSensitivity, Math.max(Math.abs(yValues[j]) * EPS, EPS));
for (int j = 0; j < nKeys; ++j) { double key = xValues[0] + interval * j; InterpolatorTestUtil.assertRelative("notClampedTest", FUNC.evaluate(resultBase, key).get(0), FUNC .evaluate(result, key) .get(0), EPS); InterpolatorTestUtil.assertArrayRelative("notClampedTest", FUNC.nodeSensitivity(resultBase, key).toArray(), FUNC.nodeSensitivity(result, key).toArray(), EPS);
@Override protected DoubleArray doParameterSensitivity(double xValue) { ArgChecker.isTrue(Math.abs(xValue) > SMALL, "magnitude of xValue must not be small"); DoubleArray resSense = FUNCTION.nodeSensitivity(polySens.get(), xValue); return resSense.multipliedBy(DoubleArray.of(resSense.size(), i -> xValues[i] / xValue)); }
@Override protected double doInterpolate(double xValue) { ArgChecker.isTrue(Math.abs(xValue) > SMALL, "magnitude of xValue must not be small"); double resValue = FUNCTION.evaluate(poly, xValue).get(0); return resValue / xValue; }
final double[] initialFirst = _function.differentiate(resultWithSensitivity, xValues).rowArray(0); final double[][] slopeSensitivity = _solver.slopeSensitivityCalculator(intervals); final DoubleArray[] initialFirstSense = _function.differentiateNodeSensitivity(resultWithSensitivity, xValues); final DoubleArray[] firstWithSensitivity = firstDerivativeWithSensitivityCalculator(yValuesSrt, intervals, initialFirst, initialFirstSense); final DoubleMatrix[] resMatrix = _solver.solveWithSensitivity(yValuesSrt, intervals, slopes, slopeSensitivity, firstWithSensitivity);
final double[] initialFirst = _function.differentiate(result, xValuesSrt).rowArray(0); final double[] first = firstDerivativeCalculator(yValuesSrt, intervals, slopes, initialFirst); final double[][] coefs = _solver.solve(yValuesSrt, intervals, slopes, first);
final PiecewisePolynomialResultsWithSensitivity resultDw = interps[k].interpolateWithSensitivity(xValues, yValuesDw); final double[] diffUp = FUNCTION.differentiateTwice(resultUp, xKeys).toArray()[0]; final double[] diffDw = FUNCTION.differentiateTwice(resultDw, xKeys).toArray()[0]; for (int i = 0; i < 10 * nData; ++i) { final double xKeyUp = xKeys[i] * (1. + EPS); final double resdiffNodeSensitivityXkeyUp = FUNCTION.differentiateNodeSensitivity(result, xKeyUp).get(j); final double resdiffNodeSensitivityXkeyDw = FUNCTION.differentiateNodeSensitivity(result, xKeyDw).get(j); final double senseFiniteXkey = 0.5 * (resdiffNodeSensitivityXkeyUp - resdiffNodeSensitivityXkeyDw) / EPS / xKeys[i]; final double resDiffTwiceNodeSensitivity = FUNCTION.differentiateTwiceNodeSensitivity(result, xKeys[i]).get(j);
if (Math.abs(xValues[xValues.length - 1] - nodes[nIntervalsAll]) < EPS) { double lastNodeX = nodes[nIntervalsAll]; double lastNodeY = FUNC.evaluate(result, lastNodeX).get(0); double extraNode = 2.0 * nodes[nIntervalsAll] - nodes[nIntervalsAll - 1]; double extraDerivative = FUNC.differentiate(result, lastNodeX).get(0); double[] newKnots = new double[nIntervalsAll + 2]; System.arraycopy(nodes, 0, newKnots, 0, nIntervalsAll + 1); if (result instanceof PiecewisePolynomialResultsWithSensitivity) { PiecewisePolynomialResultsWithSensitivity resultCast = (PiecewisePolynomialResultsWithSensitivity) result; double[] extraSense = FUNC.nodeSensitivity(resultCast, lastNodeX).toArray(); double[] extraSenseDer = FUNC.differentiateNodeSensitivity(resultCast, lastNodeX).toArray(); DoubleMatrix[] newCoefSense = new DoubleMatrix[nIntervalsAll + 1]; for (int i = 0; i < nIntervalsAll; ++i) {
keys[j] = xValues[nData - 1] + j * interval; double[] values = FUNC.evaluate(result, keys).row(0).toArray(); for (int j = 2; j < nKeys; ++j) { InterpolatorTestUtil.assertRelative("linearExtrapolationTest", values[j - 1] - values[j - 2], values[j - 1] - keys[j] = xValues[nData - 1] + j * interval; double[] values = FUNC.evaluate(result, keys).row(0).toArray(); for (int j = 2; j < nKeys; ++j) { InterpolatorTestUtil.assertRelative("linearExtrapolationTest", values[j - 1] - values[j - 2], values[j - 1] - values[j - 2], EPS); DoubleArray[] sense = FUNC.nodeSensitivity(result, keys); for (int k = 0; k < nData; ++k) { double[] yValuesUp = Arrays.copyOf(yValues, nData); PiecewisePolynomialResultsWithSensitivity resultUp = interp.interpolateWithSensitivity(xValues, yValuesUp); PiecewisePolynomialResultsWithSensitivity resultDw = interp.interpolateWithSensitivity(xValues, yValuesDw); double[] tmpUp = FUNC.evaluate(resultUp, keys).rowArray(0); double[] tmpDw = FUNC.evaluate(resultDw, keys).rowArray(0); for (int l = 0; l < nKeys; ++l) { double res = 0.5 * (tmpUp[l] - tmpDw[l]) / DELTA; // lk
@Override protected DoubleArray doParameterSensitivity(double xValue) { ArgChecker.isTrue(Math.abs(xValue) > SMALL, "magnitude of xValue must not be small"); DoubleArray resSense = FUNCTION.nodeSensitivity(polySens.get(), xValue); return resSense.multipliedBy(DoubleArray.of(resSense.size(), i -> xValues[i] / xValue)); }
@Override protected double doInterpolate(double xValue) { ArgChecker.isTrue(Math.abs(xValue) > SMALL, "magnitude of xValue must not be small"); double resValue = FUNCTION.evaluate(poly, xValue).get(0); return resValue / xValue; }