From 27866408afbe440695c9f218d422f16175696c5a Mon Sep 17 00:00:00 2001 From: Christoph Ruegg Date: Sat, 21 Dec 2013 01:39:08 +0100 Subject: [PATCH] Drop old redundant spline interpolation classes --- .../Interpolation/AkimaSplineInterpolation.cs | 261 ------------ .../CubicHermiteSplineInterpolation.cs | 203 ---------- .../Interpolation/CubicSplineInterpolation.cs | 373 ------------------ .../Interpolation/SplineInterpolation.cs | 242 ------------ src/Numerics/Numerics.csproj | 4 - 5 files changed, 1083 deletions(-) delete mode 100644 src/Numerics/Interpolation/AkimaSplineInterpolation.cs delete mode 100644 src/Numerics/Interpolation/CubicHermiteSplineInterpolation.cs delete mode 100644 src/Numerics/Interpolation/CubicSplineInterpolation.cs delete mode 100644 src/Numerics/Interpolation/SplineInterpolation.cs diff --git a/src/Numerics/Interpolation/AkimaSplineInterpolation.cs b/src/Numerics/Interpolation/AkimaSplineInterpolation.cs deleted file mode 100644 index 40ebae5d..00000000 --- a/src/Numerics/Interpolation/AkimaSplineInterpolation.cs +++ /dev/null @@ -1,261 +0,0 @@ -// -// Math.NET Numerics, part of the Math.NET Project -// http://numerics.mathdotnet.com -// http://github.com/mathnet/mathnet-numerics -// http://mathnetnumerics.codeplex.com -// -// Copyright (c) 2009-2013 Math.NET -// -// Permission is hereby granted, free of charge, to any person -// obtaining a copy of this software and associated documentation -// files (the "Software"), to deal in the Software without -// restriction, including without limitation the rights to use, -// copy, modify, merge, publish, distribute, sublicense, and/or sell -// copies of the Software, and to permit persons to whom the -// Software is furnished to do so, subject to the following -// conditions: -// -// The above copyright notice and this permission notice shall be -// included in all copies or substantial portions of the Software. -// -// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, -// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES -// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND -// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT -// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, -// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING -// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR -// OTHER DEALINGS IN THE SOFTWARE. -// - -using System; -using System.Collections.Generic; -using MathNet.Numerics.Properties; - -namespace MathNet.Numerics.Interpolation -{ - /// - /// Akima Spline Interpolation Algorithm. - /// - /// - /// This algorithm supports both differentiation and integration. - /// - public class AkimaSplineInterpolation : IInterpolation - { - /// - /// Internal Spline Interpolation - /// - readonly CubicHermiteSplineInterpolation _spline; - - /// - /// Initializes a new instance of the AkimaSplineInterpolation class. - /// - public AkimaSplineInterpolation() - { - _spline = new CubicHermiteSplineInterpolation(); - } - - /// - /// Initializes a new instance of the AkimaSplineInterpolation class. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - public AkimaSplineInterpolation(IList samplePoints, IList sampleValues) - { - _spline = new CubicHermiteSplineInterpolation(); - Initialize(samplePoints, sampleValues); - } - - /// - /// Gets a value indicating whether the algorithm supports differentiation (interpolated derivative). - /// - bool IInterpolation.SupportsDifferentiation - { - get { return true; } - } - - /// - /// Gets a value indicating whether the algorithm supports integration (interpolated quadrature). - /// - bool IInterpolation.SupportsIntegration - { - get { return true; } - } - - /// - /// Initialize the interpolation method with the given spline coefficients (sorted by the sample points t). - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - public void Initialize(IList samplePoints, IList sampleValues) - { - double[] derivatives = EvaluateSplineDerivatives(samplePoints, sampleValues); - _spline.Initialize(samplePoints, sampleValues, derivatives); - } - - /// - /// Evaluate the spline derivatives as used - /// internally by this interpolation algorithm. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Spline Derivative Vector - public static double[] EvaluateSplineDerivatives(IList samplePoints, IList sampleValues) - { - if (null == samplePoints) - { - throw new ArgumentNullException("samplePoints"); - } - - if (null == sampleValues) - { - throw new ArgumentNullException("sampleValues"); - } - - if (samplePoints.Count < 5) - { - throw new ArgumentOutOfRangeException("samplePoints"); - } - - if (samplePoints.Count != sampleValues.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength); - } - - for (var i = 1; i < samplePoints.Count; ++i) - if (samplePoints[i] <= samplePoints[i - 1]) - throw new ArgumentException(Resources.Interpolation_Initialize_SamplePointsNotStrictlyAscendingOrder, "samplePoints"); - - int n = samplePoints.Count; - - /* Prepare divided differences (diff) and weights (w) */ - - var differences = new double[n - 1]; - var weights = new double[n - 1]; - - for (int i = 0; i < differences.Length; i++) - { - differences[i] = (sampleValues[i + 1] - sampleValues[i])/(samplePoints[i + 1] - samplePoints[i]); - } - - for (int i = 1; i < weights.Length; i++) - { - weights[i] = Math.Abs(differences[i] - differences[i - 1]); - } - - /* Prepare Hermite interpolation scheme */ - - var derivatives = new double[n]; - - for (int i = 2; i < derivatives.Length - 2; i++) - { - derivatives[i] = - weights[i - 1].AlmostEqual(0.0) && weights[i + 1].AlmostEqual(0.0) - ? (((samplePoints[i + 1] - samplePoints[i])*differences[i - 1]) - + ((samplePoints[i] - samplePoints[i - 1])*differences[i])) - /(samplePoints[i + 1] - samplePoints[i - 1]) - : ((weights[i + 1]*differences[i - 1]) - + (weights[i - 1]*differences[i])) - /(weights[i + 1] + weights[i - 1]); - } - - derivatives[0] = DifferentiateThreePoint(samplePoints, sampleValues, 0, 0, 1, 2); - derivatives[1] = DifferentiateThreePoint(samplePoints, sampleValues, 1, 0, 1, 2); - derivatives[n - 2] = DifferentiateThreePoint(samplePoints, sampleValues, n - 2, n - 3, n - 2, n - 1); - derivatives[n - 1] = DifferentiateThreePoint(samplePoints, sampleValues, n - 1, n - 3, n - 2, n - 1); - - /* Build Akima spline using Hermite interpolation scheme */ - - return derivatives; - } - - /// - /// Evaluate the spline coefficients as used - /// internally by this interpolation algorithm. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Spline Coefficient Vector - public static double[] EvaluateSplineCoefficients(IList samplePoints, IList sampleValues) - { - double[] derivatives = EvaluateSplineDerivatives(samplePoints, sampleValues); - return CubicHermiteSplineInterpolation.EvaluateSplineCoefficients(samplePoints, sampleValues, derivatives); - } - - /// - /// Three-Point Differentiation Helper. - /// - /// Sample Points t. - /// Sample Values x(t). - /// Index of the point of the differentiation. - /// Index of the first sample. - /// Index of the second sample. - /// Index of the third sample. - /// The derivative approximation. - static double DifferentiateThreePoint( - IList samplePoints, IList sampleValues, - int indexT, int index0, int index1, int index2) - { - double x0 = sampleValues[index0]; - double x1 = sampleValues[index1]; - double x2 = sampleValues[index2]; - - double t = samplePoints[indexT] - samplePoints[index0]; - double t1 = samplePoints[index1] - samplePoints[index0]; - double t2 = samplePoints[index2] - samplePoints[index0]; - - double a = (x2 - x0 - (t2/t1*(x1 - x0)))/((t2*t2) - (t1*t2)); - double b = (x1 - x0 - (a*t1*t1))/t1; - return (2*a*t) + b; - } - - /// - /// Interpolate at point t. - /// - /// Point t to interpolate at. - /// Interpolated value x(t). - public double Interpolate(double t) - { - return _spline.Interpolate(t); - } - - /// - /// Differentiate at point t. - /// - /// Point t to interpolate at. - /// Interpolated first derivative at point t. - public double Differentiate(double t) - { - return _spline.Differentiate(t); - } - - /// - /// Differentiate twice at point t. - /// - /// Point t to interpolate at. - /// Interpolated second derivative at point t. - public double Differentiate2(double t) - { - return _spline.Differentiate2(t); - } - - /// - /// Indefinite integral at point t. - /// - /// Point t to integrate at. - public double Integrate(double t) - { - return _spline.Integrate(t); - } - - /// - /// Definite integral between points a and b. - /// - /// Left bound of the integration interval [a,b]. - /// Right bound of the integration interval [a,b]. - public double Integrate(double a, double b) - { - return _spline.Integrate(a, b); - } - } -} diff --git a/src/Numerics/Interpolation/CubicHermiteSplineInterpolation.cs b/src/Numerics/Interpolation/CubicHermiteSplineInterpolation.cs deleted file mode 100644 index ee54263f..00000000 --- a/src/Numerics/Interpolation/CubicHermiteSplineInterpolation.cs +++ /dev/null @@ -1,203 +0,0 @@ -// -// Math.NET Numerics, part of the Math.NET Project -// http://numerics.mathdotnet.com -// http://github.com/mathnet/mathnet-numerics -// http://mathnetnumerics.codeplex.com -// -// Copyright (c) 2009-2013 Math.NET -// -// Permission is hereby granted, free of charge, to any person -// obtaining a copy of this software and associated documentation -// files (the "Software"), to deal in the Software without -// restriction, including without limitation the rights to use, -// copy, modify, merge, publish, distribute, sublicense, and/or sell -// copies of the Software, and to permit persons to whom the -// Software is furnished to do so, subject to the following -// conditions: -// -// The above copyright notice and this permission notice shall be -// included in all copies or substantial portions of the Software. -// -// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, -// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES -// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND -// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT -// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, -// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING -// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR -// OTHER DEALINGS IN THE SOFTWARE. -// - -using System; -using System.Collections.Generic; -using MathNet.Numerics.Properties; - -namespace MathNet.Numerics.Interpolation -{ - /// - /// Cubic Hermite Spline Interpolation Algorithm. - /// - /// - /// This algorithm supports both differentiation and integration. - /// - public class CubicHermiteSplineInterpolation : IInterpolation - { - /// - /// Internal Spline Interpolation - /// - readonly SplineInterpolation _spline; - - /// - /// Initializes a new instance of the CubicHermiteSplineInterpolation class. - /// - public CubicHermiteSplineInterpolation() - { - _spline = new SplineInterpolation(); - } - - /// - /// Initializes a new instance of the CubicHermiteSplineInterpolation class. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Sample Derivatives x'(t) - public CubicHermiteSplineInterpolation(IList samplePoints, IList sampleValues, IList sampleDerivatives) - { - _spline = new SplineInterpolation(); - Initialize(samplePoints, sampleValues, sampleDerivatives); - } - - /// - /// Gets a value indicating whether the algorithm supports differentiation (interpolated derivative). - /// - bool IInterpolation.SupportsDifferentiation - { - get { return true; } - } - - /// - /// Gets a value indicating whether the algorithm supports integration (interpolated quadrature). - /// - bool IInterpolation.SupportsIntegration - { - get { return true; } - } - - /// - /// Initialize the interpolation method with the given spline coefficients (sorted by the sample points t). - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Sample Derivatives x'(t) - public void Initialize(IList samplePoints, IList sampleValues, IList sampleDerivatives) - { - double[] coefficients = EvaluateSplineCoefficients(samplePoints, sampleValues, sampleDerivatives); - _spline.Initialize(samplePoints, coefficients); - } - - /// - /// Evaluate the spline coefficients as used - /// internally by this interpolation algorithm. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Sample Derivatives x'(t) - /// Spline Coefficient Vector - public static double[] EvaluateSplineCoefficients(IList samplePoints, IList sampleValues, IList sampleDerivatives) - { - if (null == samplePoints) - { - throw new ArgumentNullException("samplePoints"); - } - - if (null == sampleValues) - { - throw new ArgumentNullException("sampleValues"); - } - - if (null == sampleDerivatives) - { - throw new ArgumentNullException("sampleDerivatives"); - } - - if (samplePoints.Count < 2) - { - throw new ArgumentOutOfRangeException("samplePoints"); - } - - if (samplePoints.Count != sampleValues.Count - || samplePoints.Count != sampleDerivatives.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength); - } - - for (var i = 1; i < samplePoints.Count; ++i) - if (samplePoints[i] <= samplePoints[i - 1]) - throw new ArgumentException(Resources.Interpolation_Initialize_SamplePointsNotStrictlyAscendingOrder, "samplePoints"); - - var coefficients = new double[4*(samplePoints.Count - 1)]; - - for (int i = 0, j = 0; i < samplePoints.Count - 1; i++, j += 4) - { - double delta = samplePoints[i + 1] - samplePoints[i]; - double delta2 = delta*delta; - double delta3 = delta*delta2; - coefficients[j] = sampleValues[i]; - coefficients[j + 1] = sampleDerivatives[i]; - coefficients[j + 2] = ((3*(sampleValues[i + 1] - sampleValues[i])) - (2*sampleDerivatives[i]*delta) - (sampleDerivatives[i + 1]*delta))/delta2; - coefficients[j + 3] = ((2*(sampleValues[i] - sampleValues[i + 1])) + (sampleDerivatives[i]*delta) + (sampleDerivatives[i + 1]*delta))/delta3; - } - - return coefficients; - } - - /// - /// Interpolate at point t. - /// - /// Point t to interpolate at. - /// Interpolated value x(t). - public double Interpolate(double t) - { - return _spline.Interpolate(t); - } - - /// - /// Differentiate at point t. - /// - /// Point t to interpolate at. - /// Interpolated first derivative at point t. - public double Differentiate(double t) - { - return _spline.Differentiate(t); - } - - /// - /// Differentiate twice at point t. - /// - /// Point t to interpolate at. - /// Interpolated second derivative at point t. - public double Differentiate2(double t) - { - return _spline.Differentiate2(t); - } - - /// - /// Indefinite integral at point t. - /// - /// Point t to integrate at. - public double Integrate(double t) - { - return _spline.Integrate(t); - } - - /// - /// Definite integral between points a and b. - /// - /// Left bound of the integration interval [a,b]. - /// Right bound of the integration interval [a,b]. - public double Integrate(double a, double b) - { - return _spline.Integrate(a, b); - } - } -} diff --git a/src/Numerics/Interpolation/CubicSplineInterpolation.cs b/src/Numerics/Interpolation/CubicSplineInterpolation.cs deleted file mode 100644 index 220d748f..00000000 --- a/src/Numerics/Interpolation/CubicSplineInterpolation.cs +++ /dev/null @@ -1,373 +0,0 @@ -// -// Math.NET Numerics, part of the Math.NET Project -// http://numerics.mathdotnet.com -// http://github.com/mathnet/mathnet-numerics -// http://mathnetnumerics.codeplex.com -// -// Copyright (c) 2009-2013 Math.NET -// -// Permission is hereby granted, free of charge, to any person -// obtaining a copy of this software and associated documentation -// files (the "Software"), to deal in the Software without -// restriction, including without limitation the rights to use, -// copy, modify, merge, publish, distribute, sublicense, and/or sell -// copies of the Software, and to permit persons to whom the -// Software is furnished to do so, subject to the following -// conditions: -// -// The above copyright notice and this permission notice shall be -// included in all copies or substantial portions of the Software. -// -// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, -// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES -// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND -// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT -// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, -// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING -// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR -// OTHER DEALINGS IN THE SOFTWARE. -// - -using System; -using System.Collections.Generic; -using MathNet.Numerics.Properties; - -namespace MathNet.Numerics.Interpolation -{ - /// - /// Cubic Spline Interpolation Algorithm with continuous first and second derivatives. - /// - /// - /// This algorithm supports both differentiation and integration. - /// - public class CubicSplineInterpolation : IInterpolation - { - /// - /// Internal Spline Interpolation - /// - readonly CubicHermiteSplineInterpolation _spline; - - /// - /// Initializes a new instance of the CubicSplineInterpolation class. - /// - public CubicSplineInterpolation() - { - _spline = new CubicHermiteSplineInterpolation(); - } - - /// - /// Initializes a new instance of the CubicSplineInterpolation class. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - public CubicSplineInterpolation(IList samplePoints, IList sampleValues) - { - _spline = new CubicHermiteSplineInterpolation(); - Initialize(samplePoints, sampleValues); - } - - /// - /// Initializes a new instance of the CubicSplineInterpolation class. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Condition of the left boundary. - /// Left boundary value. Ignored in the parabolic case. - /// Condition of the right boundary. - /// Right boundary value. Ignored in the parabolic case. - public CubicSplineInterpolation( - IList samplePoints, IList sampleValues, - SplineBoundaryCondition leftBoundaryCondition, double leftBoundary, - SplineBoundaryCondition rightBoundaryCondition, double rightBoundary) - { - _spline = new CubicHermiteSplineInterpolation(); - - Initialize( - samplePoints, sampleValues, - leftBoundaryCondition, leftBoundary, - rightBoundaryCondition, rightBoundary); - } - - /// - /// Gets a value indicating whether the algorithm supports differentiation (interpolated derivative). - /// - bool IInterpolation.SupportsDifferentiation - { - get { return true; } - } - - /// - /// Gets a value indicating whether the algorithm supports integration (interpolated quadrature). - /// - bool IInterpolation.SupportsIntegration - { - get { return true; } - } - - /// - /// Initialize the interpolation method with the given spline coefficients (sorted by the sample points t). - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - public void Initialize(IList samplePoints, IList sampleValues) - { - double[] derivatives = EvaluateSplineDerivatives(samplePoints, sampleValues, - SplineBoundaryCondition.SecondDerivative, 0.0, - SplineBoundaryCondition.SecondDerivative, 0.0); - _spline.Initialize(samplePoints, sampleValues, derivatives); - } - - /// - /// Initialize the interpolation method with the given spline coefficients (sorted by the sample points t). - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Condition of the left boundary. - /// Left boundary value. Ignored in the parabolic case. - /// Condition of the right boundary. - /// Right boundary value. Ignored in the parabolic case. - public void Initialize( - IList samplePoints, IList sampleValues, - SplineBoundaryCondition leftBoundaryCondition, double leftBoundary, - SplineBoundaryCondition rightBoundaryCondition, double rightBoundary) - { - double[] derivatives = EvaluateSplineDerivatives( - samplePoints, sampleValues, - leftBoundaryCondition, leftBoundary, - rightBoundaryCondition, rightBoundary); - _spline.Initialize(samplePoints, sampleValues, derivatives); - } - - /// - /// Evaluate the spline derivatives as used - /// internally by this interpolation algorithm. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Condition of the left boundary. - /// Left boundary value. Ignored in the parabolic case. - /// Condition of the right boundary. - /// Right boundary value. Ignored in the parabolic case. - /// Spline Derivative Vector - public static double[] EvaluateSplineDerivatives( - IList samplePoints, IList sampleValues, - SplineBoundaryCondition leftBoundaryCondition, double leftBoundary, - SplineBoundaryCondition rightBoundaryCondition, double rightBoundary) - { - if (null == samplePoints) - { - throw new ArgumentNullException("samplePoints"); - } - - if (null == sampleValues) - { - throw new ArgumentNullException("sampleValues"); - } - - if (samplePoints.Count < 2) - { - throw new ArgumentOutOfRangeException("samplePoints"); - } - - if (samplePoints.Count != sampleValues.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength); - } - - for (var i = 1; i < samplePoints.Count; ++i) - if (samplePoints[i] <= samplePoints[i - 1]) - throw new ArgumentException(Resources.Interpolation_Initialize_SamplePointsNotStrictlyAscendingOrder, "samplePoints"); - - int n = samplePoints.Count; - - // normalize special cases - if ((n == 2) - && (leftBoundaryCondition == SplineBoundaryCondition.ParabolicallyTerminated) - && (rightBoundaryCondition == SplineBoundaryCondition.ParabolicallyTerminated)) - { - leftBoundaryCondition = SplineBoundaryCondition.SecondDerivative; - leftBoundary = 0d; - rightBoundaryCondition = SplineBoundaryCondition.SecondDerivative; - rightBoundary = 0d; - } - - if (leftBoundaryCondition == SplineBoundaryCondition.Natural) - { - leftBoundaryCondition = SplineBoundaryCondition.SecondDerivative; - leftBoundary = 0d; - } - - if (rightBoundaryCondition == SplineBoundaryCondition.Natural) - { - rightBoundaryCondition = SplineBoundaryCondition.SecondDerivative; - rightBoundary = 0d; - } - - var a1 = new double[n]; - var a2 = new double[n]; - var a3 = new double[n]; - var b = new double[n]; - - // Left Boundary - switch (leftBoundaryCondition) - { - case SplineBoundaryCondition.ParabolicallyTerminated: - a1[0] = 0; - a2[0] = 1; - a3[0] = 1; - b[0] = 2*(sampleValues[1] - sampleValues[0])/(samplePoints[1] - samplePoints[0]); - break; - case SplineBoundaryCondition.FirstDerivative: - a1[0] = 0; - a2[0] = 1; - a3[0] = 0; - b[0] = leftBoundary; - break; - case SplineBoundaryCondition.SecondDerivative: - a1[0] = 0; - a2[0] = 2; - a3[0] = 1; - b[0] = (3*((sampleValues[1] - sampleValues[0])/(samplePoints[1] - samplePoints[0]))) - (0.5*leftBoundary*(samplePoints[1] - samplePoints[0])); - break; - default: - throw new NotSupportedException(Resources.InvalidLeftBoundaryCondition); - } - - // Central Conditions - for (int i = 1; i < samplePoints.Count - 1; i++) - { - a1[i] = samplePoints[i + 1] - samplePoints[i]; - a2[i] = 2*(samplePoints[i + 1] - samplePoints[i - 1]); - a3[i] = samplePoints[i] - samplePoints[i - 1]; - b[i] = (3*(sampleValues[i] - sampleValues[i - 1])/(samplePoints[i] - samplePoints[i - 1])*(samplePoints[i + 1] - samplePoints[i])) + (3*(sampleValues[i + 1] - sampleValues[i])/(samplePoints[i + 1] - samplePoints[i])*(samplePoints[i] - samplePoints[i - 1])); - } - - // Right Boundary - switch (rightBoundaryCondition) - { - case SplineBoundaryCondition.ParabolicallyTerminated: - a1[n - 1] = 1; - a2[n - 1] = 1; - a3[n - 1] = 0; - b[n - 1] = 2*(sampleValues[n - 1] - sampleValues[n - 2])/(samplePoints[n - 1] - samplePoints[n - 2]); - break; - case SplineBoundaryCondition.FirstDerivative: - a1[n - 1] = 0; - a2[n - 1] = 1; - a3[n - 1] = 0; - b[n - 1] = rightBoundary; - break; - case SplineBoundaryCondition.SecondDerivative: - a1[n - 1] = 1; - a2[n - 1] = 2; - a3[n - 1] = 0; - b[n - 1] = (3*(sampleValues[n - 1] - sampleValues[n - 2])/(samplePoints[n - 1] - samplePoints[n - 2])) + (0.5*rightBoundary*(samplePoints[n - 1] - samplePoints[n - 2])); - break; - default: - throw new NotSupportedException(Resources.InvalidRightBoundaryCondition); - } - - // Build Spline - return SolveTridiagonal(a1, a2, a3, b); - } - - /// - /// Evaluate the spline coefficients as used - /// internally by this interpolation algorithm. - /// - /// Sample Points t, sorted ascending. - /// Sample Values x(t) - /// Condition of the left boundary. - /// Left boundary value. Ignored in the parabolic case. - /// Condition of the right boundary. - /// Right boundary value. Ignored in the parabolic case. - /// Spline Coefficient Vector - public static double[] EvaluateSplineCoefficients( - IList samplePoints, IList sampleValues, - SplineBoundaryCondition leftBoundaryCondition, double leftBoundary, - SplineBoundaryCondition rightBoundaryCondition, double rightBoundary) - { - double[] derivatives = EvaluateSplineDerivatives( - samplePoints, sampleValues, - leftBoundaryCondition, leftBoundary, - rightBoundaryCondition, rightBoundary); - return CubicHermiteSplineInterpolation.EvaluateSplineCoefficients(samplePoints, sampleValues, derivatives); - } - - /// - /// Tridiagonal Solve Helper. - /// - /// The a-vector[n]. - /// The b-vector[n], will be modified by this function. - /// The c-vector[n]. - /// The d-vector[n], will be modified by this function. - /// The x-vector[n] - static double[] SolveTridiagonal(double[] a, double[] b, double[] c, double[] d) - { - for (int k = 1; k < a.Length; k++) - { - double t = a[k]/b[k - 1]; - b[k] = b[k] - (t*c[k - 1]); - d[k] = d[k] - (t*d[k - 1]); - } - - var x = new double[a.Length]; - x[x.Length - 1] = d[d.Length - 1]/b[b.Length - 1]; - for (int k = x.Length - 2; k >= 0; k--) - { - x[k] = (d[k] - (c[k]*x[k + 1]))/b[k]; - } - - return x; - } - - /// - /// Interpolate at point t. - /// - /// Point t to interpolate at. - /// Interpolated value x(t). - public double Interpolate(double t) - { - return _spline.Interpolate(t); - } - - /// - /// Differentiate at point t. - /// - /// Point t to interpolate at. - /// Interpolated first derivative at point t. - public double Differentiate(double t) - { - return _spline.Differentiate(t); - } - - /// - /// Differentiate twice at point t. - /// - /// Point t to interpolate at. - /// Interpolated second derivative at point t. - public double Differentiate2(double t) - { - return _spline.Differentiate2(t); - } - - /// - /// Indefinite integral at point t. - /// - /// Point t to integrate at. - public double Integrate(double t) - { - return _spline.Integrate(t); - } - - /// - /// Definite integral between points a and b. - /// - /// Left bound of the integration interval [a,b]. - /// Right bound of the integration interval [a,b]. - public double Integrate(double a, double b) - { - return _spline.Integrate(a, b); - } - } -} diff --git a/src/Numerics/Interpolation/SplineInterpolation.cs b/src/Numerics/Interpolation/SplineInterpolation.cs deleted file mode 100644 index 537d1117..00000000 --- a/src/Numerics/Interpolation/SplineInterpolation.cs +++ /dev/null @@ -1,242 +0,0 @@ -// -// Math.NET Numerics, part of the Math.NET Project -// http://numerics.mathdotnet.com -// http://github.com/mathnet/mathnet-numerics -// http://mathnetnumerics.codeplex.com -// -// Copyright (c) 2009-2013 Math.NET -// -// Permission is hereby granted, free of charge, to any person -// obtaining a copy of this software and associated documentation -// files (the "Software"), to deal in the Software without -// restriction, including without limitation the rights to use, -// copy, modify, merge, publish, distribute, sublicense, and/or sell -// copies of the Software, and to permit persons to whom the -// Software is furnished to do so, subject to the following -// conditions: -// -// The above copyright notice and this permission notice shall be -// included in all copies or substantial portions of the Software. -// -// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, -// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES -// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND -// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT -// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, -// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING -// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR -// OTHER DEALINGS IN THE SOFTWARE. -// - -using System; -using System.Collections.Generic; -using MathNet.Numerics.Properties; - -namespace MathNet.Numerics.Interpolation -{ - /// - /// Third-Degree Spline Interpolation Algorithm. - /// - /// - /// This algorithm supports both differentiation and integration. - /// - public class SplineInterpolation : IInterpolation - { - /// - /// Sample Points t. - /// - IList _points; - - /// - /// Spline Coefficients c(t). - /// - IList _coefficients; - - /// - /// Number of samples. - /// - int _sampleCount; - - /// - /// Initializes a new instance of the SplineInterpolation class. - /// - public SplineInterpolation() - { - } - - /// - /// Initializes a new instance of the SplineInterpolation class. - /// - /// Sample Points t (length: N), sorted ascending. - /// Spline Coefficients (length: 4*(N-1)). - public SplineInterpolation(IList samplePoints, IList splineCoefficients) - { - Initialize(samplePoints, splineCoefficients); - } - - /// - /// Gets a value indicating whether the algorithm supports differentiation (interpolated derivative). - /// - bool IInterpolation.SupportsDifferentiation - { - get { return true; } - } - - /// - /// Gets a value indicating whether the algorithm supports integration (interpolated quadrature). - /// - bool IInterpolation.SupportsIntegration - { - get { return true; } - } - - /// - /// Initialize the interpolation method with the given spline coefficients (sorted by the sample points t). - /// - /// Sample Points t (length: N), sorted ascending. - /// Spline Coefficients (length: 4*(N-1)). - public void Initialize(IList samplePoints, IList splineCoefficients) - { - if (null == samplePoints) - { - throw new ArgumentNullException("samplePoints"); - } - - if (null == splineCoefficients) - { - throw new ArgumentNullException("splineCoefficients"); - } - - if (samplePoints.Count < 2) - { - throw new ArgumentOutOfRangeException("samplePoints"); - } - - if (splineCoefficients.Count != 4*(samplePoints.Count - 1)) - { - throw new ArgumentOutOfRangeException("splineCoefficients"); - } - - for (var i = 1; i < samplePoints.Count; ++i) - if (samplePoints[i] <= samplePoints[i - 1]) - throw new ArgumentException(Resources.Interpolation_Initialize_SamplePointsNotStrictlyAscendingOrder, "samplePoints"); - - _points = samplePoints; - _coefficients = splineCoefficients; - _sampleCount = samplePoints.Count; - } - - /// - /// Interpolate at point t. - /// - /// Point t to interpolate at. - /// Interpolated value x(t). - public double Interpolate(double t) - { - int closestLeftIndex = LeftBracketIndex(t); - - // Interpolation - double offset = t - _points[closestLeftIndex]; - int k = closestLeftIndex << 2; - - return _coefficients[k] - + (offset*(_coefficients[k + 1] - + (offset*(_coefficients[k + 2] - + (offset*_coefficients[k + 3]))))); - } - - /// - /// Differentiate at point t. - /// - /// Point t to interpolate at. - /// Interpolated first derivative at point t. - public double Differentiate(double t) - { - int closestLeftIndex = LeftBracketIndex(t); - double offset = t - _points[closestLeftIndex]; - int k = closestLeftIndex << 2; - - return _coefficients[k + 1] - + (2*offset*_coefficients[k + 2]) - + (3*offset*offset*_coefficients[k + 3]); - } - - /// - /// Differentiate twice at point t. - /// - /// Point t to interpolate at. - /// Interpolated second derivative at point t. - public double Differentiate2(double t) - { - int closestLeftIndex = LeftBracketIndex(t); - double offset = t - _points[closestLeftIndex]; - int k = closestLeftIndex << 2; - - return (2*_coefficients[k + 2]) + (6*offset*_coefficients[k + 3]); - } - - /// - /// Indefinite integral at point t. - /// - /// Point t to integrate at. - public double Integrate(double t) - { - int closestLeftIndex = LeftBracketIndex(t); - - // Integration - double result = 0; - for (int i = 0, j = 0; i < closestLeftIndex; i++, j += 4) - { - double w = _points[i + 1] - _points[i]; - result += w*(_coefficients[j] - + ((w*_coefficients[j + 1]*0.5) - + (w*((_coefficients[j + 2]/3) - + (w*_coefficients[j + 3]*0.25))))); - } - - double offset = t - _points[closestLeftIndex]; - int k = closestLeftIndex << 2; - - return result + (offset*(_coefficients[k] - + (offset*_coefficients[k + 1]*0.5) - + (offset*_coefficients[k + 2]/3) - + (offset*_coefficients[k + 3]*0.25))); - } - - /// - /// Definite integral between points a and b. - /// - /// Left bound of the integration interval [a,b]. - /// Right bound of the integration interval [a,b]. - public double Integrate(double a, double b) - { - return Integrate(b) - Integrate(a); - } - - /// - /// Find the index of the greatest sample point smaller than t. - /// - /// The value to look for. - /// The sample point index. - int LeftBracketIndex(double t) - { - // Binary search in the [ t[0], ..., t[n-2] ] (t[n-1] is not included) - int low = 0; - int high = _sampleCount - 1; - while (low != high - 1) - { - int middle = (low + high)/2; - if (_points[middle] > t) - { - high = middle; - } - else - { - low = middle; - } - } - - return low; - } - } -} diff --git a/src/Numerics/Numerics.csproj b/src/Numerics/Numerics.csproj index 70003765..0057f245 100644 --- a/src/Numerics/Numerics.csproj +++ b/src/Numerics/Numerics.csproj @@ -379,15 +379,11 @@ - - - -