Algorithm 8xx: UMFPACK V3.2, an
unsymmetricpattern multifrontal method
with a column preordering strategy *
Timothy A. Davis t
January 2, 2002
Technical report TR02002. Department of Computer and Information
Science and Engineering, University of Florida.
Abstract
An ANSI C code for sparse LU factorization is presented that com
bines a leftlooking column preordering strategy with a rightlooking
unsymmetricpattern multifrontal numerical factorization. '11I. pre
ordering and symbolic analysis phase computes an upper bound on
fillin, work, and memory usage during the subsequent numerical fac
torization. Usercallable routines are provided for ordering and ana
lyzing a sparse matrix, computing the numerical factorization, solving
a system with the LU factors, transposing and permuting a sparse
matrix, and converting between sparse matrix representations. Ir.
simple user interface shields the user from the details of the com
plex sparse factorization data structures by returning simple handles
to opaque objects. Additional usercallable routines are provided for
printing and extracting the contents of these opaque objects. An even
simpler way to use the package is through its MATLAB interface. A
future version of MATLAB will incorporate UMFPACK as its default
sparse matrix factorization method.
*This work was supported by the National Science Foundation, under grants D IS
9504974 and D_ IS9803599.
tDept. of Computer and Information Science and Fii_;.. i i, Univ. of Florida,
Gainesville, FL, i S.\ email: 1,[ it, . [ti.1 ,L1 http://www.cise.ufl.edu/~davis.
Categories and Subject Descriptors: G.1.3 :'ni;til, _i .! Analysis]: .;tiil i
ical Linear Algebra linear systems (direct methods), sparse and very large
systems G.4 .l.,iI!, ii .,iii of C(i ,fi ,] .i1 111 in.l i, .i! Software algo
rithm analysis, f;. 'M .
General terms: Algorithms, Experimentation, Performance.
Keywords: sparse nonsymmetric matrices, linear equations, multifrontal
method, ordering methods.
1 Overview
U '.iFPACK Version 3.2 is a code written in ANSI C for the direct solution
of systems of linear equations, Ax = b, where A is sparse and unsymmetric.
Ti, matrix PAQ is factorized into the product LU. Ti, column ordering
Q is selected to give a good a priori upper bound on fillin and then refined
during numerical factorization (while preserving the upper bound on fillin).
Tir, row ordering P is selected during numerical factorization to maintain
numerical stability and to preserve sparsity. Thl method is a combination
of a column preordering and symbolic analysis phase based on COLA '.il)
[5, 6, 9], and a numerical factorization phase based on a modification of the
unsymmetricpattern multifrontal method (U':.iF PACK Version 2, or .1. A38
[3, 4], which does not have a preordering and symbolic analysis phase). TIi,
methods used by U':.iFPACK V3.2 are discussed in a companion paper [2].
Thi sparse matrix A factorized by U':.iFPACK must be real, square, and
nonsingular.
2 MATLAB Interface
Ti,, .iATLAB interface to U:'.iFPACK provides a replacement for '.IAT
LAB's LU routine, and the forward slash and backlash matrix operators
(when A is sparse, square, real, and nonsingular). It is typically much
faster than the builtin routines in M.iATLAB Version 6.0, uses less memory,
and returns sparser LU factors. A future version of '.IATLAB will include
U':.iFPACK as a builtin routine. Access is also provided to U'.iFPACK's
symbolic analysis phase (a version of the COLA '.il) ordering routine that
also computes information about each frontal matrix and returns the su
pernodal column elimination tree). Tih COLA '.il) ordering is suitable for
Table 1: Using U .iFPACK's .1 ATLAB interface
Function Using UMFPACK MATLAB 6.0 equivalent
Solve Ax = b. x = umfpack (A,'\',b) ; x = A \ b ;
Solve Ax = b us spparms ('autommd',0) ;
ing a different col S = spones (A) ; S = spones (A) ;
umn preordering. Q = symamd (S+S') ; Q = symamd (S+S') ;
x = umfpack (A,Q,'\',b) ; x = A (:,Q) \ b ;
x (Q) = x ;
spparms ('autommd',1)
Solve ATXT = bT. x = umfpack (b,'/',A) ; x = b / A ;
Factorize A, solve Q = colamd (A) ;
Ax = b. [L,U,P,Q] = umfpack (A) ; [L,U,P] = lu (A (:,Q))
x = U \ (L \ (b (P))) ; x = U \ (L \ (P*b)) ;
x (Q) = x ; x (Q) = x ;
most matrices, but not all, so a different input ordering can be given. Table 1
summarizes the most common uses of the U'. I iPACK mexFunction, and its
approximate equivalent in builtin M.iATLAB 6.0 routines. U':.iFPACK re
turns P and Q as permutation vectors, so that A (P, Q) is equal to L*U. In the
M.1ATLAB equivalent, P is a permutation matrix, so that P*A(: ,Q) is equal
to L*U.
3 ANSI C Interface
T1i, ANSI C U:. IiPACK library consists of 24 usercallable routines and one
include file. Twentythree of the routines come in dual versions, one that
uses into's for integers, and another that uses long's instead. Ti1, latter is
required if you wish to solve matrices requiring more than 4GB of mem
ory, assuming you can make use of the LP64 model of execution (in which
long's and pointers are both 64 bits in length). Only the int version is
described here; the long version is analogous. Both versions use double pre
cision floatingpoint arithmetic. Tli, BLAS [7] is optional, but provides the
best performance. Five UI .iPACK routines are required to solve Ax = b:
umfpack_symbolic: Preorders the columns of A to reduce fillin based
on its sparsity pattern only, finds the supernodal column elimination
tree, and postorders the tree. Returns an opaque Symbolic object as
a void pointer that can be used by umfpack_numeric to factorize A
or any other matrix with the same nonzero pattern as A. Computes
upper bounds on the nonzeros in L and U, the floatingpoint operations
required, and the memory usage of umfpacknumeric.
umfpack_numeric: iiii i i, .illy factories a sparse matrix into PAQ =
LU, using the Symbolic object. Returns an opaque Iumreric object as
a void pointer.
umfpack_solve: Solves a sparse linear system (Ax = b, ATx = b, or
systems involving just L or U), using the Iumreric object. Includes
iterative refinement with sparse backward error [1].
umfpack_free_symbolic: Frees the Symbolic object.
umfpack_free_numeric: Frees the lnumeric object.
Ti, matrix A is represented in compressed column form:
int Ap [n+1] ;
int Ai [nz] ;
double Ax [nz]
Ti, row indices and numerical values of entries in column j are stored in
Ai[Ap[j] ... Ap[j+l]l] and Ax[Ap[j] ... Ap[j+1]1], respectively. This
simple program illustrates the basic usage of U I IPACK:
#include
#include "umfpack.h"
int n = 5 ;
int Ap [ ] = {0, 2, 5, 9, 10, 12} ;
int Ai [ ] = { 0, 1, 0, 2, 4, 1, 2, 3, 4, 2, 1, 4} ;
double Ax [ ] = {2., 3., 3., 1., 4., 4., 3., 1., 2., 2., 6., 1.} ;
double b [ ] = {8., 45., 3., 3., 19.} ;
double x [5] ;
int main (int argc, char **argv)
{
double *Control = (double *) IlijiL, *Info = (double *) IlijiL ;
int i ;
void *Symbolic, *Numeric ;
umfpack_symbolic (n, Ap, Ai, &Symbolic, Control, Info)
umfpack_numeric (Ap, Ai, Ax, Symbolic, 'ilIuiit rc, Control, Info)
umfpack_free_symbolic (&Symbolic) ;
umfpack_solve ("Ax=b", Ap, Ai, Ax, x, b, Numeric, Control, Info)
umfpack_free_numeric (i lu, iric) ;
for (i = 0 ; i < n ; i++) printf ("x [%d] = 7g\n", i, x [i]) ;
Ti, Ap, Ai, and Ax arrays represent the matrix
2 3 000
3 0 406
A 0 1 3 2 0
0 0 100
0 4 201
and the solution is x = [12345]T.
Additional routines are provided for:
Cl1,.1inin.g default parameter settings, and for providing a different col
umn preordering.
Converting a matrix in triplet form to compressed column form, and
visa versa. Ti, triplet form is a simpler data structure for the user to
manipulate.
Transposing and optionally permuting a compressed column form ma
trix [8].
Getting the contents of the opaque Symbolic and Inumreric objects.
Printing and verifying the control parameters, return status, statistics,
sparse matrices, Symbolic and Iumeric objects, permutation vectors,
and dense vectors.
U:. IiPACK is available at www.cise.ufl.edu/research/sparse. T1~, pack
age includes a user guide that provides full details on how to install the
package, how to use the .iATLAB interface, and how to use the ANSI C
interface.
References
[1] : Arioli, J. W. Demmel, and I. S. Duff. Solving sparse linear systems
with sparse backward error. SIAM J. Matrix Anal. Applic., 10:165190,
1989.
[2] T. A. Davis. A column preordering strategy for the unsymmetricpattern
multifrontal method. Technical Report TR02001, Univ. of Florida,
CISE Dept., Gainesville, FL, January 2002. (www.cise.ufl.edu/tech
reports. Submitted to ACM Trans. Math. Softw.).
[3] T. A. Davis and I. S. Duff. An unsymmetricpattern multifrontal method
for sparse LU factorization. SIAM J. Matrix Anal. Applic., 18(1):140
158, 1997.
[4] T. A. Davis and I. S. Duff. A combined unifrontal/multifrontal method
for unsymmetric sparse matrices. ACM Trans. Math. Softw., 25(1):119,
1999.
[5] T. A. Davis, J. R. Gilbert, S. I. Larimore, and E. G. :... Algorithm 8xx:
COLA .i ), a column approximate minimum degree ordering algorithm.
Technical Report TR00006, Univ. of Florida, CISE Dept., Gainesville,
FL, October 2000. (www.cise.ufl.edu/techreports. Submitted to ACM
Trans. Math. Softw.).
[6] T. A. Davis, J. R. Gilbert, S. I. Larimore, and E. G. :.... A column
approximate minimum degree ordering algorithm. Technical Report TR
00005, Univ. of Florida, CISE Dept., Gainesville, FL, October 2000.
(www.cise.ufl.edu/techreports. Submitted to A CM Trans. Math. Softw.).
[7] J. J. Dongarra, J. J. Du Croz, I. S. Duff, and S. Hammarling. A set of
Level 3 Basic Linear Algebra Subprograms. ACM Trans. Math. Softw.,
16:117, 1990.
[8] F. G. Gustavson. Two fast algorithms for sparse matrices: ill; ilil ,11 i1 1i
and permuted transposition. ACM Trans. Math. Softw., 4(3):250269,
Sept. 1978.
[9] S. I. Larimore. An approximate minimum degree column ordering algo
rithm. Technical Report TR98016, Univ. of Florida, Gainesville, FL,
1998. www.cise.ufl.edu/techreports.
