A square matrix of statically-known size.

Cholesky decomposition using the Cholesky-Crout algorithm.

>>> e = F.array @[3,3] @Double [4,12,-16,12,37,-43,-16,-43,98]
>>> pretty (cholM e)
[[2.0,0.0,0.0],
 [6.0,1.0,0.0],
 [-8.0,5.0,3.0]]
>>> F.mult (cholM e) (F.transpose (cholM e)) == e
True

Inversion of a triangular matrix.

>>> t = F.array @[3,3] @Double [1,0,1,0,1,2,0,0,1]
>>> pretty (invtriM t)
[[1.0,0.0,-1.0],
 [0.0,1.0,-2.0],
 [0.0,0.0,1.0]]
>>> F.ident @[3,3] == F.mult t (invtriM t)
True

Inverse of a square matrix via Cholesky decomposition.

>>> e = F.array @[3,3] @Double [4,12,-16,12,37,-43,-16,-43,98]
>>> pretty (inverseM e)
[[49.36111111111111,-13.555555555555554,2.1111111111111107],
 [-13.555555555555554,3.7777777777777772,-0.5555555555555555],
 [2.1111111111111107,-0.5555555555555555,0.1111111111111111]]

LU decomposition with partial pivoting.

Returns (p, l, u) such that a = p^T . l . u, where p is the accumulated row-permutation matrix, l is unit lower triangular, and u is upper triangular. The algorithm is the standard iterative Schur-complement update: at each pivot k, swap the pivot row into place, then eliminate the rows below by a rank-1 update.

The Schur step a' = a - c . r is an outer product followed by subtraction, which is the matrix-level reading of the trace / feedback pattern used in circuit constructions.

LU decomposition via iterated rank-1 updates.

This is the same factorisation as luM, but the Schur-complement step is written explicitly as an outer-product contraction:

A' = A - c ⊗ r

where c is the column of multipliers below the pivot and r is the pivot row. The outer product uses expand and the subtraction is elementwise; together they are the matrix-level reading of the trace / feedback pattern.

The result agrees with luM on any invertible square matrix.

Apply a Householder reflection to zero the subdiagonal of column k.

For column x = a[:,k], compute a reflector H = I - 2 v v^T / (v^T v) such that H x = α e_k with α = -sign(x[k]) · ||x||. The action on the whole matrix is the rank-1 update

a' = a - (2 / (v^T v)) · v ⊗ (v^T a)

where v^T a is a contraction over the row axis and v ⊗ (...) is an outer product. This is the circuit-native QR step: a scheduled reflection implemented as expand contract subtract.

Newtype wrapper that gives square matrices a multiplicative (and, over a field, divisible) structure.

The unit is the identity matrix; multiplication is the usual matrix product. Division uses Cholesky-based inversion.