Group Title: Department of Computer and Information Science and Engineering Technical Reports
Title: Algorithm 8xx : UMFPACK V4.1, an unsymmetric-pattern multifrontal method
Full Citation
Permanent Link:
 Material Information
Title: Algorithm 8xx : UMFPACK V4.1, an unsymmetric-pattern multifrontal method
Series Title: Department of Computer and Information Science and Engineering Technical Reports
Physical Description: Book
Language: English
Creator: Davis, Timothy A.
Publisher: Department of Computer and Information Science and Engineering, University of Florida
Place of Publication: Gainesville, Fla.
Publication Date: May 6, 2003
Copyright Date: 2003
 Record Information
Bibliographic ID: UF00095591
Volume ID: VID00001
Source Institution: University of Florida
Holding Location: University of Florida
Rights Management: All rights reserved by the source institution and holding location.


This item has the following downloads:

2003338 ( PDF )

Full Text

Algorithm 8xx: UMFPACK V4.1, an

unsymmetric-pattern multifrontal method

Timothy A. Davis*

May 6, 2003

An ANSI C code for sparse LU factorization is presented that combines a column
pre-ordering strategy with a right-looking unsymmetric-pattern multifrontal numerical
factorization. The pre-ordering and symbolic analysis phase computes an upper bound
on fill-in, work, and memory usage during the subsequent numerical factorization.
User-callable routines are provided for ordering and analyzing a sparse matrix, com-
puting the numerical factorization, solving a system with the LU factors, transposing
and permuting a sparse matrix, and converting between sparse matrix representations.
The simple user interface shields the user from the details of the complex sparse fac-
torization data structures by returning simple handles to opaque objects. Additional
user-callable 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 in-
terface. [ FPPACK is incorporated as a built-in operator in MATLAB 6.5 as x = A\b
when A is sparse and unsymmetric.

Categories and Subject Descriptors: G.1.3 \ii'-i idil AnIi1.-i.-]: niin.!i iil Linear Al-
gebra -linear systems (direct methods), sparse and very large systems G.4 [i'-I, !i.!e ii, of
C'! 11l,,1 iil-]: Mathematical Software ,- d',i .:I/,, analysis, ti .:' '. ..
General terms: Algorithms, Experimentation, Performance.
Keywords: sparse !iii. !.-: !ii!r1 lie matrices, linear equations, multifrontal method, order-
ing methods.

1 Overview

UMIFPACK Version 4.1 is a code written in ANSI C for the direct solution of systems of
linear equations, Ax = b, where A is sparse and tii.-: iiii l!i The matrix PAQ, PRAQ,
or PR- AQ is factorized into the product LU. The column ordering Q is selected to give a
*Dept. of Computer and Information Science and En11li1.. h 1i, Univ. of Florida, Gainesville, FL, USA.
e-mail: This work was supported by the National Science
Foundation, under grants ASC-9111263, DMS-9223088, and DMS-(i2-'i..12-'7 Portions of the work were done
while on sabbatical at Stanford University and Lawrence Berkeley National Laboratory (with funding from
Stanford University and the SciDAC program).

good a priori upper bound on fill-in and then refined during numerical factorization (while
preserving the upper bound on fill-in). The row ordering P is selected during numerical
factorization to maintain numerical stability and to preserve sparsity. The diagonal matrix
R scales the rows of the matrix A. The method is a combination of a column pre-ordering
and symbolic ,i!ll.-i., phase based on COLAII) or AIlI) and a numerical factorization
phase based on a modification of the uiii.- inliii1 i:-pattern multifrontal method (UMFPACK
Version 2, or 1. 38, which does not have a pre-ordering and symbolic iii-i1,.-i., phase). The
methods used by UMFPACK V4.1 are discussed in a companion paper [1]. The sparse matrix
A factorized by UMIFPACK can be real or complex, square or rectangular, and singular or
non-singular (or any combination).
U;I FPACK ,in 1 -.1. the matrix and automatically selects between three ordering strate-
gies, described below.

unsymmetric: A column pre-ordering is computed by a modified version of COLA ,II).
The method finds a permutation Q that limits the fill-in for any subsequent choice of
P via partial pivoting. This modified version of COLA-MII) also computes the column
elimination tree and a depth-first post-ordering of the tree. During factorization, the
column pre-ordering can be modified. Columns within a single supercolumn can be
reshuffled, to reduce fill-in. Threshold partial pivoting is used with no preference given
to the diagonal entry. Within a given pivot column j, an entry aij can be chosen if
aij > 0.1 max a,j Among those numerically acceptable entries, the sparsest row i is
chosen as the pivot row.

symmetric: The column ordering is computed from AMII) applied to the pattern of
A+AT, followed by a post-ordering of the supernodal elimination tree. No modification
of the column pre-ordering is made during numerical factorization. Threshold partial
pivoting is used, with a strong preference given to the diagonal entry. The diagonal
entry is chosen if ajj > 0.001 max aj Otherwise, a sparse row is selected, using the
same method used by the tii!.-: iiiii-l1 ii: strategy (with a relative threshold of 0.1).

2-by-2: A row permutation P2 is found which attempts to reduce the number of small
diagonal entries of P2A. If aii is numerically small, the method attempts to swap
two rows i and j, such that both aij and aji are large. Once these rows are swapped
they remain in place. This does not guarantee a zero-free diagonal, but it does tend
to preserve the symmetry of the nonzero pattern. \_'-:i, the symmetric strategy (see
above) is applied to the matrix P2A.

2 MATLAB Interface

The MIATLAB interface to UMFPACK provides a replacement for ,IATLAB's LU routine,
and the forward slash and backlash matrix operators. It is typically much faster than the
built-in routines in MAI -TLAB Version 6.0, uses less memory, and returns sparser LU factors.
II ATLAB 6.5 includes UMIFPACK Version 4.0 as a built-in routine. The interface provided
here also allows access to U-,IFPACK's pre-ordering and symbolic !il-1:.-i., phase.

3 ANSI C Interface

The ANSI C U, IFPACK library consists of 31 user-callable routines and one include file.
Twc!-l -.-i:, !-i of the routines come in four versions: real or complex (both double precision),
and int or long integers. Only the double / int version is described here; the other versions
are analogous. UMIFPACK requires the BLAS (which perform dense matrix operations) for
best performance, but can be used without the BLAS. Fi e- primary UMIFPACK routines
are required to solve Ax = b:

umfpack_di_symbolic: Pre-orders the columns of A, finds the supernodal column elim-
ination tree, and post-orders the tree. Returns an opaque Symbolic object as a void *
pointer that can be used by umfpacknumeric 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 floating-point operations required, and the memory usage of umfpack_di_numeric.

umfpack_dinumeric: \niir-i i, illy factorizes a sparse matrix PAQ, PRAQ, or PR- AQ
into the product LU, using the Symbolic object. Returns an opaque Numeric object
as a void pointer.

umfpack_disolve: Solves a sparse linear system (Ax = b, ATx = b, or systems
involving just L or U), using the Numeric object. Performs iterative refinement with
sparse backward error.

umfpack_free_symbolic: Frees the Symbolic object.

umfpack_freenumeric: Frees the Numeric object.

The matrix A is represented in compressed column form:

int Ap [n+1] ;
int Ai [nz] ;
double Ax [nz]

The row indices and numerical values of entries in column j are stored in Ai[Ap[j] ...
Ap[j+1]-1] and Ax[Ap[j] ... Ap[j+1]-11, respectively. This simple program illustrates the
basic usage of UMIFPACK:

#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 (void)
double *null = (double *) NULL ;

int i ;
void *Symbolic, *Numeric ;
(void) umfpack_di_eimbJclic (n, n, Ap, Ai, Ax, rSimbJclic, null, null)
(void) umfpack_di_numeric (Ap, Ai, Ax, Symbolic, &Numeric, null, null)
umfpack_di_free_symbolic (&Symbolic) ;
(void) umfpack_di_solve (UMFPACK_A, Ap, Ai, Ax, x, b, Numeric, null, null)
umfpack_di_free_numeric (&Numeric) ;
for (i = 0 ; i < n ; i++) printf ("x [Ed] = %g\n", i, x [i]) ;
return (0) ;

The Ap, Ai, and Ax I! i,., represent the matrix

2 3 0 0 0
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:

Changing default parameter settings, and for providing a different column pre-ordering.

Converting a matrix in triplet form to compressed column form, and visa versa. The
triplet form is a simpler data structure for the user to manipulate. It consists of three
~!Ii,,.I that hold the row index, column index, and numerical value of each entry in

Transposing and optionally permuting a compressed column form matrix.

Getting the contents of the opaque Symbolic and Numeric objects, saving them to a
file, and loading them from a file.

Printing and verifying control parameters, statistics, sparse matrices, Symbolic and
Numeric objects, permutation vectors, and dense vectors.

In addition to appearing as a Collected Algorithm of the ACM, U-IFPACK is available
at The package includes a user guide that provides
full details on how to install the package, how to use the ATLAB interface, and how to
use the ANSI C interface. A basic Fortran interface is also provided.


[1] T. A. Davis. A column pre-ordering strategy for the uii:.- iiiiirl! ii-pattern multifrontal
method. ACM Trans. Math. Sof-,.. 2003 (under submission). Also TR-03-006 at

University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - - mvs