SciPy

This is documentation for an old release of SciPy (version 1.6.3). Read this page in the documentation of the latest stable release (version 1.15.1).

scipy.sparse.linalg.SuperLU

class scipy.sparse.linalg.SuperLU

LU factorization of a sparse matrix.

Factorization is represented as:

Pr * A * Pc = L * U

To construct these SuperLU objects, call the splu and spilu functions.

Notes

New in version 0.14.0.

Examples

The LU decomposition can be used to solve matrix equations. Consider:

>>>
>>> import numpy as np
>>> from scipy.sparse import csc_matrix, linalg as sla
>>> A = csc_matrix([[1,2,0,4],[1,0,0,1],[1,0,2,1],[2,2,1,0.]])

This can be solved for a given right-hand side:

>>>
>>> lu = sla.splu(A)
>>> b = np.array([1, 2, 3, 4])
>>> x = lu.solve(b)
>>> A.dot(x)
array([ 1.,  2.,  3.,  4.])

The lu object also contains an explicit representation of the decomposition. The permutations are represented as mappings of indices:

>>>
>>> lu.perm_r
array([0, 2, 1, 3], dtype=int32)
>>> lu.perm_c
array([2, 0, 1, 3], dtype=int32)

The L and U factors are sparse matrices in CSC format:

>>>
>>> lu.L.A
array([[ 1. ,  0. ,  0. ,  0. ],
       [ 0. ,  1. ,  0. ,  0. ],
       [ 0. ,  0. ,  1. ,  0. ],
       [ 1. ,  0.5,  0.5,  1. ]])
>>> lu.U.A
array([[ 2.,  0.,  1.,  4.],
       [ 0.,  2.,  1.,  1.],
       [ 0.,  0.,  1.,  1.],
       [ 0.,  0.,  0., -5.]])

The permutation matrices can be constructed:

>>>
>>> Pr = csc_matrix((np.ones(4), (lu.perm_r, np.arange(4))))
>>> Pc = csc_matrix((np.ones(4), (np.arange(4), lu.perm_c)))

We can reassemble the original matrix:

>>>
>>> (Pr.T * (lu.L * lu.U) * Pc.T).A
array([[ 1.,  2.,  0.,  4.],
       [ 1.,  0.,  0.,  1.],
       [ 1.,  0.,  2.,  1.],
       [ 2.,  2.,  1.,  0.]])
Attributes
shape

Shape of the original matrix as a tuple of ints.

nnz

Number of nonzero elements in the matrix.

perm_c

Permutation Pc represented as an array of indices.

perm_r

Permutation Pr represented as an array of indices.

L

Lower triangular factor with unit diagonal as a scipy.sparse.csc_matrix.

U

Upper triangular factor as a scipy.sparse.csc_matrix.

Methods

solve(rhs[, trans])

Solves linear system of equations with one or several right-hand sides.

Previous topic

scipy.sparse.linalg.spilu

Next topic

scipy.sparse.linalg.SuperLU.solve