Renamed GetHorizontalPoints to GetHorizontalDerivativeTValues and renamed and extended GetMaxPoint to GetMinMaxTValues. Aslo added unit tests for both functions.

Author: NetPro69

src/Numerics.Tests/InterpolationTests/CubicSplineTest.cs

@@ -39,7 +39,10 @@ public class CubicSplineTest
     {
         readonly double[] _t = { -2.0, -1.0, 0.0, 1.0, 2.0 };
         readonly double[] _y = { 1.0, 2.0, -1.0, 0.0, 1.0 };
-
+        //test data for min max values
+        readonly double[] _x = { -4, -3, -2, -1, 0, 1, 2, 3, 4 };
+        readonly double[] _z = { -7, 2, 5, 0, -3, -1, -4, 0, 6 }; 
+                                                                  
         /// <summary>
         /// Verifies that the interpolation matches the given value at all the provided sample points.
         /// </summary>
@@ -202,5 +205,37 @@ public void InterpolateAkimaSorted_MustBeThreadSafe_GitHub219([Values(8, 32, 256
 
             CollectionAssert.DoesNotContain(yipol, Double.NaN);
         }
+
+        /// <summary>
+        /// Tests that the cubic spline returns the correct t values where the derivative is 0
+        /// </summary>
+        [Test]
+        public void NaturalSplineGetHorizontalDerivativeTValues()
+        {
+            CubicSpline it = CubicSpline.InterpolateBoundaries(_x, _z, SplineBoundaryCondition.SecondDerivative, -4.0, SplineBoundaryCondition.SecondDerivative, 4.0);
+            var horizontalDerivatives = it.GetHorizontalDerivativeTValues();
+            //readonly double[] _x = { -4, -3, -2, -1,  0,  1,  2, 3, 4 };
+            //readonly double[] _z = { -7,  2,  5,  0, -3, -1, -4, 0, 6 };
+            Assert.AreEqual(4, horizontalDerivatives.Length, "Incorrect number of points with derivative value equal to 0");
+            Assert.IsTrue(horizontalDerivatives[0]>=-3 && horizontalDerivatives[0] <= -2,"Spline returns wrong t value: "+horizontalDerivatives[0]+" for first point");
+            Assert.IsTrue(horizontalDerivatives[1] >= -1 && horizontalDerivatives[1] <= 0, "Spline returns wrong t value: " + horizontalDerivatives[1] + " for second point");
+            Assert.IsTrue(horizontalDerivatives[2] >= 0 && horizontalDerivatives[2] <= 1, "Spline returns wrong t value: " + horizontalDerivatives[2] + " for third point");
+            Assert.IsTrue(horizontalDerivatives[3] >= 2 && horizontalDerivatives[3] <= 3, "Spline returns wrong t value: " + horizontalDerivatives[3] + " for fourth point");
+            Console.WriteLine("GetHorizontalDerivativeTValues checked out ok for cubic spline.");
+         }
+
+        /// <summary>
+        /// Tests that the min and max values for the natural spline are correct
+        /// </summary>
+        [Test]
+        public void NaturalSplineGetMinMaxTvalues()
+        {
+            CubicSpline it = CubicSpline.InterpolateBoundaries(_x, _z, SplineBoundaryCondition.SecondDerivative, -4.0, SplineBoundaryCondition.SecondDerivative, 4.0);
+            var minMax = it.GetMinMaxTValues();
+            Assert.AreEqual(-4, minMax.Item1, "Spline returns wrong t value for global minimum");
+            Assert.AreEqual(4, minMax.Item2, "Spline returns wrong t value for global maximum");
+            Console.WriteLine("GetMinMaxTValues checked out ok for cubic spline.");
+
+        }
     }
 }

src/Numerics/Interpolation/CubicSpline.cs

@@ -1,4 +1,4 @@
-// <copyright file="CubicSpline.cs" company="Math.NET">
+// <copyright file="CubicSpline.cs" company="Math.NET">
 // Math.NET Numerics, part of the Math.NET Project
 // http://numerics.mathdotnet.com
 // http://github.com/mathnet/mathnet-numerics
@@ -625,13 +625,13 @@ int LeftSegmentIndex(double t)
         }
 
         /// <summary>
-        /// Gets all the points (x values) where the derivative is 0
+        /// Gets all the t values where the derivative is 0
         /// </summary>
-        /// <returns>An array of X values where the derivative of the spline is 0</returns>
-        public double[] GetHorizontalPoints()
+        /// <returns>An array of t values (in the domain of teh function) where the derivative of the spline is 0</returns>
+        public double[] GetHorizontalDerivativeTValues()
         {
             List<double> points = new List<double>();
-            for (int index = 0; index < _x.Length-1; index++)
+            for (int index = 0; index < _x.Length - 1; index++)
             {
                 double a = 6 * _c3[index]; //derive ax^3 and multiply by 2
                 double b = 2 * _c2[index]; //derive bx^2
@@ -640,58 +640,70 @@ public double[] GetHorizontalPoints()
                 //first check if a is 0, if so its a linear function, this happens with quadratic condition
                 if (a.AlmostEqual(0))
                 {
-                double x = _x[index] - c / b;
-                //check if the result is in the domain
-                if (_x[index] <= x && x <= _x[index + 1]) points.Add(x);
+                    double x = _x[index] - c / b;
+                    //check if the result is in the domain
+                    if (_x[index] <= x && x <= _x[index + 1]) points.Add(x);
                 }
                 else if (d.AlmostEqual(0))//its a quadratic with a single solution
                 {
-                double x = _x[index] - b / a;
-                if (_x[index] <= x && x <= _x[index + 1]) points.Add(x);
+                    double x = _x[index] - b / a;
+                    if (_x[index] <= x && x <= _x[index + 1]) points.Add(x);
                 }
                 else if (d > 0)//only has a solution if d is greater than 0
                 {
-                d = (double)System.Math.Sqrt(d);
-                //apply quadratic equation
-                double x1 = _x[index]+ (-b + d) / a;
-                double x2 = _x[index]+ (-b - d) / a;
-                //Add any solution points that fall within the domain to the list
-                if ((_x[index] <= x1) && (x1 <= _x[index + 1])) points.Add(x1);
-                if ((_x[index] <= x2) && (x2 <= _x[index + 1])) points.Add(x2);
+                    d = (double)System.Math.Sqrt(d);
+                    //apply quadratic equation
+                    double x1 = _x[index] + (-b + d) / a;
+                    double x2 = _x[index] + (-b - d) / a;
+                    //Add any solution points that fall within the domain to the list
+                    if ((_x[index] <= x1) && (x1 <= _x[index + 1])) points.Add(x1);
+                    if ((_x[index] <= x2) && (x2 <= _x[index + 1])) points.Add(x2);
                 }
             }
             return points.ToArray();
         }
 
         /// <summary>
-        /// Returns the maximum value for the spline within its domain.
+        /// Returns the t values in the domain of the spline for which it takes the minimum and maximum value.
         /// </summary>
-        /// <returns></returns>
-        public double GetMaxPoint()
+        /// <returns>A tuple containing the t value for which the spline is minimum in the first component and maximum in the second component </returns>
+        public Tuple<double, double> GetMinMaxTValues()
         {
             double max = double.MinValue;
-            double x=0;
+            double min = double.MaxValue;
+            double minT = 0;
+            double maxT = 0;
             //go through the functions points to check if one of them has a higher value
             foreach (double p in _x)
             {
                 double y = Interpolate(p);
                 if (y > max)
                 {
-                max = y;
-                x = p;
+                    max = y;
+                    maxT = p;
+                }
+                if (y < min)
+                {
+                    min = y;
+                    minT = p;
                 }
             }
             //go through the inflexion, local minimums and local maximums
-            foreach (double p in GetHorizontalPoints())
+            foreach (double p in GetHorizontalDerivativeTValues())
             {
                 double y = Interpolate(p);
                 if (y > max)
                 {
-                max = y;
-                x = p;
+                    max = y;
+                    maxT = p;
+                }
+                if (y < min)
+                {
+                    min = y;
+                    minT = p;
                 }
             }
-            return x;
+            return new Tuple<double, double>(minT, maxT);
         }
     }
 }