110 lines
3.5 KiB
Go
110 lines
3.5 KiB
Go
// Copyright (c) 2026 Petr Balvín <opensource@petrbalvin.org> (https://petrbalvin.org)
|
|
// SPDX-License-Identifier: MIT
|
|
|
|
package linalg
|
|
|
|
import (
|
|
"math"
|
|
"testing"
|
|
)
|
|
|
|
// TestCubicSplineExactness checks the defining property of the
|
|
// spline: it passes exactly through every knot.
|
|
func TestCubicSplineExactness(t *testing.T) {
|
|
xs := mustFloats(t, []float64{0, 1, 3, 5, 8})
|
|
ys := mustFloats(t, []float64{1, 3, -2, 4, 0})
|
|
s, err := NewCubicSpline(xs, ys)
|
|
if err != nil {
|
|
t.Fatalf("NewCubicSpline: %v", err)
|
|
}
|
|
for i := range xs.Len() {
|
|
v := s.At(xs.FloatAt(i))
|
|
if math.Abs(v-ys.FloatAt(i)) > 1e-14 {
|
|
t.Fatalf("S(%v) = %v, want %v", xs.FloatAt(i), v, ys.FloatAt(i))
|
|
}
|
|
}
|
|
}
|
|
|
|
// TestCubicSplineAgainstSin checks the interpolated values against
|
|
// sin on a dense grid between knots: with a few knots the natural
|
|
// spline tracks the smooth function to the 1e-5 level.
|
|
func TestCubicSplineAgainstSin(t *testing.T) {
|
|
knots := []float64{0, 0.5, 1.0, 1.5, 2.0, 2.5, 3.14159265}
|
|
xs := mustFloats(t, knots, len(knots))
|
|
ys := make([]float64, len(knots))
|
|
for i := range knots {
|
|
ys[i] = math.Sin(knots[i])
|
|
}
|
|
ysArr := mustFloats(t, ys, len(knots))
|
|
s, err := NewCubicSpline(xs, ysArr)
|
|
if err != nil {
|
|
t.Fatalf("NewCubicSpline: %v", err)
|
|
}
|
|
for _, q := range []float64{0.25, 0.75, 1.25, 1.75, 2.25, 3.0} {
|
|
got := s.At(q)
|
|
want := math.Sin(q)
|
|
if math.Abs(got-want) > 1e-3 {
|
|
t.Fatalf("S(%v) = %v, want ≈ %v", q, got, want)
|
|
}
|
|
}
|
|
}
|
|
|
|
// TestCubicSplineNatural checks the natural boundary: the second
|
|
// derivative at both endpoints must vanish.
|
|
func TestCubicSplineNatural(t *testing.T) {
|
|
xs := mustFloats(t, []float64{0, 1, 2, 3})
|
|
ys := mustFloats(t, []float64{0, 1, 0, 1})
|
|
s, err := NewCubicSpline(xs, ys)
|
|
if err != nil {
|
|
t.Fatalf("NewCubicSpline: %v", err)
|
|
}
|
|
if s.m[0] != 0 || s.m[s.n-1] != 0 {
|
|
t.Fatalf("natural boundary: M₀ = %v, M₃ = %v, both want 0", s.m[0], s.m[s.n-1])
|
|
}
|
|
}
|
|
|
|
// TestCubicSplineEvaluate checks the array-valued entry point and the
|
|
// NaN contract outside the knot range.
|
|
func TestCubicSplineEvaluate(t *testing.T) {
|
|
xs := mustFloats(t, []float64{0, 1, 2, 3, 4, 5})
|
|
ys := mustFloats(t, []float64{0, 1, 0, 1, 0, 1})
|
|
s, err := NewCubicSpline(xs, ys)
|
|
if err != nil {
|
|
t.Fatalf("NewCubicSpline: %v", err)
|
|
}
|
|
// Inside: exact at knots.
|
|
queries := mustFloats(t, []float64{0, 2, 5}, 3)
|
|
got, err := s.Evaluate(queries)
|
|
if err != nil {
|
|
t.Fatalf("Evaluate: %v", err)
|
|
}
|
|
for i, xi := range []float64{0, 2, 5} {
|
|
if math.Abs(got.FloatAt(i)-ys.FloatAt(int(xi))) > 1e-14 {
|
|
t.Fatalf("S(%v) = %v, want %v", xi, got.FloatAt(i), ys.FloatAt(int(xi)))
|
|
}
|
|
}
|
|
// Outside: NaN.
|
|
out, _ := s.Evaluate(mustFloats(t, []float64{-1, 10}, 2))
|
|
if !math.IsNaN(out.FloatAt(0)) || !math.IsNaN(out.FloatAt(1)) {
|
|
t.Fatalf("outside the knots expected NaN, got %v and %v", out.FloatAt(0), out.FloatAt(1))
|
|
}
|
|
}
|
|
|
|
// TestCubicSplineRejectsInvalid pins the input contract: non-increasing
|
|
// abscissae and a too-short point set are refused.
|
|
func TestCubicSplineRejectsInvalid(t *testing.T) {
|
|
small := mustFloats(t, []float64{0, 1})
|
|
ys := mustFloats(t, []float64{0, 1})
|
|
if _, err := NewCubicSpline(small, ys); err == nil {
|
|
t.Fatal("expected an error for fewer than 3 points")
|
|
}
|
|
dup := mustFloats(t, []float64{0, 1, 1, 3})
|
|
if _, err := NewCubicSpline(dup, mustFloats(t, []float64{0, 1, 2, 3})); err == nil {
|
|
t.Fatal("expected an error for duplicate abscissae")
|
|
}
|
|
dec := mustFloats(t, []float64{0, 2, 1})
|
|
if _, err := NewCubicSpline(dec, mustFloats(t, []float64{0, 1, 2})); err == nil {
|
|
t.Fatal("expected an error for decreasing abscissae")
|
|
}
|
|
}
|