tensor/internal/base/scalar_abs_pin_test.go

// Copyright (c) 2026 Petr Balvín <opensource@petrbalvin.org> (https://petrbalvin.org)
// SPDX-License-Identifier: MIT

package base

import (
	"math"
	"testing"
)

// TestAbsOfEveryScalarType pins absOf for every element type the Scalar
// constraint admits. The float32 case used to fall through to the 0 the
// switch returns for an unhandled type, so a Factor[float32] pivot
// search compared every element against that zero and lost its partial
// pivoting silently.
func TestAbsOfEveryScalarType(t *testing.T) {
	// float64: the plain magnitude, including the negative zero and the
	// signed extremes.
	for _, c := range []struct {
		in   float64
		want float64
	}{
		{3.5, 3.5},
		{-3.5, 3.5},
		{0, 0},
		{math.Copysign(0, -1), 0},
		{-math.MaxFloat64, math.MaxFloat64},
		{math.Inf(-1), math.Inf(1)},
	} {
		if got := absOf(c.in); got != c.want {
			t.Errorf("absOf(float64(%v)) = %v, want %v", c.in, got, c.want)
		}
	}
	// float32: what Factor[float32] reads. The value must be the exact
	// magnitude widened to float64, not a zero.
	for _, c := range []struct {
		in   float32
		want float64
	}{
		{3.5, 3.5},
		{-3.5, 3.5},
		{0, 0},
		{float32(math.Copysign(0, -1)), 0},
		{-math.MaxFloat32, float64(math.MaxFloat32)},
		{1e-30, float64(float32(1e-30))},
	} {
		got := absOf(c.in)
		if got != c.want {
			t.Errorf("absOf(float32(%v)) = %v, want %v", c.in, got, c.want)
		}
		if got == 0 && c.in != 0 {
			t.Errorf("absOf(float32(%v)) = 0: the float32 case is missing again", c.in)
		}
	}
	// complex128: the modulus.
	for _, c := range []struct {
		in   complex128
		want float64
	}{
		{complex(3, 4), 5},
		{complex(-3, -4), 5},
		{complex(0, 0), 0},
		{complex(-0.0, 0), 0},
		{complex(1e200, 1e200), math.Sqrt2 * 1e200},
	} {
		if got := absOf(c.in); got != c.want {
			t.Errorf("absOf(complex128(%v)) = %v, want %v", c.in, got, c.want)
		}
	}
	// The pivot search is the caller that matters: a float32 matrix
	// whose largest element is negative must pick that pivot.
	m := [][]float32{
		{0.05, 0},
		{-0.1, 0},
	}
	perm, _ := Factor(m)
	if perm[0] != 1 {
		t.Errorf("Factor[float32] picked row %d as the first pivot, want the larger magnitude at row 1", perm[0])
	}
}

// TestAbsComplexExtremeScale pins AbsComplex beyond the range of the
// sqrt(r² + i²) it used to be. The squared form overflows above about
// 1.34e154 and underflows below about 1.5e-162 while the magnitude is an
// ordinary number in both directions.
func TestAbsComplexExtremeScale(t *testing.T) {
	const big = 1e200
	z := complex(big, big)
	naive := math.Sqrt(real(z)*real(z) + imag(z)*imag(z))
	if !math.IsInf(naive, 1) {
		t.Fatalf("sqrt(r² + i²) at 1e200 = %v, the overflow this test relies on is gone", naive)
	}
	if got, want := AbsComplex(z), math.Sqrt2*big; math.Abs(got-want) > 1e-15*want {
		t.Errorf("AbsComplex(1e200 + 1e200i) = %v, want %v", got, want)
	}
	// The same pair with one component zero: the magnitude is finite and
	// the squared form would still overflow.
	if got := AbsComplex(complex(big, 0)); got != big {
		t.Errorf("AbsComplex(1e200) = %v, want %v", got, big)
	}
	// The underflow side.
	const small = 1e-200
	if naive := math.Sqrt(small*small + small*small); naive != 0 {
		t.Fatalf("sqrt(r² + i²) at 1e-200 = %v, the underflow this test relies on is gone", naive)
	}
	if got, want := AbsComplex(complex(small, small)), math.Sqrt2*small; math.Abs(got-want) > 1e-15*want {
		t.Errorf("AbsComplex(1e-200 + 1e-200i) = %v, want %v", got, want)
	}
	// Ordinary magnitudes keep their exact values: the classics and a
	// Pythagorean pair whose components are exactly representable.
	for _, c := range []struct {
		in   complex128
		want float64
	}{
		{complex(3, 4), 5},
		{complex(5, 12), 13},
		{complex(1, 0), 1},
		{complex(0, -1), 1},
		{complex(0, 0), 0},
	} {
		if got := AbsComplex(c.in); got != c.want {
			t.Errorf("AbsComplex(%v) = %v, want %v", c.in, got, c.want)
		}
	}
	// A mixed pair: the large component does not swallow the small one.
	if got := AbsComplex(complex(1e200, 1e-200)); got != 1e200 {
		t.Errorf("AbsComplex(1e200 + 1e-200i) = %v, want 1e200", got)
	}
}