Linear Algebra—Wolfram Language Documentation [proxy]

Constructing Matrices	Basic Matrix Operations
Getting and Setting Pieces of Matrices	Solving Linear Systems
Scalars, Vectors, and Matrices	Eigenvalues and Eigenvectors
Operations on Scalars, Vectors, and Matrices	Advanced Matrix Operations
Multiplying Vectors and Matrices	Tensors
Vector Operations	Sparse Arrays: Linear Algebra
Matrix Inversion

Table[f,{i,m},{j,n}]	build an m×n matrix where f is a function of i and j that gives the value of the i,j ^th entry
Array[f,{m,n}]	build an m×n matrix whose i,j ^th entry is f[i,j]
ConstantArray[a,{m,n}]	build an m×n matrix with all entries equal to a
DiagonalMatrix[list]	generate a diagonal matrix with the elements of list on the diagonal
IdentityMatrix[n]	generate an n×n identity matrix
Normal[SparseArray[{{i₁,j₁}->v₁,{i₂,j₂}->v₂,…},{m,n}]]	make a matrix with nonzero values v_k at positions {i_k,j_k}

This generates a 2×2 matrix whose i,j ^th entry is a[i,j]:

This creates a 3×2 matrix of zeros:

DiagonalMatrix makes a matrix with zeros everywhere except on the leading diagonal:

This makes a 3×4 matrix with two nonzero values filled in:

MatrixForm prints the matrix in a two‐dimensional form:

Table[0,{m},{n}]	a matrix of zeros
Table[If[i>=j,1,0],{i,m},{j,n}]	a lower‐triangular matrix
RandomReal[{0,1},{m,n}]	a matrix with random numerical entries

Table evaluates If[i≥j,a++,0] separately for each element, to give a matrix with sequentially increasing entries in the lower-triangular part:

SparseArray[{},{n,n}]	a zero matrix
SparseArray[{i_,i_}->1,{n,n}]	an n×n identity matrix
SparseArray[{i_,j_}/;i>=j->1,{n,n}]	a lower‐triangular matrix

Constructing special types of matrices with SparseArray.

This sets up a general lower‐triangular matrix:

m[[i,j]]	the i,j ^th entry
m[[i]]	the i ^th row
m[[All,i]]	the i ^th column
Take[m,{i₀,i₁},{j₀,j₁}]	the submatrix with rows i₀ through i₁ and columns j₀ through j₁
m[[i₀;;i₁,j₀;;j₁]]	the submatrix with rows i₀ through i₁ and columns j₀ through j₁
m[[{i₁,…,i_r }, { j₁ , … , j_s }]]	the r×s submatrix with elements having row indices i_k and column indices j_k
Tr[m,List]	elements on the diagonal
ArrayRules[m]	positions of nonzero elements

Matrices in the Wolfram Language are represented as lists of lists. You can use all the standard Wolfram Language list‐manipulation operations on matrices.

Here is a sample 3×3 matrix:

m={{a₁₁,a₁₂,…},{a₂₁,a₂₂,…},…}	assign m to be a matrix
m[[i,j]]=a	reset element {i,j} to be a
m[[i]]=a	reset all elements in row i to be a
m[[i]]={a₁,a₂,…}	reset elements in row i to be {a₁,a₂,…}
m[[i₀;;i₁]]={v₁,v₂,…}	reset rows i₀ through i₁ to be vectors {v₁,v₂,…}
m[[All,j]]=a	reset all elements in column j to be a
m[[All,j]]={a₁,a₂,…}	reset elements in column j to be {a₁,a₂,…}
m[[i₀;;i₁,j₀;;j₁]]={{a₁₁,a₁₂,…},{a₂₁,a₂₂,…},…}	reset the submatrix with rows i₀ through i₁ and columns j₀ through j₁ to new values

Here is a 3×3 matrix:

This resets the 2, 2 element to be x, then shows the whole matrix:

This resets all elements in the second column to be z:

A range of indices can be specified by using ;; (Span).

VectorQ[expr]	give True if expr has the form of a vector, and False otherwise
MatrixQ[expr]	give True if expr has the form of a matrix, and False otherwise
Dimensions[expr]	a list of the dimensions of a vector or matrix

The list {a,b,c} has the form of a vector:

Anything that is not manifestly a list is treated as a scalar, so applying VectorQ gives False:

This is a 2×3 matrix:

For a vector, Dimensions gives a list with a single element equal to the result from Length:

Most mathematical functions in the Wolfram Language are set up to apply themselves separately to each element in a list. This is true in particular of all functions that carry the attribute Listable.

The Log applies itself separately to each element in the vector:

The differentiation function D also applies separately to each element in a list:

This adds the scalar 1 to each element of the vector:

Any object that is not manifestly a list is treated as a scalar. Here c is treated as a scalar, and added separately to each element in the vector:

This multiplies each element in the vector by the scalar k:

The object p is treated as a scalar, and added separately to each element in the vector:

This is what happens if you now replace p by the list {c,d}:

You would have gotten a different result if you had replaced p by {c,d} before you did the first operation:

cv , cm , etc.	multiply each element by a scalar
u.v , v.m , m.v , m₁.m₂ , etc.	vector and matrix multiplication
Cross[u,v]	vector cross product (also input as u×v )
Outer[Times,t,u]	outer product
KroneckerProduct[m₁,m₂,…]	Kronecker product

This multiplies each element of the vector by the scalar k:

Here are definitions for a matrix m and a vector v:

This left‐multiplies the vector v by m. The object v is effectively treated as a column vector in this case:

You can also use dot to right‐multiply v by m. Now v is effectively treated as a row vector:

You can multiply m by v on both sides to get a scalar:

For some purposes, you may need to represent vectors and matrices symbolically without explicitly giving their elements. You can use Dot to represent multiplication of such symbolic objects.

The "dot" operator gives "inner products" of vectors, matrices, and so on. In more advanced calculations, you may also need to construct outer or Kronecker products of vectors and matrices. You can use the general function Outer or KroneckerProduct to do this.

The Kronecker product of a pair of 2×2 matrices is a 4×4 matrix:

v[[i]] or Part[v,i]	give the i ^th element in the vector v
c v	scalar multiplication of c times the vector v
u.v	dot product of two vectors
Norm[v]	give the norm of v
Normalize[v]	give a unit vector in the direction of v
Standardize[v]	shift v to have zero mean and unit sample variance
Standardize[v,f₁]	shift v by f₁[v] and scale to have unit sample variance

This gives a vector u in the direction opposite to v with twice the magnitude:

This reassigns the first component of u to be its negative:

This gives the dot product of u and v:

This is the norm of v:

This is the unit vector in the same direction as v:

Transform v to have zero mean and unit sample variance:

Projection[u,v]	give the orthogonal projection of u onto v
Orthogonalize[{v₁,v_2,…}]	generate an orthonormal set from the given list of vectors

This gives the projection of u onto v:

p is a scalar multiple of v:

u-p is orthogonal to v:

Starting from the set of vectors {u,v}, this finds an orthonormal set of two vectors:

Inverse[m]

find the inverse of a square matrix

Here is a simple 2×2 matrix:

This gives the inverse of m. In producing this formula, the Wolfram Language implicitly assumes that the determinant ad-bc is nonzero:

You have to use Together to clear the denominators, and get back a standard identity matrix:

You can get rid of small off‐diagonal terms using Chop:

If you work with high‐precision approximate numbers, the Wolfram Language will keep track of the precision of matrix inverses that you generate.

This generates a 6×6 numerical matrix with entries of 20‐digit precision:

This generates a 20‐digit numerical approximation to a 6×6 Hilbert matrix. Hilbert matrices are notoriously hard to invert numerically:

Inverse works only on square matrices. "Advanced Matrix Operations" discusses the function PseudoInverse, which can also be used with nonsquare matrices.

Transpose[m]	transpose m^
ConjugateTranspose[m]	conjugate transpose m^ (Hermitian conjugate)
Inverse[m]	matrix inverse
Det[m]	determinant
Minors[m]	matrix of minors
Minors[m,k]	k ^th minors
Tr[m]	trace
MatrixRank[m]	rank of matrix

Transposing a matrix interchanges the rows and columns in the matrix. If you transpose an m×n matrix, you get an n×m matrix as the result.

Transposing a 2×3 matrix gives a 3×2 result:

Minors[m,k] gives the determinants of the k×k submatrices obtained by picking each possible set of k rows and k columns from m. Note that you can apply Minors to rectangular, as well as square, matrices.

Here is the determinant of a simple 2×2 matrix:

This generates a 3×3 matrix, whose ^th entry is a[i,j]:

Here is the determinant of m:

The trace or spur of a matrix Tr[m] is the sum of the terms on the leading diagonal.

This finds the trace of a simple 2×2 matrix:

MatrixPower[m,n]	n ^th matrix power
MatrixExp[m]	matrix exponential

Here is a 2×2 matrix:

This gives the third matrix power of m:

This gives the matrix exponential of m:

Here is an approximation to the exponential of m, based on a power series approximation:

Many calculations involve solving systems of linear equations. In many cases, you will find it convenient to write down the equations explicitly, and then solve them using Solve.

Note that if your system of equations is sparse, so that most of the entries in the matrix are zero, then it is best to represent the matrix as a SparseArray object. As discussed in "Sparse Arrays: Linear Algebra", you can convert from symbolic equations to SparseArray objects using CoefficientArrays. All the functions described here work on SparseArray objects as well as ordinary matrices.

LinearSolve[m,b]	a vector that solves the matrix equation
NullSpace[m]	a list of linearly independent vectors whose linear combinations span all solutions to the matrix equation
MatrixRank[m]	the number of linearly independent rows or columns of
RowReduce[m]	a simplified form of obtained by making linear combinations of rows

Here is a 2×2 matrix:

You can use Solve directly to solve these equations:

You can also get the vector of solutions by calling LinearSolve. The result is equivalent to the one you get from Solve:

Another way to solve the equations is to invert the matrix m, and then multiply {a,b} by the inverse. This is not as efficient as using LinearSolve:

RowReduce performs a version of Gaussian elimination and can also be used to solve the equations:

LinearSolve cannot find a solution to the equation in this case:

There is a single basis vector for the null space of m:

Multiplying the basis vector for the null space by m gives the zero vector:

There is only 1 linearly independent row in m:

NullSpace and MatrixRank have to determine whether particular combinations of matrix elements are zero. For approximate numerical matrices, the Tolerance option can be used to specify how close to zero is considered good enough. For exact symbolic matrices, you may sometimes need to specify something like ZeroTest->(FullSimplify[#]==0&) to force more to be done to test whether symbolic expressions are zero.

The basis for the null space of m contains two vectors:

Multiplying m by any linear combination of these vectors gives zero:

An important feature of functions like LinearSolve and NullSpace is that they work with rectangular, as well as square, matrices.

Underdetermined	number of equations less than the number of variables; no solutions or many solutions may exist
Overdetermined	number of independent equations more than the number of variables; solutions may or may not exist
Nonsingular	number of independent equations equal to the number of variables, and determinant nonzero; a unique solution exists
Consistent	at least one solution exists
Inconsistent	no solutions exist

LinearSolve gives one of the possible solutions to this underdetermined set of equations:

If you take the solution you get from LinearSolve, and add any linear combination of the basis vectors for the null space, you still get a solution:

The number of independent equations is the rank of the matrix MatrixRank[m]. The number of redundant equations is Length[NullSpace[m]]. Note that the sum of these quantities is always equal to the number of columns in m.

LinearSolve[m]

generate a function for solving equations of the form

You get the same result by giving the vector as an explicit second argument to LinearSolve:

But you can apply f to any vector you want:

LeastSquares[m,b]

give a vector that solves the least-squares problem

Eigenvalues[m]	a list of the eigenvalues of m
Eigenvectors[m]	a list of the eigenvectors of m
Eigensystem[m]	a list of the form {eigenvalues,eigenvectors}
Eigenvalues[N[m]], etc.	numerical eigenvalues
Eigenvalues[N[m,p]], etc.	numerical eigenvalues, starting with p‐digit precision
CharacteristicPolynomial[m,x]	the characteristic polynomial of m

Here is a 2×2 numerical matrix:

Here are the two eigenvectors of m:

Eigensystem computes the eigenvalues and eigenvectors at the same time. The assignment sets vals to the list of eigenvalues, and vecs to the list of eigenvectors:

This finds the eigenvalues of a random 4×4 matrix. For nonsymmetric matrices, the eigenvalues can have imaginary parts:

Here is a 3×3 matrix:

There is, however, only one independent eigenvector for the matrix. Eigenvectors appends two zero vectors to give a total of three vectors in this case:

Eigenvalues[m,k]	the largest k eigenvalues of m
Eigenvectors[m,k]	the corresponding eigenvectors of m
Eigensystem[m,k]	the largest k eigenvalues with corresponding eigenvectors
Eigenvalues[m,-k]	the smallest k eigenvalues of m
Eigenvectors[m,-k]	the corresponding eigenvectors of m
Eigensystem[m,-k]	the smallest k eigenvalues with corresponding eigenvectors

Eigenvalues sorts numeric eigenvalues so that the ones with large absolute value come first. In many situations, you may be interested only in the largest or smallest eigenvalues of a matrix. You can get these efficiently using Eigenvalues[m,k] and Eigenvalues[m,-k].

Eigenvalues[{m,a}]	the generalized eigenvalues of m with respect to a
Eigenvectors[{m,a}]	the generalized eigenvectors of m with respect to a
Eigensystem[{m,a}]	the generalized eigensystem of m with respect to a
CharacteristicPolynomial[{m,a},x]	the generalized characteristic polynomial of m with respect to a

The generalized eigenvalues correspond to zeros of the generalized characteristic polynomial Det[m-x a].

These two matrices share a one‐dimensional null space, so one generalized eigenvalue is Indeterminate:

SingularValueList[m]	the list of nonzero singular values of m
SingularValueList[m,k]	the k largest singular values of m
SingularValueList[{m,a}]	the generalized singular values of m with respect to a
Norm[m,p]	the p‐norm of m
Norm[m,"Frobenius"]	the Frobenius norm of m

LUDecomposition[m]	the LU decomposition
CholeskyDecomposition[m]	the Cholesky decomposition

When you create a LinearSolveFunction using LinearSolve[m], this often works by decomposing the matrix into triangular forms, and sometimes it is useful to be able to get such forms explicitly.

LU decomposition effectively factors any square matrix into a product of lower‐ and upper‐triangular matrices. Cholesky decomposition effectively factors any Hermitian positive‐definite matrix into a product of a lower‐triangular matrix and its Hermitian conjugate, which can be viewed as the analog of finding a square root of a matrix.

PseudoInverse[m]	the pseudoinverse
QRDecomposition[m]	the QR decomposition
SingularValueDecomposition[m]	the singular value decomposition
SingularValueDecomposition[{m,a}]	the generalized singular value decomposition

JordanDecomposition[m]	the Jordan decomposition
SchurDecomposition[m]	the Schur decomposition
SchurDecomposition[{m,a}]	the generalized Schur decomposition
HessenbergDecomposition[m]	the Hessenberg decomposition

Most square matrices can be reduced to a diagonal matrix of eigenvalues by applying a matrix of their eigenvectors as a similarity transformation. But even when there are not enough eigenvectors to do this, one can still reduce a matrix to a Jordan form in which there are both eigenvalues and Jordan blocks on the diagonal. Jordan decomposition in general writes any square matrix in the form as .

rank 0	scalar
rank 1	vector
rank 2	matrix
rank k	rank k tensor

A tensor of rank k is essentially a k‐dimensional table of values. To be a true rank k tensor, it must be possible to arrange the elements in the table in a k‐dimensional cuboidal array. There can be no holes or protrusions in the cuboid.

One simple way that a rank k tensor can arise is in giving a table of values for a function of k variables. In physics, the tensors that occur typically have indices which run over the possible directions in space or spacetime. Notice, however, that there is no built‐in notion of covariant and contravariant tensor indices in the Wolfram System: you have to set these up explicitly using metric tensors.

Table[f,{i₁,n₁},{i₂,n₂},…,{i_k,n_k}]
	create an n₁×n₂×…×n_k tensor whose elements are the values of f
Array[a,{n₁,n₂,…,n_k}]	create an n₁×n₂×…×n_k tensor with elements given by applying a to each set of indices
ArrayQ[t,n]	test whether t is a tensor of rank n
Dimensions[t]	give a list of the dimensions of a tensor
ArrayDepth[t]	find the rank of a tensor
MatrixForm[t]	print with the elements of t arranged in a two‐dimensional array

Here is a 2×3×2 tensor:

MatrixForm displays the elements of the tensor in a two‐dimensional array. You can think of the array as being a 2×3 matrix of column vectors:

Here is the element of the tensor:

Transpose[t]	transpose the first two indices in a tensor
Transpose[t,{p₁,p₂,…}]	transpose the indices in a tensor so that the k ^th becomes the p_k ^th
Tr[t,f]	form the generalized trace of the tensor t
Outer[f,t₁,t₂]	form the generalized outer product of the tensors t₁ and t₂ with "multiplication operator" f
t₁.t₂	form the dot product of t₁ and t₂ (last index of t₁ contracted with first index of t₂)
Inner[f,t₁,t₂,g]	form the generalized inner product, with "multiplication operator" f and "addition operator" g

You can think of a rank k tensor as having k "slots" into which you insert indices. Applying Transpose is effectively a way of reordering these slots. If you think of the elements of a tensor as forming a k‐dimensional cuboid, you can view Transpose as effectively rotating (and possibly reflecting) the cuboid.

In the most general case, Transpose allows you to specify an arbitrary reordering to apply to the indices of a tensor. The function Transpose[T,{p₁,p₂,…,p_k}] gives you a new tensor T^′ such that the value of T^′_{i₁ i₂ … i_k} is given by T_{i_p₁ i_p₂ … i_{p_k}}.

If you originally had an n_p₁×n_p₂×…×n_{p_k} tensor, then by applying Transpose, you will get an n₁×n₂×…×n_k tensor.

Here is a matrix that you can also think of as a 2×3 tensor:

Applying Transpose gives you a 3×2 tensor. Transpose effectively interchanges the two "slots" for tensor indices:

The element m[[2,3]] in the original tensor becomes the element m[[3,2]] in the transposed tensor:

This produces a 2×3×1×2 tensor:

This transposes the first two levels of t:

The result is a 3×2×1×2 tensor:

If you have a tensor that contains lists of the same length at different levels, then you can use Transpose to effectively collapse different levels.

You can also use Tr to extract diagonal elements of a tensor.

Outer products, and their generalizations, are a way of building higher‐rank tensors from lower‐rank ones. Outer products are also sometimes known as direct, tensor, or Kronecker products.

From a structural point of view, the tensor you get from Outer[f,t,u] has a copy of the structure of u inserted at the "position" of each element in t. The elements in the resulting structure are obtained by combining elements of t and u using the function f.

This gives the "outer f" of two vectors. The result is a matrix:

If you take the "outer f" of a length 3 vector with a length 2 vector, you get a 3×2 matrix:

The result of taking the "outer f" of a 2×2 matrix and a length 3 vector is a 2×2×3 tensor:

If you take the generalized outer product of an m₁×m₂×…×m_r tensor and an n₁×n₂×…×n_s tensor, you get an m₁×…×m_r×n₁×…×n_s tensor. If the original tensors have ranks r and s, your result will be a rank r+s tensor.

In terms of indices, the result of applying Outer to two tensors T_{i₁ i₂ … i_r} and U_{j₁ j₂ … j_s} is the tensor V_{i₁ i₂ … i_rj₁ j₂ … j_s} with elements f[T_{i₁ i₂ … i_r},U_{j₁ j₂ … j_s}].

In doing standard tensor calculations, the most common function f to use in Outer is Times, corresponding to the standard outer product.

Particularly in doing combinatorial calculations, however, it is often convenient to take f to be List. Using Outer, you can then get combinations of all possible elements in one tensor, with all possible elements in the other.

In constructing Outer[f,t,u] you effectively insert a copy of u at every point in t. To form Inner[f,t,u], you effectively combine and collapse the last dimension of t and the first dimension of u. The idea is to take an m₁×m₂×…×m_r tensor and an n₁×n₂×…×n_s tensor, with m_r=n₁, and get an m₁×m₂×…×m_r-1×n₂×…×n_s tensor as the result.

The simplest examples are with vectors. If you apply Inner to two vectors of equal length, you get a scalar. Inner[f,v₁,v₂,g] gives a generalization of the usual scalar product, with f playing the role of multiplication, and g playing the role of addition.

Here is a 3×2×2 tensor:

Here is a 2×3×1 tensor:

This gives a 3×2×3×1 tensor:

You can think of Inner as performing a "contraction" of the last index of one tensor with the first index of another. If you want to perform contractions across other pairs of indices, you can do so by first transposing the appropriate indices into the first or last position, then applying Inner, and then transposing the result back.

In many applications of tensors, you need to insert signs to implement antisymmetry. The function Signature[{i₁,i₂,…}], which gives the signature of a permutation, is often useful for this purpose.

Outer[f,t₁,t₂,…]	form a generalized outer product by combining the lowest‐level elements of t₁,t₂,…
Outer[f,t₁,t₂,…,n]	treat only sublists at level n as separate elements
Outer[f,t₁,t₂,…,n₁,n₂,…]	treat only sublists at level n_i in t_i as separate elements
Inner[f,t₁,t₂,g]	form a generalized inner product using the lowest‐level elements of t₁
Inner[f,t₁,t₂,g,n]	contract index n of the first tensor with the first index of the second tensor

ArrayFlatten[t,r]	create a flat rank r tensor from a rank r tensor of rank r tensors
ArrayFlatten[t]	flatten a matrix of matrices (equivalent to ArrayFlatten[t,2] )

Many large-scale applications of linear algebra involve matrices that have many elements, but comparatively few that are nonzero. You can represent such sparse matrices efficiently in the Wolfram System using SparseArray objects, as discussed in "Sparse Arrays: Manipulating Lists". SparseArray objects work by having lists of rules that specify where nonzero values appear.

SparseArray[list]	a SparseArray version of an ordinary list
SparseArray[{{i₁,j₁}->v₁,{i₂,j₂}->v₂,…},{m,n}]
	an m×n sparse array with element {i_k,j_k} having value v_k
SparseArray[{{i₁,j₁},{i₂,j₂},…}->{v₁,v₂,…},{m,n}]
	the same sparse array
Normal[array]	the ordinary list corresponding to a SparseArray

This makes a 50×50 sparse numerical matrix, with 148 nonzero elements:

You can apply most standard structural operations directly to SparseArray objects, just as you would to ordinary lists. When the results are sparse, they typically return SparseArray objects.

Dimensions[m]	the dimensions of an array
ArrayRules[m]	the rules for nonzero elements in an array
m[[i,j]]	element i, j
m[[i]]	the i ^th row
m[[All,j]]	the j ^th column
m[[i,j]]=v	reset element i, j

A few structural operations that can be done directly on SparseArray objects.

This gives the first column of m. It has only 2 nonzero elements:

This adds 3 to each element in the first column of m:

SparseArray[rules]	generate a sparse array from rules
CoefficientArrays[{eqns₁,eqns₂,…},{x₁,x₂,…}]
	get arrays of coefficients from equations
Import["file.mtx"]	import a sparse array from a file

For machine-precision numerical sparse matrices, the Wolfram System supports standard file formats such as Matrix Market (.mtx) and Harwell–Boeing. You can import and export matrices in these formats using Import and Export.