Low-level interface to UMFPACK

Routines for symbolic and numeric LU factorization of sparse matrices and for solving systems of linear equations with sparse matrices.

Tested with UMFPACK V4.4 (Jan. 28, 2005), V5.0 (May 5, 2006) Copyright (c) 2005 by Timothy A. Davis. All Rights Reserved. UMFPACK homepage: http://www.cise.ufl.edu/research/sparse/umfpack

Use ‘print UmfpackContext().funs’ to see all UMFPACK library functions the module exposes, if you need something not covered by the examples below.

References

[1] T. A. Davis, Algorithm 832: UMFPACK - an unsymmetric-pattern

multifrontal method with a column pre-ordering strategy, ACM Trans. on Mathematical Software, 30(2), 2004, pp. 196–199. https://dl.acm.org/doi/abs/10.1145/992200.992206

[2] P. Amestoy, T. A. Davis, and I. S. Duff, Algorithm 837: An approximate

minimum degree ordering algorithm, ACM Trans. on Mathematical Software, 30(3), 2004, pp. 381–388. https://dl.acm.org/doi/abs/10.1145/1024074.1024081

[3] T. A. Davis, J. R. Gilbert, S. Larimore, E. Ng, Algorithm 836: COLAMD,

an approximate column minimum degree ordering algorithm, ACM Trans. on Mathematical Software, 30(3), 2004, pp. 377–380. https://doi.org/10.1145/1024074.1024080

Module contents

UmfpackContext([family])

UMFPACK solver context

Examples

Assuming:

  • Sparse matrix in CSR or CSC format: mtx

  • Right hand side: rhs

  • Solution: sol

import scikits.umfpack as um

# Contruct the solver.
umfpack = um.UmfpackContext() # Use default 'di' family of UMFPACK routines.

# One-shot solution.
sol = umfpack( um.UMFPACK_A, mtx, rhs, autoTranspose = True )
# same as:
sol = umfpack.linsolve( um.UMFPACK_A, mtx, rhs, autoTranspose = True )

-or-

# Make LU decomposition.
umfpack.numeric( mtx )
...
# Use already LU-decomposed matrix.
sol1 = umfpack( um.UMFPACK_A, mtx, rhs1, autoTranspose = True )
sol2 = umfpack( um.UMFPACK_A, mtx, rhs2, autoTranspose = True )
# same as:
sol1 = umfpack.solve( um.UMFPACK_A, mtx, rhs1, autoTranspose = True )
sol2 = umfpack.solve( um.UMFPACK_A, mtx, rhs2, autoTranspose = True )

-or-

# Make symbolic decomposition.
umfpack.symbolic( mtx0 )
# Print statistics.
umfpack.report_symbolic()

# ...

# Make LU decomposition of mtx1 which has same structure as mtx0.
umfpack.numeric( mtx1 )
# Print statistics.
umfpack.report_numeric()

# Use already LU-decomposed matrix.
sol1 = umfpack( um.UMFPACK_A, mtx1, rhs1, autoTranspose = True )

# ...

# Make LU decomposition of mtx2 which has same structure as mtx0.
umfpack.numeric( mtx2 )
sol2 = umfpack.solve( um.UMFPACK_A, mtx2, rhs2, autoTranspose = True )

# Print all statistics.
umfpack.report_info()

-or-

# Get LU factors and permutation matrices of a matrix.
L, U, P, Q, R, do_recip = umfpack.lu( mtx )

For a given matrix A, the decomposition satisfies:

LU = PRAQ           when do_recip is true,
LU = P(R^{-1})AQ    when do_recip is false

Setting control parameters

List of control parameter names is accessible as ‘um.umfControls’ - their meaning and possible values are described in the UMFPACK documentation. To each name corresponds an attribute of the ‘um’ module, such as, for example ‘um.UMFPACK_PRL’ (controlling the verbosity of umfpack report functions). These attributes are in fact indices into the control array - to set the corresponding control array value, just do the following:

umfpack = um.UmfpackContext()
umfpack.control[um.UMFPACK_PRL] = 4 # Let's be more verbose.