Lapack

Apply Q (or Qᵀ when transpose) from qr to a length-m y, into a fresh result (LAPACK dormqr restricted to a single column). Allocates; applyQInto does not.

Apply Q (or Qᵀ when transpose) from qr to y into out, which is returned. out may be y.

Lower-triangular Cholesky decomposition A = L * LT, returned as a fresh matrix.

Throws IllegalArgumentException at the first non-positive pivot unless policy says otherwise; see CholeskyPolicy for why that is the default and what asking for CholeskyPolicy.Regularize means.

LU factorization with partial pivoting of a square a (LAPACK dgetrf); a is not modified. Allocates the factor buffers; factorInto refactorizes into existing ones.

Refactorize a into out's existing buffers, returning out. A periodic refactorization — the simplex rebuilding its basis, a filter re-decomposing a covariance — otherwise allocates an n² copy and a pivot array each time; this reuses both. out must have the same dimension as a, and its previous contents are discarded.

Invert a general matrix from its LU factorization: returns A⁻¹ given P·A = L·U (LAPACK dgetri).

The counterpart of invertSpd, which koblas had while this was missing — you could invert an SPD matrix but not a general one, for no reason beyond nobody having asked.

A⁻¹ is A · X = I, so this is the multi-RHS solveInto against an identity right-hand side. Written that way deliberately rather than as a per-column loop: it reuses the blocked solve, which means it also picks up whatever a host backend does for that routine instead of needing its own binding.

Unlike invertSpd there is no symmetry to exploit and no leading-zeros shortcut, because P scatters the unit right-hand sides.

Prefer solve when you want to apply A⁻¹ to something. An explicit inverse costs more and is less accurate than a solve against the factors — this is for the cases that genuinely need the entries, such as reading a covariance off a normal-equations matrix.

Throws

Invert an SPD matrix from its Cholesky factor: returns A^-1 given L = chol(A).

Solves A * x = e_j column by column, exploiting the unit-vector right-hand side: forward substitution starts at row j (the leading entries are provably zero) and back substitution only produces rows >= j — the strictly-upper entries of the symmetric A^-1 come from mirroring the lower triangle.

Symmetric indefinite factorization A = L·D·Lᵀ with Bunch–Kaufman partial pivoting (LAPACK dsytrf, lower). As with Blas.symv, only the lower triangle of a is read — the strictly upper triangle may hold anything — and a is not modified. Use this where the matrix is symmetric but not positive definite (KKT systems); for SPD matrices cholesky is cheaper.

QR factorization A = Q·R of an m×n a via Householder reflections (LAPACK dgeqrf); a is not modified and any shape is accepted. Rank deficiency is not detected — zero diagonal entries of R surface in solveLeastSquares as infinities/NaNs, following the triangular-solve convention.

Reciprocal condition number estimate 1 / (anorm · est(‖A⁻¹‖₁)) from a factorization (LAPACK dgecon). anorm is the 1-norm of the original, unfactored matrix (see norm1), which the caller computes before factoring. Returns 1.0 for the empty factorization and exactly 0.0 when lu is singular or anorm is zero.

This is an order-of-magnitude estimate, not an exact condition number: the estimator never exceeds the true ‖A⁻¹‖₁, so the returned value never understates the conditioning. It is the cheap signal a revised simplex uses to decide when to refactorize. The default implementation is a Hager-style 1-norm estimator over solveInto; backends may substitute a native estimator, so the exact value can differ between backends while agreeing in magnitude.

The estimator needs four vectors and runs several sweeps, so it is the allocation-heaviest routine here — which matters because a simplex calls it to decide when to refactorize. Passing a workspace reuses those buffers across calls and makes it allocation-free.

Solve A · x = b (or Aᵀ · x = b when transpose) for the factorization lu (LAPACK dgetrs); returns a fresh x. transpose serves the simplex's BTRAN (Bᵀ y = c) against a factored basis. Allocates the result; solveInto writes into a caller-owned destination instead.

Solve A · X = B (or Aᵀ · X = B when transpose) for the b.cols right-hand-side columns of b at once (LAPACK dgetrs with nrhs); returns a fresh X. The default runs the permutation and the two triangular block solves directly on the shared packed format; backends may substitute a native block solve.

Solve A · x = b for a symmetric indefinite factorization ldl (LAPACK dsytrs); returns a fresh x. Symmetry makes the transposed solve identical, so there is no transpose flag. Allocates the result; solveInto writes into a caller-owned destination instead.

Solve A · X = B for the b.cols right-hand-side columns of b at once against a symmetric indefinite factorization (LAPACK dsytrs with nrhs); returns a fresh X. Backends may substitute a native block solve.

Solve A · x = b (or Aᵀ · x = b when transpose) into out, which is returned. Nothing is allocated, so a loop that owns its destination — a simplex FTRAN/BTRAN against a factored basis, a filter update — runs without touching the collector. out may be the same array as b.

The transposed direction has to stage the solved vector before scattering it through the permutation, so pass a workspace to make that staging buffer reusable as well; without one it is the single allocation this routine still makes.

Solve A · X = B (or Aᵀ · X = B when transpose) into out, which is returned. out may be b. The transposed direction stages a block before scattering its rows through the permutation, so pass a workspace to lend that n·nrhs buffer.

Solve A · x = b into out, which is returned; allocates nothing. out may be b.

Solve A · X = B into out, which is returned. out may be b.

One dsytrs replay per right-hand side, over the contiguous column it occupies. Under row-major storage this was a single widened pass, because extracting a column meant a strided gather; here a column is already a contiguous run, so reusing the tuned vector solve costs nothing but re-reading the factor and keeps the pivot-block bookkeeping in one place instead of two.

Least-squares solve min ‖A·x − b‖₂ from the factorization (the dgels shape): requires m ≥ n and full column rank, returns the length-n solution x = R⁻¹·(Qᵀb)[0..n). This is the kernel square-root/array filters build on; it composes with cholesky rank-one updates for sliding-window problems.

Least-squares solve into out, which has length n and is returned. The Qᵀb product needs a length-m intermediate, so pass a workspace to make the call allocation-free.

Minimum-norm solution of the underdetermined consistent system A · x = b for a wide m×n A with m <= n and full row rank (LAPACK dgels's underdetermined shape, via QR of the transpose instead of LQ): pass the factorization qr(Aᵀ). With Aᵀ = Q·R we have A = Rᵀ·Qᵀ, so a forward solve Rᵀ·w = b followed by x = Q·(w padded with zeros) gives the solution of smallest 2-norm. b has length m; the result has length n. Rank deficiency is not detected and surfaces as infinities/NaNs, following the triangular-solve convention.

Minimum-norm solve into out, which has length m (the wide system's column count) and is returned. A length-n intermediate holds the forward solve, so pass a workspace to make the call allocation-free.

Solve A * x = b for x, given L = chol(A) (lower-triangular, A = L * LT). Allocates a fresh result vector; b is not modified.

Forward substitution L * y = b is column-oriented: once y[j] is final, its contribution is subtracted from the remaining right-hand side down contiguous column j via VectorKernels.axpy. Back substitution LT * x = y reads row i of LT, which is column i of L, so it uses VectorKernels.dot on the same contiguous runs.

Invert a triangular matrix in place of a fresh result (LAPACK dtrtri): returns T⁻¹ for the lower or upper triangle of the square a, taking the diagonal as 1 when unitDiag.

koblas could solve against a triangle (trsv, trsm) but not invert one. The result is triangular with the same orientation, and the opposite strict triangle of the output is zero — the input's is never read, so it may hold anything.

Unlike the triangular solves, this validates the diagonal. Triangular.kt states that the cores deliberately do not, leaving a singular triangle to produce infinities, on the grounds that the caller knows what it passed; an inverse has no such caller-supplied right-hand side to blame, and returning a matrix of infinities is worse than saying which entry was zero.

Throws