2853 lines
90 KiB
Go
2853 lines
90 KiB
Go
// Copyright (c) 2026 Petr Balvín <opensource@petrbalvin.org> (https://petrbalvin.org)
|
||
// SPDX-License-Identifier: MIT
|
||
|
||
package tensor
|
||
|
||
import (
|
||
"crypto/sha256"
|
||
"encoding/binary"
|
||
"encoding/hex"
|
||
"fmt"
|
||
"hash"
|
||
"maps"
|
||
"math"
|
||
"os"
|
||
"path/filepath"
|
||
"runtime"
|
||
"runtime/debug"
|
||
"slices"
|
||
"testing"
|
||
|
||
"sourcedock.dev/petrbalvin/tensor/internal/core"
|
||
"sourcedock.dev/petrbalvin/tensor/spmd"
|
||
)
|
||
|
||
// The oracle harness. Every domain runs one fixed workload through the
|
||
// re-exported facade only, the raw bits of every output are hashed, and
|
||
// the per-case digest is pinned in the table below. A silent behaviour
|
||
// change, a swapped reduction order, a regression in a kernel: any of
|
||
// them moves the digest and fails the test. A deliberate change to an
|
||
// algorithm re-records the table in the same commit, so the diff shows
|
||
// exactly which case moved.
|
||
//
|
||
// Inputs are fixed literals or the seeded generator, which is
|
||
// bit-stable across Go releases by construction, and the parallel
|
||
// kernels promise a fixed reduction order, so the digest is
|
||
// deterministic on a given architecture. Run with
|
||
// TENSOR_ORACLE_RECORD=1 to print the current digests instead of
|
||
// comparing them.
|
||
|
||
// oracleDigest accumulates the raw bits of the given values. Supported
|
||
// carriers: *Array, Scalar, float64, []float64, string, bool and int;
|
||
// anything else is a harness bug and panics.
|
||
type oracleDigest struct {
|
||
h hash.Hash
|
||
}
|
||
|
||
func newOracleDigest() *oracleDigest { return &oracleDigest{h: sha256.New()} }
|
||
|
||
func (d *oracleDigest) array(t *testing.T, a *core.Array) {
|
||
d.h.Write([]byte{byte(a.Dtype())})
|
||
for _, s := range a.Shape() {
|
||
var buf [8]byte
|
||
binary.LittleEndian.PutUint64(buf[:], uint64(s))
|
||
d.h.Write(buf[:])
|
||
}
|
||
switch a.Dtype() {
|
||
case core.Int:
|
||
for _, v := range a.RawInts()[:a.Len()] {
|
||
var buf [8]byte
|
||
binary.LittleEndian.PutUint64(buf[:], uint64(v))
|
||
d.h.Write(buf[:])
|
||
}
|
||
case core.Float32:
|
||
for _, v := range a.RawFloat32s()[:a.Len()] {
|
||
var buf [4]byte
|
||
binary.LittleEndian.PutUint32(buf[:], math.Float32bits(v))
|
||
d.h.Write(buf[:])
|
||
}
|
||
case core.Float16:
|
||
for _, v := range a.RawHalves()[:a.Len()] {
|
||
var buf [2]byte
|
||
binary.LittleEndian.PutUint16(buf[:], v)
|
||
d.h.Write(buf[:])
|
||
}
|
||
case core.Float:
|
||
for _, v := range a.RawFloats()[:a.Len()] {
|
||
d.f64(v)
|
||
}
|
||
case core.Complex:
|
||
for _, v := range a.RawComplexes()[:a.Len()] {
|
||
d.f64(real(v))
|
||
d.f64(imag(v))
|
||
}
|
||
case core.Bool:
|
||
// The boolean payload hashes as its own bytes, one 0 or 1 per
|
||
// element, the layout every mask round-trip is read back in.
|
||
for _, v := range a.RawBools()[:a.Len()] {
|
||
if v {
|
||
d.h.Write([]byte{1})
|
||
} else {
|
||
d.h.Write([]byte{0})
|
||
}
|
||
}
|
||
case core.Int8:
|
||
for _, v := range a.RawInt8s()[:a.Len()] {
|
||
d.h.Write([]byte{byte(v)})
|
||
}
|
||
case core.Uint8:
|
||
for _, v := range a.RawUint8s()[:a.Len()] {
|
||
d.h.Write([]byte{byte(v)})
|
||
}
|
||
case core.Int16:
|
||
for _, v := range a.RawInt16s()[:a.Len()] {
|
||
var buf [2]byte
|
||
binary.LittleEndian.PutUint16(buf[:], uint16(v))
|
||
d.h.Write(buf[:])
|
||
}
|
||
case core.Uint16:
|
||
for _, v := range a.RawUint16s()[:a.Len()] {
|
||
var buf [2]byte
|
||
binary.LittleEndian.PutUint16(buf[:], uint16(v))
|
||
d.h.Write(buf[:])
|
||
}
|
||
case core.Int32:
|
||
for _, v := range a.RawInt32s()[:a.Len()] {
|
||
var buf [4]byte
|
||
binary.LittleEndian.PutUint32(buf[:], uint32(v))
|
||
d.h.Write(buf[:])
|
||
}
|
||
case core.Uint32:
|
||
for _, v := range a.RawUint32s()[:a.Len()] {
|
||
var buf [4]byte
|
||
binary.LittleEndian.PutUint32(buf[:], uint32(v))
|
||
d.h.Write(buf[:])
|
||
}
|
||
default:
|
||
// An unhashed payload would pin nothing: fail loudly in the
|
||
// test rather than panicking or hashing the wrong bytes.
|
||
t.Fatalf("oracle harness: array payload of dtype %s (%d) has no hash rule", a.Dtype(), a.Dtype())
|
||
}
|
||
}
|
||
|
||
func (d *oracleDigest) f64(v float64) {
|
||
var buf [8]byte
|
||
binary.LittleEndian.PutUint64(buf[:], math.Float64bits(v))
|
||
d.h.Write(buf[:])
|
||
}
|
||
|
||
func (d *oracleDigest) floats(vals []float64) {
|
||
var n [8]byte
|
||
binary.LittleEndian.PutUint64(n[:], uint64(len(vals)))
|
||
d.h.Write(n[:])
|
||
for _, v := range vals {
|
||
d.f64(v)
|
||
}
|
||
}
|
||
|
||
func (d *oracleDigest) int(v int) {
|
||
var buf [8]byte
|
||
binary.LittleEndian.PutUint64(buf[:], uint64(v))
|
||
d.h.Write(buf[:])
|
||
}
|
||
|
||
// result carries one output of an oracle case into the digest.
|
||
func (d *oracleDigest) result(t *testing.T, v any) {
|
||
switch x := v.(type) {
|
||
case *Array:
|
||
d.array(t, x)
|
||
case Scalar:
|
||
d.f64(x.Float())
|
||
case float64:
|
||
d.f64(x)
|
||
case []float64:
|
||
d.floats(x)
|
||
case string:
|
||
d.h.Write([]byte(x))
|
||
case bool:
|
||
// Error paths carry the fact of the refusal, not its wording:
|
||
// the text may be reworded, the refusal must not disappear.
|
||
if x {
|
||
d.h.Write([]byte{1})
|
||
} else {
|
||
d.h.Write([]byte{0})
|
||
}
|
||
case int:
|
||
d.int(x)
|
||
default:
|
||
panic("oracle harness: unsupported carrier type")
|
||
}
|
||
}
|
||
|
||
func (d *oracleDigest) sum() string { return hex.EncodeToString(d.h.Sum(nil)) }
|
||
|
||
// mustA builds an array from fixed values.
|
||
func mustA(t *testing.T, vals []float64, shape ...int) *Array {
|
||
t.Helper()
|
||
a, err := FromFloats(vals, shape...)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return a
|
||
}
|
||
|
||
// randA draws n deterministic samples through the seeded generator.
|
||
func randA(t *testing.T, seed int64, n int) *Array {
|
||
t.Helper()
|
||
g := NewGenerator(seed)
|
||
a, err := Floats(g, n)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return a
|
||
}
|
||
|
||
// wobble is the bounded deterministic wiggle the formula-driven
|
||
// workloads use in place of a generator: a fixed function of the
|
||
// index, so a case's inputs are its own, not the generator state's.
|
||
func wobble(i int) float64 {
|
||
return 0.3*math.Sin(7.3*float64(i)+1.1)*math.Cos(2.1*float64(i)) +
|
||
0.1*math.Sin(0.7*float64(i))
|
||
}
|
||
|
||
// oracleCases is the harness workload: one entry per domain feature,
|
||
// every output it produces fed into the digest.
|
||
var oracleCases = []struct {
|
||
name string
|
||
run func(t *testing.T) []any
|
||
}{
|
||
{"spmd-shards", func(t *testing.T) []any {
|
||
// The distributed reduction answers the single-array
|
||
// reduction's exact bits at any world size: the digest pins
|
||
// the shards' answers beside the single-array ones they must
|
||
// equal, over a length that cuts the fold's partition into
|
||
// several blocks.
|
||
const gn = 200001
|
||
vals := make([]float64, gn)
|
||
for i := range vals {
|
||
vals[i] = float64((i*6559)%2001-1000) / 7.0
|
||
if i%401 == 3 {
|
||
vals[i] = math.NaN()
|
||
}
|
||
}
|
||
yVals := make([]float64, gn)
|
||
for i := range yVals {
|
||
yVals[i] = float64((i*7919)%901-450) / 5.0
|
||
}
|
||
whole, err := FromFloats(vals, gn)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
wantSum := Sum(whole)
|
||
wantMax, err := Max(whole)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out := []any{wantSum, wantMax}
|
||
for _, size := range []int{1, 3, 5} {
|
||
// Each rank writes its own slots; the digest takes them in
|
||
// rank order after the world has joined up again.
|
||
answers := make([]core.Scalar, size*5)
|
||
err := spmd.Launch(size, func(w *spmd.World) error {
|
||
span, err := spmd.Partition(gn, size, w.Rank())
|
||
if err != nil {
|
||
return err
|
||
}
|
||
local, err := Slice(whole, 0, span.Lo, span.Hi)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
gotSum, err := w.AllReduceShards(local, span, spmd.Sum)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
gotMax, err := w.AllReduceShards(local, span, spmd.Max)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
gotProd, err := w.AllReduceShards(local, span, spmd.Prod)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
gotNorm, err := w.AllReduceNormShards(local, span, 2)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
second, err := FromFloats(yVals, gn)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
secondLocal, err := Slice(second, 0, span.Lo, span.Hi)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
gotDot, err := w.AllReduceDotShards(local, secondLocal, span)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
wantProd, err := Prod(whole, 0, false)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
wantNorm, err := Norm(whole, 2, 0, false)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
wantDot, err := Dot(whole, second)
|
||
if err != nil {
|
||
return err
|
||
}
|
||
if math.Float64bits(gotSum.Float()) != math.Float64bits(wantSum.Float()) ||
|
||
math.Float64bits(gotMax.Float()) != math.Float64bits(wantMax.Float()) ||
|
||
math.Float64bits(gotProd.Float()) != math.Float64bits(wantProd.FloatAt(0)) ||
|
||
math.Float64bits(gotNorm.Float()) != math.Float64bits(wantNorm.FloatAt(0)) ||
|
||
math.Float64bits(gotDot.Float()) != math.Float64bits(wantDot.Float()) {
|
||
t.Fatalf("size %d: the sharded answers moved against the single-array reduction", size)
|
||
}
|
||
answers[w.Rank()*5] = gotSum
|
||
answers[w.Rank()*5+1] = gotMax
|
||
answers[w.Rank()*5+2] = gotProd
|
||
answers[w.Rank()*5+3] = gotNorm
|
||
answers[w.Rank()*5+4] = gotDot
|
||
return nil
|
||
})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
for _, s := range answers {
|
||
out = append(out, s)
|
||
}
|
||
out = append(out, fmt.Sprintf("size %d agreed", size))
|
||
}
|
||
return out
|
||
}},
|
||
{"core-elementwise", func(t *testing.T) []any {
|
||
a := randA(t, 1, 64)
|
||
b := randA(t, 2, 64)
|
||
sum, _ := Add(a, b)
|
||
prod, _ := Mul(a, b)
|
||
return []any{sum, prod, MulF(a, 1.5)}
|
||
}},
|
||
{"core-matmul-einsum", func(t *testing.T) []any {
|
||
a := randA(t, 3, 64)
|
||
b := randA(t, 4, 64)
|
||
mm, err := MatMul2D(mustA(t, a.RawFloats(), 8, 8), mustA(t, b.RawFloats(), 8, 8))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
e, err := Einsum("ij,jk->ik", mustA(t, a.RawFloats(), 8, 8), mustA(t, b.RawFloats(), 8, 8))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{mm, e}
|
||
}},
|
||
{"core-sort-argsort", func(t *testing.T) []any {
|
||
vals := randA(t, 5, 200).RawFloats()
|
||
vals[3] = math.NaN()
|
||
vals[17] = math.Copysign(0, -1)
|
||
vals[42] = math.Inf(1)
|
||
vals[43] = math.Inf(-1)
|
||
vals[7] = vals[8] // a tie
|
||
src := mustA(t, vals, len(vals))
|
||
s, err := Sort(src)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
idx, err := ArgSort(src)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{s, idx}
|
||
}},
|
||
{"core-reductions", func(t *testing.T) []any {
|
||
a := randA(t, 6, 1000)
|
||
b := randA(t, 7, 1000)
|
||
dot, err := Dot(a, b)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
mean, err := Mean(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
mn, err := Min(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
mx, err := Max(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
cs, err := CumSum(a, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{Sum(a), mean, mn, mx, dot, cs}
|
||
}},
|
||
{"core-shape", func(t *testing.T) []any {
|
||
a := randA(t, 8, 24)
|
||
m := mustA(t, a.RawFloats(), 4, 6)
|
||
tr, err := TransposeAxes(m, 1, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
r, err := Reshape(m, 3, 8)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sl, err := Slice(m, 0, 1, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{tr, r, sl}
|
||
}},
|
||
{"core-fft", func(t *testing.T) []any {
|
||
a := randA(t, 9, 256)
|
||
spec, err := FFT(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
back, err := IFFT(spec)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{spec, back}
|
||
}},
|
||
{"core-quasirandom", func(t *testing.T) []any {
|
||
sob, err := SobolPoints(64, 4, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sob2, err := SobolPoints(32, 4, 1000)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
hal, err := HaltonPoints(64, 4, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{sob, sob2, hal}
|
||
}},
|
||
{"linalg-solve-inv-det", func(t *testing.T) []any {
|
||
a := randA(t, 10, 25)
|
||
// Diagonal dominance keeps the solve well conditioned.
|
||
m := mustA(t, a.RawFloats(), 5, 5)
|
||
md := m.RawFloats()
|
||
for i := range 5 {
|
||
md[i*5+i] += 10
|
||
}
|
||
b := randA(t, 11, 5)
|
||
x, err := Solve(m, b)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
iv, err := Inv(m)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
det, err := Det(m)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{x, iv, det}
|
||
}},
|
||
{"linalg-factorisations", func(t *testing.T) []any {
|
||
a := randA(t, 12, 36)
|
||
m := mustA(t, a.RawFloats(), 6, 6)
|
||
// Symmetrise, then shift to positive definite.
|
||
sym := make([]float64, 36)
|
||
mm := m.RawFloats()
|
||
for i := range 6 {
|
||
for j := range 6 {
|
||
sym[i*6+j] = (mm[i*6+j] + mm[j*6+i]) / 2
|
||
}
|
||
sym[i*6+i] += 8
|
||
}
|
||
pd := mustA(t, sym, 6, 6)
|
||
q, r, err := QR(pd)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
chol, err := Cholesky(pd)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
u, sigma, vt, err := SVD(pd)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
vals, vecs, err := Eigen(pd)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{q, r, chol, u, sigma, vt, vals, vecs}
|
||
}},
|
||
{"linalg-sparse-complex", func(t *testing.T) []any {
|
||
const n = 32
|
||
// Hermitian tridiagonal: real diagonal, imaginary couplings
|
||
// conjugated across the diagonal.
|
||
var hIdx []int64
|
||
var hVal []complex128
|
||
for i := range n {
|
||
hIdx = append(hIdx, int64(i), int64(i))
|
||
hVal = append(hVal, 2+0i)
|
||
if i+1 < n {
|
||
hIdx = append(hIdx, int64(i), int64(i+1))
|
||
hVal = append(hVal, 0.5i)
|
||
hIdx = append(hIdx, int64(i+1), int64(i))
|
||
hVal = append(hVal, -0.5i)
|
||
}
|
||
}
|
||
hIndices, err := FromInts(hIdx, len(hVal), 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
hValues, err := FromComplexes(hVal, len(hVal))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
h, err := NewSparseCOO(hIndices, hValues, []int{n, n})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ones := make([]complex128, n)
|
||
for i := range ones {
|
||
ones[i] = 1
|
||
}
|
||
rhs, err := FromComplexes(ones, n)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xCG, err := SpSolveComplexCG(h, rhs, 1e-12, 200)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
// Genuinely non-Hermitian: different off-diagonal values.
|
||
var gIdx []int64
|
||
var gVal []complex128
|
||
for i := range n {
|
||
gIdx = append(gIdx, int64(i), int64(i))
|
||
gVal = append(gVal, 2+1i)
|
||
if i+1 < n {
|
||
gIdx = append(gIdx, int64(i), int64(i+1))
|
||
gVal = append(gVal, 1)
|
||
gIdx = append(gIdx, int64(i+1), int64(i))
|
||
gVal = append(gVal, 0.25i)
|
||
}
|
||
}
|
||
gIndices, err := FromInts(gIdx, len(gVal), 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
gValues, err := FromComplexes(gVal, len(gVal))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
g, err := NewSparseCOO(gIndices, gValues, []int{n, n})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xBiCG, err := SpSolveComplexBiCGSTAB(g, rhs, 1e-12, 400)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
vals, vecs, err := SpEigenComplex(h, 4, NewGenerator(11))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{xCG, xBiCG, vals, vecs}
|
||
}},
|
||
{"linalg-sparse-cholesky", func(t *testing.T) []any {
|
||
// A shuffled 5-point Laplacian: the ordering does the work,
|
||
// the solve must answer A⁻¹·(A·v) = v for both orderings.
|
||
const n = 100
|
||
g := NewGenerator(31)
|
||
perm := make([]int, n)
|
||
for i := range perm {
|
||
perm[i] = i
|
||
}
|
||
for i := n - 1; i > 0; i-- {
|
||
j := int(g.Next() % uint64(i+1))
|
||
perm[i], perm[j] = perm[j], perm[i]
|
||
}
|
||
label := func(x, y, w int) int { return perm[y*w+x] }
|
||
var idx []int64
|
||
var vals []float64
|
||
add := func(r, c int, v float64) {
|
||
idx = append(idx, int64(r), int64(c))
|
||
vals = append(vals, v)
|
||
}
|
||
const w, h = 10, 10
|
||
for y := range h {
|
||
for x := range w {
|
||
add(label(x, y, w), label(x, y, w), 4)
|
||
if x+1 < w {
|
||
add(label(x, y, w), label(x+1, y, w), -1)
|
||
add(label(x+1, y, w), label(x, y, w), -1)
|
||
}
|
||
if y+1 < h {
|
||
add(label(x, y, w), label(x, y+1, w), -1)
|
||
add(label(x, y+1, w), label(x, y, w), -1)
|
||
}
|
||
}
|
||
}
|
||
indices, err := FromInts(idx, len(vals), 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
valArray, err := FloatsFromArray(vals, len(vals))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
coo, err := NewSparseCOO(indices, valArray, []int{n, n})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
csr, err := CSRFromCOO(coo)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
v := make([]float64, n)
|
||
for i := range v {
|
||
v[i] = float64(i%13) - 6 + 0.5*float64(i%7)
|
||
}
|
||
xTrue := New(Float, n)
|
||
for i := range v {
|
||
xTrue.RawFloats()[i] = v[i]
|
||
}
|
||
b, err := csr.MatVec(xTrue)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
nat, err := NewSparseCholesky(coo, SparseOrderingNatural)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
rcm, err := NewSparseCholesky(coo, SparseOrderingReverseCuthillMcKee)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xn, err := nat.Solve(b)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xr, err := rcm.Solve(b)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
pn := nat.Permutation()
|
||
pr := rcm.Permutation()
|
||
pnf := make([]float64, len(pn))
|
||
prf := make([]float64, len(pr))
|
||
for i := range pn {
|
||
pnf[i] = float64(pn[i])
|
||
prf[i] = float64(pr[i])
|
||
}
|
||
return []any{pnf, prf, float64(nat.NNZ()), float64(rcm.NNZ()), xn, xr}
|
||
}},
|
||
{"linalg-sparse-lu", func(t *testing.T) []any {
|
||
// A nonsymmetric banded system with real fill: the factor must
|
||
// invert the map the CSR MatVec builds, pivoting included.
|
||
const n = 40
|
||
g := NewGenerator(41)
|
||
var idx []int64
|
||
var vals []float64
|
||
add := func(r, c int, v float64) {
|
||
idx = append(idx, int64(r), int64(c))
|
||
vals = append(vals, v)
|
||
}
|
||
for i := range n {
|
||
add(i, i, 5+2*g.Unit())
|
||
if i+1 < n {
|
||
add(i, i+1, -1-g.Unit())
|
||
add(i+1, i, -1-g.Unit())
|
||
}
|
||
if i+3 < n {
|
||
add(i, i+3, -0.6*g.Unit())
|
||
}
|
||
if i+4 < n {
|
||
add(i+4, i, -0.4*g.Unit())
|
||
}
|
||
}
|
||
indices, err := FromInts(idx, len(vals), 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
valArray, err := FloatsFromArray(vals, len(vals))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
coo, err := NewSparseCOO(indices, valArray, []int{n, n})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
csr, err := CSRFromCOO(coo)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
v := make([]float64, n)
|
||
for i := range v {
|
||
v[i] = float64(i%9) - 4 + 0.25*float64(i%5)
|
||
}
|
||
xTrue := New(Float, n)
|
||
for i := range v {
|
||
xTrue.RawFloats()[i] = v[i]
|
||
}
|
||
b, err := csr.MatVec(xTrue)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
f, err := NewSparseLU(coo)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
x, err := f.Solve(b)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
p := f.Permutation()
|
||
pf := make([]float64, len(p))
|
||
for i := range p {
|
||
pf[i] = float64(p[i])
|
||
}
|
||
return []any{pf, float64(f.NNZ()), x}
|
||
}},
|
||
{"integrate-fem-poisson", func(t *testing.T) []any {
|
||
// The manufactured solution u = sin(πx)·sin(πy) on a 10×10
|
||
// structured mesh of the unit square: assemble, lift the
|
||
// boundary, solve through the sparse Cholesky.
|
||
const m = 10
|
||
n := (m + 1) * (m + 1)
|
||
vertices := make([]float64, 2*n)
|
||
for j := range m + 1 {
|
||
for i := range m + 1 {
|
||
vertices[2*(j*(m+1)+i)] = float64(i) / float64(m)
|
||
vertices[2*(j*(m+1)+i)+1] = float64(j) / float64(m)
|
||
}
|
||
}
|
||
at := func(i, j int) int64 { return int64(j*(m+1) + i) }
|
||
idx := make([]int64, 0, 6*m*m)
|
||
for j := range m {
|
||
for i := range m {
|
||
idx = append(idx, at(i, j), at(i+1, j), at(i+1, j+1), at(i, j), at(i+1, j+1), at(i, j+1))
|
||
}
|
||
}
|
||
vArr, err := FromFloats(vertices, n, 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
tArr, err := FromInts(idx, m*m*2, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
mesh, err := NewTriangleMesh2D(vArr, tArr)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sol := func(x, y float64) float64 { return math.Sin(math.Pi*x) * math.Sin(math.Pi*y) }
|
||
src := func(x, y float64) float64 { return 2 * math.Pi * math.Pi * sol(x, y) }
|
||
var bound []int
|
||
var vals []float64
|
||
for j := range m + 1 {
|
||
for i := range m + 1 {
|
||
if i == 0 || i == m || j == 0 || j == m {
|
||
bound = append(bound, j*(m+1)+i)
|
||
vals = append(vals, sol(float64(i)/float64(m), float64(j)/float64(m)))
|
||
}
|
||
}
|
||
}
|
||
u, err := SolvePoissonFEM2D(mesh, src, FEMPoissonOptions{Kappa: 1, DirichletNodes: bound, DirichletValues: vals, Ordering: SparseOrderingReverseCuthillMcKee})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
first := make([]float64, 0, 12)
|
||
for i := range 12 {
|
||
first = append(first, u.FloatAt(i))
|
||
}
|
||
return []any{first, u.FloatAt(n / 2), u.FloatAt(n - 1)}
|
||
}},
|
||
{"signal-welch-savgol", func(t *testing.T) []any {
|
||
x := randA(t, 13, 512)
|
||
freqs, psd, err := WelchPSD(x, 1000, 128, 64, "hann")
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sm, err := SavitzkyGolay(x, 11, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{freqs, psd, sm}
|
||
}},
|
||
{"signal-wavelets", func(t *testing.T) []any {
|
||
x := randA(t, 14, 128)
|
||
coef, err := DWT(x, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
back, err := IDWT(coef, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
cw, err := CWT(x, Morlet, []float64{1, 2, 4, 8}, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{coef, back, cw}
|
||
}},
|
||
{"signal-lombscargle-conv", func(t *testing.T) []any {
|
||
g := NewGenerator(15)
|
||
times, err := Floats(g, 100)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
acc := 0.0
|
||
for i, v := range times.RawFloats() {
|
||
acc += v + 0.1
|
||
times.RawFloats()[i] = acc // irregular but increasing
|
||
}
|
||
vals := randA(t, 16, 100)
|
||
freqs, power, err := LombScargle(times, vals, 0.1, 5, 32)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
input := mustA(t, randA(t, 17, 64).RawFloats(), 1, 1, 64)
|
||
kernel := mustA(t, randA(t, 18, 5).RawFloats(), 1, 1, 5)
|
||
conv, err := Conv1D(input, kernel, nil, 1, 0, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{freqs, power, conv}
|
||
}},
|
||
{"signal-filters-stencils", func(t *testing.T) []any {
|
||
x := randA(t, 19, 256)
|
||
b, a, err := ButterworthLowPass(4, 1000, 100)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
filt, err := FilterApply(b, a, x)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
grad, err := Gradient1D(x, 0.5)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
m := mustA(t, randA(t, 20, 64).RawFloats(), 8, 8)
|
||
lap, err := Laplacian(m, 1.0, 1.0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{filt, grad, lap}
|
||
}},
|
||
{"stats-moments-quantiles", func(t *testing.T) []any {
|
||
a := randA(t, 21, 1000)
|
||
v, err := Var(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
s, err := Std(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
med, err := Median(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
q, err := Quantile(a, []float64{0.05, 0.25, 0.5, 0.75, 0.95})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
mean, err := Mean(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{mean, v, s, med, q}
|
||
}},
|
||
{"stats-correlation-regression", func(t *testing.T) []any {
|
||
x := randA(t, 22, 200)
|
||
ac, err := Autocorrelate(x, 16)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
pac, err := PartialAutocorrelate(x, 8)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
// The design carries its own intercept column.
|
||
xd := x.RawFloats()
|
||
design := make([]float64, 400)
|
||
for i := range 200 {
|
||
design[i*2] = 1
|
||
design[i*2+1] = xd[i]
|
||
}
|
||
lr, err := LinearRegression(mustA(t, design, 200, 2), randA(t, 23, 200))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{ac, pac, lr.Coefficients, lr.StandardErrors, lr.PValues,
|
||
lr.ResidualVariance, lr.RSquared}
|
||
}},
|
||
{"stats-poisson-regression", func(t *testing.T) []any {
|
||
// Counts from the seeded generator: one draw per row at the
|
||
// mean the true coefficients produce.
|
||
g := NewGenerator(24)
|
||
xv := randA(t, 23, 120)
|
||
design := make([]float64, 240)
|
||
y := make([]float64, 120)
|
||
for i := range 120 {
|
||
design[i*2] = 1
|
||
design[i*2+1] = xv.FloatAt(i)
|
||
counts, err := PoissonDraws(g, 1, math.Exp(0.3+0.7*xv.FloatAt(i)))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
y[i] = counts.FloatAt(0)
|
||
}
|
||
res, err := PoissonRegression(mustA(t, design, 120, 2), mustA(t, y, 120))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{res.Coefficients, res.StandardErrors, res.PValues,
|
||
res.Fitted, res.LogLikelihood, res.Iterations}
|
||
}},
|
||
{"stats-distributions", func(t *testing.T) []any {
|
||
tt, df, p, err := WelchTTest(randA(t, 24, 50), randA(t, 25, 60))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
gc, err := GammaCDF(2.5, 3, 1.5)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ec, err := ExponentialCDF(1.25, 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sc, err := StudentTCDF(1.5, 7)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
pc, err := PoissonCDF(3, 2.5)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
bc, err := BinomialCDF(4, 10, 0.35)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
nq, err := NormalQuantile(0.975)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{NormalCDF(1.96), tt, df, p, gc, ec, sc, pc, bc, nq}
|
||
}},
|
||
{"integrate-quadrature-ode", func(t *testing.T) []any {
|
||
q, qerr, err := IntegrateFunction(func(x float64) (float64, error) {
|
||
return math.Exp(-x * x), nil
|
||
}, 0, 1, QuadratureOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
decay := func(t float64, y *Array) (*Array, error) { return MulF(y, -1), nil }
|
||
y0 := mustA(t, []float64{1}, 1)
|
||
end, err := IntegrateODE(decay, 0, 2, y0, ODEOptions{MaxSteps: 10000})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
nd, err := IntegrateND(func(x []float64) float64 {
|
||
return math.Exp(-x[0] - x[1]*x[1])
|
||
}, []float64{0, 0}, []float64{1, 1}, CubatureOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
heat, err := IntegrateHeat1D(mustA(t, []float64{1, 2, 3, 4, 5, 6, 7, 8}, 8),
|
||
0.1, 0.1, 0.05, 0.005, 3, 0, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{q, qerr, end, nd, heat}
|
||
}},
|
||
{"optim-roots-minima", func(t *testing.T) []any {
|
||
root, err := FindRoot(func(x float64) float64 { return math.Cos(x) - x }, 0, 1, 1e-12)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
rosen := func(p *Array) (float64, error) {
|
||
x := p.FloatAt(0)
|
||
y := p.FloatAt(1)
|
||
return (1-x)*(1-x) + 100*(y-x*x)*(y-x*x), nil
|
||
}
|
||
pt, val, err := Minimise(rosen, mustA(t, []float64{-1.2, 1}, 2),
|
||
MinimiseOptions{MaxIterations: 500, Tolerance: 1e-12})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{root, pt, val}
|
||
}},
|
||
{"optim-bounded-minima", func(t *testing.T) []any {
|
||
// The bounded solver on the wall case: with x forced past 1.5
|
||
// the Rosenbrock neck at (1, 1) is infeasible and the minimum
|
||
// sits at (1.5, 2.25) with value 0.25, the first coordinate
|
||
// pinned by the projection and released never.
|
||
rosen := func(p *Array) (float64, error) {
|
||
x := p.FloatAt(0)
|
||
y := p.FloatAt(1)
|
||
return (1-x)*(1-x) + 100*(y-x*x)*(y-x*x), nil
|
||
}
|
||
pt, val, err := MinimiseLBFGS(rosen, nil, mustA(t, []float64{-1.2, 1}, 2),
|
||
LBFGSOptions{Tolerance: 1e-10, Lower: []float64{1.5, math.Inf(-1)}})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
refused := false
|
||
_, _, err = MinimiseLBFGS(rosen, nil, mustA(t, []float64{0, 0}, 2),
|
||
LBFGSOptions{Lower: []float64{2}, Upper: []float64{1}})
|
||
if err != nil {
|
||
refused = true
|
||
}
|
||
return []any{pt, val, refused}
|
||
}},
|
||
{"optim-constrained-minima", func(t *testing.T) []any {
|
||
// The augmented Lagrangian on an active inequality: the bowl
|
||
// around (2, -1) cut by x + y ≤ 0 bottoms out on the wall at
|
||
// (1.5, -1.5) with value 0.5.
|
||
bowl := func(p *Array) (float64, error) {
|
||
dx := p.FloatAt(0) - 2
|
||
dy := p.FloatAt(1) + 1
|
||
return dx*dx + dy*dy, nil
|
||
}
|
||
A := mustA(t, []float64{1, 1}, 1, 2)
|
||
pt, val, err := MinimiseConstrained(bowl, nil, mustA(t, []float64{4, 4}, 2),
|
||
LinearConstraints{A: A, Lower: []float64{math.Inf(-1)}, Upper: []float64{0}},
|
||
LBFGSOptions{Tolerance: 1e-10})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
refused := false
|
||
_, _, err = MinimiseConstrained(bowl, nil, mustA(t, []float64{0, 0}, 2),
|
||
LinearConstraints{Lower: []float64{0}, Upper: []float64{1}}, LBFGSOptions{})
|
||
if err != nil {
|
||
refused = true
|
||
}
|
||
return []any{pt, val, refused}
|
||
}},
|
||
{"grad-backward", func(t *testing.T) []any {
|
||
x := FromArray(mustA(t, []float64{1.5, -2, 3}, 3), true)
|
||
loss, err := x.Mul(x)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
l, err := loss.Sum()
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
if err := l.Backward(); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{x.Grad()}
|
||
}},
|
||
{"io-roundtrips", func(t *testing.T) []any {
|
||
dir := t.TempDir()
|
||
a := randA(t, 26, 24)
|
||
m := mustA(t, a.RawFloats(), 4, 6)
|
||
fitsPath := filepath.Join(dir, "o.fits")
|
||
if err := SaveFITS(fitsPath, m, map[string]string{"OBJECT": "oracle"}); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
backFits, hdr, err := LoadFITS(fitsPath)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ncPath := filepath.Join(dir, "o.nc")
|
||
dims := []NetCDFDim{{Name: "row", Length: 4}, {Name: "col", Length: 6}}
|
||
vars := []NetCDFVar{{Name: "field", Dims: []string{"row", "col"}, Values: m}}
|
||
if err := SaveNetCDF(ncPath, dims, vars, map[string]string{"title": "oracle"}); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
_, backNc, ncAttrs, err := LoadNetCDF(ncPath)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{backFits, hdr["OBJECT"], backNc[0].Values, ncAttrs["title"]}
|
||
}},
|
||
{"core-views-and-int-precision", func(t *testing.T) []any {
|
||
// Views, the integer range above 2^53 and a widened query
|
||
// array: the 2026-09 review found defects in all three.
|
||
x := randA(t, 41, 12)
|
||
view, err := Slice(x, 0, 3, 9)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
pow, err := PowI(view, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
big, err := FromInts([]int64{1 << 53, 1<<53 + 1, -(1 << 53)}, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
hi, err := ArgMax(big)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
lo, err := ArgMin(big)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
// A float32 view whose payload is longer than its extent: the
|
||
// widening inside Interpolate2D used to walk the payload.
|
||
xs, err := FromFloat32s(make([]float32, x.Len()), x.Len())
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
for i := range xs.Len() {
|
||
// Stay inside the (3, 4) grid's closed domain.
|
||
xs.SetFloatAt(i, float64(i%4))
|
||
}
|
||
viewXs, err := Slice(xs, 0, 0, view.Len())
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
grid, err := FromFloats(make([]float64, 12), 3, 4)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
// The y queries need their own range: the grid's y domain is
|
||
// [0, rows−1], narrower than x's.
|
||
ys, err := FromFloat32s(make([]float32, viewXs.Len()), viewXs.Len())
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
for i := range ys.Len() {
|
||
ys.SetFloatAt(i, float64(i%3))
|
||
}
|
||
interp, err := Interpolate2D(grid, viewXs, ys, 0, 0, 1, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{pow, hi, lo, interp}
|
||
}},
|
||
{"io-hdf5-fixture", func(t *testing.T) []any {
|
||
// The HDF5 reference-library fixture: contiguous ints, a chunked
|
||
// deflate+shuffle float vector, a float32 matrix in a group,
|
||
// and attributes inherited from the groups.
|
||
sets, err := LoadHDF5(filepath.Join("io", "testdata", "h5", "fixture.h5"))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out := []any{}
|
||
for _, d := range sets {
|
||
out = append(out, d.Path, d.Values)
|
||
for _, k := range slices.Sorted(maps.Keys(d.Attrs)) {
|
||
out = append(out, k, d.Attrs[k])
|
||
}
|
||
}
|
||
return out
|
||
}},
|
||
{"io-hostile-inputs", func(t *testing.T) []any {
|
||
// A header that lies about its sizes must be refused, never
|
||
// allocated and never panicked over: both files here are a few
|
||
// dozen bytes and claim data no file could hold.
|
||
dir := t.TempDir()
|
||
ncPath := filepath.Join(dir, "hostile.nc")
|
||
if err := os.WriteFile(ncPath, hostileOracleHeader(), 0o644); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
_, _, _, ncErr := LoadNetCDF(ncPath)
|
||
fitsPath := filepath.Join(dir, "hostile.fits")
|
||
if err := os.WriteFile(fitsPath, hostileOracleTable(), 0o644); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
_, fitsErr := LoadFITSTable(fitsPath)
|
||
return []any{ncErr != nil, fitsErr != nil}
|
||
}},
|
||
{"core-float16", func(t *testing.T) []any {
|
||
// The half dtype end to end: narrowing of values that straddle
|
||
// the format's corners, arithmetic through the promotion
|
||
// ladder, ordering that must widen before it compares.
|
||
vals := []float64{1, -2.5, 0.25, 65504, -65504, 1.0 / 16384, 1.0 / 16777216, 0.1, 2048.5}
|
||
h, err := FromFloat16s(vals, len(vals))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
s := Sum(h)
|
||
sorted, err := Sort(h)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
prod := MulF(h, 2)
|
||
widened, err := Add(h, mustA(t, make([]float64, len(vals)), len(vals)))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ints, err := Astype(h, Int)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{h, s, sorted, prod, widened, ints,
|
||
HalfToFloat64(HalfFromFloat64(0.1)), HalfToFloat64(0x7BFF),
|
||
HalfToFloat64(0x0001), HalfFromFloat64(65520) == 0x7C00,
|
||
HalfFromFloat64(-1.0/33554432.0) == 0x8000}
|
||
}},
|
||
{"signal-kalman-arma", func(t *testing.T) []any {
|
||
// The linear filter on a scalar random walk with noisy
|
||
// measurements from the seeded generator, the unscented filter
|
||
// over the same exactly linear model, and the AR fit with its
|
||
// theoretical spectrum.
|
||
g := NewGenerator(50)
|
||
z, err := Floats(g, 40)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
one := mustA(t, []float64{1}, 1, 1)
|
||
qv := mustA(t, []float64{0.01}, 1, 1)
|
||
rv := mustA(t, []float64{1}, 1, 1)
|
||
zero := mustA(t, []float64{0}, 1)
|
||
pv := mustA(t, []float64{1}, 1, 1)
|
||
res, err := KalmanFilter(z, one, one, KalmanOptions{
|
||
InitialState: zero, InitialCovariance: pv,
|
||
ProcessNoise: qv, MeasurementNoise: rv})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
uk, err := UnscentedKalmanFilter(z,
|
||
func(x *Array) (*Array, error) { return MulF(x, 1), nil },
|
||
func(x *Array) (*Array, error) { return MulF(x, 1), nil },
|
||
KalmanOptions{
|
||
InitialState: zero, InitialCovariance: pv,
|
||
ProcessNoise: qv, MeasurementNoise: rv})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
x, err := Floats(NewGenerator(51), 300)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ar, err := EstimateAR(x, 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
freqs, psd, err := ARMASpectrum(ar, 32)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{res.States, res.Covariances, res.Innovations, res.LogLikelihood,
|
||
uk.LogLikelihood, uk.States, ar.AR, ar.InnovationVariance, ar.AIC, freqs, psd}
|
||
}},
|
||
{"signal-windows-filtfilt-dwt", func(t *testing.T) []any {
|
||
// The window catalogue, the zero-phase filter, the 2-D median
|
||
// and the Daubechies transform over seeded inputs.
|
||
hann, err := WindowHann(16, false)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
kaiser, err := WindowKaiser(16, 4.5, true)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
x := randA(t, 52, 128)
|
||
b, a, err := ButterworthLowPass(3, 100, 20)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ff, err := Filtfilt(b, a, x)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
img := mustA(t, randA(t, 53, 64).RawFloats(), 8, 8)
|
||
med, err := MedianFilter2D(img, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
coef, err := DaubechiesDWT(x, DB4, 2, DWTPeriodic)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
back, err := DaubechiesIDWT(coef, DB4, 2, DWTPeriodic)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{hann, kaiser, ff, med, coef, back}
|
||
}},
|
||
{"stats-distributions2-multipletest", func(t *testing.T) []any {
|
||
// The second-league distributions and the corrections, on
|
||
// fixed parameters and the classic Benjamini-Hochberg vector.
|
||
wc, err := WeibullCDF(1.5, 2, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
wq, err := WeibullQuantile(0.4, 2, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
lc, err := LognormalCDF(1, 0, 0.5)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
pq, err := ParetoQuantile(0.5, 1, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
nb, err := NegativeBinomialCDF(5, 3, 0.4)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ncx, err := NoncentralChiSquareCDF(9, 4, 2.5)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
nct, err := NoncentralTCDF(2, 8, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ncf, err := NoncentralFQuantile(0.5, 4, 10, 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
draws, err := DirichletDraws(NewGenerator(54), 12, []float64{2, 3, 4})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
px := randA(t, 55, 100)
|
||
py := randA(t, 56, 100)
|
||
rho, err := SpearmanRho(px, py)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
tau, err := KendallTau(px, py)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
p := []float64{0.001, 0.008, 0.039, 0.041, 0.042, 0.06, 0.074, 0.205,
|
||
0.212, 0.216, 0.222, 0.251, 0.269, 0.275, 0.34, 0.341, 0.384,
|
||
0.456, 0.657, 0.876}
|
||
bonf, err := Bonferroni(p)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
holm, err := Holm(p)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
bh, err := BenjaminiHochberg(p)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{wc, wq, lc, pq, nb, ncx, nct, ncf, draws, rho, tau, bonf, holm, bh}
|
||
}},
|
||
{"stats-regression2", func(t *testing.T) []any {
|
||
// The regularised, robust and quantile fits on one seeded
|
||
// design, with one gross outlier the robust fit must survive.
|
||
g := NewGenerator(57)
|
||
xv, err := Floats(g, 60)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
n := xv.Len()
|
||
design := make([]float64, 0, 2*n)
|
||
y := make([]float64, n)
|
||
for i := range n {
|
||
design = append(design, 1, xv.FloatAt(i))
|
||
y[i] = 1 + 2*xv.FloatAt(i) + 0.3*(float64(i%7)-3)
|
||
}
|
||
y[7] += 100
|
||
dArr := mustA(t, design, n, 2)
|
||
yArr := mustA(t, y, n)
|
||
xOnly := make([]float64, n)
|
||
for i := range n {
|
||
xOnly[i] = xv.FloatAt(i)
|
||
}
|
||
xArr := mustA(t, xOnly, n, 1)
|
||
lasso, err := Lasso(xArr, yArr, 0.05)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
en, err := ElasticNet(xArr, yArr, 0.05, 0.5)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
huber, err := HuberRegression(dArr, yArr)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
si, sl, err := TheilSenRegression(xv, yArr)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
qr, err := QuantileRegression(dArr, yArr, 0.75)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{lasso.Intercept, lasso.Coefficients, lasso.Iterations,
|
||
en.Coefficients, huber.Coefficients, huber.Scale, huber.Iterations,
|
||
si, sl, qr.Coefficients, qr.Objective}
|
||
}},
|
||
{"stats-unsupervised", func(t *testing.T) []any {
|
||
// PCA over a seeded cloud, k-means over two blobs, the mixture
|
||
// over a univariate draw and the Gaussian-process posterior.
|
||
g := NewGenerator(58)
|
||
a1, err := Normal(g, 48, 0, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
a2, err := Normal(g, 48, 0, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
a3, err := Normal(g, 48, 0, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
cloud := make([]float64, 144)
|
||
for i := range 48 {
|
||
cloud[3*i] = 3 * a1.FloatAt(i)
|
||
cloud[3*i+1] = a2.FloatAt(i)
|
||
cloud[3*i+2] = 0.1 * a3.FloatAt(i)
|
||
}
|
||
pca, err := PCA(mustA(t, cloud, 48, 3))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
bg := NewGenerator(59)
|
||
b1, err := Normal(bg, 40, 0, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
b2, err := Normal(bg, 40, 0, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
b3, err := Normal(bg, 40, 8, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
b4, err := Normal(bg, 40, 8, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
blobs := make([]float64, 160)
|
||
for i := range 40 {
|
||
blobs[2*i] = b1.FloatAt(i)
|
||
blobs[2*i+1] = b2.FloatAt(i)
|
||
blobs[40+2*i] = b3.FloatAt(i)
|
||
blobs[40+2*i+1] = b4.FloatAt(i)
|
||
}
|
||
km, err := KMeans(NewGenerator(60), mustA(t, blobs, 80, 2), 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
centres := make([]float64, 0, 4)
|
||
for _, c := range km.Centres {
|
||
centres = append(centres, c...)
|
||
}
|
||
mg := NewGenerator(61)
|
||
hi, err := Normal(mg, 40, 6, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
lo, err := Normal(mg, 40, 0, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
mix := make([]float64, 80)
|
||
for i := range 40 {
|
||
mix[2*i] = hi.FloatAt(i)
|
||
mix[2*i+1] = lo.FloatAt(i)
|
||
}
|
||
gmm, err := GaussianMixture(NewGenerator(62), mustA(t, mix, 80, 1), 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
means := make([]float64, 0, len(gmm.Means))
|
||
for _, m := range gmm.Means {
|
||
means = append(means, m...)
|
||
}
|
||
trainX := mustA(t, randA(t, 63, 16).RawFloats(), 16, 1)
|
||
sinVals, err := Sin(trainX)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
trainY, err := Reshape(sinVals, 16)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
testX := mustA(t, randA(t, 64, 8).RawFloats(), 8, 1)
|
||
kern, err := SquaredExponentialKernel(1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
gp, err := GaussianProcessRegression(kern, trainX, trainY, 0.05, testX)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
mll, err := MarginalLogLikelihood(kern, trainX, trainY, 0.05)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{pca.Loadings, pca.ExplainedVarianceRatio, float64(km.Iterations),
|
||
centres, km.Inertia, gmm.Weights, means, gmm.BIC, gp.Mean, gp.Variance, mll}
|
||
}},
|
||
{"stats-contingency", func(t *testing.T) []any {
|
||
// The exact and the asymptotic table tests over fixed tables,
|
||
// every alternative, and the effect sizes beside them.
|
||
out := []any{}
|
||
tables := [][]float64{
|
||
{12, 5, 7, 10},
|
||
{8, 2, 1, 5},
|
||
{17, 0, 0, 23},
|
||
{21, 8, 3, 9, 15, 11, 4, 12, 30},
|
||
}
|
||
for ti, tab := range tables[:3] {
|
||
a, _, err := FisherExactTest(mustA(t, tab, 2, 2), TwoSided)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
less, _, err := FisherExactTest(mustA(t, tab, 2, 2), Less)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
greater, or, err := FisherExactTest(mustA(t, tab, 2, 2), Greater)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, a, less, greater, or)
|
||
p, err := McNemarTest(mustA(t, tab, 2, 2))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, p)
|
||
if ti < 2 {
|
||
v, err := CramersV(mustA(t, tab, 2, 2))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, v)
|
||
}
|
||
}
|
||
chi2, df, p, err := ChiSquareIndependence(mustA(t, tables[3], 3, 3))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
v, err := CramersV(mustA(t, tables[3], 3, 3))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, chi2, df, p, v)
|
||
return out
|
||
}},
|
||
{"stats-mixedmodel", func(t *testing.T) []any {
|
||
// A balanced random-intercept fit and a random-slope fit, over
|
||
// formula-driven data no generator touches.
|
||
const groups8 = 8
|
||
const per = 6
|
||
y := make([]float64, 0, groups8*per)
|
||
labels := make([]int, 0, groups8*per)
|
||
for g := range groups8 {
|
||
effect := 2 * math.Sin(1.7*float64(g)+0.4)
|
||
for i := range per {
|
||
y = append(y, 5+effect+wobble(g*per+i))
|
||
labels = append(labels, g)
|
||
}
|
||
}
|
||
ones := make([]float64, len(y))
|
||
for i := range ones {
|
||
ones[i] = 1
|
||
}
|
||
res, err := LinearMixedModel(mustA(t, y, len(y)), mustA(t, ones, len(y), 1),
|
||
mustA(t, ones, len(y), 1), labels)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out := []any{res.Coefficients, res.StandardErrors,
|
||
res.RandomEffects[0], res.RandomEffects[7], res.RandomCovariance,
|
||
res.ResidualVariance, res.LogLikelihood, res.Converged,
|
||
fmt.Sprintf("%v", res.GroupLabels)}
|
||
// The slope fit: the random design carries the covariate
|
||
// alone, the identified configuration.
|
||
slopes := []float64{0.9, -1.1, 1.9, -0.3, -1.8, 0.4}
|
||
sy := make([]float64, 0, len(slopes)*per)
|
||
sx := make([]float64, 0, 2*len(slopes)*per)
|
||
sz := make([]float64, 0, len(slopes)*per)
|
||
slabels := make([]int, 0, len(slopes)*per)
|
||
jitter := make([]float64, 0, len(slopes)*per)
|
||
for g := range slopes {
|
||
for _, x := range []float64{-1, -0.6, -0.2, 0.2, 0.6, 1} {
|
||
jitter = append(jitter, wobble(len(jitter)+13)*0.1)
|
||
slabels = append(slabels, g)
|
||
sz = append(sz, x)
|
||
}
|
||
}
|
||
mean := 0.0
|
||
for _, j := range jitter {
|
||
mean += j
|
||
}
|
||
mean /= float64(len(jitter))
|
||
for g, s := range slopes {
|
||
for k := range per {
|
||
i := g*per + k
|
||
x := sz[i]
|
||
sy = append(sy, 1+2*x+0.5*s*x+0.1*(jitter[i]-mean))
|
||
sx = append(sx, 1, x)
|
||
}
|
||
}
|
||
slope, err := LinearMixedModel(mustA(t, sy, len(sy)), mustA(t, sx, len(sy), 2),
|
||
mustA(t, sz, len(sz), 1), slabels)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return append(out, slope.Coefficients, slope.RandomCovariance,
|
||
slope.ResidualVariance, slope.Converged)
|
||
}},
|
||
{"stats-hmm", func(t *testing.T) []any {
|
||
// The scaled recursions, the decode and the Baum-Welch fit on
|
||
// a fixed sequence under a fixed model.
|
||
model, err := NewHiddenMarkovModel([]float64{0.6, 0.4},
|
||
[]float64{0.7, 0.3, 0.2, 0.8}, []float64{0.9, 0.1, 0.25, 0.75})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
observations := []int{0, 0, 1, 0, 1, 1, 0, 1, 1, 1, 0, 0, 1, 0}
|
||
filtered, ll, err := model.Forward(observations)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
smoothed, sll, err := model.Smooth(observations)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
path, pathProb, err := model.Viterbi(observations)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out := []any{filtered[len(observations)-1], smoothed[0], ll, sll,
|
||
fmt.Sprintf("%v", path), pathProb}
|
||
seq := make([]int, 240)
|
||
for i := range seq {
|
||
// The wiggle stays inside ±0.4, so the shifted scale is
|
||
// positive and the modulo lands inside the symbol set.
|
||
seq[i] = int(wobble(i)*100+200) % 3
|
||
}
|
||
fit, err := FitHiddenMarkovModel(NewGenerator(77), seq, 2, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return append(out, fit.Model.Initial, fit.Model.Transition, fit.Model.Emission,
|
||
fit.LogLikelihood, fit.Converged)
|
||
}},
|
||
{"stats-hierarchy", func(t *testing.T) []any {
|
||
// The five linkages over one fixed two-column sample, with the
|
||
// cuts the dendrogram answers.
|
||
sample := make([]float64, 24)
|
||
for i := range 12 {
|
||
sample[2*i] = wobble(2*i+5)*3 + float64(i%4)
|
||
sample[2*i+1] = wobble(2*i+31) * 2
|
||
}
|
||
out := []any{}
|
||
for _, method := range []Linkage{SingleLinkage, CompleteLinkage, AverageLinkage, CentroidLinkage, WardLinkage} {
|
||
d, err := HierarchicalClustering(mustA(t, sample, 12, 2), method)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
labels, err := d.Cut(3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, d.Heights, fmt.Sprintf("%v", labels), d.Sizes[10])
|
||
}
|
||
first, err := HierarchicalClustering(mustA(t, sample, 12, 2), WardLinkage)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
byHeight, err := first.CutHeight(first.Heights[8])
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return append(out, fmt.Sprintf("%v", byHeight))
|
||
}},
|
||
{"integrate-stiff-solvers", func(t *testing.T) []any {
|
||
// The stiff batch on one decay problem, the index-1 DAE on its
|
||
// circuit, the collocation on a linear two-point problem, the
|
||
// symplectic pair and the advection pair.
|
||
decay := func(t float64, y *Array) (*Array, error) { return MulF(y, -1), nil }
|
||
y0 := mustA(t, []float64{1}, 1)
|
||
stats := &BDFVarStats{}
|
||
bdf, err := IntegrateBDFVar(decay, 0, 2, y0, BDFVarOptions{MaxSteps: 10000, Stats: stats})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ros, err := IntegrateROS4(decay, 0, 2, y0, ODEOptions{MaxSteps: 10000})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
daeRHS := func(t float64, y *Array) (*Array, error) {
|
||
out, err := Zeros(Float, 2)
|
||
if err != nil {
|
||
return nil, err
|
||
}
|
||
out.SetFloatAt(0, 1-y.FloatAt(0))
|
||
out.SetFloatAt(1, y.FloatAt(1)-y.FloatAt(0))
|
||
return out, nil
|
||
}
|
||
mass, err := FromFloats([]float64{1, 0, 0, 0}, 2, 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
dae, err := IntegrateDAE(daeRHS, mass, 0, 1, mustA(t, []float64{0, 0}, 2), 50, DAEOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
colRHS := func(t float64, y *Array) (*Array, error) {
|
||
out, err := Zeros(Float, 2)
|
||
if err != nil {
|
||
return nil, err
|
||
}
|
||
out.SetFloatAt(0, y.FloatAt(1))
|
||
return out, nil
|
||
}
|
||
col, err := SolveBoundaryCollocation(colRHS, 0, 1, mustA(t, []float64{0, 0}, 2),
|
||
BoundaryConditions{Start: []int{0}, End: []int{0}, EndValues: []float64{1}},
|
||
CollocationOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
accel := func(q *Array) (*Array, error) { return MulF(q, -1), nil }
|
||
positions, momenta, err := IntegrateYoshida4(accel, 0, 6.283185307179586,
|
||
mustA(t, []float64{1}, 1), mustA(t, []float64{0}, 1), 200)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
gradH := func(z *Array) (*Array, error) {
|
||
out, err := Zeros(Float, 2)
|
||
if err != nil {
|
||
return nil, err
|
||
}
|
||
out.SetFloatAt(0, z.FloatAt(1))
|
||
out.SetFloatAt(1, z.FloatAt(0))
|
||
return out, nil
|
||
}
|
||
mids, _, err := IntegrateMidpoint(gradH, 0, 1, mustA(t, []float64{1}, 1),
|
||
mustA(t, []float64{1}, 1), 32, MidpointOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
n := 64
|
||
dx := 1.0 / 65
|
||
pulse := make([]float64, n)
|
||
for i := range n {
|
||
c := (float64(i) + 0.5) * dx
|
||
if c >= 0.3 && c <= 0.6 {
|
||
pulse[i] = 1
|
||
}
|
||
}
|
||
upwind, err := IntegrateUpwindAdvection1D(mustA(t, pulse, n), 1, dx, 0.2, 0.9*dx, 2, 0, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
koren, err := IntegrateAdvection1D(mustA(t, pulse, n), 1, dx, 0.2, 0.9*dx, 2, 0, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{bdf, float64(stats.MaxOrder), ros, dae,
|
||
col.Values[len(col.Values)-1], float64(len(col.Mesh)),
|
||
positions[len(positions)-1], momenta[len(momenta)-1], mids[len(mids)-1],
|
||
upwind, koren}
|
||
}},
|
||
{"integrate-fem3d", func(t *testing.T) []any {
|
||
// The tetrahedral Poisson assembly on the unit box with the
|
||
// manufactured sin(πx)·sin(πy)·sin(πz) solution.
|
||
mesh, err := BoxTetraMesh3D(0, 0, 0, 1, 1, 1, 4, 4, 4)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
src := func(x, y, z float64) float64 {
|
||
return 3 * math.Pi * math.Pi * math.Sin(math.Pi*x) * math.Sin(math.Pi*y) * math.Sin(math.Pi*z)
|
||
}
|
||
var bound []int
|
||
var vals []float64
|
||
for i := range mesh.Vertices3() {
|
||
x, y, z := mesh.Vertices[3*i], mesh.Vertices[3*i+1], mesh.Vertices[3*i+2]
|
||
if x == 0 || x == 1 || y == 0 || y == 1 || z == 0 || z == 1 {
|
||
bound = append(bound, i)
|
||
vals = append(vals, 0)
|
||
}
|
||
}
|
||
u, err := SolvePoissonFEM3D(mesh, src, FEMPoisson3DOptions{
|
||
Kappa: 1, DirichletNodes: bound, DirichletValues: vals,
|
||
Ordering: SparseOrderingReverseCuthillMcKee})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{u.FloatAt(u.Len() / 2), u.FloatAt(u.Len() - 1), u.FloatAt(7)}
|
||
}},
|
||
{"linalg-sparse-lsqr-rrqr-update", func(t *testing.T) []any {
|
||
// Sparse least squares against a consistent right-hand side,
|
||
// the pivoted QR on a seeded matrix and the sparse rank-one
|
||
// round trip on a banded factor.
|
||
g := NewGenerator(65)
|
||
const m, n = 12, 8
|
||
var idx []int64
|
||
var vals []float64
|
||
for i := range m {
|
||
idx = append(idx, int64(i), int64(i%(n)))
|
||
vals = append(vals, 2+g.Unit())
|
||
if i+1 < m && i%(n)+1 < n {
|
||
idx = append(idx, int64(i), int64(i%n+1))
|
||
vals = append(vals, -0.5-0.5*g.Unit())
|
||
}
|
||
}
|
||
indices, err := FromInts(idx, len(vals), 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
valArr, err := FloatsFromArray(vals, len(vals))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
coo, err := NewSparseCOO(indices, valArr, []int{m, n})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
csr, err := CSRFromCOO(coo)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xTrue := New(Float, n)
|
||
for i := range n {
|
||
xTrue.RawFloats()[i] = float64(i%5) - 2 + 0.25*float64(i%3)
|
||
}
|
||
b, err := csr.MatVec(xTrue)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xl, infoL, err := SpLSQR(coo, b, 1e-12, 300, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xm, infoM, err := SpLSMR(coo, b, 1e-12, 300, 0)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
dm := mustA(t, randA(t, 66, 60).RawFloats(), 10, 6)
|
||
q, r, perm, rank, err := RRQR(dm)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
permF := make([]float64, len(perm))
|
||
for i, p := range perm {
|
||
permF[i] = float64(p)
|
||
}
|
||
b6 := randA(t, 67, 10)
|
||
xr, err := SolveRRQR(dm, b6)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
const bn = 30
|
||
var bIdx []int64
|
||
var bVals []float64
|
||
for i := range bn {
|
||
bIdx = append(bIdx, int64(i), int64(i))
|
||
bVals = append(bVals, 4)
|
||
if i+1 < bn {
|
||
bIdx = append(bIdx, int64(i), int64(i+1))
|
||
bVals = append(bVals, -1)
|
||
bIdx = append(bIdx, int64(i+1), int64(i))
|
||
bVals = append(bVals, -1)
|
||
}
|
||
}
|
||
bIndices, err := FromInts(bIdx, len(bVals), 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
bValArr, err := FloatsFromArray(bVals, len(bVals))
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
bCoo, err := NewSparseCOO(bIndices, bValArr, []int{bn, bn})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
factor, err := NewSparseCholesky(bCoo, SparseOrderingNatural)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
bv := New(Float, bn)
|
||
for i := range bn {
|
||
bv.RawFloats()[i] = float64(i%7) - 3
|
||
}
|
||
before, err := factor.Solve(bv)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
upd := New(Float, bn)
|
||
upd.RawFloats()[7] = 2
|
||
if err := factor.Update(upd); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
after, err := factor.Solve(bv)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
if err := factor.Downdate(upd); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
restored, err := factor.Solve(bv)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{xl, infoL.Criterion, infoL.ResidualNorm, float64(infoL.Iterations),
|
||
xm, infoM.Criterion, q, r, permF, float64(rank), xr,
|
||
before, after, restored}
|
||
}},
|
||
{"optim-lp-qp-global2", func(t *testing.T) []any {
|
||
// The simplex through the two-sided rows, the active-set QP on
|
||
// an active wall, the seeded global pair on one rotated bowl,
|
||
// the nonlinear equality and the scalar bracket.
|
||
a := mustA(t, []float64{1, 1, 1, 0}, 2, 2)
|
||
xLP, vLP, err := MinimiseLinearRows(mustA(t, []float64{1, 1}, 2),
|
||
LinearConstraints{A: a, Lower: []float64{2, math.Inf(-1)}, Upper: []float64{math.Inf(1), 3}},
|
||
LinearProgramOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
h := mustA(t, []float64{2, 0, 0, 2}, 2, 2)
|
||
xQP, vQP, multQP, err := MinimiseQP(h, mustA(t, []float64{-2, -4}, 2),
|
||
LinearConstraints{A: mustA(t, []float64{1, 0}, 1, 2),
|
||
Lower: []float64{math.Inf(-1)}, Upper: []float64{1}},
|
||
nil, QPOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
cmaF := func(p *Array) (float64, error) {
|
||
var s float64
|
||
w := []float64{4, 3, 2, 1}
|
||
for i := range 4 {
|
||
di := p.FloatAt(i) - float64(i+1)
|
||
s += w[i] * di * di
|
||
if i+1 < 4 {
|
||
s += 0.5 * di * (p.FloatAt(i+1) - float64(i+2))
|
||
}
|
||
}
|
||
return s, nil
|
||
}
|
||
xCMA, vCMA, err := MinimiseCMAES(cmaF, mustA(t, []float64{2, -2, 1, -1}, 4),
|
||
CMAESOptions{Sigma0: 0.5, Generations: 300, Tolerance: 1e-10, Seed: 7})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xSA, vSA, err := MinimiseSimulatedAnnealing(cmaF, mustA(t, []float64{2, -2, 1, -1}, 4),
|
||
SimulatedAnnealingOptions{Steps: 20000, Tolerance: 1e-6, Seed: 42, AllowBudgetExit: true})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
xNL, vNL, multNL, err := MinimiseNonlinearConstrained(
|
||
func(p *Array) (float64, error) { return p.FloatAt(0), nil },
|
||
nil, mustA(t, []float64{0.5, 0.5}, 2),
|
||
NonlinearConstraints{Equalities: []func(*Array) (float64, error){
|
||
func(p *Array) (float64, error) {
|
||
x, y := p.FloatAt(0), p.FloatAt(1)
|
||
return x*x + y*y - 1, nil
|
||
},
|
||
}}, LBFGSOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
root, err := FindRootBrent(func(x float64) float64 { return math.Cos(x) - x }, -10, 10, BrentOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
broydenRHS := func(x *Array) (*Array, error) {
|
||
out, err := Zeros(Float, 2)
|
||
if err != nil {
|
||
return nil, err
|
||
}
|
||
xf, yf := x.FloatAt(0), x.FloatAt(1)
|
||
out.SetFloatAt(0, xf*xf+yf*yf-1)
|
||
out.SetFloatAt(1, xf-yf)
|
||
return out, nil
|
||
}
|
||
xBr, resBr, err := FindRootSystem(broydenRHS, mustA(t, []float64{2, 0.5}, 2),
|
||
RootSystemOptions{UseBroyden: true})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{xLP, vLP, xQP, vQP, multQP, xCMA, vCMA, xSA, vSA,
|
||
xNL, vNL, multNL, root, xBr, resBr}
|
||
}},
|
||
{"io-hdf5-write", func(t *testing.T) []any {
|
||
// The writer round trip: chunked, shuffled, deflated floats, a
|
||
// float32 matrix, wide integers, a nested group and attributes,
|
||
// read back through the unchanged reader in both superblocks.
|
||
dir := t.TempDir()
|
||
f64 := randA(t, 68, 32)
|
||
f32, err := FromFloat32s([]float32{1.5, -2.5, 3.25, 4, 5, 6}, 2, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
wide, err := FromInts([]int64{5, -3, 1 << 40}, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
nested, err := FromFloats([]float64{1, 2, 3, 4}, 2, 2)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sets := []HDF5Dataset{
|
||
{Path: "/signals/main", Shape: f64.Shape(), Values: f64, Attrs: map[string]string{"units": "m"}},
|
||
{Path: "/image", Shape: f32.Shape(), Values: f32},
|
||
{Path: "/wide", Shape: wide.Shape(), Values: wide},
|
||
{Path: "/g/sub/deep", Shape: nested.Shape(), Values: nested},
|
||
}
|
||
attrs := map[string]map[string]string{"/": {"title": "oracle"}, "/g": {"note": "middle"}}
|
||
classic := filepath.Join(dir, "classic.h5")
|
||
if err := SaveHDF5(classic, sets, attrs, HDF5WriteOptions{Gzip: 6, Shuffle: true}); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
back, err := LoadHDF5(classic)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out := []any{}
|
||
for _, d := range back {
|
||
out = append(out, d.Path, d.Values)
|
||
for _, k := range slices.Sorted(maps.Keys(d.Attrs)) {
|
||
out = append(out, k, d.Attrs[k])
|
||
}
|
||
}
|
||
latest := filepath.Join(dir, "latest.h5")
|
||
if err := SaveHDF5(latest, sets, attrs, HDF5WriteOptions{Latest: true}); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
backLatest, err := LoadHDF5(latest)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, len(backLatest) == len(back))
|
||
for _, d := range backLatest {
|
||
if d.Path == "/signals/main" {
|
||
out = append(out, d.Values)
|
||
}
|
||
}
|
||
return out
|
||
}},
|
||
{"io-hdf5-narrow-write-roundtrip", func(t *testing.T) []any {
|
||
// The seven narrow dtypes through the writer and back, the
|
||
// values at each dtype's extremes: bool through the boolean
|
||
// enumeration convention, every integer at its stored width and
|
||
// signedness. The same sets are written contiguous, in the
|
||
// latest layout and chunked with shuffle and deflate, and every
|
||
// read-back must land the written dtype and values.
|
||
dir := t.TempDir()
|
||
sets := []HDF5Dataset{}
|
||
must := func(values *Array, err error) *Array {
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return values
|
||
}
|
||
put := func(path string, values *Array) {
|
||
sets = append(sets, HDF5Dataset{Path: path, Shape: values.Shape(), Values: values})
|
||
}
|
||
put("/bool", must(FromBools([]bool{false, true, true, false}, 2, 2)))
|
||
put("/i8", must(FromInt8s([]int8{-128, 127, 0, -1}, 2, 2)))
|
||
put("/u8", must(FromUint8s([]uint8{0, 255, 7, 128}, 2, 2)))
|
||
put("/i16", must(FromInt16s([]int16{-32768, 32767, 0, -2}, 2, 2)))
|
||
put("/u16", must(FromUint16s([]uint16{0, 65535, 5, 32768}, 2, 2)))
|
||
put("/i32", must(FromInt32s([]int32{-2147483648, 2147483647, 0, -3}, 2, 2)))
|
||
put("/u32", must(FromUint32s([]uint32{0, 4294967295, 9, 2147483648}, 2, 2)))
|
||
out := []any{}
|
||
for _, c := range []struct {
|
||
file string
|
||
opts HDF5WriteOptions
|
||
}{
|
||
{"contiguous.h5", HDF5WriteOptions{}},
|
||
{"latest.h5", HDF5WriteOptions{Latest: true}},
|
||
{"chunked.h5", HDF5WriteOptions{Gzip: 1, Shuffle: true, ChunkBytes: 8}},
|
||
} {
|
||
path := filepath.Join(dir, c.file)
|
||
if err := SaveHDF5(path, sets, nil, c.opts); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
back, err := LoadHDF5(path)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, c.file)
|
||
for _, d := range back {
|
||
out = append(out, d.Path, d.Values)
|
||
}
|
||
}
|
||
return out
|
||
}},
|
||
{"io-hdf5-uint8-roundtrip", func(t *testing.T) []any {
|
||
// A uint8 dataset through the HDF5 read path: the fixture pins
|
||
// the reader's decode directly, hand-built in the classic
|
||
// layout, and must land uint8 on the values below. The writer
|
||
// path is covered by io-hdf5-narrow-write-roundtrip above,
|
||
// which writes every narrow dtype and reads it back.
|
||
path := filepath.Join(t.TempDir(), "u8.h5")
|
||
file := oracleHDF5Narrow(oracleHDF5Fixed(1, false), []uint64{5},
|
||
[]byte{0, 1, 200, 255, 42})
|
||
if err := os.WriteFile(path, file, 0o644); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sets, err := LoadHDF5(path)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{len(sets), sets[0].Path, sets[0].Values}
|
||
}},
|
||
{"io-hdf5-int16-roundtrip", func(t *testing.T) []any {
|
||
// An int16 dataset through the same hand-built fixture path:
|
||
// the read-back must land int16 on the extremes below.
|
||
path := filepath.Join(t.TempDir(), "i16.h5")
|
||
file := oracleHDF5Narrow(oracleHDF5Fixed(2, true), []uint64{4},
|
||
[]byte{0xfd, 0xff, 0x00, 0x80, 0xff, 0x7f, 0x07, 0x00})
|
||
if err := os.WriteFile(path, file, 0o644); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sets, err := LoadHDF5(path)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{len(sets), sets[0].Path, sets[0].Values}
|
||
}},
|
||
{"io-hdf5-bool-roundtrip", func(t *testing.T) []any {
|
||
// A boolean enumeration dataset, the convention HDF5 writers
|
||
// carry booleans in: the read-back must land core.Bool.
|
||
path := filepath.Join(t.TempDir(), "bool.h5")
|
||
file := oracleHDF5Narrow(oracleHDF5EnumBool(), []uint64{3}, []byte{1, 0, 1})
|
||
if err := os.WriteFile(path, file, 0o644); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
sets, err := LoadHDF5(path)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{len(sets), sets[0].Path, sets[0].Values}
|
||
}},
|
||
{"io-netcdf-classic-natives", func(t *testing.T) []any {
|
||
// A NetCDF classic file carrying byte, short, int, char, float
|
||
// and double values, read natively. The NetCDF writer stores
|
||
// float64, float32 and int64 only, so the fixture is written
|
||
// here as classic-format bytes and the oracle pins the landings.
|
||
path := filepath.Join(t.TempDir(), "natives.nc")
|
||
if err := os.WriteFile(path, oracleNetCDFNatives(), 0o644); err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
_, vars, _, err := LoadNetCDF(path)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out := []any{len(vars)}
|
||
for _, v := range vars {
|
||
out = append(out, v.Name, v.Values)
|
||
}
|
||
return out
|
||
}},
|
||
{"optim-lm-fit", func(t *testing.T) []any {
|
||
// The Levenberg-Marquardt fit through the Fit surface: an
|
||
// exact exponential decay with a requested covariance, then a
|
||
// diagonal-Sigma run whose first observation carries a quarter
|
||
// of the variance, so the whitening moves the answer. The
|
||
// digest pins the parameters, both chi2 values, the statuses
|
||
// and the covariance.
|
||
xs := []float64{0, 1, 2, 3, 4, 5, 6, 7}
|
||
ys := make([]float64, len(xs))
|
||
for i, x := range xs {
|
||
ys[i] = 2.5*math.Exp(-0.7*float64(x)) + 0.3
|
||
}
|
||
residual := func(p *Array) (*Array, error) {
|
||
out := New(Float, len(xs))
|
||
for i := range len(xs) {
|
||
out.RawFloats()[i] = ys[i] - (p.FloatAt(0)*math.Exp(-p.FloatAt(1)*float64(xs[i])) + p.FloatAt(2))
|
||
}
|
||
return out, nil
|
||
}
|
||
p0, err := FromFloats([]float64{2, 0.5, 0.1}, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
res, err := LevenbergMarquardtFit(residual, p0, LMOptions{RequestCovariance: true})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
if res.Status != FitConverged {
|
||
t.Fatalf("status %d, want FitConverged", res.Status)
|
||
}
|
||
out := []any{res.Parameters, res.Chi2, int(res.Status), res.Covariance}
|
||
sigma, err := FromFloats([]float64{0.25, 1, 1, 1, 1, 1, 1, 1}, 8)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
wres, err := LevenbergMarquardtFit(residual, p0, LMOptions{Sigma: sigma})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
if wres.Status != FitConverged {
|
||
t.Fatalf("weighted status %d, want FitConverged", wres.Status)
|
||
}
|
||
return append(out, wres.Parameters, wres.Chi2, int(wres.Status))
|
||
}},
|
||
{"core-cosm1", func(t *testing.T) []any {
|
||
// The small-argument cosine deficit across its regimes: the
|
||
// digest pins the factored series from the underflowed end to
|
||
// the crossover and the direct subtraction above it.
|
||
a, err := FromFloats([]float64{
|
||
1e-300, 1e-12, 1e-9, 1e-8, 1e-4, 0.01, 0.5,
|
||
math.Pi / 4, 0.79, 1, 2, -0.3, -2.5, 1e6,
|
||
}, 14)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
v, err := Cosm1(a)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{v}
|
||
}},
|
||
{"core-substream", func(t *testing.T) []any {
|
||
// The substream family and the scalar splitmix64 chain: the
|
||
// digest pins the first draws of three members of one seed's
|
||
// family beside the chain of states the scalar surface hands
|
||
// out, so a change to the seeding moves loudly.
|
||
out := []any{}
|
||
for _, i := range []int{0, 1, 7} {
|
||
g, err := Substream(42, i)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
a, err := Floats(g, 8)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, a)
|
||
}
|
||
state, v := Splitmix64(0)
|
||
for range 3 {
|
||
out = append(out, fmt.Sprintf("%#016x", state), fmt.Sprintf("%#016x", v))
|
||
state, v = Splitmix64(state)
|
||
}
|
||
return out
|
||
}},
|
||
{"signal-chirp", func(t *testing.T) []any {
|
||
// The chirp synthesiser: the digest pins a rising sweep, a
|
||
// falling one and the constant tone that must be the plain
|
||
// sine.
|
||
rise, err := Chirp(64, 5, 60, 128)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
fall, err := Chirp(16, 900, 100, 8000)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
tone, err := Chirp(32, 11, 11, 256)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{rise, fall, tone}
|
||
}},
|
||
{"core-besseljreal", func(t *testing.T) []any {
|
||
// The real-order Bessel J across its regimes: the series at a
|
||
// small argument, the climb seeded from the expansion, the
|
||
// fractional Miller walk past the argument, and the integer
|
||
// delegation, pinned beside the points themselves.
|
||
type point struct {
|
||
nu, x float64
|
||
}
|
||
pts := []point{
|
||
{0.5, 2}, {2.5, 0.05}, {1.5, 12.1}, {2.3, 15},
|
||
{7.7, 40}, {20.5, 15}, {1.7, 1000}, {3, 20},
|
||
}
|
||
out := []any{len(pts)}
|
||
for _, p := range pts {
|
||
v, err := BesselJRealOrder(p.nu, p.x)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
out = append(out, v)
|
||
}
|
||
return out
|
||
}},
|
||
{"integrate-filon", func(t *testing.T) []any {
|
||
// The Filon quadrature across its regimes: polynomial
|
||
// amplitudes, where the construction is exact whatever the
|
||
// carrier does, an exponential amplitude under a carrier of
|
||
// five hundred, and the zero-frequency degeneration beside
|
||
// them.
|
||
poly := func(x float64) (float64, error) { return x, nil }
|
||
c1, s1, err := IntegrateFilon(poly, 2, 7, 500, FilonOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
exp := func(x float64) (float64, error) { return math.Exp(0.5 * x), nil }
|
||
c2, s2, err := IntegrateFilon(exp, 2, 7, 500, FilonOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
c3, s3, err := IntegrateFilon(exp, 0, 3, 0, FilonOptions{})
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{c1, s1, c2, s2, c3, s3}
|
||
}},
|
||
{"core-narrow", func(t *testing.T) []any {
|
||
// The narrow integer dtypes end to end: every width over the
|
||
// corners of its own range, promotion across widths and across
|
||
// the class boundary, the exact int64 folds, the bool masks the
|
||
// comparisons answer and the selection and logic they feed.
|
||
i8, err := Int8sFromArray([]int8{-128, -1, 0, 1, 127, -128, 127, 0, 3, -7}, 10)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
u8, err := Uint8sFromArray([]uint8{0, 1, 127, 128, 255, 254, 3, 0, 9, 200}, 10)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
i16, err := Int16sFromArray([]int16{-32768, -1, 0, 1, 32767, -32768, 32767, 0, 500, -900}, 10)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
u16, err := Uint16sFromArray([]uint16{0, 1, 32767, 32768, 65535, 65534, 7, 0, 40, 60000}, 10)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
i32, err := Int32sFromArray([]int32{-2147483648, -1, 0, 1, 2147483647, -2147483648, 2147483647, 0, 70000, -300}, 10)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
u32, err := Uint32sFromArray([]uint32{0, 1, 2147483647, 2147483648, 4294967295, 4294967294, 11, 0, 5, 4000000000}, 10)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
cross1, err := Add(i8, u16)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
cross2, err := Sub(i32, u8)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
cross3, err := Mul(i16, u8)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
cross4, err := Add(u32, i16)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ltMask, err := Lt(i8, u16)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
eqMask, err := EqI(u8, 1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
leMask, err := LeI(i32, -1)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
geMask, err := GeF(i8, 0.5)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
ltfMask, err := LtF(i16, -0.5)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
joined, err := And(ltMask, eqMask)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
flipped, err := Not(eqMask)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
picked, err := Select(i8, joined)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
chosen, err := Where(eqMask, u8, i8)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
widened, err := Astype(u32, Float)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
wide, err := FromInts([]int64{2147483648, -2147483649, -1}, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
_, narrowErr := Astype(wide, Int32)
|
||
if narrowErr == nil {
|
||
t.Fatal("Astype narrowed 2^31 into int32")
|
||
}
|
||
inRange, err := FromInts([]int64{2147483647, -2147483648, -1}, 3)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
narrowed, err := Astype(inRange, Int32)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
count, err := CountNonzero(leMask)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
minU8, err := Min(u8)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
maxI32, err := Max(i32)
|
||
if err != nil {
|
||
t.Fatal(err)
|
||
}
|
||
return []any{i8, u8, i16, u16, i32, u32,
|
||
cross1, cross2, cross3, cross4,
|
||
Sum(i8), Sum(i16), Sum(u32), minU8, maxI32,
|
||
ltMask, eqMask, leMask, geMask, ltfMask, joined, flipped,
|
||
picked, chosen, widened, narrowed, narrowErr != nil, count}
|
||
}},
|
||
}
|
||
|
||
// hostileNetCDF builds a CDF-1 file whose single variable spans three
|
||
// dimensions of 2^21, 2^20 and 2^20 elements: their product times the
|
||
// 8-byte element width wraps to zero in 64-bit arithmetic, so a reader
|
||
// that multiplies before it bounds would allocate what the header
|
||
// claims. The file is 112 bytes.
|
||
func hostileOracleHeader() []byte {
|
||
u32 := func(v uint32) []byte {
|
||
var b [4]byte
|
||
binary.BigEndian.PutUint32(b[:], v)
|
||
return b[:]
|
||
}
|
||
name := func(s string) []byte {
|
||
b := u32(uint32(len(s)))
|
||
b = append(b, s...)
|
||
for len(s)%4 != 0 {
|
||
b = append(b, 0)
|
||
s += " "
|
||
}
|
||
return b
|
||
}
|
||
var b []byte
|
||
b = append(b, 'C', 'D', 'F', 1)
|
||
b = append(b, u32(0)...) // numrecs
|
||
b = append(b, u32(10)...) // NC_DIMENSION
|
||
b = append(b, u32(3)...)
|
||
for i, l := range []uint32{1 << 21, 1 << 20, 1 << 20} {
|
||
b = append(b, name(string(rune('a'+i)))...)
|
||
b = append(b, u32(l)...)
|
||
}
|
||
b = append(b, u32(0)...) // absent attribute list
|
||
b = append(b, u32(0)...)
|
||
b = append(b, u32(11)...) // NC_VARIABLE
|
||
b = append(b, u32(1)...)
|
||
b = append(b, name("v")...)
|
||
b = append(b, u32(3)...)
|
||
b = append(b, u32(0)...)
|
||
b = append(b, u32(1)...)
|
||
b = append(b, u32(2)...)
|
||
b = append(b, u32(0)...) // absent variable attributes
|
||
b = append(b, u32(0)...)
|
||
b = append(b, u32(6)...) // NC_DOUBLE
|
||
b = append(b, u32(0)...) // vsize
|
||
b = append(b, u32(0)...) // begin
|
||
return b
|
||
}
|
||
|
||
// hostileFITSTable builds a BINTABLE header with a zero-width TFORM and
|
||
// no data block: the row arithmetic must not be able to divide by the
|
||
// row size, and the truncated file must be refused.
|
||
func hostileOracleTable() []byte {
|
||
var b []byte
|
||
for _, body := range []string{
|
||
"XTENSION= 'BINTABLE'", "BITPIX = 8", "NAXIS = 2",
|
||
"NAXIS1 = 4", "NAXIS2 = 2", "TFIELDS = 1",
|
||
"TFORM1 = 'X '", "END",
|
||
} {
|
||
b = append(b, []byte(body)...)
|
||
for len(b)%80 != 0 {
|
||
b = append(b, ' ')
|
||
}
|
||
}
|
||
return b
|
||
}
|
||
|
||
// oracleHDF5Fixed renders a version 1 fixed-point datatype message of
|
||
// the given element width and signedness, with the bit offset and bit
|
||
// precision the HDF5 file format specification's fixed-point property
|
||
// table defines behind the header.
|
||
func oracleHDF5Fixed(size uint32, signed bool) []byte {
|
||
m := []byte{0x10, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0} // version 1, class 0
|
||
if signed {
|
||
m[1] = 0x08 // class bit field: bit 3 marks two's complement
|
||
}
|
||
binary.LittleEndian.PutUint32(m[4:], size)
|
||
binary.LittleEndian.PutUint16(m[10:], uint16(8*size)) // bit precision
|
||
return m
|
||
}
|
||
|
||
// oracleHDF5EnumBool renders the boolean enumeration datatype message:
|
||
// class 8 over a one-byte unsigned base, the member names padded in
|
||
// their own fields to multiples of eight bytes, the member values
|
||
// packed behind them, exactly as the HDF5 file format specification's
|
||
// enumeration class defines.
|
||
func oracleHDF5EnumBool() []byte {
|
||
m := []byte{0x18, 2, 0, 0, 1, 0, 0, 0} // class 8, two members, size 1
|
||
m = append(m, oracleHDF5Fixed(1, false)...) // the base type
|
||
for _, name := range []string{"TRUE", "FALSE"} {
|
||
start := len(m)
|
||
m = append(m, name...)
|
||
m = append(m, 0)
|
||
for (len(m)-start)%8 != 0 {
|
||
m = append(m, 0)
|
||
}
|
||
}
|
||
return append(m, 1, 0) // TRUE = 1, FALSE = 0
|
||
}
|
||
|
||
// oracleHDF5Narrow builds a minimal classic-layout HDF5 file carrying
|
||
// one dataset: a version 0 superblock of eight-byte addresses and
|
||
// lengths, a version 1 object header with the given datatype message,
|
||
// contiguous storage and the payload bytes. The io tests build their
|
||
// hostile fixtures the same way; these bytes are fixed, so the digest
|
||
// of what the reader lands is stable across runs.
|
||
func oracleHDF5Narrow(dtypeMsg []byte, dims []uint64, payload []byte) []byte {
|
||
const headerAt, dataAt = 96, 448
|
||
align8 := func(n int) int { return (n + 7) / 8 * 8 }
|
||
n := dataAt + len(payload)
|
||
f := make([]byte, n)
|
||
copy(f, []byte{0x89, 'H', 'D', 'F', '\r', '\n', 0x1a, '\n'})
|
||
f[8] = 0 // superblock version 0, the classic layout
|
||
f[13], f[14] = 8, 8 // address and length sizes
|
||
put := func(at int, v uint64) { binary.LittleEndian.PutUint64(f[at:], v) }
|
||
put(32, math.MaxUint64) // free space undefined
|
||
put(40, uint64(n)) // end of file
|
||
put(48, math.MaxUint64) // driver information undefined
|
||
put(64, headerAt) // the root object header
|
||
space := make([]byte, 8+8*len(dims))
|
||
space[0] = 1 // dataspace version 1
|
||
space[1] = byte(len(dims))
|
||
for i, d := range dims {
|
||
binary.LittleEndian.PutUint64(space[8+8*i:], d)
|
||
}
|
||
lay := make([]byte, 18)
|
||
lay[0], lay[1] = 3, 1 // layout version 3, contiguous
|
||
binary.LittleEndian.PutUint64(lay[2:], dataAt)
|
||
binary.LittleEndian.PutUint64(lay[10:], uint64(len(payload)))
|
||
msgs := []struct {
|
||
typ uint16
|
||
body []byte
|
||
}{
|
||
{1, space}, // dataspace
|
||
{3, dtypeMsg},
|
||
{8, lay}, // data layout
|
||
}
|
||
off := headerAt
|
||
f[off] = 1 // object header version 1
|
||
binary.LittleEndian.PutUint16(f[off+2:], uint16(len(msgs)))
|
||
binary.LittleEndian.PutUint32(f[off+4:], 1) // reference count
|
||
region := off + 16
|
||
size := 0
|
||
for _, m := range msgs {
|
||
size = align8(size + 8 + len(m.body))
|
||
}
|
||
binary.LittleEndian.PutUint32(f[off+8:], uint32(size))
|
||
for _, m := range msgs {
|
||
binary.LittleEndian.PutUint16(f[region:], m.typ)
|
||
binary.LittleEndian.PutUint16(f[region+2:], uint16(len(m.body)))
|
||
copy(f[region+8:], m.body)
|
||
region += align8(8 + len(m.body))
|
||
}
|
||
copy(f[dataAt:], payload)
|
||
return f
|
||
}
|
||
|
||
// oracleNetCDFNatives builds a CDF-1 file by hand carrying one
|
||
// variable per classic type code: NC_BYTE, NC_SHORT, NC_INT, NC_CHAR,
|
||
// NC_FLOAT and NC_DOUBLE over one shared four-element dimension, the
|
||
// integers at the extremes of every width and the floats at
|
||
// recognisable values. The NetCDF writer stores float64, float32 and
|
||
// int64 only, so the fixture bytes live here exactly as the foreign
|
||
// fixtures the io package carries.
|
||
func oracleNetCDFNatives() []byte {
|
||
u32 := func(v uint32) []byte {
|
||
var b [4]byte
|
||
binary.BigEndian.PutUint32(b[:], v)
|
||
return b[:]
|
||
}
|
||
name := func(s string) []byte {
|
||
b := u32(uint32(len(s)))
|
||
b = append(b, s...)
|
||
for len(b)%4 != 0 {
|
||
b = append(b, 0)
|
||
}
|
||
return b
|
||
}
|
||
type classicVar struct {
|
||
name string
|
||
code uint32
|
||
data []byte
|
||
}
|
||
vars := []classicVar{
|
||
{"b", 1, []byte{0x80, 0xff, 0x00, 0x7f}}, // NC_BYTE: -128, -1, 0, 127
|
||
{"s", 3, []byte{0x80, 0x00, 0xff, 0xff, 0x00, 0x01, 0x7f, 0xff}}, // NC_SHORT extremes
|
||
{"i", 4, []byte{0x80, 0, 0, 0, 0xff, 0xff, 0xff, 0xff, 0, 0, 0, 1, 0x7f, 0xff, 0xff, 0xff}}, // NC_INT
|
||
{"c", 2, []byte{0x00, 0x41, 0xc8, 0xff}}, // NC_CHAR: raw bytes
|
||
{"f", 5, []byte{ // NC_FLOAT: 1.5, -2.5, 0, 3.25
|
||
0x3f, 0xc0, 0x00, 0x00, 0xc0, 0x20, 0x00, 0x00,
|
||
0x00, 0x00, 0x00, 0x00, 0x40, 0x50, 0x00, 0x00}},
|
||
{"d", 6, []byte{ // NC_DOUBLE: 1.5, -2.5, 0, 4
|
||
0x3f, 0xf8, 0, 0, 0, 0, 0, 0, 0xc0, 0x04, 0, 0, 0, 0, 0, 0,
|
||
0, 0, 0, 0, 0, 0, 0, 0, 0x40, 0x10, 0, 0, 0, 0, 0, 0}},
|
||
}
|
||
var b []byte
|
||
b = append(b, 'C', 'D', 'F', 1)
|
||
b = append(b, u32(0)...) // numrecs
|
||
b = append(b, u32(10)...) // NC_DIMENSION
|
||
b = append(b, u32(1)...) // one dimension
|
||
b = append(b, name("n")...)
|
||
b = append(b, u32(4)...) // length 4
|
||
b = append(b, u32(0)...) // gatt_list ABSENT
|
||
b = append(b, u32(0)...) // gatt_list ABSENT
|
||
b = append(b, u32(11)...) // NC_VARIABLE
|
||
b = append(b, u32(uint32(len(vars)))...) // six variables
|
||
// The header ends where the first payload begins: one 36-byte
|
||
// variable record per variable behind the fixed prefix.
|
||
begin := len(b) + len(vars)*36
|
||
off := begin
|
||
for _, v := range vars {
|
||
b = append(b, name(v.name)...)
|
||
b = append(b, u32(1)...) // rank 1
|
||
b = append(b, u32(0)...) // dimid 0
|
||
b = append(b, u32(0)...) // variable attributes ABSENT
|
||
b = append(b, u32(0)...) // variable attributes ABSENT
|
||
b = append(b, u32(v.code)...)
|
||
b = append(b, u32(uint32(len(v.data)))...) // vsize
|
||
b = append(b, u32(uint32(off))...)
|
||
off += len(v.data)
|
||
}
|
||
for _, v := range vars {
|
||
b = append(b, v.data...)
|
||
}
|
||
return b
|
||
}
|
||
|
||
// pinnedOracleDigests is the recorded behaviour of the harness, keyed
|
||
// by runtime.GOOS/GOARCH plus the code generation the binary was built
|
||
// with: floating-point kernels may differ between architectures (FMA
|
||
// availability, libm corners), and GOAMD64 decides whether the compiler
|
||
// may contract a multiply and an add into one fused operation, which
|
||
// changes the last bits without changing the algorithm. A deliberate
|
||
// algorithm change re-records the affected entries on the configuration
|
||
// it was made on; an accidental one fails the test there. A build whose
|
||
// combination has no block skips loudly until it is recorded, so
|
||
// determinism is measured per configuration, never claimed in general.
|
||
//
|
||
// The portable build is the product and the only build, and its
|
||
// digests are recorded at the toolchain default level (v1): the
|
||
// compiler has no auto-vectoriser, so a pinned higher level would buy
|
||
// only scalar FMA contraction, which the bit-pinned kernels already
|
||
// suppress by spelling. A build pinned to another GOAMD64 level has no
|
||
// block and skips with the recording instructions rather than failing.
|
||
//
|
||
// A re-recorded digest is sometimes the point itself: an old digest
|
||
// can have pinned a bug, and three entries exist because their
|
||
// predecessors did. signal-welch-savgol's magnitude goes through
|
||
// math.Hypot, so it no longer rounds through a squared sum.
|
||
// optim-roots-minima and optim-constrained-minima carry the L-BFGS
|
||
// line search that no longer returns its start point when a stiff
|
||
// objective rejects every trial step, the Nelder-Mead that no longer
|
||
// stops on a simplex lying on one level set, and MinimiseConstrained
|
||
// returning f at the answer instead of the augmented Lagrangian value.
|
||
var pinnedOracleDigests = map[string]map[string]string{
|
||
"linux/amd64+v1": {
|
||
"spmd-shards": "0bec5217a9a79b182241ac0937f32e9d20c89a5b392fc797ab3787a03fc8c9ed",
|
||
"io-hdf5-narrow-write-roundtrip": "f2e87eaadcbab846e4acf5473b3023faca3015cd3b3f0919c4c3792114589ee9",
|
||
"io-hdf5-uint8-roundtrip": "54550518525945c92d970a38d2d81a5005fb8a899b597819c6060610ae77dd54",
|
||
"io-hdf5-int16-roundtrip": "b4163c6e9ed7a04ddb4a2328d394afe074237533738d19aa5ecaa68bff9e2685",
|
||
"io-hdf5-bool-roundtrip": "fa56fd9a421d2e78f17711e717f851c4b19dc1f6946d5c4e1dddad85cc42d347",
|
||
"io-netcdf-classic-natives": "c0219b60aa963176ce509f5633c22b96bc070529db1035d53a643b67d614d01b",
|
||
"core-elementwise": "0812ac01c62b084326d0e590941a46327136e9022a101fad68998781b68f7014",
|
||
"core-matmul-einsum": "a62b07d724f06603fe4c5ebf3d448b61d7638b93d1b58df605d94096d30ba4b4",
|
||
"core-sort-argsort": "24fc532b4140bc6633b5fabb8a19c708c1919932c75ae8fce048c413c94aff0f",
|
||
"core-reductions": "764eac75a3bcc1a38f7975ed4de1e5df39d70889dfb845507ef1cfe1ff456eac",
|
||
"core-shape": "c84940b2aff4e4a25e7a4292fa6d76c2a0306f328ba334e3a6085843d94fcc3b",
|
||
"core-fft": "1848e23763a0a72f282eff45c1cc54a3d2486f4ff53b12b4a8f462bcdc773061",
|
||
"core-quasirandom": "000620721b8ce6fff07287da6301b89dbcdaa7293d68a0783aae60fa4370aa6d",
|
||
"linalg-solve-inv-det": "5b73c4a57a6221f2b5abd3b6769d6aa13723db0ada0373e870d551f795c61bf3",
|
||
"linalg-factorisations": "f70c03b3759703a6fcd1ff1ef3cbab024acd000781e2065331c18f95c45b4d15",
|
||
"linalg-sparse-cholesky": "da49d6554341ffbf40500a99f241e0133528c2e6fadfc385cfa6d7d59685da28",
|
||
"linalg-sparse-lu": "3ccf6e7064b410c3c6df0a6b2697b4c66f3722aec99491ed9298db6212d65fbb",
|
||
"integrate-fem-poisson": "065cb65bb7fba9838ed6f84b703d7ffb69d877564bfbe006905e81fdea6b9175",
|
||
"linalg-sparse-complex": "d81a8dd242a23ce7026ded53c44a4fea3a58d9b20c8be4ed2c03147c740e2130",
|
||
"signal-welch-savgol": "36d6581ac2b06191b2d6eb743e042283d8cd0cbc187c54b6d40f490eef5984a4",
|
||
"signal-wavelets": "f324710c1f67c8e5bb099731c8647c944744a62e9219ef8fdbb070412df434ab",
|
||
"signal-lombscargle-conv": "4c3f3eb0f292b978d701dcea384ccbfc26be06acff348756ec0bda19c4cca20e",
|
||
"signal-filters-stencils": "00bccd2f3b5d61f07c3383280dc8111897e1eaf57a6e16a0d70081a703a8b158",
|
||
"stats-moments-quantiles": "61d196193641395bbb855685f2727cfe29521848bba111d89de9f480ef7264cc",
|
||
"stats-correlation-regression": "f398fea69ae522ffaafdf0d307edd02d657aaf7a452c6c75d174307becbf7218",
|
||
"stats-poisson-regression": "e385d9fda18acef4669ac1d1f1e8b12a5abeb682eff37a06f6ec140dd7b9f91a",
|
||
"stats-distributions": "8986f622708f3f1f28c1c00fed9225014c17a90a8e0714e2ce6e2f181ee03321",
|
||
"integrate-quadrature-ode": "66dd47bbc30de692902937078808da5322c76e6de1facaf5099fd72344f248a5",
|
||
"optim-roots-minima": "9af044a7e71416f0ab1f8dcbd4073f199163791ad28779b59c7faabdf982e6f4",
|
||
"optim-bounded-minima": "368dbb82d31c5e53c7dad9f944921e741ccccedafa7fda8d8b180f9daa58628e",
|
||
"optim-constrained-minima": "884020c56542ac84e5460f0cbd171f64904f86254ebbae6479e3bd6662cc6662",
|
||
"grad-backward": "0f3a1be94160019b796b457aa42531817685b7bc478ed041cc6c73ec5810b744",
|
||
"io-roundtrips": "de1a36e2b2e9235a40af37b4ee5823e526f6a46583c34af6245dbbd322e3f1e3",
|
||
"core-views-and-int-precision": "ca3c789ead2894ed2c53416f128b8aa665aed9dc3d4d67e38ae9d540a1ac5ccd",
|
||
"io-hdf5-fixture": "01e96137069f7d9dd7d4cd2e236de8f5fbcf4ba6121a1f121f21c041197cf117",
|
||
"io-hostile-inputs": "9dcf97a184f32623d11a73124ceb99a5709b083721e878a16d78f596718ba7b2",
|
||
"core-float16": "6af37847cc6187b64f2bb20b7fefdd3b28a039e42988708d942c41bf62576949",
|
||
"signal-kalman-arma": "874b9c7bfd3dae3111cb49576a51913972ad572963a2ab7ad1eaa3187f15196b",
|
||
"signal-windows-filtfilt-dwt": "b435d694350e7a7f06303ded85c2b63e0e562795cf69686f2ae854f2c6dc81c5",
|
||
"stats-distributions2-multipletest": "6d2f2661fa3a2bdaefbdf16f3d26b41ef69b29c93347e16e6df8542f30b14078",
|
||
"stats-regression2": "cabb52d66ccb5d92b8d0885f5a5a166c1814d364dc75a4e8cac89819820edba3",
|
||
"stats-unsupervised": "c8a17031b7133d58a6b6b6ec21a5624281d8f291717749869e41251b3e8ef68f",
|
||
"integrate-stiff-solvers": "48276e65576eeefa8d7c3da8da42ee782441462ce34cf6bb814e494c0e00668f",
|
||
"integrate-fem3d": "8323b5c89ac9187441cfac694193ef6394aff0508e71b53052d1b76fd59a535b",
|
||
"linalg-sparse-lsqr-rrqr-update": "806455d7e53765a272afca8cfef90aa233ab78e62e4ad682c085ff7e117231fc",
|
||
"optim-lp-qp-global2": "510ce65940755a164784ed4577502db2b60604aa9cbe184f22b2ab2a4c03bd57",
|
||
"io-hdf5-write": "21482614860f885f277b05666924ab3b63a2c55b8b24734cbadc4577f2c7ac90",
|
||
"optim-lm-fit": "cdbf15c3a4c5615e1a364a9fe4e5ec431d0cb65d9744cab52cbca66af51c58cf",
|
||
"core-cosm1": "f995296d5c1f285ec6f2e719ee62f8eaf7c9a8e751cea8702cbe16b3df66ddb2",
|
||
"core-substream": "ce49fae36f5e15b37aacc25d62b97c9f800d51aa08dee553cee1a5b6e1c01798",
|
||
"signal-chirp": "ff9af5783d32eed88bf2658bde5660c247a0bd30a699c403f1755da290830a4b",
|
||
"core-besseljreal": "49698a21228f091079cd65ae94cda134b16c3b0bec58a47642f6ad3fec7461d7",
|
||
"integrate-filon": "4d5242520a5c144e4bdf111df83c2f2faf4fadce0859f87738acab91e3ed93e0",
|
||
"core-narrow": "80913507636c4aabcdc57e5f90a75cd18499c86e244d84006650de3e1273da92",
|
||
"stats-contingency": "9ac61a7a29227ddae02861076c85a669fd467f466a8c1118c52965ad4ac20a1b",
|
||
"stats-mixedmodel": "6db27ea7908863a1f2d0725325d8a72a93a5ec83ad3882bca28a5d3d2b614227",
|
||
"stats-hmm": "779e2c9cb8ee1e1871b17a04a42b2d574ff3de6132430bae8c58be502d450b75",
|
||
"stats-hierarchy": "c4171ae9c61b76c676d016356ae319c5e3c15c0e32f7d0c26f0b1b38700cb96e",
|
||
},
|
||
"linux/amd64+v1+race": {
|
||
"spmd-shards": "0bec5217a9a79b182241ac0937f32e9d20c89a5b392fc797ab3787a03fc8c9ed",
|
||
"io-hdf5-narrow-write-roundtrip": "f2e87eaadcbab846e4acf5473b3023faca3015cd3b3f0919c4c3792114589ee9",
|
||
"io-hdf5-uint8-roundtrip": "54550518525945c92d970a38d2d81a5005fb8a899b597819c6060610ae77dd54",
|
||
"io-hdf5-int16-roundtrip": "b4163c6e9ed7a04ddb4a2328d394afe074237533738d19aa5ecaa68bff9e2685",
|
||
"io-hdf5-bool-roundtrip": "fa56fd9a421d2e78f17711e717f851c4b19dc1f6946d5c4e1dddad85cc42d347",
|
||
"io-netcdf-classic-natives": "c0219b60aa963176ce509f5633c22b96bc070529db1035d53a643b67d614d01b",
|
||
"core-elementwise": "0812ac01c62b084326d0e590941a46327136e9022a101fad68998781b68f7014",
|
||
"core-matmul-einsum": "a62b07d724f06603fe4c5ebf3d448b61d7638b93d1b58df605d94096d30ba4b4",
|
||
"core-sort-argsort": "24fc532b4140bc6633b5fabb8a19c708c1919932c75ae8fce048c413c94aff0f",
|
||
"core-reductions": "764eac75a3bcc1a38f7975ed4de1e5df39d70889dfb845507ef1cfe1ff456eac",
|
||
"core-shape": "c84940b2aff4e4a25e7a4292fa6d76c2a0306f328ba334e3a6085843d94fcc3b",
|
||
"core-fft": "1848e23763a0a72f282eff45c1cc54a3d2486f4ff53b12b4a8f462bcdc773061",
|
||
"core-quasirandom": "000620721b8ce6fff07287da6301b89dbcdaa7293d68a0783aae60fa4370aa6d",
|
||
"linalg-solve-inv-det": "5b73c4a57a6221f2b5abd3b6769d6aa13723db0ada0373e870d551f795c61bf3",
|
||
"linalg-factorisations": "f70c03b3759703a6fcd1ff1ef3cbab024acd000781e2065331c18f95c45b4d15",
|
||
"linalg-sparse-cholesky": "da49d6554341ffbf40500a99f241e0133528c2e6fadfc385cfa6d7d59685da28",
|
||
"linalg-sparse-lu": "3ccf6e7064b410c3c6df0a6b2697b4c66f3722aec99491ed9298db6212d65fbb",
|
||
"integrate-fem-poisson": "065cb65bb7fba9838ed6f84b703d7ffb69d877564bfbe006905e81fdea6b9175",
|
||
"linalg-sparse-complex": "d81a8dd242a23ce7026ded53c44a4fea3a58d9b20c8be4ed2c03147c740e2130",
|
||
"signal-welch-savgol": "36d6581ac2b06191b2d6eb743e042283d8cd0cbc187c54b6d40f490eef5984a4",
|
||
"signal-wavelets": "f324710c1f67c8e5bb099731c8647c944744a62e9219ef8fdbb070412df434ab",
|
||
"signal-lombscargle-conv": "4c3f3eb0f292b978d701dcea384ccbfc26be06acff348756ec0bda19c4cca20e",
|
||
"signal-filters-stencils": "00bccd2f3b5d61f07c3383280dc8111897e1eaf57a6e16a0d70081a703a8b158",
|
||
"stats-moments-quantiles": "61d196193641395bbb855685f2727cfe29521848bba111d89de9f480ef7264cc",
|
||
"stats-correlation-regression": "f398fea69ae522ffaafdf0d307edd02d657aaf7a452c6c75d174307becbf7218",
|
||
"stats-poisson-regression": "e385d9fda18acef4669ac1d1f1e8b12a5abeb682eff37a06f6ec140dd7b9f91a",
|
||
"stats-distributions": "8986f622708f3f1f28c1c00fed9225014c17a90a8e0714e2ce6e2f181ee03321",
|
||
"integrate-quadrature-ode": "66dd47bbc30de692902937078808da5322c76e6de1facaf5099fd72344f248a5",
|
||
"optim-roots-minima": "9af044a7e71416f0ab1f8dcbd4073f199163791ad28779b59c7faabdf982e6f4",
|
||
"optim-bounded-minima": "368dbb82d31c5e53c7dad9f944921e741ccccedafa7fda8d8b180f9daa58628e",
|
||
"optim-constrained-minima": "884020c56542ac84e5460f0cbd171f64904f86254ebbae6479e3bd6662cc6662",
|
||
"grad-backward": "0f3a1be94160019b796b457aa42531817685b7bc478ed041cc6c73ec5810b744",
|
||
"io-roundtrips": "de1a36e2b2e9235a40af37b4ee5823e526f6a46583c34af6245dbbd322e3f1e3",
|
||
"core-views-and-int-precision": "ca3c789ead2894ed2c53416f128b8aa665aed9dc3d4d67e38ae9d540a1ac5ccd",
|
||
"io-hdf5-fixture": "01e96137069f7d9dd7d4cd2e236de8f5fbcf4ba6121a1f121f21c041197cf117",
|
||
"io-hostile-inputs": "9dcf97a184f32623d11a73124ceb99a5709b083721e878a16d78f596718ba7b2",
|
||
"core-float16": "6af37847cc6187b64f2bb20b7fefdd3b28a039e42988708d942c41bf62576949",
|
||
"signal-kalman-arma": "874b9c7bfd3dae3111cb49576a51913972ad572963a2ab7ad1eaa3187f15196b",
|
||
"signal-windows-filtfilt-dwt": "b435d694350e7a7f06303ded85c2b63e0e562795cf69686f2ae854f2c6dc81c5",
|
||
"stats-distributions2-multipletest": "6d2f2661fa3a2bdaefbdf16f3d26b41ef69b29c93347e16e6df8542f30b14078",
|
||
"stats-regression2": "cabb52d66ccb5d92b8d0885f5a5a166c1814d364dc75a4e8cac89819820edba3",
|
||
"stats-unsupervised": "c8a17031b7133d58a6b6b6ec21a5624281d8f291717749869e41251b3e8ef68f",
|
||
"integrate-stiff-solvers": "48276e65576eeefa8d7c3da8da42ee782441462ce34cf6bb814e494c0e00668f",
|
||
"integrate-fem3d": "8323b5c89ac9187441cfac694193ef6394aff0508e71b53052d1b76fd59a535b",
|
||
"linalg-sparse-lsqr-rrqr-update": "806455d7e53765a272afca8cfef90aa233ab78e62e4ad682c085ff7e117231fc",
|
||
"optim-lp-qp-global2": "510ce65940755a164784ed4577502db2b60604aa9cbe184f22b2ab2a4c03bd57",
|
||
"io-hdf5-write": "21482614860f885f277b05666924ab3b63a2c55b8b24734cbadc4577f2c7ac90",
|
||
"optim-lm-fit": "cdbf15c3a4c5615e1a364a9fe4e5ec431d0cb65d9744cab52cbca66af51c58cf",
|
||
"core-cosm1": "f995296d5c1f285ec6f2e719ee62f8eaf7c9a8e751cea8702cbe16b3df66ddb2",
|
||
"core-substream": "ce49fae36f5e15b37aacc25d62b97c9f800d51aa08dee553cee1a5b6e1c01798",
|
||
"signal-chirp": "ff9af5783d32eed88bf2658bde5660c247a0bd30a699c403f1755da290830a4b",
|
||
"core-besseljreal": "49698a21228f091079cd65ae94cda134b16c3b0bec58a47642f6ad3fec7461d7",
|
||
"integrate-filon": "4d5242520a5c144e4bdf111df83c2f2faf4fadce0859f87738acab91e3ed93e0",
|
||
"core-narrow": "80913507636c4aabcdc57e5f90a75cd18499c86e244d84006650de3e1273da92",
|
||
"stats-contingency": "9ac61a7a29227ddae02861076c85a669fd467f466a8c1118c52965ad4ac20a1b",
|
||
"stats-mixedmodel": "6db27ea7908863a1f2d0725325d8a72a93a5ec83ad3882bca28a5d3d2b614227",
|
||
"stats-hmm": "779e2c9cb8ee1e1871b17a04a42b2d574ff3de6132430bae8c58be502d450b75",
|
||
"stats-hierarchy": "c4171ae9c61b76c676d016356ae319c5e3c15c0e32f7d0c26f0b1b38700cb96e",
|
||
},
|
||
}
|
||
|
||
// oracleBuild names the code generation this binary was built with: the
|
||
// GOAMD64 level decides whether the compiler may contract a multiply and
|
||
// an add into one fused operation, and the race detector's
|
||
// instrumentation changes which loops the compiler still contracts.
|
||
// Each changes the last bits of ordinary arithmetic without changing
|
||
// the algorithm, so the digests are keyed by both alongside the
|
||
// architecture: a build whose combination has no recorded block skips
|
||
// loudly instead of failing, exactly as an unseen architecture does.
|
||
func oracleBuild() string {
|
||
level := "v1"
|
||
if bi, ok := debug.ReadBuildInfo(); ok {
|
||
for _, s := range bi.Settings {
|
||
if s.Key == "GOAMD64" && s.Value != "" {
|
||
level = s.Value
|
||
}
|
||
}
|
||
}
|
||
return "+" + level + oracleRaceSuffix
|
||
}
|
||
|
||
// oracleKey names this architecture and code-generation block.
|
||
func oracleKey() string {
|
||
return runtime.GOOS + "/" + runtime.GOARCH + oracleBuild()
|
||
}
|
||
|
||
func TestOracle(t *testing.T) {
|
||
record := os.Getenv("TENSOR_ORACLE_RECORD") == "1"
|
||
current := make(map[string]string, len(oracleCases))
|
||
for _, c := range oracleCases {
|
||
d := newOracleDigest()
|
||
for _, v := range c.run(t) {
|
||
d.result(t, v)
|
||
}
|
||
current[c.name] = d.sum()
|
||
}
|
||
if record {
|
||
for _, c := range oracleCases {
|
||
t.Logf("oracle %s %s", c.name, current[c.name])
|
||
}
|
||
t.Logf("record these into pinnedOracleDigests[%q] and commit", oracleKey())
|
||
return
|
||
}
|
||
pinned, known := pinnedOracleDigests[oracleKey()]
|
||
if !known {
|
||
t.Skipf("no oracle digests recorded for %s; run TENSOR_ORACLE_RECORD=1 go test -run TestOracle -v . on it, "+
|
||
"then commit the pinnedOracleDigests[%q] block", oracleKey(), oracleKey())
|
||
}
|
||
if len(pinned) != len(oracleCases) {
|
||
t.Fatalf("pinnedOracleDigests[%q] holds %d entries, the harness runs %d; record with TENSOR_ORACLE_RECORD=1",
|
||
oracleKey(), len(pinned), len(oracleCases))
|
||
}
|
||
for _, c := range oracleCases {
|
||
want, ok := pinned[c.name]
|
||
if !ok {
|
||
t.Errorf("oracle case %q is not pinned", c.name)
|
||
continue
|
||
}
|
||
if current[c.name] != want {
|
||
t.Errorf("oracle case %q moved:\n got %s\nwant %s\n"+
|
||
"If the change is deliberate, re-record with TENSOR_ORACLE_RECORD=1 in the same commit.",
|
||
c.name, current[c.name], want)
|
||
}
|
||
}
|
||
}
|