190 lines
6.8 KiB
Go
190 lines
6.8 KiB
Go
// Copyright (c) 2026 Petr Balvín <opensource@petrbalvin.org> (https://petrbalvin.org)
|
|
// SPDX-License-Identifier: MIT
|
|
|
|
package integrate
|
|
|
|
import (
|
|
"math"
|
|
"testing"
|
|
|
|
"sourcedock.dev/petrbalvin/tensor/internal/base"
|
|
)
|
|
|
|
// intFilon closed forms: the antiderivatives the referents come from.
|
|
// intFilon1 integrates 1·cos(kx) and 1·sin(kx); intFilonX integrates x
|
|
// against the same kernels; intFilonExp integrates e^{ax} against
|
|
// them. All are exact calculus, evaluated independently of the code
|
|
// under test.
|
|
func intFilon1(a, b, k float64) (c, s float64) {
|
|
return (math.Sin(k*b) - math.Sin(k*a)) / k,
|
|
(math.Cos(k*a) - math.Cos(k*b)) / k
|
|
}
|
|
|
|
func intFilonX(a, b, k float64) (c, s float64) {
|
|
cb, sb := math.Cos(k*b), math.Sin(k*b)
|
|
ca, sa := math.Cos(k*a), math.Sin(k*a)
|
|
c = (cb+k*b*sb)/k/k - (ca+k*a*sa)/k/k
|
|
s = (sb-k*b*cb)/k/k - (sa-k*a*ca)/k/k
|
|
return c, s
|
|
}
|
|
|
|
func intFilonExp(a, b, amp, k float64) (c, s float64) {
|
|
cb, sb := math.Cos(k*b), math.Sin(k*b)
|
|
ca, sa := math.Cos(k*a), math.Sin(k*a)
|
|
eb, ea := math.Exp(amp*b), math.Exp(amp*a)
|
|
c = eb*(amp*cb+k*sb)/(amp*amp+k*k) - ea*(amp*ca+k*sa)/(amp*amp+k*k)
|
|
s = eb*(amp*sb-k*cb)/(amp*amp+k*k) - ea*(amp*sa-k*ca)/(amp*amp+k*k)
|
|
return c, s
|
|
}
|
|
|
|
func wantFilon(t *testing.T, label string, gotC, gotS, wantC, wantS, tol float64) {
|
|
t.Helper()
|
|
if d := math.Abs(gotC - wantC); d > tol*math.Max(1, math.Abs(wantC)) {
|
|
t.Fatalf("%s: cos part = %.17g, want %.17g (absolute %.3g)", label, gotC, wantC, d)
|
|
}
|
|
if d := math.Abs(gotS - wantS); d > tol*math.Max(1, math.Abs(wantS)) {
|
|
t.Fatalf("%s: sin part = %.17g, want %.17g (absolute %.3g)", label, gotS, wantS, d)
|
|
}
|
|
}
|
|
|
|
// TestIntegrateFilonExactAmplitudes pins the exactness the method
|
|
// promises: unit and linear amplitudes are polynomials below the
|
|
// default degree, so every frequency from one to a thousand must land
|
|
// on the closed form at the rounding floor, whatever the carrier does
|
|
// between the samples.
|
|
func TestIntegrateFilonExactAmplitudes(t *testing.T) {
|
|
one := func(float64) (float64, error) { return 1, nil }
|
|
identity := func(x float64) (float64, error) { return x, nil }
|
|
for _, c := range []struct {
|
|
label string
|
|
a, b float64
|
|
k float64
|
|
f func(float64) (float64, error)
|
|
ref func(a, b, k float64) (c, s float64)
|
|
}{
|
|
{"unit on [0, π]", 0, math.Pi, 1, one, intFilon1},
|
|
{"unit on [2, 7]", 2, 7, 100, one, intFilon1},
|
|
{"unit on [0, 1]", 0, 1, 1000, one, intFilon1},
|
|
{"x on [0, π]", 0, math.Pi, 1, identity, intFilonX},
|
|
{"x on [2, 7]", 2, 7, 500, identity, intFilonX},
|
|
} {
|
|
gotC, gotS, err := IntegrateFilon(c.f, c.a, c.b, c.k, FilonOptions{})
|
|
if err != nil {
|
|
t.Fatalf("%s: %v", c.label, err)
|
|
}
|
|
wantC, wantS := c.ref(c.a, c.b, c.k)
|
|
wantFilon(t, c.label, gotC, gotS, wantC, wantS, 1e-12)
|
|
}
|
|
}
|
|
|
|
// TestIntegrateFilonZeroFrequency pins the degeneration at k = 0: the
|
|
// sine part is exactly zero and the cosine part is the plain integral
|
|
// of the amplitude.
|
|
func TestIntegrateFilonZeroFrequency(t *testing.T) {
|
|
f := func(x float64) (float64, error) { return x * x, nil }
|
|
gotC, gotS, err := IntegrateFilon(f, 0, 3, 0, FilonOptions{})
|
|
if err != nil {
|
|
t.Fatalf("IntegrateFilon: %v", err)
|
|
}
|
|
if gotS != 0 {
|
|
t.Fatalf("the sine part at k = 0 is %g, want 0", gotS)
|
|
}
|
|
if d := math.Abs(gotC - 9); d > 1e-12 {
|
|
t.Fatalf("the cosine part at k = 0 is %.17g, want 9", gotC)
|
|
}
|
|
}
|
|
|
|
// TestIntegrateFilonExponential holds a non-polynomial amplitude
|
|
// against the exact antiderivative at a frequency whose carrier the
|
|
// automatic panel count must respect: two hundred and fifty
|
|
// oscillations over the interval, answered from a few thousand
|
|
// amplitude samples.
|
|
func TestIntegrateFilonExponential(t *testing.T) {
|
|
const amp = 0.5
|
|
f := func(x float64) (float64, error) { return math.Exp(amp * x), nil }
|
|
gotC, gotS, err := IntegrateFilon(f, 2, 7, 500, FilonOptions{})
|
|
if err != nil {
|
|
t.Fatalf("IntegrateFilon: %v", err)
|
|
}
|
|
wantC, wantS := intFilonExp(2, 7, amp, 500)
|
|
wantFilon(t, "exp amplitude", gotC, gotS, wantC, wantS, 1e-11)
|
|
}
|
|
|
|
// TestIntegrateFilonOrientation pins the reversed interval and the
|
|
// empty one: reversing negates both parts and an empty interval
|
|
// integrates to nothing.
|
|
func TestIntegrateFilonOrientation(t *testing.T) {
|
|
f := func(x float64) (float64, error) { return math.Exp(0.2 * x), nil }
|
|
fc, fs, err := IntegrateFilon(f, 0, 3, 40, FilonOptions{})
|
|
if err != nil {
|
|
t.Fatal(err)
|
|
}
|
|
rc, rs, err := IntegrateFilon(f, 3, 0, 40, FilonOptions{})
|
|
if err != nil {
|
|
t.Fatal(err)
|
|
}
|
|
if rc != -fc || rs != -fs {
|
|
t.Fatalf("the reversed interval gave (%.17g, %.17g), want the negation of (%.17g, %.17g)", rc, rs, fc, fs)
|
|
}
|
|
if ec, es, err := IntegrateFilon(f, 2, 2, 40, FilonOptions{}); err != nil || ec != 0 || es != 0 {
|
|
t.Fatalf("the empty interval gave (%g, %g, %v)", ec, es, err)
|
|
}
|
|
}
|
|
|
|
// TestIntegrateFilonDeterministic redraws one integral and requires
|
|
// the same bits, the contract every entry point here carries.
|
|
func TestIntegrateFilonDeterministic(t *testing.T) {
|
|
f := func(x float64) (float64, error) { return math.Exp(0.1 * x), nil }
|
|
one := func() (float64, float64) {
|
|
c, s, err := IntegrateFilon(f, 0, 5, 300, FilonOptions{})
|
|
if err != nil {
|
|
t.Fatal(err)
|
|
}
|
|
return c, s
|
|
}
|
|
c1, s1 := one()
|
|
c2, s2 := one()
|
|
if c1 != c2 || s1 != s2 {
|
|
t.Fatalf("the same call moved: (%.17g, %.17g) against (%.17g, %.17g)", c1, s1, c2, s2)
|
|
}
|
|
}
|
|
|
|
// TestIntegrateFilonErrors pins the contract: NaN bounds or frequency,
|
|
// a node count out of range, a forced panel count whose panels carry
|
|
// more carrier than the weights can be built within, and a failing or
|
|
// non-finite amplitude all surface as errors naming themselves.
|
|
func TestIntegrateFilonErrors(t *testing.T) {
|
|
f := func(x float64) (float64, error) { return 1, nil }
|
|
if _, _, err := IntegrateFilon(f, math.NaN(), 1, 10, FilonOptions{}); err == nil {
|
|
t.Fatal("a NaN bound: want an error")
|
|
}
|
|
if _, _, err := IntegrateFilon(f, 0, 1, math.NaN(), FilonOptions{}); err == nil {
|
|
t.Fatal("a NaN frequency: want an error")
|
|
}
|
|
if _, _, err := IntegrateFilon(f, 0, 1, math.Inf(1), FilonOptions{}); err == nil {
|
|
t.Fatal("an infinite frequency: want an error")
|
|
}
|
|
if _, _, err := IntegrateFilon(f, 0, 1, 10, FilonOptions{Nodes: 1}); err == nil {
|
|
t.Fatal("one node: want an error")
|
|
}
|
|
if _, _, err := IntegrateFilon(f, 0, 1, 10, FilonOptions{Nodes: 33}); err == nil {
|
|
t.Fatal("33 nodes: want an error")
|
|
}
|
|
if _, _, err := IntegrateFilon(f, 0, 1, 1e6, FilonOptions{Panels: 2}); err == nil {
|
|
t.Fatal("two panels under a carrier of 1e6: want an error")
|
|
}
|
|
boom := func(float64) (float64, error) { return 0, base.Errf("amplitude failed") }
|
|
if _, _, err := IntegrateFilon(boom, 0, 1, 10, FilonOptions{}); err == nil {
|
|
t.Fatal("a failing amplitude: want the error to propagate")
|
|
}
|
|
bad := func(x float64) (float64, error) {
|
|
if x > 0.5 {
|
|
return math.NaN(), nil
|
|
}
|
|
return 1, nil
|
|
}
|
|
if _, _, err := IntegrateFilon(bad, 0, 1, 10, FilonOptions{}); err == nil {
|
|
t.Fatal("a non-finite amplitude value: want an error")
|
|
}
|
|
}
|