scikits.umfpack.UmfpackLU

class scikits.umfpack.UmfpackLU(A)

LU factorization of a sparse matrix.

Factorization is represented as:

Pr * (R^-1) * A * Pc = L * U

Parameters:

Acsc_matrix or csr_matrix

Matrix to decompose

Examples

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

>>> import numpy as np
>>> from scipy.sparse import csc_matrix
>>> from scikits import umfpack
>>> 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 = umfpack.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((4, 4))
>>> Pr[lu.perm_r, np.arange(4)] = 1
>>> Pc = csc_matrix((4, 4))
>>> Pc[np.arange(4), lu.perm_c] = 1

Similarly for the row scalings:

>>> R = csc_matrix((4, 4))
>>> R.setdiag(lu.R)

We can reassemble the original matrix:

>>> (Pr.T * R * (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

Combined number of nonzeros in L and U: L.nnz + U.nnz

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.

R

Row scaling factors, as a 1D array.

Methods

solve(b)	Solve linear equation A x = b for x
solve_sparse(B)	Solve linear equation of the form A X = B.