uncomplicate.neanderthal.linalg

Contains type-agnostic linear algebraic functions roughly corresponding to the functionality usually defined in LAPACK (factorizations, solvers, etc.). This namespace works similarly to the uncomplicate.neanderthal.core namespace; see there for more details about the intended use.

Cheat Sheet

Also see:

con

(con lu nrm nrm1?)(con a nrm1-?)(con a)

Computes the reciprocal of the condition number of a triangularized matrix lu.

If nrm1? is true, computes the reciprocal based on 1 norm, if not, on infinity norm. nrm and matching nrm1? have to be supplied.

If the LU has been stale, or norm is not possible to compute, throws ex-info.

See related info about gecon.

det

(det a)

Computes the determinant of a triangularized matrix a.

If the matrix is not square, throws ex-info.

es

(es a)

Computes eigenvalues using Schur factorization without destroying a’s data.

es!

(es! a w vs)(es! a w)(es! a)

Computes the eigenvalues and Schur factorization of a matrix a.

On exit, a is overwritten with the Schur form T. The first 2 columns of a column-oriented GE matrix w are overwritten with eigenvalues of a. If v GE matrix is provided, it will be overwritten by the orthogonal matrix Z of Schur vectors. If vs is nil, only eigenvalues are computed.

If the QR algorithm failed to compute all the eigenvalues, throws ex-info, with the information on the index of the first eigenvalue that converged. If w is not column-oriented, or does not have 2 columns, throws ex-info. If v’s dimensions do not fit with a’s dimensions, throws ex-info. If some value in the native call is illegal, throws ex-info.

See related info about gees.

ev

(ev a)

Computes eigenvalues without destroying a’s data.

ev!

(ev! a w vl vr)(ev! a w)(ev! a)

Computes the eigenvalues and left and right eigenvectors of a matrix a.

a - source matrix and computed orthogonal factorization (m x m). w - computed eigenvalues (m x 2 or k x 1 for symmetric matrices). vl and/or vr - left and right eigenvectors (m x m or m x k for symmetric matrices).

On exit, a is overwritten with QR factors. The first 2 columns of a column-oriented GE matrix w are overwritten with eigenvalues of a. If vl and vr GE matrices are provided, they will be overwritten with left and right eigenvectors. For symmetric matrices, computes the first k eigenvalues (k <= m). If vl and/or vr are nil, eigenvectors are not computed.

If the QR algorithm failed to compute all the eigenvalues, throws ex-info, with the information on the index of the first eigenvalue that converged. If w is not column-oriented, or does not have 2 columns, throws ex-info. If vl or vr dimensions do not fit with a’s dimensions, throws ex-info. If some value in the native call is illegal, throws ex-info.

See related info about geev.

gls

(gls a b d)

Solves the generalized linear least squares problem using a generalized QR factorization.

Minimizes the 2-norm ||y|| subject to constraints d=ax+by.

See gls!.

gls!

(gls! a b d x y)

Destructively solves the generalized linear least squares problem using a generalized QR factorization.

Minimizes the 2-norm ||y|| subject to constraints d=ax+by.

Overwrites a with the R (rqf), x with the solution, b with T (rqf), d with garbage, x ad y with solution to the GLS problem.

If the dimensions of arguments do not fit, throws ex-info. If a and b do not have the same layout (column or row oriented), throws ex-info. If the least squares cannot be computed the function throws ex-info. If some value in the native call is illegal, throws ex-info.

See related info about gglse.

lqf

(lqf a)

Computes the LQ factorization of a GE m x n matrix.

See qrf!.

See related info about gelqf.

lqf!

(lqf! a)

Destructively computes the LQ factorization of a GE m x n matrix.

See qrf!.

See related info about gelqf.

ls

(ls a b)

Solves an overdetermined or underdetermined linear system ax=b with full rank matrix using QR or LQ factorization.

Uses a temporary copy ofa for the factorization data: - QR if m >= n; - LQ if m < n.

Returns a new instance containing: - the least squares solution vectors if m >= n - minimum norm solution vectors if m < n.

If a and b do not have the same layout (column or row oriented), throws ex-info. If the least squares cannot be computed the function throws ex-info. If some value in the native call is illegal, throws ex-info.

See ls!. See related info about gels.

ls!

(ls! a b)

Destructively solves an overdetermined or underdetermined linear system AX = B with full rank matrix using QR or LQ factorization.

Overwrites a with the factorization data: - QR if m >= n; - LQ if m < n.

Overwrites b with: - the least squares solution vectors if m >= n - minimum norm solution vectors if m < n.

If a and b do not have the same layout (column or row oriented), throws ex-info. If the least squares cannot be computed the function throws ex-info. If some value in the native call is illegal, throws ex-info.

See related info about gels.

lse

(lse a b c d)

Solves the generalized linear least squares problem with equality constraints using RQ factorization.

See lse! See related info about gglse.

lse!

(lse! a b c d x)

Destructively solves the generalized linear least squares problem with equality constraints using RQ factorization.

Minimizes the 2-norm ||ax-c|| subject to constraints bx=d.

Overwrites a with the T (rqf), x with the solution, b with R (rqf), d with garbage, c with the residual sum of squares.

If the dimensions of arguments do not fit, throws ex-info. If a and b do not have the same layout (column or row oriented), throws ex-info. If the least squares cannot be computed the function throws ex-info. If some value in the native call is illegal, throws ex-info.

See related info about gglse.

org

(org or)

Generates the real orthogonal matrix Q of the QR, RQ, QL or LQ factorization formed by qrf! or qrfp! (or rqf! etc.).

The input is a structure containing orthogonalized matrix :lu and vector :tau that were previously processed by qrf! (and friends).

See related info about orgqr.

org!

(org! or)

Destructively generates the real orthogonal matrix Q of the QR, RQ, QL or LQ factorization formed by qrf! or qrfp! (or rqf! etc.).

The input is a structure containing orthogonalized matrix :lu and vector :tau that were previously processed by qrf! (and friends). Overwrites the input with the appropriate portion of the resulting matrix Q.

See related info about orgqr.

psv

(psv a b)

Solves a system of linear equations with a positive definite symmetric coefficient matrix a, and multiple right sides matrix b. Returns the solution matrix in a new instance.

If G is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about psv!.

psv!

(psv! a b)

Destructively solves a system of linear equations with a positive definite symmetric coefficient matrix a, and multiple right sides matrix b. Overwrites b by the solution matrix.

Overwrites a with G, and b with the solution.

If G is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about gesv.

ptrf!

(ptrf! a)

Triangularizes a positive definite symmetric matrix a. Destructively computes the Cholesky factorization of a nxn matrix a, and places it in a record with the key gg.

Overwrites a with G in the lower triangle, or G transposed in the upper triangle, depending on whether a is :lower or :upper. Cholesky does not need pivoting.

If a is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about potrf.

qlf

(qlf a)

Computes the QL factorization of a GE m x n matrix.

See qrf!.

See related info about geqlf.

qlf!

(qlf! a)

Destructively computes the QL factorization of a GE m x n matrix.

See qrf!.

See related info about geqlf.

qpf

(qpf a)

Purely computes the QR factorization of a GE m x n matrix and places it in record that contains :or and :tau.

It is similar to qrf! and can replace it, with the caveat that the results have swapped columns (because of pivoting) so the solutions based on it have to be permuted back.

See qpf! and qrf!.

qpf!

(qpf! a)

Destructively computes the QR factorization with pivoting of a m x n matrix and places it in record that contains :or and :tau.

It is similar to qrf! and can replace it, with the caveat that the results have swapped columns (because of pivoting) so the solutions based on it have to be permuted back.

See related info about geqpf.

qrf

(qrf a)

Purely computes the QR factorization of a GE m x n matrix and places it in record that contains :or and :tau.

See qrf!

qrf!

(qrf! a)

Destructively computes the QR factorization of a m x n matrix and places it in record that contains :or and :tau.

The input is a matrix a. The output overwrites the contents of a. Output QR is laid out in a in the following way: The elements in the upper triangle (or trapezoid) contain the (min m n) x n upper triangular (or trapezoidal) matrix R. The elements in the lower triangle (or trapezoid) below the diagonal, with the vector tau contain Q as a product of (min m n) elementary reflectors. Other routines work with Q in this representation. If you need to compute q, call org! or org.

See related info about geqrf.

qrfp

(qrfp a)

Purely computes the QR factorization of a GE m x n matrix, with non-negative diagonal elements.

See qrfp!. See related info about geqrfp.

qrfp!

(qrfp! a)

Destructively computes the QR factorization of a GE m x n matrix, with non-negative diagonal elements.

See qrf!. See related info about geqrfp.

rqf

(rqf a)

Computes the RQ factorization of a GE m x n.

See qrf!.

See related info about gerqf.

rqf!

(rqf! a)

Destructively computes the RQ factorization of a GE m x n.

See qrf!.

See related info about gerqf.

sv

(sv a b)

Solves a system of linear equations with a square coefficient matrix a, and multiple right hand sides matrix b. Returns the solution in a new matrix instance.

If a is symmetric, tries to do Cholesky factorization first, and only does LDLt if it turns out not to be positive definite.

If U is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about sv!.

sv!

(sv! a b)

Destructively solves a system of linear equations with a square coefficient matrix a, and multiple right hand sides matrix b. Overwrites b by the solution matrix.

Overwrites a with L and U, and b with the solution. L is stored as a lower unit triangle, and U as an upper triangle. Pivot is not retained. If you need to reuse LU, use trf, and trs, or their destructive versions.

If U is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about gesv.

svd

(svd a u? vt? sdd?)(svd a u? vt?)(svd a sdd?)(svd a)

Computes the singular value decomposition of a matrix a, and returns a SVD record containing :sigma (the singular values), u (if u? is true) and vt (if vt? is true).

If sdd? is true, uses the faster divide and conquer SDD method. If sdd? is not provided, uses ordinary SVD if both u? and vt? are false.

If the reduction to bidiagonal form failed to converge, throws ex-info, with the information on the number of converged superdiagonals. If some value in the native call is illegal, throws ex-info.

See related info about svd!

svd!

(svd! a sigma u vt superb)(svd! a sigma u vt)(svd! a sigma superb)(svd! a sigma)

Computes the singular value decomposition of a matrix a.

If superb is not provided, uses faster but slightly less precise divide and conquer method (SDD).

On exit, a’s contents is destroyed, or, if u or vt are nil, overwritten with U or transposed V singular vectors of a. s is a diagonal matrix populated with sorted singular values. If the factorization does not converge, a diagonal matrix superb is populated with the unconverged superdiagonal elements (see LAPACK documentation). If called without u and vt, U and transposed V are not computed.

If the reduction to bidiagonal form failed to converge, throws ex-info, with the information on the number of converged superdiagonals. If some value in the native call is illegal, throws ex-info.

See related info about gesvd.

trf

(trf a)

Triangularizes a non-TR matrix a. Computes the LU (or LDLt, or UDUt) factorization of a mxn matrix a, and places it in a record that contains :lu, :a and :ipiv. If a is positive definite symmetric, compute the Cholesky factorization, that contains only :gg and no ipiv is needed.

Pivot is placed into the :ipiv, a vector of integers or longs (if applicable).

If U is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about trf! and ptrf!.

trf!

(trf! a)

Triangularizes a non-triangular matrix a. Destructively computes the LU (or LDLt, or UDUt, or GGt) factorization of a mxn matrix a, and places it in a record that contains :lu and :ipiv.

Overwrites a with L and U. L is stored as a lower unit triangle, and U as an upper triangle. Pivot is placed into the :ipiv, a vector of integers or longs.

If U is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about getrf.

tri

(tri a)

Computes the inverse of a triangularized matrix a. Returns the results in a new matrix instance.

If U is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info. If the matrix is not square, throws ex-info.

See related info about tri!.

tri!

(tri! a)

Destructively computes the inverse of a triangularized matrix a. Overwrites a (or :lu) with a nxn inverse.

If U is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info. If the matrix is not square, throws ex-info.

See related info about getri.

trs

(trs a b)

Solves a system of linear equations with a triangularized matrix a, with multiple right hand sides matrix b. Returns the results in a new matrix instance.

If U is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about trs!.

trs!

(trs! a b)

Destructively solves a system of linear equations with a triangularized matrix a, with multiple right hand sides matrix b. Overwrites b by the solution matrix.

If U is exactly singular (it can’t be used for solving a system of linear equations), throws ex-info.

See related info about getrs.