Constructing Matrices

Functions for constructing matrices.

This generates a 2×2 matrix whosei,j^th entry isa[i,j]:

Here is another way to produce the same matrix:

This creates a 3×2 matrix of zeros:

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

IdentityMatrix[n] produces ann×n identity matrix:

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

MatrixForm prints the matrix in a two‐dimensional form:

Constructing special types of matrices.

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

Constructing special types of matrices withSparseArray.

This sets up a general lower‐triangular matrix:

Getting and Setting Pieces of Matrices

Ways to get pieces of matrices.

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:

This picks out the second row of the matrix:

Here is the second column of the matrix:

This picks out a submatrix:

Resetting parts of matrices.

Here is a 3×3 matrix:

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

This resets all elements in the second column to bez:

This separately resets the three elements in the second column:

This increments all the values in the second column:

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

This resets the first two rows to be new vectors:

This resets elements in the first and third columns of each row:

This resets elements in the first and third columns of rows 2 through 3:

Scalars, Vectors, and Matrices

The Wolfram Language represents matrices and vectors using lists. Anything that is not a list the Wolfram Language considers as a scalar.

A vector in the Wolfram Language consists of a list of scalars. A matrix consists of a list of vectors, representing each of its rows. In order to be a valid matrix, all the rows must be the same length, so that the elements of the matrix effectively form a rectangular array.

Functions for testing the structure of vectors and matrices.

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

Anything that is not manifestly a list is treated as a scalar, so applyingVectorQ givesFalse:

This is a 2×3 matrix:

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

This object does not count as a matrix because its rows are of different lengths:

Operations on Scalars, Vectors, and Matrices

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 attributeListable.

A consequence is that most mathematical functions are applied element by element to matrices and vectors.

TheLog applies itself separately to each element in the vector:

The same is true for a matrix, or, for that matter, for any nested list:

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

The sum of two vectors is carried out element by element:

If you try to add two vectors with different lengths, you get an error:

This adds the scalar1 to each element of the vector:

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

This multiplies each element in the vector by the scalark:

It is important to realize that the Wolfram Language treats an object as a vector in a particular operation only if the object is explicitly a list at the time when the operation is done. If the object is not explicitly a list, the Wolfram Language always treats it as a scalar. This means that you can get different results, depending on whether you assign a particular object to be a list before or after you do a particular operation.

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

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

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

Multiplying Vectors and Matrices

Different kinds of vector and matrix multiplication.

This multiplies each element of the vector by the scalark:

The "dot" operator gives the scalar product of two vectors:

You can also use dot to multiply a matrix by a vector:

Dot is also the notation for matrix multiplication in the Wolfram Language:

It is important to realize that you can use "dot" for both left‐ and right‐multiplication of vectors by matrices. The Wolfram Language makes no distinction between "row" and "column" vectors.Dot carries out whatever operation is possible. (In formal terms, contracts the last index of the tensor with the first index of.)

Here are definitions for a matrixm and a vectorv:

This left‐multiplies the vectorv bym. The objectv is effectively treated as a column vector in this case:

You can also use dot to right‐multiplyv bym. Nowv is effectively treated as a row vector:

You can multiplym byv 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 useDot to represent multiplication of such symbolic objects.

Dot effectively acts here as a noncommutative form of multiplication:

It is, nevertheless, associative:

Dot products of sums are not automatically expanded out:

You can apply the distributive law in this case using the functionDistribute, as discussed in"Structural Operations":

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 functionOuter orKroneckerProduct to do this.

The outer product of two vectors is a matrix:

The outer product of a matrix and a vector is a rank three tensor:

Outer products are discussed in more detail in"Tensors".

The Kronecker product of a matrix and a vector is a matrix:

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

Vector Operations

Basic vector operations.

This is a vector in three dimensions:

This gives a vectoru in the direction opposite tov with twice the magnitude:

This reassigns the first component ofu to be its negative:

This gives the dot product ofu andv:

This is the norm ofv:

This is the unit vector in the same direction asv:

This verifies that the norm is 1:

Transformv to have zero mean and unit sample variance:

This shows the transformed values have mean 0 and variance 1:

Two vectors are orthogonal if their dot product is zero. A set of vectors is orthonormal if they are all unit vectors and are pairwise orthogonal.

Orthogonal vector operations.

This gives the projection ofu ontov:

p is a scalar multiple ofv:

u-p is orthogonal tov:

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

When one of the vectors is linearly dependent on the vectors preceding it, the corresponding position in the result will be a zero vector:

Matrix Inversion

Matrix inversion.

Here is a simple 2×2 matrix:

This gives the inverse ofm. In producing this formula, the Wolfram Language implicitly assumes that the determinantad-bc is nonzero:

Multiplying the inverse by the original matrix should give the identity matrix:

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

Here is a matrix of rational numbers:

The Wolfram Language finds the exact inverse of the matrix:

Multiplying by the original matrix gives the identity matrix:

If you try to invert a singular matrix, the Wolfram Language prints a warning message, and returns the input unchanged:

If you give a matrix with exact symbolic or numerical entries, the Wolfram Language gives the exact inverse. If, on the other hand, some of the entries in your matrix are approximate real numbers, then the Wolfram Language finds an approximate numerical result.

Here is a matrix containing approximate real numbers:

This finds the numerical inverse:

Multiplying by the original matrix gives you an identity matrix with small round-off errors:

You can get rid of small off‐diagonal terms usingChop:

When you try to invert a matrix with exact numerical entries, the Wolfram Language can always tell whether or not the matrix is singular. When you invert an approximate numerical matrix, The Wolfram Language can usually not tell for certain whether or not the matrix is singular: all it can tell is, for example, that the determinant is small compared to the entries of the matrix. When the Wolfram Language suspects that you are trying to invert a singular numerical matrix, it prints a warning.

The Wolfram Language prints a warning if you invert a numerical matrix that it suspects is singular:

This matrix is singular, but the warning is different, and the result is useless:

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 takes the matrix, multiplies it by its inverse, and shows the first row of the result:

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

The result is still correct, but the zeros now have lower accuracy:

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

Basic Matrix Operations

Some basic matrix operations.

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

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

Det[m] gives the determinant of a square matrixm.Minors[m] is the matrix whose^th element gives the determinant of the submatrix obtained by deleting the^th row and the^th column ofm. The^th cofactor ofm is times the^th element of the matrix of minors.

Minors[m,k] gives the determinants of thek×k submatrices obtained by picking each possible set ofk rows andk columns fromm. Note that you can applyMinors 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 isa[i,j]:

Here is the determinant ofm:

Thetrace orspur of a matrixTr[m] is the sum of the terms on the leading diagonal.

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

Therank of a matrix is the number of linearly independent rows or columns.

This finds the rank of a matrix:

Powers and exponentials of matrices.

Here is a 2×2 matrix:

This gives the third matrix power ofm:

It is equivalent to multiplying three copies of the matrix:

Here is the millionth matrix power:

The matrix exponential of a matrixm is, where indicates a matrix power.

This gives the matrix exponential ofm:

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

Solving Linear Systems

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 usingSolve.

In some cases, however, you may prefer to convert the system of linear equations into a matrix equation, and then apply matrix manipulation operations to solve it. This approach is often useful when the system of equations arises as part of a general algorithm, and you do not know in advance how many variables will be involved.

A system of linear equations can be stated in matrix form as, where is the vector of variables.

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 aSparseArray object. As discussed in"Sparse Arrays: Linear Algebra", you can convert from symbolic equations toSparseArray objects usingCoefficientArrays. All the functions described here work onSparseArray objects as well as ordinary matrices.

Solving and analyzing linear systems.

Here is a 2×2 matrix:

This gives two linear equations:

You can useSolve directly to solve these equations:

You can also get the vector of solutions by callingLinearSolve. The result is equivalent to the one you get fromSolve:

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

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

If you have a square matrix with a nonzero determinant, then you can always find a unique solution to the matrix equation for any. If, however, the matrix has determinant zero, then there may be either no vector, or an infinite number of vectors which satisfy for a particular. This occurs when the linear equations embodied in are not independent.

When has determinant zero, it is nevertheless always possible to find nonzero vectors that satisfy. The set of vectors satisfying this equation form thenull space orkernel of the matrix. Any of these vectors can be expressed as a linear combination of a particular set of basis vectors, which can be obtained usingNullSpace[m].

Here is a simple matrix, corresponding to two identical linear equations:

The matrix has determinant zero:

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

There is a single basis vector for the null space ofm:

Multiplying the basis vector for the null space bym gives the zero vector:

There is only1 linearly independent row inm:

NullSpace andMatrixRank have to determine whether particular combinations of matrix elements are zero. For approximate numerical matrices, theTolerance 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 likeZeroTest->(FullSimplify[#]==0&) to force more to be done to test whether symbolic expressions are zero.

Here is a simple symbolic matrix with determinant zero:

The basis for the null space ofm contains two vectors:

Multiplyingm by any linear combination of these vectors gives zero:

An important feature of functions likeLinearSolve andNullSpace is that they work withrectangular, as well assquare, matrices.

When you represent a system of linear equations by a matrix equation of the form, the number of columns in gives the number of variables, and the number of rows gives the number of equations. There are a number of cases.

Classes of linear systems represented by rectangular matrices.

This asks for the solution to the inconsistent set of equations and:

This matrix represents two equations, for three variables:

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

When a matrix represents an underdetermined system of equations, the matrix has a nontrivial null space. In this case, the null space is spanned by a single vector:

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

The number of independent equations is therank of the matrixMatrixRank[m]. The number of redundant equations isLength[NullSpace[m]]. Note that the sum of these quantities is always equal to the number of columns inm.

GeneratingLinearSolveFunction objects.

In some applications, you will want to solve equations of the form many times with the same, but different. You can do this efficiently in the Wolfram Language by usingLinearSolve[m] to create a singleLinearSolveFunction that you can apply to as many vectors as you want.

This creates aLinearSolveFunction:

You can apply this to a vector:

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

But you can applyf to any vector you want:

Solving least-squares problems.

This linear system is inconsistent:

LeastSquares finds a vector that minimizes in the least-squares sense:

Eigenvalues and Eigenvectors

Eigenvalues and eigenvectors.

The eigenvalues of a matrix are the values for which one can find nonzero vectors such that. The eigenvectors are the vectors.

Thecharacteristic polynomialCharacteristicPolynomial[m,x] for an× matrix is given byDet[m-xIdentityMatrix[n]]. The eigenvalues are the roots of this polynomial.

Finding the eigenvalues of an× matrix in general involves solving an^th‐degree polynomial equation. For, therefore, the results cannot in general be expressed purely in terms of explicit radicals.Root objects can nevertheless always be used, although except for fairly sparse or otherwise simple matrices the expressions obtained are often unmanageably complex.

Even for a matrix as simple as this, the explicit form of the eigenvalues is quite complicated:

If you give a matrix of approximate real numbers, the Wolfram Language will find the approximate numerical eigenvalues and eigenvectors.

Here is a 2×2 numerical matrix:

The matrix has two eigenvalues, in this case both real:

Here are the two eigenvectors ofm:

Eigensystem computes the eigenvalues and eigenvectors at the same time. The assignment setsvals to the list of eigenvalues, andvecs to the list of eigenvectors:

This verifies that the first eigenvalue and eigenvector satisfy the appropriate condition:

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

The functionEigenvalues always gives you a list of eigenvalues for an× matrix. The eigenvalues correspond to the roots of the characteristic polynomial for the matrix, and may not necessarily be distinct.Eigenvectors, on the other hand, gives a list of eigenvectors which are guaranteed to be independent. If the number of such eigenvectors is less than, thenEigenvectors supplements the list with zero vectors, so that the total length of the list is always.

Here is a 3×3 matrix:

The matrix has three eigenvalues, all equal to zero:

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:

This gives the characteristic polynomial of the matrix:

Finding largest and smallest eigenvalues.

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 usingEigenvalues[m,k] andEigenvalues[m,-k].

This computes the exact eigenvalues of an integer matrix:

The eigenvalues are sorted in decreasing order of size:

This gives the three eigenvalues with largest absolute value:

Generalized eigenvalues, eigenvectors, and characteristic polynomial.

The generalized eigenvalues for a matrix with respect to a matrix are defined to be those for which.

The generalized eigenvalues correspond to zeros of the generalized characteristic polynomialDet[m-xa].

Note that while ordinary matrix eigenvalues always have definite values, some generalized eigenvalues will always beIndeterminate if the generalized characteristic polynomial vanishes, which happens if and share a null space. Note also that generalized eigenvalues can be infinite.

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

This gives a generalized characteristic polynomial:

Advanced Matrix Operations

Finding singular values and norms of matrices.

Thesingular values of a matrix are the square roots of the eigenvalues of, where denotes Hermitian transpose. The number of such singular values is the smaller dimension of the matrix.SingularValueList sorts the singular values from largest to smallest. Very small singular values are usually numerically meaningless. With the option settingTolerance->t,SingularValueList drops singular values that are less than a fractiont of the largest singular value. For approximate numerical matrices, the tolerance is by default slightly greater than zero.

If you multiply the vector for each point in a unit sphere in‐dimensional space by an× matrix, then you get an‐dimensional ellipsoid, whose principal axes have lengths given by the singular values of.

The2‐norm of a matrixNorm[m,2] is the largest principal axis of the ellipsoid, equal to the largest singular value of the matrix. This is also the maximum 2‐norm length of for any possible unit vector.

The‐norm of a matrixNorm[m,p] is in general the maximum‐norm length of that can be attained. The cases most often considered are,, and. Also sometimes considered is the Frobenius normNorm[m,"Frobenius"], which is the square root of the trace of.

Decomposing square matrices into triangular forms.

When you create aLinearSolveFunction usingLinearSolve[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.

Orthogonal decompositions of matrices.

The standard definition for the inverse of a matrix fails if the matrix is not square or is singular. Thepseudoinverse of a matrix can however still be defined. It is set up to minimize the sum of the squares of all entries in, where is the identity matrix. The pseudoinverse is sometimes known as the generalized inverse, or the Moore–Penrose inverse. It is particularly used for problems related to least‐squares fitting.

QR decomposition writes any matrix as a product, where is an orthonormal matrix, denotes Hermitian transpose, and is a triangular matrix, in which all entries below the leading diagonal are zero.

Singular value decomposition, orSVD, is an underlying element in many numerical matrix algorithms. The basic idea is to write any matrix in the form, where is a matrix with the singular values of on its diagonal, and are orthonormal matrices, and is the Hermitian transpose of.

Functions related to eigenvalue problems.

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 aJordan 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.

Numerically more stable is theSchur decomposition, which writes any square matrix in the form, where is an orthonormal matrix, and is block upper‐triangular. Also related is theHessenberg decomposition, which writes a square matrix in the form, where is an orthonormal matrix, and can have nonzero elements down to the diagonal below the leading diagonal.

Tensors

Tensors are mathematical objects that give generalizations of vectors and matrices. In the Wolfram System, a tensor is represented as a set of lists, nested to a certain number of levels. The nesting level is therank of the tensor.

Interpretations of nested lists.

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

Theindices that specify a particular element in the tensor correspond to the coordinates in the cuboid. Thedimensions of the tensor correspond to the side lengths of the cuboid.

One simple way that a rankk tensor can arise is in giving a table of values for a function ofk 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.

Functions for creating and testing the structure of tensors.

Here is a 2×3×2 tensor:

This is another way to produce the same 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:

Dimensions gives the dimensions of the tensor:

Here is the element of the tensor:

ArrayDepth gives the rank of the tensor:

The rank of a tensor is equal to the number of indices needed to specify each element. You can pick out subtensors by using a smaller number of indices.

Tensor manipulation operations.

You can think of a rankk tensor as havingk "slots" into which you insert indices. ApplyingTranspose is effectively a way of reordering these slots. If you think of the elements of a tensor as forming ak‐dimensional cuboid, you can viewTranspose 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 functionTranspose[T,{p₁,p₂,…,p_k}] gives you a new tensorT^′ such that the value ofT^′_i1i2…ik is given byT_ip1ip2…ipk.

If you originally had ann_p1×n_p2×…×n_pk tensor, then by applyingTranspose, you will get ann₁×n₂×…×n_k tensor.

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

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

The elementm[[2,3]] in the original tensor becomes the elementm[[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 useTranspose to effectively collapse different levels.

This collapses all three levels, giving a list of the elements on the "main diagonal":

This collapses only the first two levels:

You can also useTr to extract diagonal elements of a tensor.

This forms the ordinary trace of a rank 3 tensor:

Here is a generalized trace, with elements combined into a list:

This combines diagonal elements only down to level 2:

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 fromOuter[f,t,u] has a copy of the structure ofu inserted at the "position" of each element int. The elements in the resulting structure are obtained by combining elements oft andu using the functionf.

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

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

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

Here are the dimensions of the tensor:

If you take the generalized outer product of anm₁×m₂×…×m_r tensor and ann₁×n₂×…×n_s tensor, you get anm₁×…×m_r×n₁×…×n_s tensor. If the original tensors have ranksr ands, your result will be a rankr+s tensor.

In terms of indices, the result of applyingOuter to two tensorsT_i1i2…ir andU_j1j2…js is the tensorV_{i1i2…irj1j2…js} with elementsf[T_i1i2…ir,U_j1j2…js].

In doing standard tensor calculations, the most common functionf to use inOuter isTimes, corresponding to the standard outer product.

Particularly in doing combinatorial calculations, however, it is often convenient to takef to beList. UsingOuter, you can then get combinations of all possible elements in one tensor, with all possible elements in the other.

In constructingOuter[f,t,u] you effectively insert a copy ofu at every point int. To formInner[f,t,u], you effectively combine and collapse the last dimension oft and the first dimension ofu. The idea is to take anm₁×m₂×…×m_r tensor and ann₁×n₂×…×n_s tensor, withm_r=n₁, and get anm₁×m₂×…×m_r-1×n₂×…×n_s tensor as the result.

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

This gives a generalization of the standard scalar product of two vectors:

This gives a generalization of a matrix product:

Here is a 3×2×2 tensor:

Here is a 2×3×1 tensor:

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

Here are the dimensions of the result:

You can think ofInner 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 applyingInner, and then transposing the result back.

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

Treating only certain sublists in tensors as separate elements.

Here every single symbol is treated as a separate element:

But here only sublists at level 1 are treated as separate elements:

Flattening block tensors.

Here is a block matrix (a matrix of matrices that can be viewed as blocks that fit edge to edge within a larger matrix):

Here is the matrix formed by piecing the blocks together:

Sparse Arrays: Linear Algebra

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 usingSparseArray objects, as discussed in"Sparse Arrays: Manipulating Lists".SparseArray objects work by having lists of rules that specify where nonzero values appear.

Specifying sparse arrays.

As discussed in"Sparse Arrays: Manipulating Lists", you can use patterns to specify collections of elements in sparse arrays. You can also have sparse arrays that correspond to tensors of any rank.

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

This shows a visual representation of the matrix elements:

Here are the four largest eigenvalues of the matrix:

Dot gives aSparseArray result:

You can extract parts just like in an ordinary array:

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

A few structural operations that can be done directly onSparseArray objects.

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

This adds 3 to each element in the first column ofm:

Now all the elements in the first column are nonzero:

This gives the rules for the nonzero elements on the second row:

Typical ways to get sparse arrays.

This generates a tridiagonal random matrix:

Even the tenth power of the matrix is still fairly sparse:

This extracts the coefficients as sparse arrays:

Here are the corresponding ordinary arrays:

This reproduces the original forms:

CoefficientArrays can handle general polynomial equations:

The coefficients of the quadratic part are given in a rank 3 tensor:

This reproduces the original forms:

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 usingImport andExport.