Sparse Matrices¶
SparseMatrix Class Reference¶
- class sympy.matrices.sparse.SparseMatrix[source]¶
A sparse matrix (a matrix with a large number of zero elements).
Examples
>>> from sympy.matrices import SparseMatrix, ones >>> SparseMatrix(2, 2, range(4)) Matrix([ [0, 1], [2, 3]]) >>> SparseMatrix(2, 2, {(1, 1): 2}) Matrix([ [0, 0], [0, 2]])
A SparseMatrix can be instantiated from a ragged list of lists:
>>> SparseMatrix([[1, 2, 3], [1, 2], [1]]) Matrix([ [1, 2, 3], [1, 2, 0], [1, 0, 0]])
For safety, one may include the expected size and then an error will be raised if the indices of any element are out of range or (for a flat list) if the total number of elements does not match the expected shape:
>>> SparseMatrix(2, 2, [1, 2]) Traceback (most recent call last): ... ValueError: List length (2) != rows*columns (4)
Here, an error is not raised because the list is not flat and no element is out of range:
>>> SparseMatrix(2, 2, [[1, 2]]) Matrix([ [1, 2], [0, 0]])
But adding another element to the first (and only) row will cause an error to be raised:
>>> SparseMatrix(2, 2, [[1, 2, 3]]) Traceback (most recent call last): ... ValueError: The location (0, 2) is out of designated range: (1, 1)
To autosize the matrix, pass None for rows:
>>> SparseMatrix(None, [[1, 2, 3]]) Matrix([[1, 2, 3]]) >>> SparseMatrix(None, {(1, 1): 1, (3, 3): 3}) Matrix([ [0, 0, 0, 0], [0, 1, 0, 0], [0, 0, 0, 0], [0, 0, 0, 3]])
Values that are themselves a Matrix are automatically expanded:
>>> SparseMatrix(4, 4, {(1, 1): ones(2)}) Matrix([ [0, 0, 0, 0], [0, 1, 1, 0], [0, 1, 1, 0], [0, 0, 0, 0]])
A ValueError is raised if the expanding matrix tries to overwrite a different element already present:
>>> SparseMatrix(3, 3, {(0, 0): ones(2), (1, 1): 2}) Traceback (most recent call last): ... ValueError: collision at (1, 1)
See also
- property CL¶
Alternate faster representation
- LDLdecomposition(hermitian=True)[source]¶
Returns the LDL Decomposition (matrices
L
andD
) of matrixA
, such thatL * D * L.T == A
.A
must be a square, symmetric, positive-definite and non-singular.This method eliminates the use of square root and ensures that all the diagonal entries of L are 1.
Examples
>>> from sympy.matrices import SparseMatrix >>> A = SparseMatrix(((25, 15, -5), (15, 18, 0), (-5, 0, 11))) >>> L, D = A.LDLdecomposition() >>> L Matrix([ [ 1, 0, 0], [ 3/5, 1, 0], [-1/5, 1/3, 1]]) >>> D Matrix([ [25, 0, 0], [ 0, 9, 0], [ 0, 0, 9]]) >>> L * D * L.T == A True
- property RL¶
Alternate faster representation
- applyfunc(f)[source]¶
Apply a function to each element of the matrix.
Examples
>>> from sympy.matrices import SparseMatrix >>> m = SparseMatrix(2, 2, lambda i, j: i*2+j) >>> m Matrix([ [0, 1], [2, 3]]) >>> m.applyfunc(lambda i: 2*i) Matrix([ [0, 2], [4, 6]])
- as_mutable()[source]¶
Returns a mutable version of this matrix.
Examples
>>> from sympy import ImmutableMatrix >>> X = ImmutableMatrix([[1, 2], [3, 4]]) >>> Y = X.as_mutable() >>> Y[1, 1] = 5 # Can set values in Y >>> Y Matrix([ [1, 2], [3, 5]])
- cholesky(hermitian=True)[source]¶
Returns the Cholesky decomposition L of a matrix A such that L * L.T = A
A must be a square, symmetric, positive-definite and non-singular matrix
Examples
>>> from sympy.matrices import SparseMatrix >>> A = SparseMatrix(((25,15,-5),(15,18,0),(-5,0,11))) >>> A.cholesky() Matrix([ [ 5, 0, 0], [ 3, 3, 0], [-1, 1, 3]]) >>> A.cholesky() * A.cholesky().T == A True
The matrix can have complex entries:
>>> from sympy import I >>> A = SparseMatrix(((9, 3*I), (-3*I, 5))) >>> A.cholesky() Matrix([ [ 3, 0], [-I, 2]]) >>> A.cholesky() * A.cholesky().H Matrix([ [ 9, 3*I], [-3*I, 5]])
Non-hermitian Cholesky-type decomposition may be useful when the matrix is not positive-definite.
>>> A = SparseMatrix([[1, 2], [2, 1]]) >>> L = A.cholesky(hermitian=False) >>> L Matrix([ [1, 0], [2, sqrt(3)*I]]) >>> L*L.T == A True
- col_list()[source]¶
Returns a column-sorted list of non-zero elements of the matrix.
Examples
>>> from sympy.matrices import SparseMatrix >>> a=SparseMatrix(((1, 2), (3, 4))) >>> a Matrix([ [1, 2], [3, 4]]) >>> a.CL [(0, 0, 1), (1, 0, 3), (0, 1, 2), (1, 1, 4)]
- liupc()[source]¶
Liu’s algorithm, for pre-determination of the Elimination Tree of the given matrix, used in row-based symbolic Cholesky factorization.
Examples
>>> from sympy.matrices import SparseMatrix >>> S = SparseMatrix([ ... [1, 0, 3, 2], ... [0, 0, 1, 0], ... [4, 0, 0, 5], ... [0, 6, 7, 0]]) >>> S.liupc() ([[0], [], [0], [1, 2]], [4, 3, 4, 4])
References
Symbolic Sparse Cholesky Factorization using Elimination Trees, Jeroen Van Grondelle (1999) http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.39.7582
- row_list()[source]¶
Returns a row-sorted list of non-zero elements of the matrix.
Examples
>>> from sympy.matrices import SparseMatrix >>> a = SparseMatrix(((1, 2), (3, 4))) >>> a Matrix([ [1, 2], [3, 4]]) >>> a.RL [(0, 0, 1), (0, 1, 2), (1, 0, 3), (1, 1, 4)]
- row_structure_symbolic_cholesky()[source]¶
Symbolic cholesky factorization, for pre-determination of the non-zero structure of the Cholesky factororization.
Examples
>>> from sympy.matrices import SparseMatrix >>> S = SparseMatrix([ ... [1, 0, 3, 2], ... [0, 0, 1, 0], ... [4, 0, 0, 5], ... [0, 6, 7, 0]]) >>> S.row_structure_symbolic_cholesky() [[0], [], [0], [1, 2]]
References
Symbolic Sparse Cholesky Factorization using Elimination Trees, Jeroen Van Grondelle (1999) http://citeseerx.ist.psu.edu/viewdoc/summary?doi=10.1.1.39.7582
- solve(rhs, method='LDL')[source]¶
Return solution to self*soln = rhs using given inversion method.
For a list of possible inversion methods, see the .inv() docstring.
- solve_least_squares(rhs, method='LDL')[source]¶
Return the least-square fit to the data.
By default the cholesky_solve routine is used (method=’CH’); other methods of matrix inversion can be used. To find out which are available, see the docstring of the .inv() method.
Examples
>>> from sympy.matrices import SparseMatrix, Matrix, ones >>> A = Matrix([1, 2, 3]) >>> B = Matrix([2, 3, 4]) >>> S = SparseMatrix(A.row_join(B)) >>> S Matrix([ [1, 2], [2, 3], [3, 4]])
If each line of S represent coefficients of Ax + By and x and y are [2, 3] then S*xy is:
>>> r = S*Matrix([2, 3]); r Matrix([ [ 8], [13], [18]])
But let’s add 1 to the middle value and then solve for the least-squares value of xy:
>>> xy = S.solve_least_squares(Matrix([8, 14, 18])); xy Matrix([ [ 5/3], [10/3]])
The error is given by S*xy - r:
>>> S*xy - r Matrix([ [1/3], [1/3], [1/3]]) >>> _.norm().n(2) 0.58
If a different xy is used, the norm will be higher:
>>> xy += ones(2, 1)/10 >>> (S*xy - r).norm().n(2) 1.5
ImmutableSparseMatrix Class Reference¶
- class sympy.matrices.immutable.ImmutableSparseMatrix(*args, **kwargs)[source]¶
Create an immutable version of a sparse matrix.
Examples
>>> from sympy import eye >>> from sympy.matrices.immutable import ImmutableSparseMatrix >>> ImmutableSparseMatrix(1, 1, {}) Matrix([[0]]) >>> ImmutableSparseMatrix(eye(3)) Matrix([ [1, 0, 0], [0, 1, 0], [0, 0, 1]]) >>> _[0, 0] = 42 Traceback (most recent call last): ... TypeError: Cannot set values of ImmutableSparseMatrix >>> _.shape (3, 3)