%------------------------------------------------------------------------------- % The UserGuide.m4tex file. Processed into UserGuide.tex via m4. %------------------------------------------------------------------------------- % \documentstyle[12pt,fullpage]{article} \documentstyle[12pt,fullpage,times]{article} \newcommand{\m}[1]{{\bf{#1}}} % for matrices and vectors \newcommand{\tr}{^{\sf T}} % transpose \begin{document} \author{Timothy A. Davis \\ Dept. of Computer and Information Science and Engineering \\ Univ. of Florida, Gainesville, FL} \title{UMFPACK Version 3.2 User Guide} \date{January 1, 2002} \maketitle %------------------------------------------------------------------------------- \begin{abstract} UMFPACK is a set of routines for solving unsymmetric sparse linear systems, \newline $\m{Ax}=\m{b}$, using the Unsymmetric MultiFrontal method and direct sparse LU factorization. It is written in ANSI/ISO C, with a MATLAB (Version 6.0 or 6.1) interface. UMFPACK relies on the Level-3 Basic Linear Algebra Subprograms (dense matrix multiply) for its performance. \end{abstract} %------------------------------------------------------------------------------- Copyright\copyright 2002 by Timothy A. Davis, University of Florida, davis@cise.ufl.edu. All Rights Reserved. {\bf UMFPACK License:} Your use or distribution of UMFPACK or any derivative code implies that you agree to this License. THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY EXPRESSED OR IMPLIED. ANY USE IS AT YOUR OWN RISK. Permission is hereby granted to use or copy this program, provided that the Copyright, this License, and the Availability of the original version is retained on all copies. User documentation of any code that uses this code or any derivative code must cite the Copyright, this License, the Availability note, and "Used by permission." If this code or any derivative code is accessible from within MATLAB, then typing "help umfpack" must cite the Copyright, and "type umfpack" must also cite this License and the Availability note. Permission to modify the code and to distribute modified code is granted, provided the Copyright, this License, and the Availability note are retained, and a notice that the code was modified is included. This software was developed with support from the National Science Foundation, and is provided to you free of charge. {\bf Acknowledgments:} This work was supported by the National Science Foundation, under grants DMS-9504974 and DMS-9803599. %------------------------------------------------------------------------------- \newpage %------------------------------------------------------------------------------- \tableofcontents %------------------------------------------------------------------------------- \newpage \section{Overview} %------------------------------------------------------------------------------- UMFPACK Version 3.2 is a set of routines for solving systems of linear equations, $\m{Ax}=\m{b}$, when $\m{A}$ is sparse and unsymmetric. It is based on the Unsymmetric MultiFrontal method \cite{DavisDuff97,DavisDuff99}, which factorizes $\m{PAQ}$ into the product $\m{LU}$, where $\m{L}$ and $\m{U}$ are lower and upper triangular, respectively, and $\m{P}$ are $\m{Q}$ are permutation matrices. Both $\m{P}$ and $\m{Q}$ are chosen to reduce fill-in (new nonzeros in $\m{L}$ and $\m{U}$ that are not present in $\m{A}$). The permutation $\m{P}$ has the dual role of reducing fill-in and maintaining numerical accuracy (via relaxed partial pivoting and row interchanges). UMFPACK first finds a column pre-ordering that reduces fill-in, without regard to numerical values, with a modified version of COLAMD \cite{DavisGilbertLarimoreNg00_algo,DavisGilbertLarimoreNg00,Larimore98}. The method finds a symmetric permutation $\m{Q}$ of the matrix $\m{A}\tr\m{A}$ (without forming $\m{A}\tr\m{A}$ explicitly). This is a good choice for $\m{Q}$, since the Cholesky factors of $\m{(AQ)\tr(AQ)}$ are an upper bound (in terms of nonzero pattern) of the factor $\m{U}$ for the unsymmetric LU factorization ($\m{PAQ}=\m{LU}$) regardless of the choice of $\m{P}$ \cite{GeorgeNg85,GeorgeNg87,GilbertNg93}. Next, the method breaks the factorization of the matrix $\m{A}$ down into a sequence of dense rectangular frontal matrices. The frontal matrices are related to each other by a supernodal column elimination tree, in which each node in the tree represents one frontal matrix. This analysis phase also determines upper bounds on the memory usage, the floating-point operation count, and the number of nonzeros in the LU factors. UMFPACK factorizes each {\em chain} of frontal matrices in a single working array, similar to how the unifrontal method \cite{dusc:96} factorizes the whole matrix. A chain of frontal matrices is a sequence of fronts where the parent of front $i$ is $i$+1 in the supernodal column elimination tree. UMFPACK is an outer-product based, right-looking method. At the $k$-th step of Gaussian elimination, it represents the updated submatrix $\m{A}_k$ as an implicit summation of a set of dense submatrices (referred to as {\em elements}, borrowing a phrase from finite-element methods) that arise when the frontal matrices are factorized and their pivot rows and columns eliminated. Each frontal matrix represents the elimination of one or more columns; each column of $\m{A}$ will be eliminated in a specific frontal matrix, and which frontal matrix will be used for each column is determined by the pre-analysis phase. The pre-analysis phase also determines the worst-case size of each frontal matrix so that they can hold any candidate pivot column and any candidate pivot row. From the perspective of the analysis phase, any candidate pivot column in the frontal matrix is identical (in terms of nonzero pattern), and so is any row. However, the numerical factorization phase has more information than the analysis phase. It uses this information to reorder the columns within each frontal matrix to reduce fill-in. Similarly, since the number of nonzeros in each row and column are maintained (more precisely, COLMMD-style approximate degrees \cite{GilbertMolerSchreiber}), a pivot row can be selected based on sparsity-preserving criteria (low degree) as well as numerical considerations (relaxed threshold partial pivoting). This information about row and column degrees is not available to left-looking methods such as {\tt SuperLU} \cite{SuperLU99} or MATLAB's {\tt LU} \cite{GilbertMolerSchreiber,GilbertPeierls88}. More details of the method, including experimental results, are described in \cite{Davis02,Davis02_algo}, available at www.cise.ufl.edu/tech-reports. %------------------------------------------------------------------------------- \section{Availability} %------------------------------------------------------------------------------- UMFPACK Version 3.2 is available at www.cise.ufl.edu/research/sparse, and has been submitted as a collected algorithm of the ACM \cite{Davis02,Davis02_algo}. It makes use of a modified version of COLAMD V2.0 by Timothy A.~Davis, Stefan Larimore, John Gilbert, and Esmond Ng. The original COLAMD V2.0 is available in MATLAB V6.0 (or later), and at www.cise.ufl.edu/research/sparse. These codes are also available in Netlib \cite{netlib} at www.netlib.org. Prior versions of UMFPACK, co-authored with Iain Duff, are available at www.cise.ufl.edu/research/sparse and as MA38 (functionally equivalent to Version 2.2.1) in the Harwell Subroutine Library. %------------------------------------------------------------------------------- \section{Using UMFPACK in MATLAB} %------------------------------------------------------------------------------- The easiest way to use UMFPACK is within MATLAB. This discussion assumes that you have MATLAB Version 6.0 or later (which includes the BLAS, and the {\tt colamd} ordering routine). To compile the UMFPACK mexFunction, you can either type {\tt make umfpack} in the Unix system shell, or type {\tt umfpack\_make} in MATLAB (which should work on any system, including Windows). See Section~\ref{Install} for more details on how to install UMFPACK. Once installed, the UMFPACK mexFunction can analyze, factor, and solve linear systems. Table~\ref{matlab} summarizes some of the more common uses of UMFPACK within MATLAB. \begin{table} \caption{Using UMFPACK's MATLAB interface} \label{matlab} \vspace{0.1in} {\footnotesize \begin{tabular}{l|l|l} \hline Function & Using UMFPACK & MATLAB 6.0 equivalent \\ \hline & & \\ \begin{minipage}[t]{1.5in} Solve $\m{Ax}=\m{b}$. \end{minipage} & \begin{minipage}[t]{2.2in} \begin{verbatim} x = umfpack (A,'\',b) ; \end{verbatim} \end{minipage} & \begin{minipage}[t]{2.2in} \begin{verbatim} x = A \ b ; \end{verbatim} \end{minipage} \\ & & \\ \hline & & \\ \begin{minipage}[t]{1.5in} Solve $\m{Ax}=\m{b}$ using a different column pre-ordering. \end{minipage} & \begin{minipage}[t]{2.2in} \begin{verbatim} S = spones (A) ; Q = symamd (S+S') ; x = umfpack (A,Q,'\',b) ; \end{verbatim} \end{minipage} & \begin{minipage}[t]{2.2in} \begin{verbatim} spparms ('autommd',0) ; S = spones (A) ; Q = symamd (S+S') ; x = A (:,Q) \ b ; x (Q) = x ; spparms ('autommd',1) ; \end{verbatim} \end{minipage} \\ & & \\ \hline & & \\ \begin{minipage}[t]{1.5in} Solve $\m{A}\tr\m{x}\tr = \m{b}\tr$. \end{minipage} & \begin{minipage}[t]{2.2in} \begin{verbatim} x = umfpack (b,'/',A) ; \end{verbatim} \end{minipage} & \begin{minipage}[t]{2.2in} \begin{verbatim} x = b / A ; \end{verbatim} \end{minipage} \\ & & \\ \hline & & \\ \begin{minipage}[t]{1.5in} Factorize $\m{A}$, then solve $\m{Ax}=\m{b}$. \end{minipage} & \begin{minipage}[t]{2.2in} \begin{verbatim} [L,U,P,Q] = umfpack (A) ; x = U \ (L \ (b (P))) ; x (Q) = x ; \end{verbatim} \end{minipage} & \begin{minipage}[t]{2.2in} \begin{verbatim} Q = colamd (A) ; [L,U,P] = lu (A (:,Q)) ; x = U \ (L \ (P*b)) ; x (Q) = x ; \end{verbatim} \end{minipage} \\ & & \\ \hline \end{tabular} } \end{table} UMFPACK requires {\tt A} to be square, sparse, nonsingular and not complex, and it requires {\tt b} to be a dense column vector (and not complex). This is more restrictive than what you can do with MATLAB's backslash or {\tt LU}. Future releases of UMFPACK may remove some of these restrictions. MATLAB's {\tt [L,U,P] = lu (A)} returns a lower triangular {\tt L}, an upper triangular {\tt U}, and a permutation matrix {\tt P} such that {\tt P*A} is equal to {\tt L*U}. UMFPACK behaves differently; it returns {\tt P} and {\tt Q} such that {\tt A (P,Q)} is equal to {\tt L*U}, where {\tt P} and {\tt Q} are permutation vectors. If you prefer permutation matrices, use the following MATLAB code: {\footnotesize \begin{verbatim} [L,U,P,Q] = umfpack (A) ; n = size (A,1) ; I = speye (n) ; P = I (P,:) ; Q = I (:,Q) ; \end{verbatim} } Now {\tt P*A*Q} is equal to {\tt L*U}. Note that {\tt x = umfpack (A, '}$\backslash${\tt ', b)} requires that {\tt b} be a dense column vector. If you wish to use the LU factors from UMFPACK to solve a linear system, $\m{Ax}=\m{b}$ where $\m{b}$ is a either a dense or sparse matrix with more than one column, do this: {\footnotesize \begin{verbatim} [L,U,P,Q] = umfpack (A) ; x = U \ (L \ (b (P,:))) ; x (Q,:) = x ; \end{verbatim} } The above two examples do not make use of the iterative refinement that is built into {\tt x = umfpack (A, '}$\backslash${\tt ', b)} however. Since the numeric factorization refines its column pre-ordering, the {\tt Q} in {\tt [L,U,P,Q] = umfpack (A)} and {\tt [Q,F,C] = umfpack (A, 'symbolic')} will in general be different. There are more options; you can provide your own column pre-ordering (in which case UMFPACK does not call COLAMD), you can modify other control settings (similar to the {\tt spparms} in MATLAB), and you get various statistics on the analysis, factorization, and solution of the linear system. Type {\tt help umfpack\_details} and {\tt help umfpack\_report} in MATLAB for more information. Two demo m-files are provided. Just type {\tt umfpack\_simple} and {\tt umfpack\_demo} to run them. They roughly correspond to the C programs {\tt umfpack\_simple.c} and {\tt umfpack\_demo.c}. You may want to type {\tt more on} before running the {\tt umfpack\_simple} demo since it generates lots of output. A simple M-file ({\tt umfpack\_btf}) is provided that first permutes the matrix to upper block triangular form, using MATLAB's {\tt dmperm} routine. It solves $\m{Ax}=\m{b}$; the LU factors are not returned. Its usage is simple: {\tt x = umfpack\_btf (A, b)}. Type {\tt help umfpack\_btf} for more options. One issue you may encounter is how UMFPACK allocates its memory when being used in a mexFunction. One part of its working space is of variable size. The symbolic analysis phase determines an upper bound on the size of this memory, but not all of this memory will typically be used in the numerical factorization. UMFPACK tries to allocate a decent amount of working space (70\% of the upper bound, by default), with some elbow room so that it can run more efficiently. If this fails, it reduces its request and uses less memory. However, {\tt mxMalloc} aborts the {\tt umfpack} mexFunction if it fails, so this strategy doesn't work in MATLAB. The strategy works fine when {\tt malloc} is used instead. If you run out of memory in MATLAB, try reducing {\tt Control(7)} to be less than 0.70, and try again. Alternatively, set {\tt Control(7)} to 1.0 or 1.05 to avoid all reallocations of memory. Type {\tt help umfpack\_details} and {\tt umfpack\_report} for more information, and refer to the {\tt Control [UMFPACK\_ALLOC\_INIT]} parameter described in {\tt umfpack\_numeric} in Section~\ref{Primary}, below. There is a solution to this problem, but it relies on undocumented internal routines. See the {\tt -DMATHWORKS} option in {\tt umf\_config.h} in Section~\ref{Config} for details. Memory allocation on a PC is notoriously bad, so I recommend setting {\tt Control(7)} to a non-default value of 1.0 or even 1.05. This will avoid most reallocations of memory. %------------------------------------------------------------------------------- \section{Using UMFPACK in a C program} %------------------------------------------------------------------------------- The C-callable UMFPACK library consists of 24 user-callable routines and one include file. Twenty-three of the routines come in dual versions, with different sizes of integers. All user-callable routine names begin with {\tt umfpack\_} or {\tt umfpack\_l\_}; other routine names beginning with {\tt umf\_} or {\tt umfl\_} are internal to the package, and should not be called by the user. The include file {\tt umfpack.h}, listed in Section~\ref{Include}, must be included in any C program that uses UMFPACK. %------------------------------------------------------------------------------- \subsection{The size of an integer} %------------------------------------------------------------------------------- There are two versions of each user-callable routine (except for one routine). The routine names starting with just {\tt umfpack\_} use {\tt int} integer arguments; those starting with {\tt umfpack\_l\_} use {\tt long} integer arguments. If you compile UMFPACK in the standard ILP32 mode (32-bit {\tt int}'s, {\tt long}'s, an pointers) then the versions are essentially identical. You will be able to solve problems using up to 4GB of memory. If you compile UMFPACK in the standard LP64 mode, the size of an {\tt int} remains 32-bits, but the size of a {\tt long} and a pointer both get promoted to 64-bits. In the LP64 mode, the {\tt umfpack\_l\_*} routines can solve huge problems (not limited to 4GB), limited of course by the amount of available memory. The only drawback to the 64-bit mode is that few BLAS libraries support 64-bit integers. This limits the performance you will obtain. Both versions are discussed below. Use only one version for any one problem; do not attempt to use one version to analyze the matrix and another version to factorize the matrix, for example. %------------------------------------------------------------------------------- \subsection{Primary routines, and a simple example} %------------------------------------------------------------------------------- Five primary UMFPACK routines are required to solve $\m{Ax}=\m{b}$. They are fully described in Section~\ref{Primary}: \begin{itemize} \item {\tt umfpack\_symbolic}, {\tt umfpack\_l\_symbolic}: Pre-orders the columns of $\m{A}$ to reduce fill-in, based on its sparsity pattern only, finds the supernodal column elimination tree, and post-orders the tree. Returns an opaque {\tt Symbolic} object as a {\tt void *} pointer. The object contains the symbolic analysis and is needed for the numerical factorization. This routine requires only $O(|\m{A}|)$ space, where $|\m{A}|$ is the number of nonzero entries in the matrix. It computes upper bounds on the nonzeros in $\m{L}$ and $\m{U}$, the floating-point operations required, and the memory usage of {\tt umfpack\_numeric}. The {\tt Symbolic} object is small; it contains just the column pre-ordering, the supernodal column elimination tree, and information about each frontal matrix, and is no larger than about $6n$ integers (where $\m{A}$ is $n$-by-$n$). The matrix must be structurally non-singular (more precisely, each row and column must have at least one entry). \item {\tt umfpack\_numeric}, {\tt umfpack\_l\_numeric}: Numerically factorizes a sparse matrix into $\m{PAQ}=\m{LU}$. Requires the symbolic ordering and analysis computed by {\tt umfpack\_symbolic} or {\tt umfpack\_qsymbolic}. Returns an opaque {\tt Numeric} object as a {\tt void *} pointer. The object contains the numerical factorization and is used by {\tt umfpack\_solve}. You can factorize a new matrix with a different values (but identical pattern) as the matrix analyzed by {\tt umfpack\_symbolic} or {\tt umfpack\_qsymbolic} by re-using the {\tt Symbolic} object (this feature is available when using UMFPACK in a C program, but not in MATLAB). The matrix must be non-singular. \item {\tt umfpack\_solve}, {\tt umfpack\_l\_solve}: Solves a sparse linear system ($\m{Ax}=\m{b}$, $\m{A}\tr\m{x}=\m{b}$, or systems involving just $\m{L}$ or $\m{U}$), using the numeric factorization computed by {\tt umfpack\_numeric}. Iterative refinement with sparse backward error \cite{ardd:89} is used by default. \item {\tt umfpack\_free\_symbolic}, {\tt umfpack\_l\_free\_symbolic}: Frees the {\tt Symbolic} object created by {\tt umfpack\_symbolic} or {\tt umfpack\_qsymbolic}. \item {\tt umfpack\_free\_numeric}, {\tt umfpack\_l\_free\_numeric}: Frees the {\tt Numeric} object created by {\tt umfpack\_numeric}. \end{itemize} Be careful not to free a {\tt Symbolic} object with {\tt umfpack\_free\_numeric}. Nor should you attempt to free a {\tt Numeric} object with {\tt umfpack\_free\_symbolic}. Failure to free these objects will lead to memory leaks. The matrix $\m{A}$ is represented in compressed column form, which is identical to the sparse matrix representation used by MATLAB. It consists of three arrays, where the matrix is {\tt n}-by{\tt n}, with {\tt nz} entries. For the {\tt int} version of UMFPACK: {\footnotesize \begin{verbatim} int Ap [n+1] ; int Ai [nz] ; double Ax [nz] ; \end{verbatim} } For the {\tt long} version of UMFPACK: {\footnotesize \begin{verbatim} long Ap [n+1] ; long Ai [nz] ; double Ax [nz] ; \end{verbatim} } All nonzeros are entries, but an entry may be numerically zero. The row indices of entries in column {\tt j} are stored in {\tt Ai[Ap[j]} ... {\tt Ap[j+1]-1]}. The corresponding numerical values are stored in {\tt Ax[Ap[j]} ... {\tt Ap[j+1]-1]}. No duplicate row indices may be present, and the row indices in any given column must be sorted in ascending order. The first entry {\tt Ap [0]} must be zero. The total number of entries in the matrix is thus {\tt nz = Ap [n]}. Except for the fact that extra zero entries can be included, there is thus a unique compressed column representation of any given matrix $\m{A}$. Here is a simple main program, {\tt umfpack\_simple.c}, that illustrates the basic usage of UMFPACK. {\footnotesize \begin{verbatim} include(umfpack_simple.c) \end{verbatim} } It solves the same linear system as the {\tt umfpack\_simple.m} MATLAB m-file. The {\tt Ap}, {\tt Ai}, and {\tt Ax} arrays represent the matrix \[ \m{A} = \left[ \begin{array}{rrrrr} 2 & 3 & 0 & 0 & 0 \\ 3 & 0 & 4 & 0 & 6 \\ 0 & -1 & -3 & 2 & 0 \\ 0 & 0 & 1 & 0 & 0 \\ 0 & 4 & 2 & 0 & 1 \\ \end{array} \right]. \] and the solution is $\m{x} = [1 \, 2 \, 3 \, 4 \, 5]\tr$. The program uses default control settings and does not return any statistics about the ordering, factorization, or solution ({\tt Control} and {\tt Info} are both {\tt (double *) NULL}). %------------------------------------------------------------------------------- \subsection{Alternative routines} %------------------------------------------------------------------------------- Three alternative routines are provided that modify UMFPACK's default behavior. They are fully described in Section~\ref{Alternative}: \begin{itemize} \item {\tt umfpack\_defaults}, {\tt umfpack\_l\_defaults}: Sets the default control parameters in the {\tt Control} array. These can then be modified as desired before passing the array to the other UMFPACK routines. Control parameters are fully described in Section~\ref{defaults}. One particular parameter deserves special notice. UMFPACK uses relaxed partial pivoting, where a candidate pivot entry is numerically acceptable if its magnitude is greater than or equal to a tolerance parameter times the magnitude of the largest entry in the same column. The parameter {\tt Info [UMFPACK\_PIVOT\_TOLERANCE]} has a default value of 0.1. This may be too small for some matrices, particularly for ill-conditioned or poorly scaled ones. With the default pivot tolerance and default iterative refinement, {\tt x = umfpack (A, '}$\backslash${\tt ', b)} is just as accurate as {\tt x = A}$\backslash${\tt b} in MATLAB for nearly all matrices. \item {\tt umfpack\_qsymbolic}, {\tt umfpack\_l\_qsymbolic}: An alternative to {\tt umfpack\_symbolic}. Allows the user to specify his or her own column pre-ordering, rather than using the default COLAMD pre-ordering. \item {\tt umfpack\_wsolve}, {\tt umfpack\_l\_wsolve}: An alternative to {\tt umfpack\_solve} which does not dynamically allocate any memory. Requires the user to pass several additional size-$n$ work arrays. \end{itemize} %------------------------------------------------------------------------------- \subsection{Matrix manipulation routines} %------------------------------------------------------------------------------- The compressed column data structure is compact, and simplifies the UMFPACK routines that operate on the sparse matrix $\m{A}$. However, it can be inconvenient for the user to generate. Section~\ref{Manipulate} presents the details of routines for manipulating sparse matrices in {\em triplet} form, compressed column form, and compressed row form (the transpose of the compressed column form). The triplet form of a matrix consists of three arrays. For the {\tt int} version of UMFPACK: {\footnotesize \begin{verbatim} int Ti [nz] ; int Tj [nz] ; double Tx [nz] ; \end{verbatim} } For the {\tt long} version: {\footnotesize \begin{verbatim} long Ti [nz] ; long Tj [nz] ; double Tx [nz] ; \end{verbatim} } The {\tt k}-th triplet is $(i,j,a_{ij})$, where $i =$ {\tt Ti [k]}, $j =$ {\tt Tj [k]}, and $a_{ij} =$ {\tt Tx [k]}. The triplets can be in any order in the {\tt Ti}, {\tt Tj}, and {\tt Tx} arrays, and duplicate entries may exist. Any duplicate entries are summed when the triplet form is converted to compressed column form. This is a convenient way to create a matrix arising in finite-element methods, for example. Three routines are provided for manipulating sparse matrices: \begin{itemize} \item {\tt umfpack\_triplet\_to\_col}, {\tt umfpack\_l\_triplet\_to\_col}: Converts a triplet form of a matrix to compressed column form (ready for input to {\tt umfpack\_symbolic}, {\tt umfpack\_qsymbolic}, and {\tt umfpack\_numeric}). Identical to {\tt A = spconvert (i,j,x)} in MATLAB, except that zero entries are not removed, so that the pattern of entries in the compressed column form of $\m{A}$ are fully under user control. This is important if you want to factorize a new matrix with the {\tt Symbolic} object from a prior matrix with the same pattern as the new one. MATLAB never stores explicitly zero entries. \item {\tt umfpack\_col\_to\_triplet}, {\tt umfpack\_l\_col\_to\_triplet}: The opposite of {\tt umfpack\_triplet\_to\_col}. Identical to {\tt [i,j,x] = find (A)} in MATLAB, except that numerically zero entries may be included. \item {\tt umfpack\_transpose}, {\tt umfpack\_l\_transpose}: Transposes and optionally permutes a column form matrix \cite{Gustavson78}. Identical to {\tt B = A (P,Q)'} in MATLAB, except for the presence of numerically zero entries. \end{itemize} It is quite easy to add matrices in triplet form, transpose them, and permute them. See the discussion in {\tt umfpack\_triplet\_to\_col} in Section~\ref{Manipulate} for more details. All of the matrix manipulation routines can correctly operate on singular matrices. %------------------------------------------------------------------------------- \subsection{Getting the contents of opaque objects} %------------------------------------------------------------------------------- There are cases where the user would like to do more with the LU factorization of a matrix than solve a linear system. The opaque {\tt Symbolic} and {\tt Numeric} objects are just that - opaque. In addition, the LU factors are stored in the {\tt Numeric} object in a compact way that does not store all of the row and column indices \cite{GeorgeLiu}. These objects may not be dereferenced by the user, and even if they were, it would be difficult for a typical user to understand how the LU factors are stored. Thus, three routines are provided for copying their contents into user-provided arrays using simpler data structures. They are fully described in Section~\ref{Get}: \begin{itemize} \item {\tt umfpack\_get\_lunz}, {\tt umfpack\_l\_get\_lunz}: Returns the number of nonzeros in $\m{L}$ and $\m{U}$. \item {\tt umfpack\_get\_numeric}, {\tt umfpack\_l\_get\_numeric}: Copies $\m{L}$, $\m{U}$, $\m{P}$, and $\m{Q}$ from the {\tt Numeric} object into arrays provided by the user. The matrix $\m{L}$ is returned in compressed row form (with the column indices in each row sorted in ascending order). The matrix $\m{U}$ is returned in compressed column form (also with sorted columns). There are no explicit zero entries in $\m{L}$ and $\m{U}$, but such entries may exist in the {\tt Numeric} object. The permutations $\m{P}$ and $\m{Q}$ are represented as permutation vectors, where {\tt P [k] = i} means that row {\tt i} of the original matrix is the {\tt k}-th pivot row (or the {\tt k}-th row of $\m{PAQ}$), and where {\tt Q [k] = j} means that column {\tt j} of the original matrix is the {\tt k}-th pivot column. This is identical to how MATLAB uses permutation vectors. \item {\tt umfpack\_get\_symbolic}, {\tt umfpack\_l\_get\_symbolic}: Copies the contents of the {\tt Symbolic} object (the initial column preordering, and supernodal column elimination tree, and information about each frontal matrix) into arrays provided by the user. \end{itemize} UMFPACK itself does not make use of the output of the {\tt umfpack\_get\_*} routines; they are provided solely for returning the contents of the opaque {\tt Symbolic} and {\tt Numeric} objects to the user. %------------------------------------------------------------------------------- \subsection{Reporting routines} %------------------------------------------------------------------------------- None of the UMFPACK routines discussed so far prints anything, even when an error occurs. UMFPACK provides you with nine routines for printing the input and output arguments (including the {\tt Control} settings and {\tt Info} statistics) of UMFPACK routines discussed above. They are fully described in Section~\ref{Report}: \begin{itemize} \item {\tt umfpack\_report\_status}, {\tt umfpack\_l\_report\_status}: Prints the status (return value) of other {\tt umfpack\_*} routines. \item {\tt umfpack\_report\_info}, {\tt umfpack\_l\_report\_info}: Prints the statistics returned in the {\tt Info} array by {\tt umfpack\_*symbolic}, {\tt umfpack\_numeric}, and {\tt umfpack\_*solve}. \item {\tt umfpack\_report\_control}, {\tt umfpack\_l\_report\_control}: Prints the {\tt Control} settings. \item {\tt umfpack\_report\_matrix}, {\tt umfpack\_l\_report\_matrix}: Verifies and prints a compressed column-form or compressed row-form sparse matrix. \item {\tt umfpack\_report\_triplet}, {\tt umfpack\_l\_report\_triplet}: Verifies and prints a matrix in triplet form. \item {\tt umfpack\_report\_symbolic}, {\tt umfpack\_l\_report\_symbolic}: Verifies and prints a {\tt Symbolic} object. \item {\tt umfpack\_report\_numeric}, {\tt umfpack\_l\_report\_numeric}: Verifies and prints a {\tt Numeric} object. \item {\tt umfpack\_report\_perm}, {\tt umfpack\_l\_report\_perm}: Verifies and prints a permutation vector. \item {\tt umfpack\_report\_vector}, {\tt umfpack\_l\_report\_vector}: Verifies and prints a real vector. \end{itemize} The {\tt umfpack\_report\_*} routines behave slightly differently when compiled into the C-callable UMFPACK library than when used in the MATLAB mexFunction. MATLAB stores its sparse matrices using the same compressed column data structure discussed above, where row and column indices are in the range {\tt 0} to {\tt n-1}, but it prints them as if they are in the range {\tt 1} to {\tt n}. The UMFPACK mexFunction behaves the same way. You can control how much the {\tt umfpack\_report\_*} routines print by modifying the {\tt Control [UMFPACK\_PRL]} parameter. Its default value is {\tt UMFPACK\_DEFAULT\_PRL} which is equal to 1. Here is a summary of how the routines use this print level parameter: \begin{itemize} \item {\tt umfpack\_report\_status}, {\tt umfpack\_l\_report\_status}: No output if the print level is 0 or less, even when an error occurs. If 1, then error messages are printed, and nothing is printed if the status is {\tt UMFPACK\_OK}. If 2 or more, then the status is always printed. If 4 or more, then the UMFPACK Copyright is printed. If 6 or more, then the UMFPACK License is printed. See also the first page of this User Guide for the Copyright and License. \item {\tt umfpack\_report\_control}, {\tt umfpack\_l\_report\_control}: No output if the print level is 1 or less. If 2 or more, all of {\tt Control} is printed. \item {\tt umfpack\_report\_info}, {\tt umfpack\_l\_report\_info}: No output if the print level is 1 or less. If 2 or more, all of {\tt Info} is printed. \item all other {\tt umfpack\_report\_*} routines: If the print level is 2 or less, then these routines return silently without checking their inputs. If 3 or more, the inputs are fully verified and a short status summary is printed. If 4, then the first few entries of the input arguments are printed. If 5, then all of the input arguments are printed. \end{itemize} %------------------------------------------------------------------------------- \subsection{Utility routines} %------------------------------------------------------------------------------- UMFPACK includes a routine that returns the time used by the process, {\tt umfpack\_timer}. The routine uses either {\tt getrusage} (which is preferred), or the ANSI C {\tt clock} routine if that is not available. It is fully described in Section~\ref{Utility}. It is the only routine that is identical in both {\tt int} and {\tt long} versions (there is no {\tt umfpack\_l\_timer} routine). %------------------------------------------------------------------------------- \subsection{Control parameters} %------------------------------------------------------------------------------- UMFPACK uses an optional {\tt double} array, {\tt Control}, to modify its control parameters. These may be modified by the user (see {\tt umfpack\_defaults} and {\tt umfpack\_l\_defaults}). Each user-callable routine includes a complete description of how each control setting modifies its behavior. Table~\ref{control} summarizes the entire contents of the {\tt Control} array. Future versions may make use of additional entries in the {\tt Control} array. Note that ANSI C uses 0-based indexing, while MATLAB user's 1-based indexing. Thus, {\tt Control(1)} in MATLAB is the same as {\tt Control[0]} or {\tt Control[UMFPACK\_PRL]} in ANSI C. \begin{table} \caption{UMFPACK Control parameters} \label{control} {\footnotesize \begin{tabular}{llll} \hline MATLAB & ANSI C & default & description \\ \hline \multicolumn{4}{l}{Used by reporting routines:} \\ {\tt Control(1)} & {\tt Control[UMFPACK\_PRL]} & 1 & printing level \\ \hline % /* used in UMFPACK_*symbolic only: */ \multicolumn{4}{l}{Used by {\tt umfpack\_*symbolic:}} \\ {\tt Control(2)} & {\tt Control[UMFPACK\_DENSE\_ROW]} & 0.2 & dense row threshold \\ {\tt Control(3)} & {\tt Control[UMFPACK\_DENSE\_COL]} & 0.2 & dense column threshold \\ \hline % /* used in UMFPACK_numeric only: */ \multicolumn{4}{l}{Used by {\tt umfpack\_*numeric:}} \\ {\tt Control(4)} & {\tt Control[UMFPACK\_PIVOT\_TOLERANCE]} & 0.1 & partial pivoting tolerance \\ {\tt Control(5)} & {\tt Control[UMFPACK\_BLOCK\_SIZE]} & 24 & BLAS block size \\ {\tt Control(6)} & {\tt Control[UMFPACK\_RELAXED\_AMALGAMATION]} & 0.25 & amalgamation \\ {\tt Control(7)} & {\tt Control[UMFPACK\_ALLOC\_INIT]} & 0.7 & initial memory allocation \\ {\tt Control(13)} & {\tt Control[UMFPACK\_PIVOT\_OPTION]} & 0 & symmetric pivot preference \\ {\tt Control(14)} & {\tt Control[UMFPACK\_RELAXED2\_AMALGAMATION]} & 0.1 & amalgamation \\ {\tt Control(15)} & {\tt Control[UMFPACK\_RELAXED3\_AMALGAMATION]} & 0.125 & amalgamation \\ \hline % /* used in UMFPACK_*solve only: */ \multicolumn{4}{l}{Used by {\tt umfpack\_*solve:}} \\ {\tt Control(8)} & {\tt Control[UMFPACK\_IRSTEP]} & 2 & max iter. refinement steps \\ \hline % % /* compile-time settings - Control [8..11] cannot be changed at run time: */ \multicolumn{4}{l}{Can only be changed at compile time:} \\ {\tt Control(9)} & {\tt Control[UMFPACK\_COMPILED\_WITH\_BLAS]} & - & true if BLAS is used \\ {\tt Control(10)} & {\tt Control[UMFPACK\_COMPILED\_FOR\_MATLAB]} & - & true for mexFunction \\ {\tt Control(11)} & {\tt Control[UMFPACK\_COMPILED\_WITH\_GETRUSAGE]} & - & true if {\tt getrusage} used \\ {\tt Control(12)} & {\tt Control[UMFPACK\_COMPILED\_IN\_DEBUG\_MODE]} & - & true if debug mode enabled \\ % %/* Control [15...19] unused */ \hline \end{tabular} } \end{table} %------------------------------------------------------------------------------- \subsection{Larger examples} %------------------------------------------------------------------------------- A full example of all user-callable UMFPACK routines (the {\tt int} routines) is available in the C main program, {\tt umfpack\_demo.c} listed in Section~\ref{Demo}. A nearly identical program that uses the {\tt long} integer version of UMFPACK is in {\tt umfpack\_l\_demo.c}. Another example is the UMFPACK mexFunction, {\tt umfpackmex.c}. The mexFunction accesses only the user-callable C interface to UMFPACK. The only features that it does not use are the support for the triplet form (MATLAB's sparse arrays are already in the compressed column form) and the ability to reuse the {\tt Symbolic} object to numerically factorize a matrix whose pattern is the same as a prior matrix analyzed by {\tt umfpack\_symbolic} or {\tt umfpack\_qsymbolic}. The latter is an important feature, but the mexFunction does not return its opaque {\tt Symbolic} and {\tt Numeric} objects to MATLAB. Instead, it gets the contents of these objects after extracting them via the {\tt umfpack\_get\_*} routines, and returns them as MATLAB sparse matrices. %------------------------------------------------------------------------------- \section{Synopsis of all C-callable routines ({\tt int} version)} %------------------------------------------------------------------------------- Each subsection, below, summarizes the input variables, output variables, return values, and calling sequences of the routines in one category. Variables with the same name as those already listed in a prior category have the same size and type. %------------------------------------------------------------------------------- \subsection{Primary routines} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} #include "umfpack.h" int status, n, nz, Ap [n+1], Ai [nz] ; double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO], Ax [nz], X [n], B [n] ; void *Symbolic, *Numeric ; char *sys ; status = umfpack_symbolic (n, Ap, Ai, &Symbolic, Control, Info) ; status = umfpack_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ; status = umfpack_solve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info) ; umfpack_free_symbolic (&Symbolic) ; umfpack_free_numeric (&Numeric) ; \end{verbatim} } %------------------------------------------------------------------------------- \subsection{Alternative routines} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} int Qinit [n], Wi [n] ; double W [n], Y [n], Z [n], S [n] ; umfpack_defaults (Control) ; status = umfpack_qsymbolic (n, Ap, Ai, Qinit, &Symbolic, Control, Info) ; status = umfpack_wsolve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info, Wi, W, Y, Z, S) ; \end{verbatim} } %------------------------------------------------------------------------------- \subsection{Matrix manipulation routines} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} int Ti [nz], Tj [nz], Bp [n+1], Bi [max(n,nz)], P [n], Q [n], Cp [n+1], Ci [nz] ; double Tx [nz], Cx [nz], Bx [nz] ; status = umfpack_col_to_triplet (n, Ap, Tj) ; status = umfpack_triplet_to_col (n, nz, Ti, Tj, Tx, Bp, Bi, Bx) ; status = umfpack_transpose (n, Ap, Ai, Ax, P, Q, Cp, Ci, Cx) ; \end{verbatim} } %------------------------------------------------------------------------------- \subsection{Getting the contents of opaque objects} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} int lnz, unz, Lp [n+1], Li [lnz], Up [n+1], Ui [unz] ; double Lx [lnz], Ux [unz] ; int nfr, nchains, nsparse_col, Qtree [n], Front_npivots [n], Front_parent [n], Chain_start [n], Chain_maxrows [n], Chain_maxcols [n] ; status = umfpack_get_lunz (&lnz, &unz, &n, Numeric) ; status = umfpack_get_numeric (Lp, Li, Lx, Up, Ui, Ux, P, Q, Numeric) ; status = umfpack_get_symbolic (&n, &nz, &nfr, &nchains, &nsparse_col, Qtree, Front_npivots, Front_parent, Chain_start, Chain_maxrows, Chain_maxcols, Symbolic) ; \end{verbatim} } Note: the {\tt nsparse\_col} argument is no longer relevant. It is always equal to {\tt n} in this version. %------------------------------------------------------------------------------- \subsection{Reporting routines} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} char *name, *form ; umfpack_report_status (Control, status) ; umfpack_report_control (Control) ; umfpack_report_info (Control, Info) ; status = umfpack_report_matrix (name, n, Ap, Ai, Ax, form, Control) ; status = umfpack_report_numeric (name, Numeric, Control) ; status = umfpack_report_perm (name, n, P, Control) ; status = umfpack_report_symbolic (name, Symbolic, Control) ; status = umfpack_report_triplet (name, n, nz, Ti, Tj, Tx, Control) ; status = umfpack_report_vector (name, n, X, Control) ; \end{verbatim} } %------------------------------------------------------------------------------- \section{Synopsis of all C-callable routines ({\tt long} version)} %------------------------------------------------------------------------------- Each subsection, below, summarizes the input variables, output variables, return values, and calling sequences of the routines in one category. Variables with the same name as those already listed in a prior category have the same size and type. Note that the include file, {\tt umfpack.h}, is the same for both {\tt int} and {\tt long} versions of UMFPACK. %------------------------------------------------------------------------------- \subsection{Primary routines} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} #include "umfpack.h" long status, n, nz, Ap [n+1], Ai [nz] ; double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO], Ax [nz], X [n], B [n] ; void *Symbolic, *Numeric ; char *sys ; status = umfpack_l_symbolic (n, Ap, Ai, &Symbolic, Control, Info) ; status = umfpack_l_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ; status = umfpack_l_solve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info) ; umfpack_l_free_symbolic (&Symbolic) ; umfpack_l_free_numeric (&Numeric) ; \end{verbatim} } %------------------------------------------------------------------------------- \subsection{Alternative routines} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} long Qinit [n], Wi [n] ; double W [n], Y [n], Z [n], S [n] ; umfpack_l_defaults (Control) ; status = umfpack_l_qsymbolic (n, Ap, Ai, Qinit, &Symbolic, Control, Info) ; status = umfpack_l_wsolve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info, Wi, W, Y, Z, S) ; \end{verbatim} } %------------------------------------------------------------------------------- \subsection{Matrix manipulation routines} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} long Ti [nz], Tj [nz], Bp [n+1], Bi [max(n,nz)], P [n], Q [n], Cp [n+1], Ci [nz] ; double Tx [nz], Cx [nz], Bx [nz] ; status = umfpack_l_col_to_triplet (n, Ap, Tj) ; status = umfpack_l_triplet_to_col (n, nz, Ti, Tj, Tx, Bp, Bi, Bx) ; status = umfpack_l_transpose (n, Ap, Ai, Ax, P, Q, Cp, Ci, Cx) ; \end{verbatim} } %------------------------------------------------------------------------------- \subsection{Getting the contents of opaque objects} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} long lnz, unz, Lp [n+1], Li [lnz], Up [n+1], Ui [unz] ; double Lx [lnz], Ux [unz] ; long nfr, nchains, nsparse_col, Qtree [n], Front_npivots [n], Front_parent [n], Chain_start [n], Chain_maxrows [n], Chain_maxcols [n] ; status = umfpack_l_get_lunz (&lnz, &unz, &n, Numeric) ; status = umfpack_l_get_numeric (Lp, Li, Lx, Up, Ui, Ux, P, Q, Numeric) ; status = umfpack_l_get_symbolic (&n, &nz, &nfr, &nchains, &nsparse_col, Qtree, Front_npivots, Front_parent, Chain_start, Chain_maxrows, Chain_maxcols, Symbolic) ; \end{verbatim} } Note: the {\tt nsparse\_col} argument is no longer relevant. It is always equal to {\tt n} in this version. %------------------------------------------------------------------------------- \subsection{Reporting routines} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} char *name, *form ; umfpack_l_report_status (Control, status) ; umfpack_l_report_control (Control) ; umfpack_l_report_info (Control, Info) ; status = umfpack_l_report_matrix (name, n, Ap, Ai, Ax, form, Control) ; status = umfpack_l_report_numeric (name, Numeric, Control) ; status = umfpack_l_report_perm (name, n, P, Control) ; status = umfpack_l_report_symbolic (name, Symbolic, Control) ; status = umfpack_l_report_triplet (name, n, nz, Ti, Tj, Tx, Control) ; status = umfpack_l_report_vector (name, n, X, Control) ; \end{verbatim} } %------------------------------------------------------------------------------- \section{Synopsis of utility routines} %------------------------------------------------------------------------------- This routine is the same in both {\tt int} and {\tt long} versions of UMFPACK. {\footnotesize \begin{verbatim} double t ; t = umfpack_timer ( ) ; \end{verbatim} } %------------------------------------------------------------------------------- \section{Installation} \label{Install} %------------------------------------------------------------------------------- UMFPACK comes with a {\tt Makefile} for compiling the C-callable {\tt umfpack.a} library and the {\tt umfpack} mexFunction on Unix. System-dependent configurations are controlled by the {\tt Makefile}, and defined in {\tt umf\_config.h} listed in Section~\ref{Config}. You should not have to modify {\tt umf\_config.h}. To compile {\tt umfpack.a} on most Unix systems, all you need to do is to type {\tt make}. This will use the generic configuration, in {\tt Make.generic}. The three demo programs will be executed, and the output of {\tt umfpack\_demo.c} and {\tt umfpack\_l\_demo.c} will be compared with {\tt umfpack\_demo.out} and {\tt umfpack\_l\_demo.out}. These two demo programs are identical, except that {\tt umfpack\_demo.c} uses the {\tt int} version, while {\tt umfpack\_l\_demo.out} uses the {\tt long} version of UMFPACK. Expect to see a few differences, such as residual norms, compile-time control settings, and perhaps memory usage differences. (The Compaq Alpha uses the LP64 model by default, so if you're using that computer compare your output with the 64-bit Solaris output in {\tt umfpack\_demo.out64} and {\tt umfpack\_l\_demo.out64}). The BLAS \cite{DaydeDuff99,ACM679a,ATLAS} will not be used, so the performance of UMFPACK will not be as high as possible. For better performance, edit the {\tt Makefile} and un-comment the {\tt include Make.*} statement that is specific to your computer. For example, {\footnotesize \begin{verbatim} # include Make.generic # include Make.linux # include Make.sgi include Make.solaris # include Make.alpha # include Make.rs6000 \end{verbatim} } will include the Solaris-specific configurations, which uses the Sun Performance Library BLAS ({\tt sunperf}), and compiler optimizations that are different than the generic settings. If you change the {\tt Makefile} or your system-specific {\tt Make.*} file, be sure to type {\tt make purge} before recompiling. Here are the various parameters that you can control in your {\tt Make.*} file; more details are in {\tt umf\_config.h} listed in Section~\ref{Config}: \begin{itemize} \item {\tt CC = } your C compiler, usually, {\tt cc}. If you don't modify this string at all in your {\tt Make.*}, then the {\tt make} program will use your default C compiler (if {\tt make} is installed properly). \item {\tt RANLIB = } your system's {\tt ranlib} program, if needed. \item {\tt CFLAGS = } optimization flags, such as {\tt -O}. \item {\tt CONFIG = } configuration settings. \item {\tt LIB = } your libraries, such as {\tt -lm} or {\tt -lblas}. \end{itemize} The {\tt CONFIG} string can include combinations of the following: \begin{itemize} \item {\tt -DNBLAS} if you do not have any BLAS at all. By default, {\tt umf\_config.h} assumes you have some version of the BLAS. The BLAS are de-selected in {\tt Make.generic} with the statement {\tt CONFIG = -DNBLAS}. \item {\tt -DNCBLAS} if you do not have the C-BLAS \cite{ATLAS}. The interface to the C-BLAS is identical on any system (Unix or Windows). By default, {\tt umf\_config.h} assumes you have the C-BLAS, except for Solaris (which has {\tt sunperf}) and MATLAB, which has its own BLAS for compiling the MATLAB mexFunction on any system. \item {\tt -DNSUNPERF} if you are on Solaris but do not have {\tt sunperf}. \item {\tt -DLONGBLAS} if your BLAS can take {\tt long} integer input arguments. If not defined, then the {\tt umfpack\_l\_*} version of UMFPACK that uses {\tt long} integers does not call the BLAS. \item {\tt -DGETRUSAGE} if you have the {\tt getrusage} function. This should exist on any UNIX system. \item Options for controlling how C calls the Fortran BLAS: {\tt -DBLAS\_BY\_VALUE}, \newline {\tt -DBLAS\_NO\_UNDERSCORE}, and {\tt -DBLAS\_CHAR\_ARG}. These are set automatically for Sun Solaris, SGI Irix, Red Hat Linux, Compaq Alpha, and AIX (the IBM RS 6000). \end{itemize} To compile the {\tt umfpack} mexFunction on Unix, type {\tt make umfpack}. The MATLAB {\tt mex} command will select the appropriate compiler and compiler flags for your system, and the BLAS internal to MATLAB will be used. The {\tt mexopts.sh} file in your UMFPACK directory has been modified from the MATLAB default; the unmodified version is in {\tt mexopts.sh.orig} for comparison. If you're running Windows, and all you want to do is use UMFPACK in MATLAB, then just type {\tt umfpack\_make} in MATLAB. MATLAB Version 6.0 or higher is required. You won't be able to use the BLAS when compiling with the {\tt lcc} compiler provided with MATLAB Version 6.0); you will get an error stating that {\tt \_dgemm} is undefined. There is no work-around for this problem. Either use a different C compiler, or don't use the BLAS. %------------------------------------------------------------------------------- \section{Future work} \label{Future} %------------------------------------------------------------------------------- Here are a few features that are not in UMFPACK Version 3.2, in no particular order. They may appear in a future release of UMFPACK. If you are interested, let me know and I could consider including them: \begin{enumerate} \item Future versions may have different default {\tt Control} parameters. \item a condition number estimator. You can write your own in MATLAB by making a copy of the built-in MATLAB {\tt condest.m} routine and replacing {\tt LU} with {\tt umfpack}. Be sure to do so only if your MATLAB license allows you to, and do not distribute the derivative MATLAB code without direct permission from The Mathworks, Inc. \item an estimate of the 1-norm of $\m{PAQ}-\m{LU}$. \newline See {\tt ftp://ftp.mathworks.com/pub/contrib/v4/linalg/normest1.m} for a similar algorithm that computes the 1-norm estimate of $\sigma \m{I}+\m{AA}\tr-\m{LL}\tr$. It can easily be modified to compute the 1-norm estimate of $\m{PAQ}-\m{LU}$. See also \cite{DavisHager99}. \item a complex version. \item when using iterative refinement, the residual $\m{Ax}-\m{b}$ could be returned by {\tt umfpack\_solve} ({\tt umfpack\_wsolve} already does so, but this is not documented). \item the solve routines could handle multiple right-hand sides, and sparse right-hand sides. \item an option to redirect the error and diagnostic output to something other than standard output. \item permutation to block-triangular-form \cite{Duff78a} for the C-callable interface. \item the symbolic and numeric factorization could handle singular matrices, just like MATLAB's {\tt LU}. \item the ability to use user-provided {\tt malloc}, {\tt free}, and {\tt realloc} memory allocation routines. Note that UMFPACK makes very few calls to these routines. \item the ability to use user-provided work arrays, so that {\tt malloc}, {\tt free}, and {\tt realloc} realloc are not called. The {\tt umfpack\_wsolve} routine is one example. \item future versions may return more statistics in the {\tt Info} array, and they may use more entries in the {\tt Control} array. \item use a method that takes time proportional to the number of nonzeros in $\m{A}$ to analyze $\m{A}$ when {\tt Qinit} is provided (or when {\tt Qinit} is not provided and {\tt umf\_colamd} ignores "dense" rows) \cite{GilbertNgPeyton94}. The current method in {\tt umf\_analyze.c} takes time proportional to the number of nonzeros in the upper bound of $\m{U}$. \item an option of extracting the diagonal of $\m{U}$ (or other subsets of $\m{L}$ and $\m{U}$) from the {\tt Numeric} object without having to extract the entire LU factorization. \item a Fortran interface (this would probably require modifying UMFPACK to use user-provided work arrays). \item a C++ interface. \item a parallel version using MPI. \end{enumerate} %------------------------------------------------------------------------------- \newpage \section{The primary UMFPACK routines} \label{Primary} %------------------------------------------------------------------------------- The include files are the same for both {\tt int} and {\tt long} versions of UMFPACK. The generic integer type is {\tt Int}, which is an {\tt int} or {\tt long}, depending on which version of UMFPACK you are using. \subsection{umfpack\_symbolic and umfpack\_l\_symbolic} {\footnotesize \begin{verbatim} include(umfpack_symbolic.h) \end{verbatim} } \newpage \subsection{umfpack\_numeric and umfpack\_l\_numeric} {\footnotesize \begin{verbatim} include(umfpack_numeric.h) \end{verbatim} } \newpage \subsection{umfpack\_solve and umfpack\_l\_solve} {\footnotesize \begin{verbatim} include(umfpack_solve.h) \end{verbatim} } \newpage \subsection{umfpack\_free\_symbolic and umfpack\_l\_free\_symbolic} {\footnotesize \begin{verbatim} include(umfpack_free_symbolic.h) \end{verbatim} } \newpage \subsection{umfpack\_free\_numeric and umfpack\_l\_free\_numeric} {\footnotesize \begin{verbatim} include(umfpack_free_numeric.h) \end{verbatim} } %------------------------------------------------------------------------------- \newpage \section{Alternatives routines} \label{Alternative} %------------------------------------------------------------------------------- \subsection{umfpack\_defaults and umfpack\_l\_defaults} \label{defaults} {\footnotesize \begin{verbatim} include(umfpack_defaults.h) \end{verbatim} } \newpage \subsection{umfpack\_qsymbolic and umfpack\_l\_qsymbolic} {\footnotesize \begin{verbatim} include(umfpack_qsymbolic.h) \end{verbatim} } \newpage \subsection{umfpack\_wsolve and umfpack\_l\_wsolve} {\footnotesize \begin{verbatim} include(umfpack_wsolve.h) \end{verbatim} } %------------------------------------------------------------------------------- \newpage \section{Matrix manipulation routines} \label{Manipulate} %------------------------------------------------------------------------------- \subsection{umfpack\_col\_to\_triplet and umfpack\_l\_col\_to\_triplet} {\footnotesize \begin{verbatim} include(umfpack_col_to_triplet.h) \end{verbatim} } \newpage \subsection{umfpack\_triplet\_to\_col and umfpack\_l\_triplet\_to\_col} {\footnotesize \begin{verbatim} include(umfpack_triplet_to_col.h) \end{verbatim} } \newpage \subsection{umfpack\_transpose and umfpack\_l\_transpose} {\footnotesize \begin{verbatim} include(umfpack_transpose.h) \end{verbatim} } %------------------------------------------------------------------------------- \newpage \section{Getting the contents of opaque objects} \label{Get} %------------------------------------------------------------------------------- \subsection{umfpack\_get\_lunz and umfpack\_l\_get\_lunz} {\footnotesize \begin{verbatim} include(umfpack_get_lunz.h) \end{verbatim} } \newpage \subsection{umfpack\_get\_numeric and umfpack\_l\_get\_numeric} {\footnotesize \begin{verbatim} include(umfpack_get_numeric.h) \end{verbatim} } \newpage \subsection{umfpack\_get\_symbolic and umfpack\_l\_get\_symbolic} {\footnotesize \begin{verbatim} include(umfpack_get_symbolic.h) \end{verbatim} } %------------------------------------------------------------------------------- \newpage \section{Reporting routines} \label{Report} %------------------------------------------------------------------------------- \subsection{umfpack\_report\_status and umfpack\_l\_report\_status} {\footnotesize \begin{verbatim} include(umfpack_report_status.h) \end{verbatim} } \newpage \subsection{umfpack\_report\_control and umfpack\_l\_report\_control} {\footnotesize \begin{verbatim} include(umfpack_report_control.h) \end{verbatim} } \newpage \subsection{umfpack\_report\_info and umfpack\_l\_report\_info} {\footnotesize \begin{verbatim} include(umfpack_report_info.h) \end{verbatim} } \newpage \subsection{umfpack\_report\_matrix and umfpack\_l\_report\_matrix} {\footnotesize \begin{verbatim} include(umfpack_report_matrix.h) \end{verbatim} } \newpage \subsection{umfpack\_report\_numeric and umfpack\_l\_report\_numeric} {\footnotesize \begin{verbatim} include(umfpack_report_numeric.h) \end{verbatim} } \newpage \subsection{umfpack\_report\_perm and umfpack\_l\_report\_perm} {\footnotesize \begin{verbatim} include(umfpack_report_perm.h) \end{verbatim} } \newpage \subsection{umfpack\_report\_symbolic and umfpack\_l\_report\_symbolic} {\footnotesize \begin{verbatim} include(umfpack_report_symbolic.h) \end{verbatim} } \newpage \subsection{umfpack\_report\_triplet and umfpack\_l\_report\_triplet} {\footnotesize \begin{verbatim} include(umfpack_report_triplet.h) \end{verbatim} } \newpage \subsection{umfpack\_report\_vector and umfpack\_l\_report\_vector} {\footnotesize \begin{verbatim} include(umfpack_report_vector.h) \end{verbatim} } %------------------------------------------------------------------------------- \newpage \section{Utility routines} \label{Utility} %------------------------------------------------------------------------------- \subsection{umfpack\_timer} {\footnotesize \begin{verbatim} include(umfpack_timer.h) \end{verbatim} } %------------------------------------------------------------------------------- \newpage \section{{\tt umfpack.h} include file} \label{Include} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} include(umfpack.h) \end{verbatim} } %------------------------------------------------------------------------------- \newpage \section{Demo C main program, {\tt umfpack\_demo.c}} \label{Demo} %------------------------------------------------------------------------------- The {\tt umfpack\_l\_demo.c} is identical except that all the routine names are changed, and {\tt long}'s are used instead of {\tt int}'s. {\footnotesize \begin{verbatim} include(umfpack_demo.c) \end{verbatim} } %------------------------------------------------------------------------------- \newpage \section{Configuration file {\tt umf\_config.h}} \label{Config} %------------------------------------------------------------------------------- {\footnotesize \begin{verbatim} include(umf_config.h) \end{verbatim} } %------------------------------------------------------------------------------- \newpage % References %------------------------------------------------------------------------------- \bibliographystyle{plain} \bibliography{UserGuide} \end{document}