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:

ev!

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

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

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.

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

See related info about lapacke_?geev.

glq!

(glq! a tau)

Generates the real orthogonal matrix Q of the LQ factorization formed by lqf!.

See grq!.

See related info about lapacke_?orglq.

gql!

(gql! a tau)

Generates the real orthogonal matrix Q of the QL factorization formed by qlf!.

See gqr!.

See related info about lapacke_?orgql.

gqr!

(gqr! a tau)

Generates the real orthogonal matrix Q of the QR factorization formed by qrf! or qrfp!.

The input is a GE matrix a and vector tau that were previously processed by qrf!. The matrix Q is square, m x m dimensional, but if the dimension of factorization is m x n, only the leading n columns of Q are computed. The output overwrites a with the leading n columns of Q which form an orthonormal basis in the space spanned by columns of a.

If the condition (<= (dim tau) (ncols a) (mrows a)) is not satisfied, throws ExceptionInfo.’’ If the stride of tau is not 1, or some value is illegal, throws ExceptionInfo.

See related info about lapacke_?orgqr.

grq!

(grq! a tau)

Generates the real orthogonal matrix Q of the RQ factorization formed by rqf!.

The input is a GE matrix a and vector tau that were previously processed by rqf!. The matrix Q is square, m x m dimensional, but if the dimension of factorization is m x n, only the last m rows of Q are computed. The output overwrites a with the last m rows of Q which form an orthonormal basis in the space spanned by rows of a.

If the condition (<= (dim tau) (mrows a) (ncols a)) is not satisfied, throws ExceptionInfo.’’ If the stride of tau is not 1, or some value is illegal, throws ExceptionInfo.

See related info about lapacke_?orgrq.

inv!

(inv! a)

Computes the inverse of a square matrix a, using LU factorization.

Overwrites a with its inverse.

If a is not square, or not invertible, throws ExceptionInfo.

lqf!

(lqf! a tau)(lqf! a)

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

See qrf!.

See related info about lapacke_?gelqf.

ls!

(ls! a b)

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 order (column or row oriented), throws ExceptionInfo. If the i-th element of the triangular factor of a is zero, so that a does not have full rank, the least squares cannot be computed, and the function throws ExceptionInfo. If some value in the native call is illegal, throws ExceptionInfo.

See related info about lapacke_?gels.

mlq!

(mlq! x y z)

Multiplies a real matrix by the orthogonal matrix Q formed by lqf!.

See mrq!.

See related info about lapacke_?ormlq.

mql!

(mql! x y z)

Multiplies a real matrix by the orthogonal matrix Q formed by qlf!.

See mqr!.

See related info about lapacke_?ormql.

mqr!

(mqr! x y z)

Multiplies a real matrix by the orthogonal matrix Q formed by qrf! or qrfp!.

The input is a GE matrix a and vector tau that were previously processed by qrf!, and a matrix C. The side of the multiplication (CQ or QC) is determined by the position of the vector (also the output of qrf!), which, by convenience, is after Q: the input x y z is either C Q tau, or Q tau C. The output overwrites C with the product.

The dimensions of C, Q, and tau must satisfy the following conditions, or ExceptionInfo is thrown:

  • (= (ncols a) (dim tau))
  • for Q tau C (multiplication from the left):
    • (<= (dim tau) (mrows c))
    • (= (mrows C) (mrows a))
  • for C Q tau (multiplication from the right): (<= (dim tau) (ncols c))
    • (= (ncols C) (mrows a))

If the stride of tau is not 1, or some value is illegal, throws ExceptionInfo.

See related info about lapacke_?ormqr.

mrq!

(mrq! x y z)

Multiplies a real matrix by the orthogonal matrix Q formed by rqf!.

The input is a GE matrix a and vector tau that were previously processed by rqf!, and a matrix C. The side of the multiplication (CQ or QC) is determined by the position of the vector (also the output of rqf!), which, by convenience, is after Q: the input x y z is either C Q tau, or Q tau C. The output overwrites C with the product.

The dimensions of C, Q, and tau must satisfy the following conditions, or ExceptionInfo is thrown: * (= (mrows C) (ncols a)) * (= (mrows a) (dim tau)) * for Q tau C (multiplication from the left): - (<= (dim tau) (mrows c)) * for C Q tau (multiplication from the right): - (<= (dim tau) (ncols c))

If the stride of tau is not 1, or some value is illegal, throws ExceptionInfo.

See related info about lapacke_?ormrq.

qlf!

(qlf! a tau)(qlf! a)

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

See qrf!.

See related info about lapacke_?geqlf.

qrf!

(qrf! a tau)(qrf! a)

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

The input is a GE 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 gqr!.

If the stride of tau is not 1, or some value is illegal, throws ExceptionInfo.

See related info about lapacke_?geqrf.

qrfp!

(qrfp! a tau)(qrfp! a)

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

See qrf!.

See related info about lapacke_?geqrfp.

rqf!

(rqf! a tau)(rqf! a)

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

See qrf!.

See related info about lapacke_?gerqf.

sv!

(sv! a b ipiv)(sv! a b)

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, b with the solution, and ipiv with a pivoting vector. 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 the stride of ipiv is not 1, or its dimension is not equeals to a’s number of columns, or some value is illegal, throws ExceptionInfo. If U is exactly singular (it can’t be used for solving a system of linear equations), throws ExceptionInfo.

See related info about lapacke_?gesv.

svd!

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

Computes the singular value decomposition of a matrix a.

On exit, a is destroyed, or, if u or vt are nil, overwritten with U or transposed V singular vectors of a. s is populated with sorted singular values. If the factorization does not converge, 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 ExceptionInfo, with the information on the number of converged superdiagonals. If some value in the native call is illegal, throws ExceptionInfo.

See related info about lapacke_?gesvd.

trf!

(trf! a ipiv)(trf! a)

Computes the LU factorization of a GE mxn matrix a.

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 the stride of ipiv is not 1, or its dimension is not equals to a’s number of columns, or some value is illegal, throws ExceptionInfo. If U is exactly singular (it can’t be used for solving a system of linear equations), throws ExceptionInfo.

See related info about lapacke_?getrf.

tri!

(tri! a ipiv)

Computes the inverse of a matrix a, previously LU-factorized by trf!. ipiv also originates from trf!

Overwrites a with a nxn inverse of pre-LU a. If the stride of ipiv is not 1, or its dimension is not equals to a’s number of columns, or some value is illegal, throws ExceptionInfo. If U is exactly singular (it can’t be used for solving a system of linear equations), throws ExceptionInfo.

See related info about lapacke_?getri.

trs!

(trs! lu b ipiv)(trs! lu b)

Solves a system of linear equations with an LU factored matrix lu, with multiple right hand sides matrix b. Overwrites b by the solution matrix.

If the stride of ipiv is not 1, or its dimension is not equeals to lu’s number of columns, or some value is illegal, throws ExceptionInfo. If U is exactly singular (it can’t be used for solving a system of linear equations), throws ExceptionInfo.

See related info about lapacke_?getrs.