Linear Algebra¶
Linear algebra functions in Julia are largely implemented by calling functions from LAPACK. Sparse factorizations call functions from SuiteSparse.
-
*
(A, B)¶ Matrix multiplication
-
\
(A, B)¶ Matrix division using a polyalgorithm. For input matrices
A
andB
, the resultX
is such thatA*X == B
whenA
is square. The solver that is used depends upon the structure ofA
. A direct solver is used for upper- or lower triangularA
. For HermitianA
(equivalent to symmetricA
for non-complexA
) the BunchKaufman factorization is used. Otherwise an LU factorization is used. For rectangularA
the result is the minimum-norm least squares solution computed by reducingA
to bidiagonal form and solving the bidiagonal least squares problem. For sparse, squareA
the LU factorization (from UMFPACK) is used.
-
dot
(x, y)¶ Compute the dot product
-
cross
(x, y)¶ Compute the cross product of two 3-vectors
-
norm
(a)¶ Compute the norm of a
Vector
or aMatrix
-
lu
(A) → L, U, P¶ Compute the LU factorization of
A
, such thatP*A = L*U
.
-
lufact
(A) → LU¶ Compute the LU factorization of
A
, returning anLU
object for denseA
or anUmfpackLU
object for sparseA
. The individual components of the factorizationF
can be accesed by indexing:F[:L]
,F[:U]
, andF[:P]
(permutation matrix) orF[:p]
(permutation vector). AnUmfpackLU
object has additional componentsF[:q]
(the left permutation vector) andRs
the vector of scaling factors. The following functions are available for bothLU
andUmfpackLU
objects:size
,\
anddet
. ForLU
there is also aninv
method. The sparse LU factorization is such thatL*U
is equal to``diagmm(Rs,A)[p,q]``.
-
lufact!
(A) → LU¶ lufact!
is the same aslufact
but saves space by overwriting the input A, instead of creating a copy. For sparseA
thenzval
field is not overwritten but the index fields,colptr
androwval
are decremented in place, converting from 1-based indices to 0-based indices.
-
chol
(A[, LU]) → F¶ Compute Cholesky factorization of a symmetric positive-definite matrix
A
and return the matrixF
. IfLU
isL
(Lower),A = L*L'
. IfLU
isU
(Upper),A = R'*R
.
-
cholfact
(A[, LU]) → Cholesky¶ Compute the Cholesky factorization of a dense symmetric positive-definite matrix
A
and return aCholesky
object.LU
may be ‘L’ for using the lower part or ‘U’ for the upper part. The default is to use ‘U’. The triangular matrix can be obtained from the factorizationF
with:F[:L]
andF[:U]
. The following functions are available forCholesky
objects:size
,\
,inv
,det
. ALAPACK.PosDefException
error is thrown in case the matrix is not positive definite.
-
cholfact
(A[, ll]) → CholmodFactor Compute the sparse Cholesky factorization of a sparse matrix
A
. IfA
is Hermitian its Cholesky factor is determined. IfA
is not Hermitian the Cholesky factor ofA*A'
is determined. A fill-reducing permutation is used. Methods forsize
,solve
,\
,findn_nzs
,diag
,det
andlogdet
. One of the solve methods includes an integer argument that can be used to solve systems involving parts of the factorization only. The optional boolean argument,ll
determines whether the factorization returned is of theA[p,p] = L*L'
form, whereL
is lower triangular orA[p,p] = diagmm(L,D)*L'
form whereL
is unit lower triangular andD
is a non-negative vector. The default is LDL.
-
cholfact!
(A[, LU]) → Cholesky¶ cholfact!
is the same ascholfact
but saves space by overwriting the input A, instead of creating a copy.
-
cholpfact
(A[, LU]) → CholeskyPivoted¶ Compute the pivoted Cholesky factorization of a symmetric positive semi-definite matrix
A
and return aCholeskyPivoted
object.LU
may be ‘L’ for using the lower part or ‘U’ for the upper part. The default is to use ‘U’. The triangular factors containted in the factorizationF
can be obtained withF[:L]
andF[:U]
, whereas the permutation can be obtained withF[:P]
orF[:p]
. The following functions are available forCholeskyPivoted
objects:size
,\
,inv
,det
. ALAPACK.RankDeficientException
error is thrown in case the matrix is rank deficient.
-
cholpfact!
(A[, LU]) → CholeskyPivoted¶ cholpfact!
is the same ascholpfact
but saves space by overwriting the input A, instead of creating a copy.
-
qr
(A[, thin]) → Q, R¶ Compute the QR factorization of
A
such thatA = Q*R
. Also seeqrfact
. The default is to compute a thin factorization.
-
qrfact
(A)¶ Compute the QR factorization of
A
and return aQR
object. The coomponents of the factorizationF
can be accessed as follows: the orthogonal matrixQ
can be extracted withF[:Q]
and the triangular matrixR
withF[:R]
. The following functions are available forQR
objects:size
,\
. WhenQ
is extracted, the resulting type is theQRPackedQ
object, and has the*
operator overloaded to support efficient multiplication byQ
andQ'
.
-
qrfact!
(A)¶ qrfact!
is the same asqrfact
but saves space by overwriting the input A, instead of creating a copy.
-
qrp
(A[, thin]) → Q, R, P¶ Compute the QR factorization of
A
with pivoting, such thatA*P = Q*R
, Also seeqrpfact
. The default is to compute a thin factorization.
-
qrpfact
(A) → QRPivoted¶ Compute the QR factorization of
A
with pivoting and return aQRPivoted
object. The components of the factorizationF
can be accessed as follows: the orthogonal matrixQ
can be extracted withF[:Q]
, the triangular matrixR
withF[:R]
, and the permutation withF[:P]
orF[:p]
. The following functions are available forQRPivoted
objects:size
,\
. WhenQ
is extracted, the resulting type is theQRPivotedQ
object, and has the*
operator overloaded to support efficient multiplication byQ
andQ'
. AQRPivotedQ
matrix can be converted into a regular matrix withfull
.
-
qrpfact!
(A) → QRPivoted¶ qrpfact!
is the same asqrpfact
but saves space by overwriting the input A, instead of creating a copy.
-
sqrtm
(A)¶ Compute the matrix square root of
A
. IfB = sqrtm(A)
, thenB*B == A
within roundoff error.
-
eig
(A) → D, V¶ Compute eigenvalues and eigenvectors of A
-
eigvals
(A)¶ Returns the eigenvalues of
A
.
-
eigmax
(A)¶ Returns the largest eigenvalue of
A
.
-
eigmin
(A)¶ Returns the smallest eigenvalue of
A
.
-
eigvecs
(A[, eigvals])¶ Returns the eigenvectors of
A
.For SymTridiagonal matrices, if the optional vector of eigenvalues
eigvals
is specified, returns the specific corresponding eigenvectors.
-
eigfact
(A)¶ Compute the eigenvalue decomposition of
A
and return anEigen
object. IfF
is the factorization object, the eigenvalues can be accessed withF[:values]
and the eigenvectors withF[:vectors]
. The following functions are available forEigen
objects:inv
,det
.
-
eigfact!
(A)¶ eigfact!
is the same aseigfact
but saves space by overwriting the input A, instead of creating a copy.
-
hessfact
(A)¶ Compute the Hessenberg decomposition of
A
and return aHessenberg
object. IfF
is the factorization object, the unitary matrix can be accessed withF[:Q]
and the Hessenberg matrix withF[:H]
. WhenQ
is extracted, the resulting type is theHessenbergQ
object, and may be converted to a regular matrix withfull
.
-
hessfact!
(A)¶ hessfact!
is the same ashessfact
but saves space by overwriting the input A, instead of creating a copy.
-
schurfact
(A) → Schur¶ Computes the Schur factorization of the matrix
A
. The (quasi) triangular Schur factor can be obtained from theSchur
objectF
with eitherF[:Schur]
orF[:T]
and the unitary/orthogonal Schur vectors can be obtained withF[:vectors]
orF[:Z]
such thatA=F[:vectors]*F[:Schur]*F[:vectors]'
. The eigenvalues ofA
can be obtained withF[:values]
.
-
schur
(A) → Schur[:T], Schur[:Z], Schur[:values]¶ See schurfact
-
schurfact
(A, B) → GeneralizedSchur Computes the Generalized Schur (or QZ) factorization of the matrices
A
andB
. The (quasi) triangular Schur factors can be obtained from theSchur
objectF
withF[:S]
andF[:T]
, the left unitary/orthogonal Schur vectors can be obtained withF[:left]
orF[:Q]
and the right unitary/orthogonal Schur vectors can be obtained withF[:right]
orF[:Z]
such thatA=F[:left]*F[:S]*F[:right]'
andB=F[:left]*F[:T]*F[:right]'
. The generalized eigenvalues ofA
andB
can be obtained withF[:alpha]./F[:beta]
.
-
schur
(A, B) → GeneralizedSchur[:S], GeneralizedSchur[:T], GeneralizedSchur[:Q], GeneralizedSchur[:Z] See schurfact
-
svdfact
(A[, thin]) → SVD¶ Compute the Singular Value Decomposition (SVD) of
A
and return anSVD
object.U
,S
,V
andVt
can be obtained from the factorizationF
withF[:U]
,F[:S]
,F[:V]
andF[:Vt]
, such thatA = U*diagm(S)*Vt
. Ifthin
istrue
, an economy mode decomposition is returned. The algorithm producesVt
and henceVt
is more efficient to extract thanV
. The default is to produce a thin decomposition.
-
svdfact!
(A[, thin]) → SVD¶ svdfact!
is the same assvdfact
but saves space by overwriting the input A, instead of creating a copy. Ifthin
istrue
, an economy mode decomposition is returned. The default is to produce a thin decomposition.
-
svd
(A[, thin]) → U, S, V¶ Compute the SVD of A, returning
U
, vectorS
, andV
such thatA == U*diagm(S)*V'
. Ifthin
istrue
, an economy mode decomposition is returned.
-
svdvals
(A)¶ Returns the singular values of
A
.
-
svdvals!
(A)¶ Returns the singular values of
A
, while saving space by overwriting the input.
-
svdfact
(A, B) → GeneralizedSVD Compute the generalized SVD of
A
andB
, returning aGeneralizedSVD
Factorization object, such thatA = U*D1*R0*Q'
andB = V*D2*R0*Q'
.
-
svd
(A, B) → U, V, Q, D1, D2, R0 Compute the generalized SVD of
A
andB
, returningU
,V
,Q
,D1
,D2
, andR0
such thatA = U*D1*R0*Q'
andB = V*D2*R0*Q'
.
-
svdvals
(A, B) Return only the singular values from the generalized singular value decomposition of
A
andB
.
-
triu
(M)¶ Upper triangle of a matrix
-
tril
(M)¶ Lower triangle of a matrix
-
diag
(M[, k])¶ The
k
-th diagonal of a matrix, as a vector
-
diagm
(v[, k])¶ Construct a diagonal matrix and place
v
on thek
-th diagonal
-
diagmm
(matrix, vector)¶ Multiply matrices, interpreting the vector argument as a diagonal matrix. The arguments may occur in the other order to multiply with the diagonal matrix on the left.
-
Tridiagonal
(dl, d, du)¶ Construct a tridiagonal matrix from the lower diagonal, diagonal, and upper diagonal
-
Bidiagonal
(dv, ev, isupper)¶ Constructs an upper (isupper=true) or lower (isupper=false) bidiagonal matrix using the given diagonal (dv) and off-diagonal (ev) vectors
-
Woodbury
(A, U, C, V)¶ Construct a matrix in a form suitable for applying the Woodbury matrix identity
-
rank
(M)¶ Compute the rank of a matrix
-
norm
(A[, p]) Compute the
p
-norm of a vector or a matrix.p
is2
by default, if not provided. IfA
is a vector,norm(A, p)
computes thep
-norm.norm(A, Inf)
returns the largest value inabs(A)
, whereasnorm(A, -Inf)
returns the smallest. IfA
is a matrix, valid values forp
are1
,2
, orInf
. In order to compute the Frobenius norm, usenormfro
.
-
normfro
(A)¶ Compute the Frobenius norm of a matrix
A
.
-
cond
(M[, p])¶ Matrix condition number, computed using the p-norm.
p
is 2 by default, if not provided. Valid values forp
are1
,2
, orInf
.
-
trace
(M)¶ Matrix trace
-
det
(M)¶ Matrix determinant
-
inv
(M)¶ Matrix inverse
-
pinv
(M)¶ Moore-Penrose inverse
-
null
(M)¶ Basis for null space of M.
-
repmat
(A, n, m)¶ Construct a matrix by repeating the given matrix
n
times in dimension 1 andm
times in dimension 2.
-
kron
(A, B)¶ Kronecker tensor product of two vectors or two matrices.
-
linreg
(x, y)¶ Determine parameters
[a, b]
that minimize the squared error betweeny
anda+b*x
.
-
linreg
(x, y, w) Weighted least-squares linear regression.
-
expm
(A)¶ Matrix exponential.
-
issym
(A)¶ Test whether a matrix is symmetric.
-
isposdef
(A)¶ Test whether a matrix is positive-definite.
-
istril
(A)¶ Test whether a matrix is lower-triangular.
-
istriu
(A)¶ Test whether a matrix is upper-triangular.
-
ishermitian
(A)¶ Test whether a matrix is hermitian.
-
transpose
(A)¶ The transpose operator (.’).
-
ctranspose
(A)¶ The conjugate transpose operator (‘).
BLAS Functions¶
This module provides wrappers for some of the BLAS functions for
linear algebra. Those BLAS functions that overwrite one of the input
arrays have names ending in '!'
.
Usually a function has 4 methods defined, one each for Float64
,
Float32
, Complex128
and Complex64
arrays.
-
copy!
(n, X, incx, Y, incy)¶ Copy
n
elements of arrayX
with strideincx
to arrayY
with strideincy
. ReturnsY
.
-
dot
(n, X, incx, Y, incy) Dot product of two vectors consisting of
n
elements of arrayX
with strideincx
andn
elements of arrayY
with strideincy
. There are nodot
methods forComplex
arrays.
-
nrm2
(n, X, incx)¶ 2-norm of a vector consisting of
n
elements of arrayX
with strideincx
.
-
axpy!
(n, a, X, incx, Y, incy)¶ Overwrite
Y
witha*X + Y
. ReturnsY
.
-
syrk!
(uplo, trans, alpha, A, beta, C)¶ Rank-k update of the symmetric matrix
C
asalpha*A*A.' + beta*C
oralpha*A.'*A + beta*C
according to whethertrans
is ‘N’ or ‘T’. Whenuplo
is ‘U’ the upper triangle ofC
is updated (‘L’ for lower triangle). ReturnsC
.
-
syrk
(uplo, trans, alpha, A)¶ Returns either the upper triangle or the lower triangle, according to
uplo
(‘U’ or ‘L’), ofalpha*A*A.'
oralpha*A.'*A
, according totrans
(‘N’ or ‘T’).
-
herk!
(uplo, trans, alpha, A, beta, C)¶ Methods for complex arrays only. Rank-k update of the Hermitian matrix
C
asalpha*A*A' + beta*C
oralpha*A'*A + beta*C
according to whethertrans
is ‘N’ or ‘T’. Whenuplo
is ‘U’ the upper triangle ofC
is updated (‘L’ for lower triangle). ReturnsC
.
-
herk
(uplo, trans, alpha, A)¶ Methods for complex arrays only. Returns either the upper triangle or the lower triangle, according to
uplo
(‘U’ or ‘L’), ofalpha*A*A'
oralpha*A'*A
, according totrans
(‘N’ or ‘T’).
-
gbmv!
(trans, m, kl, ku, alpha, A, x, beta, y)¶ Update vector
y
asalpha*A*x + beta*y
oralpha*A'*x + beta*y
according totrans
(‘N’ or ‘T’). The matrixA
is a general band matrix of dimensionm
bysize(A,2)
withkl
sub-diagonals andku
super-diagonals. Returns the updatedy
.
-
gbmv
(trans, m, kl, ku, alpha, A, x, beta, y)¶ Returns
alpha*A*x
oralpha*A'*x
according totrans
(‘N’ or ‘T’). The matrixA
is a general band matrix of dimensionm
bysize(A,2)
withkl
sub-diagonals andku
super-diagonals.
-
sbmv!
(uplo, k, alpha, A, x, beta, y)¶ Update vector
y
asalpha*A*x + beta*y
whereA
is a a symmetric band matrix of ordersize(A,2)
withk
super-diagonals stored in the argumentA
. The storage layout forA
is described the reference BLAS module, level-2 BLAS at <http://www.netlib.org/lapack/explore-html/>.Returns the updated
y
.
-
sbmv
(uplo, k, alpha, A, x)¶ Returns
alpha*A*x
whereA
is a symmetric band matrix of ordersize(A,2)
withk
super-diagonals stored in the argumentA
.
-
gemm!
(tA, tB, alpha, A, B, beta, C)¶ Update
C
asalpha*A*B + beta*C
or the other three variants according totA
(transposeA
) andtB
. Returns the updatedC
.
-
gemm
(tA, tB, alpha, A, B)¶ Returns
alpha*A*B
or the other three variants according totA
(transposeA
) andtB
.