diff --git a/src/Numerics/Interpolation/Barycentric.cs b/src/Numerics/Interpolation/Barycentric.cs index 90b3418d7..2bb50763e 100644 --- a/src/Numerics/Interpolation/Barycentric.cs +++ b/src/Numerics/Interpolation/Barycentric.cs @@ -334,6 +334,13 @@ public double Interpolate(double t) /// Interpolated second derivative at point t. double IInterpolation.Differentiate2(double t) => throw new NotSupportedException(); + /// + /// Differentiate three times at point t. NOT SUPPORTED. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + double IInterpolation.Differentiate3(double t) => throw new NotSupportedException(); + /// /// Indefinite integral at point t. NOT SUPPORTED. /// diff --git a/src/Numerics/Interpolation/BulirschStoerRationalInterpolation.cs b/src/Numerics/Interpolation/BulirschStoerRationalInterpolation.cs index 9dfeb767a..ec84aca90 100644 --- a/src/Numerics/Interpolation/BulirschStoerRationalInterpolation.cs +++ b/src/Numerics/Interpolation/BulirschStoerRationalInterpolation.cs @@ -182,6 +182,13 @@ public double Interpolate(double t) /// Interpolated second derivative at point t. double IInterpolation.Differentiate2(double t) => throw new NotSupportedException(); + /// + /// Differentiate three times at point t. NOT SUPPORTED. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + double IInterpolation.Differentiate3(double t) => throw new NotSupportedException(); + /// /// Indefinite integral at point t. NOT SUPPORTED. /// diff --git a/src/Numerics/Interpolation/CubicSpline.cs b/src/Numerics/Interpolation/CubicSpline.cs index d1f87e46c..f09b3c48a 100644 --- a/src/Numerics/Interpolation/CubicSpline.cs +++ b/src/Numerics/Interpolation/CubicSpline.cs @@ -576,7 +576,18 @@ public double Differentiate2(double t) { int k = LeftSegmentIndex(t); var x = t - _x[k]; - return 2*_c2[k] + x*6*_c3[k]; + return 2 * _c2[k] + x * 6 * _c3[k]; + } + + /// + /// Differentiate three times at point t. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + public double Differentiate3(double t) + { + int k = LeftSegmentIndex(t); + return 6 * _c3[k]; } /// diff --git a/src/Numerics/Interpolation/IInterpolation.cs b/src/Numerics/Interpolation/IInterpolation.cs index 03c4e66c3..60659dbf6 100644 --- a/src/Numerics/Interpolation/IInterpolation.cs +++ b/src/Numerics/Interpolation/IInterpolation.cs @@ -65,6 +65,13 @@ public interface IInterpolation /// Interpolated second derivative at point t. double Differentiate2(double t); + /// + /// Differentiate three times at point t. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + double Differentiate3(double t); + /// /// Indefinite integral at point t. /// diff --git a/src/Numerics/Interpolation/LinearSpline.cs b/src/Numerics/Interpolation/LinearSpline.cs index 0f192c3e3..eecddb43f 100644 --- a/src/Numerics/Interpolation/LinearSpline.cs +++ b/src/Numerics/Interpolation/LinearSpline.cs @@ -152,6 +152,13 @@ public double Differentiate(double t) /// Interpolated second derivative at point t. public double Differentiate2(double t) => 0d; + /// + /// Differentiate three times at point t. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + public double Differentiate3(double t) => 0d; + /// /// Indefinite integral at point t. /// diff --git a/src/Numerics/Interpolation/LogLinear.cs b/src/Numerics/Interpolation/LogLinear.cs index e40702e08..156923b61 100644 --- a/src/Numerics/Interpolation/LogLinear.cs +++ b/src/Numerics/Interpolation/LogLinear.cs @@ -152,6 +152,24 @@ public double Differentiate2(double t) return secondDerivative; } + /// + /// Differentiate three times at point t. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + public double Differentiate3(double t) + { + var linearFirstDerivative = _spline.Differentiate(t); + var linearSecondDerivative = _spline.Differentiate2(t); + var linearThirdDerivative = _spline.Differentiate3(t); + + var thirdDerivative = Differentiate2(t) * linearFirstDerivative + + 2 * Differentiate(t) * linearSecondDerivative + + Interpolate(t) * linearThirdDerivative; + + return thirdDerivative; + } + /// /// Indefinite integral at point t. /// diff --git a/src/Numerics/Interpolation/NevillePolynomialInterpolation.cs b/src/Numerics/Interpolation/NevillePolynomialInterpolation.cs index 389ea8456..3be6a10a9 100644 --- a/src/Numerics/Interpolation/NevillePolynomialInterpolation.cs +++ b/src/Numerics/Interpolation/NevillePolynomialInterpolation.cs @@ -197,6 +197,36 @@ public double Differentiate2(double t) return ddx[0]; } + /// + /// Differentiate three times at point t. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + public double Differentiate3(double t) + { + var x = new double[_y.Length]; + var dx = new double[_y.Length]; + var ddx = new double[_y.Length]; + var dddx = new double[_y.Length]; + _y.CopyTo(x, 0); + + for (int level = 1; level < x.Length; level++) + { + for (int i = 0; i < x.Length - level; i++) + { + double hp = t - _x[i + level]; + double ho = _x[i] - t; + double den = _x[i] - _x[i + level]; + dddx[i] = ((hp * dddx[i]) + (3 * ddx[i]) + (ho * dddx[i + 1]) - (3 * ddx[i + 1])) / den; + ddx[i] = ((hp * ddx[i]) + (2 * dx[i]) + (ho * ddx[i + 1]) - (2 * dx[i + 1])) / den; + dx[i] = ((hp * dx[i]) + x[i] + (ho * dx[i + 1]) - x[i + 1]) / den; + x[i] = ((hp * x[i]) + (ho * x[i + 1])) / den; + } + } + + return dddx[0]; + } + /// /// Indefinite integral at point t. NOT SUPPORTED. /// diff --git a/src/Numerics/Interpolation/QuadraticSpline.cs b/src/Numerics/Interpolation/QuadraticSpline.cs index f25edf55b..491931915 100644 --- a/src/Numerics/Interpolation/QuadraticSpline.cs +++ b/src/Numerics/Interpolation/QuadraticSpline.cs @@ -110,6 +110,13 @@ public double Differentiate2(double t) return 2*_c2[k]; } + /// + /// Differentiate three times at point t. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + public double Differentiate3(double t) => 0d; + /// /// Indefinite integral at point t. /// diff --git a/src/Numerics/Interpolation/StepInterpolation.cs b/src/Numerics/Interpolation/StepInterpolation.cs index c439ae87a..0d33bfbec 100644 --- a/src/Numerics/Interpolation/StepInterpolation.cs +++ b/src/Numerics/Interpolation/StepInterpolation.cs @@ -139,6 +139,13 @@ public double Differentiate(double t) /// Interpolated second derivative at point t. public double Differentiate2(double t) => Differentiate(t); + /// + /// Differentiate three times at point t. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + public double Differentiate3(double t) => Differentiate2(t); + /// /// Indefinite integral at point t. /// diff --git a/src/Numerics/Interpolation/TransformedInterpolation.cs b/src/Numerics/Interpolation/TransformedInterpolation.cs index 36b58a6f8..d0099c889 100644 --- a/src/Numerics/Interpolation/TransformedInterpolation.cs +++ b/src/Numerics/Interpolation/TransformedInterpolation.cs @@ -146,6 +146,13 @@ public static TransformedInterpolation Interpolate( /// Interpolated second derivative at point t. double IInterpolation.Differentiate2(double t) => throw new NotSupportedException(); + /// + /// Differentiate three times at point t. NOT SUPPORTED. + /// + /// Point t to interpolate at. + /// Interpolated third derivative at point t. + double IInterpolation.Differentiate3(double t) => throw new NotSupportedException(); + /// /// Indefinite integral at point t. NOT SUPPORTED. ///