octave-nkf: liboctave/UMFPACK/UMFPACK/Doc/UserGuide.stex comparison

comparison liboctave/UMFPACK/UMFPACK/Doc/UserGuide.stex @ 5164:57077d0ddc8e

[project @ 2005-02-25 19:55:24 by jwe]

author	jwe
date	Fri, 25 Feb 2005 19:55:28 +0000
parents
children

comparison

equal deleted inserted replaced

-:9f3299378193
+:57077d0ddc8e
+%-------------------------------------------------------------------------------
+% The UserGuide.stex file.  Processed into UserGuide.tex via sed.
+%-------------------------------------------------------------------------------
+\documentclass[11pt]{article}
+\newcommand{\m}[1]{{\bf{#1}}}       % for matrices and vectors
+\newcommand{\tr}{^{\sf T}}          % transpose
+\newcommand{\he}{^{\sf H}}          % complex conjugate transpose
+\newcommand{\implies}{\rightarrow}
+\topmargin 0in
+\textheight 9in
+\oddsidemargin 0pt
+\evensidemargin 0pt
+\textwidth 6.5in
+\begin{document}
+\author{Timothy A. Davis \\
+Dept. of Computer and Information Science and Engineering \\
+Univ. of Florida, Gainesville, FL}
+\title{UMFPACK Version 4.4 User Guide}
+\date{Jan. 28, 2005}
+\maketitle
+%-------------------------------------------------------------------------------
+\begin{abstract}
+UMFPACK is a set of routines for solving unsymmetric sparse linear
+systems, $\m{Ax}=\m{b}$, using the Unsymmetric MultiFrontal method
+and direct sparse LU factorization.  It is written in ANSI/ISO C, with a
+MATLAB interface.  UMFPACK relies on the Level-3 Basic
+Linear Algebra Subprograms (dense matrix multiply) for its performance.
+This code works on Windows and many versions of Unix (Sun Solaris,
+Red Hat Linux, IBM AIX, SGI IRIX, and Compaq Alpha).
+\end{abstract}
+%-------------------------------------------------------------------------------
+Technical Report TR-04-003 (revised)
+UMFPACK Version 4.4 (Jan. 28, 2005), Copyright\copyright 2005 by Timothy A.
+Davis.  All Rights Reserved.
+{\bf UMFPACK License:}
+Your use or distribution of UMFPACK or any modified version of
+UMFPACK 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 UMFPACK or any modified version of UMFPACK code must cite the
+Copyright, this License, the Availability note, and ``Used by permission.''
+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 Availability:}
+http://www.cise.ufl.edu/research/sparse/umfpack
+{\bf Acknowledgments:}
+This work was supported by the National Science Foundation, under
+grants DMS-9504974, DMS-9803599, and CCR-0203270.
+The upgrade to Version 4.1 and the inclusion of the
+symmetric and 2-by-2 pivoting strategies
+were done while the author was on sabbatical at
+Stanford University and Lawrence Berkeley National Laboratory.
+%-------------------------------------------------------------------------------
+\newpage
+%-------------------------------------------------------------------------------
+\tableofcontents
+%-------------------------------------------------------------------------------
+\newpage
+\section{Overview}
+%-------------------------------------------------------------------------------
+UMFPACK\footnote{Pronounced with two syllables: umph-pack}
+Version 4.4 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-pattern MultiFrontal method \cite{DavisDuff97,DavisDuff99}.
+UMFPACK factorizes
+$\m{PAQ}$, $\m{PRAQ}$, or $\m{PR}^{-1}\m{AQ}$ into the product $\m{LU}$,
+where $\m{L}$ and $\m{U}$
+are lower and upper triangular, respectively, $\m{P}$ and $\m{Q}$ are
+permutation matrices, and $\m{R}$ is a diagonal matrix of row scaling factors
+(or $\m{R}=\m{I}$ if row-scaling is not used).  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).
+The sparse matrix $\m{A}$ can be square or rectangular, singular
+or non-singular, and real or complex (or any combination).  Only square
+matrices $\m{A}$ can be used to solve $\m{Ax}=\m{b}$ or related systems.
+Rectangular matrices can only be factorized.
+UMFPACK first finds a column pre-ordering that reduces fill-in, without regard
+to numerical values.  It scales and analyzes the matrix, and then automatically
+selects one of three strategies for pre-ordering the rows and columns:
+{\em unsymmetric},
+{\em 2-by-2}, and
+{\em symmetric}.  These strategies are described below.
+First, all pivots with zero Markowitz cost are eliminated and placed in the
+LU factors.  The remaining submatrix $\m{S}$ is then analyzed.
+The following rules are applied, and the first one that matches defines
+the strategy.
+\begin{itemize}
+\item Rule 1: $\m{A}$ rectangular $\implies$ unsymmetric.
+\item Rule 2:
+If the zero-Markowitz elimination results in a rectangular $\m{S}$,
+or an $\m{S}$ whose diagonal has not been preserved, the
+unsymmetric strategy is used.
+\item The symmetry $\sigma_1$ of $\m{S}$ is computed.  It is defined as
+the number of {\em matched} off-diagonal entries, divided by the
+total number of off-diagonal entries.  An entry $s_{ij}$ is matched
+if $s_{ji}$ is also an entry.  They need not be numerically equal.
+An {\em entry} is a value in $\m{A}$ which is present
+in the input data structure.  All nonzeros are entries, but some entries
+may be numerically zero.
+Rule 3: $\sigma_1 < 0.1 \implies$ unsymmetric.
+The matrix is very unsymmetric.
+\item Let $d$ be the number of nonzero entries on the diagonal of $\m{S}$.
+Let $\m{S}$ be $\nu$-by-$\nu$.
+Rule 4: $(\sigma_1 \ge 0.7) \:\wedge\: (d = \nu) \implies$ symmetric.
+The matrix has a nearly symmetric nonzero pattern, and a zero-free
+diagonal.
+\end{itemize}
+If the strategy has not yet been determined,
+the 2-by-2 strategy is attempted.  A row permutation $\m{P}_2$
+is found which attempts to reduce the number of small
+diagonal entries of $\m{P}_2 \m{S}$.
+An entry $s_{ij}$ is determined to be small if
+$|s_{ij}| < 0.01 \max |s_{*j}|$, or large otherwise.
+If $s_{ii}$ is numerically small, the method attempts to swap
+two rows $i$ and $j$, such that both $s_{ij}$ and $s_{ji}$ are large.
+Once these rows are swapped,
+they remain in place.  Let $\sigma_2$ be the symmetry of $\m{P}_2 \m{S}$,
+and let $d_2$ be the number of nonzero entries (either small or large)
+on the diagonal of $\m{P}_2 \m{S}$.
+\begin{itemize}
+\item Rule 5:
+($\sigma_2 > 1.1 \sigma_1) \:\wedge\: (d_2 > 0.9 \nu) \implies$ 2-by-2.
+The 2-by-2 permutation has made the matrix significantly more symmetric.
+\item Rule 6: $\sigma_2 < 0.7 \sigma_1 \implies$ unsymmetric.
+The 2-by-2 strategy has significantly deteriorated the symmetry,
+\item Rule 7: $\sigma_2 < 0.25 \implies$ unsymmetric.
+The matrix is still very unsymmetric.
+\item Rule 8: $\sigma_2 \ge 0.51 \implies$ 2-by-2.
+The matrix is roughly symmetric.
+\item Rule 9: $\sigma_2 \ge 0.999 \sigma_1 \implies$ 2-by-2.
+The 2-by-2 permutation has preserved symmetry, or made it only
+slightly worse.
+\item Rule 10: if no rule has yet triggered, use the unsymmetric strategy.
+\end{itemize}
+Each strategy is described below:
+\begin{itemize}
+\item {\em unsymmetric}:
+The column pre-ordering of $\m{S}$ is computed by a modified version of COLAMD
+\cite{DavisGilbertLarimoreNg00_algo,DavisGilbertLarimoreNg00,Larimore98}.
+The method finds a symmetric permutation $\m{Q}$ of the matrix $\m{S}\tr\m{S}$
+(without forming $\m{S}\tr\m{S}$ explicitly).  This is a good choice for
+$\m{Q}$, since the Cholesky factors of $\m{(SQ)\tr(SQ)}$ are an upper bound (in
+terms of nonzero pattern) of the factor $\m{U}$ for the unsymmetric LU
+factorization ($\m{PSQ}=\m{LU}$) regardless of the choice of $\m{P}$
+\cite{GeorgeNg85,GeorgeNg87,GilbertNg93}.  This modified version of
+COLAMD also computes the column elimination tree and post-orders the
+tree.  It finds the upper bound on the number of nonzeros in L and U.
+It also has a different threshold for determining dense rows and columns.
+During factorization, the column pre-ordering can be modified.
+Columns within a single super-column 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 $a_{ij}$ can be chosen if
+$|a_{ij}| \ge 0.1 \max |a_{*j}|$.  Among those numerically acceptable
+entries, the sparsest row $i$ is chosen as the pivot row.
+\item {\em 2-by-2}:
+The symmetric strategy (see below) is applied to the matrix $\m{P}_2 \m{S}$,
+rather than $\m{S}$.
+\item {\em symmetric}:
+The column ordering is computed from AMD
+\cite{AmestoyDavisDuff96,AmestoyDavisDuff03},
+applied to the pattern of $\m{S}+\m{S}\tr$
+followed by a post-ordering of the supernodal elimination
+tree of $\m{S}+\m{S}\tr$.
+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
+$a_{jj} \ge 0.001 \max |a_{*j}|$.  Otherwise, a sparse row is selected,
+using the same method used by the unsymmetric strategy.
+\end{itemize}
+The symmetric and 2-by-2 strategies, and their automatic selection,
+are new to Version 4.1.  Version 4.0 only used the unsymmetric strategy.
+Once the strategy is selected,
+the factorization of the matrix $\m{A}$ is broken down into the factorization
+of 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.  For the
+nonsingular matrices factorized with the unsymmetric strategy, there are
+exactly the same number of chains as there are leaves 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 sub-matrices (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 which 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 numeric 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).
+When the symmetric or 2-by-2 strategies are used,
+the column preordering is not refined during numeric factorization.
+Row pivoting for sparsity and numerical accuracy is performed if the
+diagonal entry is too small.
+More details of the method, including experimental results, are
+described in \cite{Davis03,Davis03_algo}, available at
+http://www.cise.ufl.edu/tech-reports.
+%-------------------------------------------------------------------------------
+\section{Availability}
+%-------------------------------------------------------------------------------
+In addition to appearing as a Collected Algorithm of the ACM,
+UMFPACK Version 4.4 is available at http://www.cise.ufl.edu/research/sparse.
+Version 4.3 is included as a built-in routine in MATLAB
+7.1.  Version 4.0 (in MATLAB 6.5)
+does not have the symmetric or 2-by-2 strategies and it takes
+less advantage of the level-3
+BLAS \cite{DaydeDuff99,ACM679a,ATLAS,GotoVandeGeijn02}.
+Version 4.4 (and v4.3 through v4.1) tend to be much faster than Version 4.0,
+particularly on unsymmetric matrices with mostly symmetric
+nonzero pattern (such as finite element and circuit simulation matrices).
+Version 3.0 and following make
+use of a modified version of COLAMD V2.0 by Timothy A.~Davis, Stefan
+Larimore, John Gilbert, and Esmond Ng.  The original COLAMD V2.1 is available in
+as a built-in routine in MATLAB V6.0 (or later), and at
+http://www.cise.ufl.edu/research/sparse.
+These codes are also available in Netlib \cite{netlib} at
+http://www.netlib.org.
+UMFPACK Versions 2.2.1 and earlier, co-authored with Iain Duff,
+are available at http://www.cise.ufl.edu/research/sparse and as
+MA38 (functionally equivalent to Version 2.2.1) in the Harwell
+Subroutine Library.
+%-------------------------------------------------------------------------------
+\section{Primary changes from prior versions}
+%-------------------------------------------------------------------------------
+A detailed list of changes is in the {\tt ChangeLog} file.
+%-------------------------------------------------------------------------------
+\subsection{Version 4.4}
+%-------------------------------------------------------------------------------
+Bug fix in strategy selection in {\tt umfpack\_*\_qsymbolic}.
+Added packed complex case for all complex input/output arguments.
+Added {\tt umfpack\_get\_determinant}.
+Added minimal support for Microsoft Visual Studio
+(the {\tt umf\_multicompile.c} file).
+%-------------------------------------------------------------------------------
+\subsection{Version 4.3.1}
+%-------------------------------------------------------------------------------
+Minor bug fix in the forward/backsolve.  This bug had the effect of turning
+off iterative refinement when solving $\m{A}\tr\m{x}=\m{b}$ after factorizing
+$\m{A}$.  UMFPACK mexFunction now factorizes $\m{A}\tr$ in its forward-slash
+operation.
+%-------------------------------------------------------------------------------
+\subsection{Version 4.3}
+%-------------------------------------------------------------------------------
+No changes are visible to the C or MATLAB user, except the presence of
+one new control parameter in the {\tt Control} array,
+and three new statistics in the {\tt Info} array.
+The primary change is the addition of an (optional) drop tolerance.
+%-------------------------------------------------------------------------------
+\subsection{Version 4.1}
+%-------------------------------------------------------------------------------
+The following is a summary of the main changes that are visible to the C
+or MATLAB user:
+\begin{enumerate}
+\item New ordering strategies added.  No changes are required in user code
+(either C or MATLAB) to use the new default strategy, which is an automatic
+selection of the unsymmetric, symmetric, or 2-by-2 strategies.
+\item Row scaling added.  This is only visible to the MATLAB caller when using
+the form {\tt [L,U,P,Q,R] = umfpack (A)}, to retrieve the LU factors.
+Likewise, it is only visible to the C caller when the LU factors are
+retrieved, or when solving systems with just $\m{L}$ or $\m{U}$.
+New C-callable and MATLAB-callable routines are included to get and to
+apply the scale factors computed by UMFPACK.  Row scaling is enabled by
+default, but can be disabled.  Row scaling usually leads to a better
+factorization, particularly when the symmetric strategy is used.
+\item Error code {\tt UMFPACK\_ERROR\_problem\_to\_large} removed.
+Version 4.0 would generate this error when the upper bound memory usage
+exceeded 2GB (for the {\tt int} version), even when the actual memory
+usage was less than this.  The new version properly handles this case,
+and can successfully factorize the matrix if sufficient memory is
+available.
+\item New control parameters and statistics provided.
+\item The AMD symmetric approximate minimum degree ordering routine added
+\cite{AmestoyDavisDuff96,AmestoyDavisDuff03}.
+It is used by UMFPACK, and can also be called independently from C or
+MATLAB.
+\item The {\tt umfpack} mexFunction now returns permutation matrices, not
+permutation vectors, when using the form {\tt [L,U,P,Q] = umfpack (A)}
+or the new form {\tt [L,U,P,Q,R] = umfpack (A)}.
+\item New arguments added to the user-callable routines
+{\tt umfpack\_*\_symbolic},
+{\tt umfpack\_*\_qsymbolic},
+{\tt umfpack\_*\_get\_numeric}, and
+{\tt umfpack\_*\_get\_symbolic}.
+The symbolic analysis now makes use of the numerical values of the matrix
+$\m{A}$, to guide the 2-by-2 strategy.  The subsequent matrix passed to
+the numeric factorization step does not have to have the same numerical
+values.  All of the new arguments are optional.  If you do not wish to
+include them, simply pass {\tt NULL} pointers instead.  The 2-by-2 strategy
+will assume all entries are numerically large, for example.
+\item New routines added to save and load the {\tt Numeric} and {\tt Symbolic}
+objects to and from a binary file.
+\item A Fortran interface added.  It provides access to a subset of
+UMFPACK's features.
+\item You can compute an incomplete LU factorization, by dropping small
+entries from $\m{L}$ and $\m{U}$.  By default, no nonzero entry is
+dropped, no matter how small in absolute value.  This feature is new
+to Version 4.3.
+\end{enumerate}
+%-------------------------------------------------------------------------------
+\section{Using UMFPACK in MATLAB}
+%-------------------------------------------------------------------------------
+The easiest way to use UMFPACK is within MATLAB.  Version 4.3 is a built-in
+routine in MATLAB 7.1, and is used in {\tt x = A}$\backslash${\tt b} when
+{\tt A} is sparse, square, unsymmetric (or symmetric but not positive definite),
+and with nonzero entries that are not confined in a narrow band.
+It is also used for the {\tt [L,U,P,Q] = lu (A)} usage of {\tt lu}.
+Type {\tt help lu} in MATLAB 6.5 or later for more details.
+To use the UMFPACK mexFunction, you must download and compile it,
+since the mexFunction itself is not part of MATLAB.
+The following discussion assumes that
+you have MATLAB Version 6.0 or later (which includes the BLAS, and the
+{\tt colamd} ordering routine).  To compile both the UMFPACK and AMD
+mexFunctions, just type {\tt make} in the Unix system shell,
+while in the {\tt UMFPACK} directory.
+You can also type {\tt umfpack\_make} in MATLAB, if you are in the
+{\tt UMFPACK/MATLAB} directory, or if that directory is in your MATLAB path.
+This works on any system with MATLAB, 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 the UMFPACK mexFunction within MATLAB.
+An optional input argument can be used to modify the control parameters for
+UMFPACK, and an optional output argument provides statistics on the
+factorization.
+Refer to the AMD User Guide for more details about the AMD mexFunction.
+\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 row and column pre-ordering
+(symmetric ordering).
+\end{minipage}
+&
+\begin{minipage}[t]{2.2in}
+\begin{verbatim}
+S = spones (A) ;
+Q = symamd (S+S') ;
+Control = umfpack ;
+Control (6) = 3 ;
+x = umfpack (A,Q,'\',b,Control) ;
+\end{verbatim}
+\end{minipage}
+&
+\begin{minipage}[t]{2.2in}
+\begin{verbatim}
+spparms ('autommd',0) ;
+S = spones (A) ;
+Q = symamd (S+S') ;
+x = A (Q,Q) \ b (Q) ;
+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}
+Note: $\m{A}$ is factorized.
+\end{minipage}
+&
+\begin{minipage}[t]{2.2in}
+\begin{verbatim}
+x = b / A ;
+\end{verbatim}
+Note: $\m{A}\tr$ is factorized.
+\end{minipage}
+\\
+& & \\
+\hline
+& & \\
+\begin{minipage}[t]{1.5in}
+Scale and factorize $\m{A}$, then solve $\m{Ax}=\m{b}$.
+\end{minipage}
+&
+\begin{minipage}[t]{2.2in}
+\begin{verbatim}
+[L,U,P,Q,R] = umfpack (A) ;
+c = P * (R \ b) ;
+x = Q * (U \ (L \ c)) ;
+\end{verbatim}
+\end{minipage}
+&
+\begin{minipage}[t]{2.2in}
+\begin{verbatim}
+[m n] = size (A) ;
+r = full (sum (abs (A), 2)) ;
+r (find (r == 0)) = 1 ;
+R = spdiags (r, 0, m, m) ;
+I = speye (n) ;
+Q = I (:, colamd (A)) ;
+[L,U,P] = lu ((R\A)*Q) ;
+c = P * (R \ b) ;
+x = Q * (U \ (L \ c)) ;
+\end{verbatim}
+\end{minipage}
+\\
+& & \\
+\hline
+\end{tabular}
+}
+\end{table}
+Note: in MATLAB 6.5 or later, use {\tt spparms ('autoamd',0)} in addition to
+{\tt spparms ('autommd',0)}, in Table~\ref{matlab}, to turn off MATLAB's
+default reordering.
+UMFPACK requires
+{\tt b} to be a dense vector (real or complex) of the appropriate dimension.
+This is more restrictive than what you can do with MATLAB's
+backslash or forward slash.  See {\tt umfpack\_solve} for an M-file that
+removes this restriction.
+This restriction does not apply to the built-in backslash operator
+in MATLAB 6.5 or later, which uses UMFPACK to factorize the matrix.
+You can do this yourself in MATLAB:
+{\footnotesize
+\begin{verbatim}
+[L,U,P,Q,R] = umfpack (A) ;
+x = Q * (U \ (L \ (P * (R \ b)))) ;
+\end{verbatim}
+}
+or, with no row scaling:
+{\footnotesize
+\begin{verbatim}
+[L,U,P,Q] = umfpack (A) ;
+x = Q * (U \ (L \ (P * b))) ;
+\end{verbatim}
+}
+The above examples do not make use of the iterative refinement
+that is built into
+{\tt x = }{\tt umfpack (A,'}$\backslash${\tt ',b)}
+however.
+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.  By default, it scales
+the rows of {\tt A} and reorders the columns of {\tt A} prior to
+factorization, so that {\tt L*U} is equal to {\tt P*(R}$\backslash${\tt A)*Q},
+where {\tt R} is a diagonal sparse matrix of scale factors for the rows
+of {\tt A}.  The scale factors {\tt R} are applied to {\tt A} via the MATLAB
+expression {\tt R}$\backslash${\tt A} to avoid multiplying by
+the reciprocal, which can be numerically inaccurate.
+There are more options; you can provide your own column pre-ordering (in which
+case UMFPACK does not call COLAMD or AMD), you can modify other control settings
+(similar to the {\tt spparms} in MATLAB), and you can get various statistics on
+the analysis, factorization, and solution of the linear system.  Type
+{\tt umfpack\_details} and {\tt 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.
+The output of these two programs should be about the same
+as the files {\tt umfpack\_simple.m.out} and {\tt umfpack\_demo.m.out}
+that are provided.
+Factorizing {\tt A'} (or {\tt A.'}) and using the transposed factors can
+sometimes be faster than factorizing {\tt A}.  It can also be preferable to
+factorize {\tt A'} if {\tt A} is rectangular.  UMFPACK pre-orders the columns
+to maintain sparsity; the row ordering is not determined until the matrix
+is factorized.  Thus, if {\tt A} is {\tt m} by {\tt n} with rank {\tt m}
+and {\tt m} $<$ {\tt n}, then {\tt umfpack} might not find a factor
+{\tt U} with a zero-free diagonal.  Unless the matrix ill-conditioned or
+poorly scaled, factorizing {\tt A'} in this case will guarantee that both
+factors will have zero-free diagonals.  Here's how you can factorize {\tt A'}
+and get the factors of {\tt A} instead:
+\begin{verbatim}
+[l,u,p,q] = umfpack (A') ;
+L = u' ;
+U = l' ;
+P = q ;
+Q = p ;
+clear l u p q
+\end{verbatim}
+This is an alternative to {\tt [L,U,P,Q]=umfpack(A)}.
+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, and then
+solves each block.  The LU factors are not returned.  Its usage is simple:
+{\tt x = umfpack\_btf(A,b)}.  Type {\tt help umfpack\_btf} for more options.
+An estimate of the 1-norm of {\tt L*U-P*A*Q} can be computed in MATLAB
+as {\tt lu\_normest(P*A*Q,L,U)}, using the {\tt lu\_normest.m} M-file
+by Hager and Davis \cite{DavisHager99} that is included with the
+UMFPACK distribution.  With row scaling enabled, use
+{\tt lu\_normest(P*(R}$\backslash${\tt A)*Q,L,U)} instead.
+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.
+This is 70\% of the upper bound, by default, for the unsymmetric strategy.
+For the symmetric strategy, the fraction of the upper bound is computed
+automatically (assuming a best-case scenario with no numerical pivoting
+required during numeric factorization).
+If this initial allocation fails, it reduces its request
+and uses less memory.   If the space is not large enough during factorization,
+it is increased via {\tt mxRealloc}.
+However, {\tt mxMalloc} and {\tt mxRealloc} abort the {\tt umfpack} mexFunction
+if they fail, so this strategy does not work in MATLAB.  The strategy works fine
+when {\tt malloc} or the internal memory allocator {\tt utMalloc} are used
+instead, since those routines return {\tt NULL} on failure, and do not terminate
+the mexFunction.  The {\tt umfpack} mexFunction can be compiled to use
+{\tt utMalloc}, but this is an internal undocumented utility routine in MATLAB,
+and thus using {\tt utMalloc} might not always be successful.
+To use the documented {\tt mxMalloc} routine instead, compile the
+mexFunction with the {\tt -DNUTIL} flag enabled.
+To compute the determinant with UMFPACK:
+\begin{verbatim}
+d = umfpack (A, 'det') ;
+[d e] = umfpack (A, 'det') ;
+\end{verbatim}
+The first case is identical to MATLAB's {\tt det}.
+The second case returns the determinant in the form
+$d \times 10^e$, which avoids overflow if $e$ is large.
+%-------------------------------------------------------------------------------
+\section{Using UMFPACK in a C program}
+\label{C}
+%-------------------------------------------------------------------------------
+The C-callable UMFPACK library consists of 32 user-callable routines and one
+include file.  All but three of the routines come in four versions, with
+different sizes of integers and for real or complex floating-point numbers:
+\begin{enumerate}
+\item {\tt umfpack\_di\_*}: real double precision, {\tt int} integers.
+\item {\tt umfpack\_dl\_*}: real double precision, {\tt long} integers.
+\item {\tt umfpack\_zi\_*}: complex double precision, {\tt int} integers.
+\item {\tt umfpack\_zl\_*}: complex double precision, {\tt long} integers.
+\end{enumerate}
+where {\tt *} denotes the specific name of one of the routines.
+Routine names beginning with {\tt umf\_} are internal to the package,
+and should not be called by the user.  The include file {\tt umfpack.h}
+must be included in any C program that uses UMFPACK.
+The other three routines are the same for all four versions.
+In addition, the C-callable AMD library distributed with UMFPACK
+includes 4 user-callable routines (in two versions with {\tt int} and
+{\tt long} integers) and one include file.  Refer to the AMD documentation
+for more details.
+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.
+The notation {\tt umfpack\_di\_*} refers to all user-callable routines
+for the real double precision and {\tt int} integer case.  The notation
+{\tt umfpack\_*\_numeric}, for example, refers all four versions
+(real/complex, int/long) of a single operation
+(in this case numeric factorization).
+%-------------------------------------------------------------------------------
+\subsection{The size of an integer}
+%-------------------------------------------------------------------------------
+The {\tt umfpack\_di\_*} and {\tt umfpack\_zi\_*} routines use {\tt int} integer
+arguments; those starting with {\tt umfpack\_dl\_} or {\tt umfpack\_zl\_}
+use {\tt long} integer arguments.  If you compile UMFPACK in the standard
+ILP32 mode (32-bit {\tt int}'s, {\tt long}'s, and pointers) then the versions
+are essentially identical.  You will be able to solve problems using up to 2GB
+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\_dl\_*}
+and {\tt umfpack\_zl\_*} routines can solve huge
+problems (not limited to 2GB), limited of course by the amount of available
+memory.  The only drawback to the 64-bit mode is that not all BLAS libraries
+support 64-bit integers.  This limits the performance you will obtain.
+Those that do support 64-bit integers are specific to particular
+architectures, and are not portable.  UMFPACK and AMD should be compiled
+in the same mode.
+If you compile UMFPACK and AMD in the LP64 mode,
+be sure to add {\tt -DLP64} to the compilation command.  See the examples in
+{\tt Make.alpha}, {\tt Make.sgi}, and {\tt Make.solaris}.
+%-------------------------------------------------------------------------------
+\subsection{Real and complex floating-point}
+%-------------------------------------------------------------------------------
+The {\tt umfpack\_di\_*} and {\tt umfpack\_dl\_*} routines take (real) double
+precision arguments, and return double precision arguments.  In the
+{\tt umfpack\_zi\_*} and {\tt umfpack\_zl\_*} routines, these same arguments
+hold the real part of the matrices; and second double precision arrays hold
+the imaginary part of the input and output matrices.  Internally, complex
+numbers are stored in arrays with their real and imaginary parts interleaved,
+as required by the BLAS (``packed'' complex form).
+New to Version 4.4 is the option of providing input/output arguments
+in packed complex form.
+%-------------------------------------------------------------------------------
+\subsection{Primary routines, and a simple example}
+%-------------------------------------------------------------------------------
+Five primary UMFPACK routines are required to factorize $\m{A}$ or
+solve $\m{Ax}=\m{b}$.  They are fully described in Section~\ref{Primary}:
+\begin{itemize}
+\item {\tt umfpack\_*\_symbolic}:
+Pre-orders the columns of $\m{A}$ to reduce fill-in.
+Returns an opaque {\tt Symbolic} object as a {\tt void *}
+pointer.  The object contains the symbolic analysis and is needed for the
+numeric 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. It is no larger than about $13n$ integers if $\m{A}$ is
+$n$-by-$n$.
+\item {\tt umfpack\_*\_numeric}:
+Numerically scales and then factorizes a sparse matrix into
+$\m{PAQ}$, $\m{PRAQ}$, or $\m{PR}^{-1}\m{AQ}$ into the product $\m{LU}$,
+where
+$\m{P}$ and $\m{Q}$ are permutation matrices, $\m{R}$ is a diagonal
+matrix of scale factors, $\m{L}$ is lower triangular with unit diagonal,
+and $\m{U}$ is upper triangular.  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 or Fortran program, but not in MATLAB).
+The matrix
+$\m{U}$ will have zeros on the diagonal if $\m{A}$ is singular; this
+produces a warning, but the factorization is still valid.
+\item {\tt umfpack\_*\_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.  The matrix $\m{A}$ must
+be square.  If it is singular, then a divide-by-zero will occur, and your
+solution with contain IEEE Inf's or NaN's in the appropriate places.
+\item {\tt umfpack\_*\_free\_symbolic}:
+Frees the {\tt Symbolic} object created by {\tt umfpack\_*\_symbolic}
+or {\tt umfpack\_*\_qsymbolic}.
+\item {\tt umfpack\_*\_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 or four arrays, where the matrix is {\tt m}-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}
+}
+The complex versions add another array for the imaginary part:
+{\footnotesize
+\begin{verbatim}
+double Az [nz] ;
+\end{verbatim}
+}
+Alternatively, if {\tt Az} is {\tt NULL},
+the real part of the $k$th entry is located in
+{\tt Ax[2*k]} and the imaginary part is located in
+{\tt Ax[2*k+1]}, and the {\tt Ax} array is of size {\tt 2*nz}.
+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]} \ldots {\tt Ap[j+1]-1]}.
+The corresponding numerical values are stored in
+{\tt Ax[Ap[j]} \ldots {\tt Ap[j+1]-1]}.
+The imaginary part, for the complex versions, is stored in
+{\tt Az[Ap[j]} \ldots {\tt Ap[j+1]-1]}
+(see above for the packed complex case).
+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}$.
+For a more flexible method for providing an input matrix to UMFPACK,
+see Section~\ref{triplet}.
+Here is a simple main program, {\tt umfpack\_simple.c}, that illustrates the
+basic usage of UMFPACK.  See Section~\ref{Synopsis} for a short description
+of each calling sequence, including a list of options for the first
+argument of {\tt umfpack\_di\_solve}.
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_simple.c via sed
+\end{verbatim}
+}
+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 to $\m{Ax}=\m{b}$ 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}).  It also ignores the status value returned by
+most user-callable UMFPACK routines.
+%-------------------------------------------------------------------------------
+\subsection{A note about zero-sized arrays}
+%-------------------------------------------------------------------------------
+UMFPACK uses many user-provided arrays of
+size {\tt m} or {\tt n} (the order of the matrix), and of size
+{\tt nz} (the number of nonzeros in a matrix).  UMFPACK does not handle
+zero-dimensioned arrays;
+it returns an error code if {\tt m} or {\tt n}
+are zero.  However, {\tt nz} can be zero, since all singular matrices are
+handled correctly.  If you attempt to {\tt malloc} an array of size {\tt nz}
+= 0, however, {\tt malloc} will return a null pointer which UMFPACK will report
+as a missing argument.  If you {\tt malloc} an array of
+size {\tt nz} to pass to UMFPACK, make sure that you handle the {\tt nz} = 0
+case correctly (use a size equal to the maximum of {\tt nz} and 1, or use a
+size of {\tt nz+1}).
+%-------------------------------------------------------------------------------
+\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}:
+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 summarized in Section~\ref{control_param}.
+Three particular parameters deserve 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 Control [UMFPACK\_PIVOT\_TOLERANCE]} has a
+default value of 0.1, and is used for the unsymmetric strategy.
+For complex matrices, a cheap approximation of the absolute value is
+used for the threshold pivoting test
+($|a| \approx |a_{\mbox{real}}|+|a_{\mbox{imag}}|$).
+For the symmetric strategy, a second tolerance is used for diagonal
+entries: \newline {\tt Control [UMFPACK\_SYM\_PIVOT\_TOLERANCE]}, with
+a default value of 0.001.  The first parameter (with a default of 0.1)
+is used for any off-diagonal candidate pivot entries.
+These two parameters may be too small for some matrices, particularly for
+ill-conditioned or poorly scaled ones.  With the default pivot tolerances
+and default iterative refinement,
+{\tt x = umfpack (A,'}$\backslash${\tt ',b)}
+is just as accurate as (or more accurate) than
+{\tt x = A}$\backslash${\tt b}
+in MATLAB 6.1 for nearly all matrices.
+If {\tt Control [UMFPACK\_PIVOT\_TOLERANCE]} is zero, than any
+nonzero entry is acceptable as a pivot (this is changed from Version 4.0,
+which treated a value of 0.0 the same as 1.0).  If the symmetric strategy is
+used, and {\tt Control [UMFPACK\_SYM\_PIVOT\_TOLERANCE]} is zero, then any
+nonzero entry on the diagonal is accepted as a pivot.  Off-diagonal pivoting
+will still occur if the diagonal entry is exactly zero.  The
+{\tt Control [UMFPACK\_SYM\_PIVOT\_TOLERANCE]} parameter is new to Version
+4.1.  It is similar in function to the pivot tolerance for left-looking
+methods (the MATLAB {\tt THRESH} option in {\tt [L,U,P] = lu (A, THRESH)},
+and the pivot tolerance parameter in SuperLU).
+The parameter {\tt Control [UMFPACK\_STRATEGY]} can be used to bypass
+UMFPACK's automatic strategy selection.  The automatic strategy nearly
+always selects the best method.  When it does not, the different methods
+nearly always give about the same quality of results.  There may be
+cases where the automatic strategy fails to pick a good strategy. Also,
+you can save some computing time if you know the right strategy for your
+set of matrix problems.
+\item {\tt umfpack\_*\_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
+or AMD pre-orderings.  For example, a graph partitioning-based order
+of $\m{A}\tr\m{A}$ would be suitable for UMFPACK's unsymmetric strategy.
+A partitioning of $\m{A}+\m{A}\tr$ would be suitable for UMFPACK's
+symmetric or 2-by-2 strategies.
+\item {\tt umfpack\_*\_wsolve}:
+An alternative to {\tt umfpack\_*\_solve} which does not dynamically
+allocate any memory.  Requires the user to pass two additional work
+arrays.
+\end{itemize}
+%-------------------------------------------------------------------------------
+\subsection{Matrix manipulation routines}
+\label{triplet}
+%-------------------------------------------------------------------------------
+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 or
+four 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 complex versions use another array to hold the imaginary part:
+{\footnotesize
+\begin{verbatim}
+double Tz [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]}.  For the complex versions,
+{\tt Tx[k]} is the real part of $a_{ij}$ and
+{\tt Tz[k]} is the imaginary part.
+The triplets can be in any
+order in the {\tt Ti}, {\tt Tj}, and {\tt Tx} arrays (and {\tt Tz} for
+the complex versions), and duplicate entries may
+exist.
+If {\tt Tz} is NULL, then the array {\tt Tx} becomes of size {\tt 2*nz},
+and the real and imaginary parts of the
+{\tt k}-th triplet are located in {\tt Tx[2*k]} and {\tt Tx[2*k+1]},
+respectively.
+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.
+Four routines are provided for manipulating sparse matrices:
+\begin{itemize}
+\item {\tt umfpack\_*\_triplet\_to\_col}:
+Converts a triplet form of a matrix to compressed column form (ready for
+input to \newline
+{\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.
+\item {\tt umfpack\_*\_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}:
+Transposes and optionally permutes a column form matrix \cite{Gustavson78}.
+Identical to
+{\tt R = A(P,Q)'} (linear algebraic transpose, using the complex conjugate)
+or {\tt R = A(P,Q).'} (the array transpose)
+in MATLAB, except for the presence of numerically zero entries.
+Factorizing $\m{A}\tr$ and then solving $\m{Ax}=\m{b}$ with the transposed
+factors can sometimes be much faster or much slower than factorizing
+$\m{A}$.  It is highly dependent on your particular matrix.
+\item {\tt umfpack\_*\_scale}:
+Applies the row scale factors to a user-provided vector.  This is not
+required to solve the sparse linear system $\m{Ax}=\m{b}$ or
+$\m{A}\tr\m{x}=\m{b}$, since {\tt umfpack\_*\_solve} applies the scale
+factors for those systems.
+\end{itemize}
+It is quite easy to add matrices in triplet form, subtract them, transpose
+them, permute them, construct a submatrix, and multiply a triplet-form matrix
+times a vector.  UMFPACK does not provide code for these basic operations,
+however.  Refer to the discussion of
+{\tt umfpack\_*\_triplet\_to\_col} in Section~\ref{Manipulate} for more details
+on how to compute these operations in your own code.
+The only primary matrix operation not provided by UMFPACK is the
+multiplication of two sparse matrices \cite{Gustavson78}.
+A future package under development (as of Jan. 2005), CHOLMOD,
+will provide many of these matrix operations, which
+can then be used in conjunction with UMFPACK.
+Watch my web page for details.
+%-------------------------------------------------------------------------------
+\subsection{Getting the contents of opaque objects}
+%-------------------------------------------------------------------------------
+There are cases where you may wish 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.  You cannot do anything with them
+except to pass them back to subsequent calls to UMFPACK.  Three routines
+are provided for copying their contents into user-provided arrays using simpler
+data structures.  Four routines are provided for saving and loading the
+{\tt Numeric} and {\tt Symbolic} objects to/from binary files.
+An additional routine is provided that computes the determinant.
+They are fully described in Section~\ref{Get}:
+\begin{itemize}
+\item {\tt umfpack\_*\_get\_lunz}:
+Returns the number of nonzeros in $\m{L}$ and $\m{U}$.
+\item {\tt umfpack\_*\_get\_numeric}:
+Copies $\m{L}$, $\m{U}$, $\m{P}$, $\m{Q}$, and $\m{R}$
+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 (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
+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 column of $\m{PAQ}$.  This is identical to how MATLAB uses
+permutation vectors (type {\tt help colamd} in MATLAB 6.1 or later).
+\item {\tt umfpack\_*\_get\_symbolic}:
+Copies the contents of the {\tt Symbolic} object (the initial row and column
+preordering, supernodal column elimination tree, and information
+about each frontal matrix) into arrays provided by the user.
+\item {\tt umfpack\_*\_get\_determinant}:
+Computes the determinant from the diagonal of $\m{U}$ and the permutations
+$\m{P}$ and $\m{Q}$.  This is mostly of theoretical interest.
+It is not a good test to determine if your matrix is singular or not.
+\item {\tt umfpack\_*\_save\_numeric}:
+Saves a copy of the {\tt Numeric} object to a file, in binary format.
+\item {\tt umfpack\_*\_load\_numeric}:
+Creates a {\tt Numeric} object by loading it from a file created
+by {\tt umfpack\_*\_save\_numeric}.
+\item {\tt umfpack\_*\_save\_symbolic}:
+Saves a copy of the {\tt Symbolic} object to a file, in binary format.
+\item {\tt umfpack\_*\_load\_symbolic}:
+Creates a {\tt Symbolic} object by loading it from a file created
+by {\tt umfpack\_*\_save\_symbolic}.
+\end{itemize}
+UMFPACK itself does not make use of these routines;
+they are provided solely for returning the contents of the opaque
+{\tt Symbolic} and {\tt Numeric} objects to the user, and saving/loading
+them to/from a binary file.  None of them do any computation, except for
+{\tt umfpack\_*\_get\_determinant}.
+%-------------------------------------------------------------------------------
+\subsection{Reporting routines}
+\label{Reporting}
+%-------------------------------------------------------------------------------
+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}:
+Prints the status (return value) of other {\tt umfpack\_*} routines.
+\item {\tt umfpack\_*\_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}:
+Prints the {\tt Control} settings.
+\item {\tt umfpack\_*\_report\_matrix}:
+Verifies and prints a compressed column-form or compressed row-form sparse
+matrix.
+\item {\tt umfpack\_*\_report\_triplet}:
+Verifies and prints a matrix in triplet form.
+\item {\tt umfpack\_*\_report\_symbolic}:
+Verifies and prints a {\tt Symbolic} object.
+\item {\tt umfpack\_*\_report\_numeric}:
+Verifies and prints a {\tt Numeric} object.
+\item {\tt umfpack\_*\_report\_perm}:
+Verifies and prints a permutation vector.
+\item {\tt umfpack\_*\_report\_vector}:
+Verifies and prints a real or complex 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 of an $m$-by-$n$
+matrix are in the range 0 to $m-1$ or $n-1$, respectively\footnote{Complex
+matrices in MATLAB use the split array form, with one {\tt double} array
+for the real part and another array for the imaginary part.  UMFPACK
+supports that format, as well as the packed complex format (new to Version 4.4).}
+It prints them as if they are in the range 1 to $m$ or $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 1.
+Here is a summary of how the routines use this print level parameter:
+\begin{itemize}
+\item {\tt umfpack\_*\_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}.  A warning message is printed if
+the matrix is singular.  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}:
+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}:
+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}
+This print level parameter has an additional effect on the MATLAB mexFunction.
+If zero, then no warnings of singular or nearly singular matrices are
+printed (similar to the MATLAB commands
+{\tt warning off MATLAB:singularMatrix} and
+{\tt warning off MATLAB:nearlySingularMatrix}).
+%-------------------------------------------------------------------------------
+\subsection{Utility routines}
+%-------------------------------------------------------------------------------
+UMFPACK v4.0 included 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 still available in
+UMFPACK v4.1 and following, but not used internally.
+Two new timing routines are provided in UMFPACK Version 4.1 and following,
+{\tt umfpack\_tic} and {\tt umfpack\_toc}.  They use POSIX-compliant
+{\tt sysconf} and {\tt times} routines to find both the CPU time
+and wallclock time.
+These three routines are the only user-callable
+routine that is identical in all four {\tt int}/{\tt long}, real/complex
+versions (there is no {\tt umfpack\_di\_timer} routine, for example).
+%-------------------------------------------------------------------------------
+\subsection{Control parameters}
+\label{control_param}
+%-------------------------------------------------------------------------------
+UMFPACK uses an optional {\tt double} array (currently of size 20)
+to modify its control parameters.  If you pass {\tt (double *) NULL} instead
+of a {\tt Control} array, then defaults are used.  These defaults provide
+nearly optimal performance (both speed, memory usage, and numerical accuracy)
+for a wide range of matrices from real applications.
+This array will almost certainly grow in size in future releases,
+so be sure to dimension your {\tt Control} array to be of size
+{\tt UMFPACK\_CONTROL}.  That constant is currently defined to be 20,
+but may increase in future versions, since all 20 entries are in use.
+The contents of this array may be modified by the user
+(see {\tt umfpack\_*\_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.
+Note that ANSI C uses 0-based indexing, while MATLAB uses 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
+{\tt Control(1)}  & {\tt Control[UMFPACK\_PRL]} & 1 & printing level \\
+{\tt Control(2)}  & {\tt Control[UMFPACK\_DENSE\_ROW]} & 0.2 & dense row parameter \\
+{\tt Control(3)}  & {\tt Control[UMFPACK\_DENSE\_COL]} & 0.2 & dense column parameter \\
+{\tt Control(4)}  & {\tt Control[UMFPACK\_PIVOT\_TOLERANCE]} & 0.1 & partial pivoting tolerance \\
+{\tt Control(5)}  & {\tt Control[UMFPACK\_BLOCK\_SIZE]} & 32 & BLAS block size \\
+{\tt Control(6)}  & {\tt Control[UMFPACK\_STRATEGY]} & 0 (auto) & select strategy \\
+{\tt Control(7)}  & {\tt Control[UMFPACK\_ALLOC\_INIT]} & 0.7 & initial memory allocation  \\
+{\tt Control(8)}  & {\tt Control[UMFPACK\_IRSTEP]} & 2 & max iter. refinement steps \\
+{\tt Control(13)} & {\tt Control[UMFPACK\_2BY2\_TOLERANCE]} & 0.01 & defines ``large'' entries \\
+{\tt Control(14)} & {\tt Control[UMFPACK\_FIXQ]} & 0 (auto) & fix or modify Q \\
+{\tt Control(15)} & {\tt Control[UMFPACK\_AMD\_DENSE]} & 10 & AMD dense row/column parameter \\
+{\tt Control(16)} & {\tt Control[UMFPACK\_SYM\_PIVOT\_TOLERANCE]} & 0.001 & for diagonal entries \\
+{\tt Control(17)} & {\tt Control[UMFPACK\_SCALE]} & 1 (sum) & row scaling (none, sum, or max) \\
+{\tt Control(18)} & {\tt Control[UMFPACK\_FRONT\_ALLOC\_INIT]} & 0.5 & frontal matrix allocation ratio \\
+{\tt Control(19)} & {\tt Control[UMFPACK\_DROPTOL]} & 0 & drop tolerance \\
+{\tt Control(20)} & {\tt Control[UMFPACK\_AGGRESSIVE]} & 1 (yes) & aggressive absorption \\
+& & & in AMD and COLAMD \\
+%
+\hline
+\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]} & - & 1 if {\tt getrusage} used \\
+{\tt Control(12)} & {\tt Control[UMFPACK\_COMPILED\_IN\_DEBUG\_MODE]} & - & true if debug mode enabled \\
+\hline
+\end{tabular}
+}
+\end{table}
+Let $\alpha_r = ${\tt Control [UMFPACK\_DENSE\_ROW]},
+$\alpha_c = ${\tt Control [UMFPACK\_DENSE\_COL]}, and
+$\alpha = ${\tt Control [UMFPACK\_AMD\_DENSE]}.
+Suppose the submatrix $\m{S}$, obtained after eliminating pivots with
+zero Markowitz cost, is $m$-by-$n$.
+Then a row is considered ``dense'' if it has more than
+$\max (16, 16 \alpha_r \sqrt{n})$ entries.
+A column is considered ``dense'' if it has more than
+$\max (16, 16 \alpha_c \sqrt{m})$ entries.
+These rows and columns are treated different in COLAMD and during numerical
+factorization.   In COLAMD, dense columns are placed last in their natural
+order, and dense rows are ignored.  During numerical factorization, dense
+rows are stored differently.
+In AMD, a row/column of the square matrix $\m{S}+\m{S}\tr$ is
+considered ``dense'' if it has more than $\max (16, \alpha \sqrt{n})$ entries.
+These rows/columns are placed last in AMD's output ordering.
+For more details on the control parameters, refer to the documentation of
+{\tt umfpack\_*\_qsymbolic}, {\tt umfpack\_*\_numeric}, {\tt umfpack\_*\_solve},
+and the {\tt umfpack\_*\_report\_*} routines,
+in Sections~\ref{Primary}~through~\ref{Report}, below.
+%-------------------------------------------------------------------------------
+\subsection{Error codes}
+\label{error_codes}
+%-------------------------------------------------------------------------------
+Many of the routines return a {\tt status} value.
+This is also returned as the first entry in the {\tt Info} array, for
+those routines with that argument.  The following list summarizes
+all of the error codes in UMFPACK.  Each error code is given a
+specific name in the {\tt umfpack.h} include file, so you can use
+those constants instead of hard-coded values in your program.
+Future versions may report additional error codes.
+A value of zero means everything was successful, and the matrix is
+non-singular.  A value greater than zero means the routine was successful,
+but a warning occurred.
+A negative value means the routine was not successful.
+In this case, no {\tt Symbolic} or {\tt Numeric} object was created.
+\begin{itemize}
+\item {\tt UMFPACK\_OK},  (0):  UMFPACK was successful.
+\item {\tt UMFPACK\_WARNING\_singular\_matrix},  (1):  Matrix is singular.
+There are exact zeros on the diagonal of $\m{U}$.
+\item {\tt UMFPACK\_WARNING\_determinant\_underflow}, (2):
+The determinant is nonzero, but smaller in magnitude than
+the smallest positive floating-point number.
+\item {\tt UMFPACK\_WARNING\_determinant\_overflow}, (3):
+The determinant is larger in magnitude than
+the largest positive floating-point number (IEEE Inf).
+\item {\tt UMFPACK\_ERROR\_out\_of\_memory},  (-1):  Not enough memory.
+The ANSI C {\tt malloc} or {\tt realloc} routine failed.
+\item {\tt UMFPACK\_ERROR\_invalid\_Numeric\_object},  (-3):
+Routines that take a {\tt Numeric} object as input (or load it
+from a file) check this object and return this error code if it is
+invalid.  This can be caused by a memory leak or overrun in your
+program, which can overwrite part of the Numeric object.  It can also
+be caused by passing a Symbolic object by mistake, or some other pointer.
+If you try to factorize a matrix using one version of UMFPACK and
+then use the factors in another version, this error code will trigger as
+well.  You cannot factor your matrix using
+version 4.0 and then solve with version 4.1, for example.\footnote{
+Exception: v4.3, v4.3.1, and v4.4 use identical data structures
+for the {\tt Numeric} and {\tt Symbolic} objects}.
+You cannot use different precisions of the same version
+(real and complex, for example).
+It is possible for the {\tt Numeric} object to be corrupted by your
+program in subtle ways that are not detectable by this quick check.
+In this case, you may see an
+{\tt UMFPACK\_ERROR\_different\_pattern} error code, or even an
+{\tt UMFPACK\_ERROR\_internal\_error}.
+\item {\tt UMFPACK\_ERROR\_invalid\_Symbolic\_object},  (-4):
+Routines that take a {\tt Symbolic} object as input (or load it
+from a file) check this object and return this error code if it is
+invalid.  The causes of this error are analogous to the
+{\tt UMFPACK\_ERROR\_invalid\_Numeric\_object} error described above.
+\item {\tt UMFPACK\_ERROR\_argument\_missing},  (-5):
+Some arguments of some are optional (you can pass a {\tt NULL} pointer
+instead of an array).  This error code occurs if you pass a {\tt NULL}
+pointer when that argument is required to be present.
+\item {\tt UMFPACK\_ERROR\_n\_nonpositive}  (-6):
+The number of rows or columns of the matrix must be greater than zero.
+\item {\tt UMFPACK\_ERROR\_invalid\_matrix}  (-8):
+The matrix is invalid.  For the column-oriented input, this error
+code will occur if the contents of {\tt Ap} and/or {\tt Ai} are invalid.
+{\tt Ap} is an integer array of size {\tt n\_col+1}.
+On input, it holds the
+``pointers'' for the column form of the sparse matrix $\m{A}$.
+Column {\tt j} of
+the matrix A is held in {\tt Ai [(Ap [j])} \ldots {\tt (Ap [j+1]-1)]}.
+The first entry, {\tt Ap [0]}, must be zero,
+and {\tt Ap [j]} $\le$ {\tt Ap [j+1]} must hold for all
+{\tt j} in the range 0 to {\tt n\_col-1}.
+The value {\tt nz = Ap [n\_col]} is thus the
+total number of entries in the pattern of the matrix A.
+{\tt nz} must be greater than or equal to zero.
+The nonzero pattern (row indices) for column {\tt j} is stored in
+{\tt Ai [(Ap [j])} \ldots {\tt (Ap [j+1]-1)]}.  The row indices in a given
+column {\tt j}
+must be in ascending order, and no duplicate row indices may be present.
+Row indices must be in the range 0 to {\tt n\_row-1}
+(the matrix is 0-based).
+Some routines take a triplet-form input, with arguments
+{\tt nz}, {\tt Ti}, and {\tt Tj}.  This error code is returned
+if {\tt nz} is less than zero,
+if any row    index in {\tt Ti} is outside the range 0 to {\tt n\_col-1}, or
+if any column index in {\tt Tj} is outside the range 0 to {\tt n\_row-1}.
+\item {\tt UMFPACK\_ERROR\_different\_pattern},  (-11):
+The most common cause of this error is that the pattern of the
+matrix has changed between the symbolic and numeric factorization.
+It can also occur if the {\tt Numeric} or {\tt Symbolic} object has
+been subtly corrupted by your program.
+\item {\tt UMFPACK\_ERROR\_invalid\_system},  (-13):
+The {\tt sys} argument provided to one of the solve routines is invalid.
+\item {\tt UMFPACK\_ERROR\_invalid\_permutation},  (-15):
+The permutation vector provided as input is invalid.
+\item {\tt UMFPACK\_ERROR\_file\_IO},  (-17):
+This error code is returned by the routines that save and load
+the {\tt Numeric} or {\tt Symbolic} objects to/from a file, if a
+file I/O error has occurred.  The file may not exist or may not be readable,
+you may be trying to create a file that you don't have permission to create,
+or you may be out of disk space.  The file you are trying to read might
+be the wrong one, and an earlier end-of-file condition would then result
+in this error.
+\item {\tt UMFPACK\_ERROR\_internal\_error},  (-911):
+An internal error has occurred, of unknown cause.  This is either a bug
+in UMFPACK, or the result of a memory overrun from your program.
+Try modifying the file {\tt AMD/Source/amd\_internal.h} and adding
+the statement {\tt \#undef NDEBUG}, to enable the debugging mode.
+Recompile UMFPACK and rerun your program.
+A failed assertion might occur which
+can give you a better indication as to what is going wrong.  Be aware that
+UMFPACK will be extraordinarily slow when running in debug mode.
+If all else fails, contact the developer (davis@cise.ufl.edu) with
+as many details as possible.
+\end{itemize}
+%-------------------------------------------------------------------------------
+\subsection{Larger examples}
+%-------------------------------------------------------------------------------
+Full examples of all user-callable UMFPACK routines
+are available in four stand-alone C main programs, {\tt umfpack\_*\_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.
+The {\tt umf4.c} program for reading matrices in Harwell/Boeing format
+\cite{DuffGrimesLewis87b} is provided.  It requires three Fortran 77 programs
+({\tt readhb.f}, {\tt readhb\_nozeros.f}, and {\tt readhb\_size.f})
+for reading in the sample Harwell/Boeing files in the {\tt UMFPACK/Demo/HB}
+directory.  More matrices are available at
+http://www.cise.ufl.edu/research/sparse/matrices.
+Type {\tt make hb} in the {\tt UMFPACK/Demo/HB} directory
+to compile and run this demo.  This program was used for the experimental
+results in \cite{Davis03}.
+%-------------------------------------------------------------------------------
+\section{Synopsis of C-callable routines}
+\label{Synopsis}
+%-------------------------------------------------------------------------------
+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.
+The real, {\tt long} integer {\tt umfpack\_dl\_*} routines are
+identical to the real, {\tt int} routines, except that {\tt \_di\_} is replaced
+with {\tt \_dl\_} in the name, and all {\tt int} arguments become {\tt long}.
+Similarly, the complex, {\tt long} integer {\tt umfpack\_zl\_*} routines are
+identical to the complex, {\tt int} routines, except that {\tt \_zi\_} is
+replaced
+with {\tt \_zl\_} in the name, and all {\tt int} arguments become {\tt long}.
+Only the real and complex {\tt int} versions are listed in the synopsis below.
+The matrix $\m{A}$ is {\tt m}-by-{\tt n} with {\tt nz} entries.
+The {\tt sys} argument of {\tt umfpack\_*\_solve}
+is an integer in the range 0 to 14 which defines which linear system is
+to be solved.
+\footnote{Integer values for {\tt sys} are used instead of strings (as in LINPACK
+and LAPACK) to avoid C-to-Fortran portability issues.}
+Valid values are listed in Table~\ref{sys}.
+The notation $\m{A}\he$ refers to the matrix transpose, which is the
+complex conjugate transpose for complex matrices ({\tt A'} in MATLAB).
+The array transpose is $\m{A}\tr$, which is {\tt A.'} in MATLAB.
+\begin{table}
+\begin{center}
+\caption{UMFPACK {\tt sys} parameter}
+\label{sys}
+{\footnotesize
+\begin{tabular}{ll|l}
+\hline
+Value & & system \\
+\hline
+& & \\
+{\tt UMFPACK\_A}      &  (0) & $\m{Ax}=\m{b}$ \\
+{\tt UMFPACK\_At}     &  (1) & $\m{A}\he\m{x}=\m{b}$ \\
+{\tt UMFPACK\_Aat}    &  (2) & $\m{A}\tr\m{x}=\m{b}$ \\
+& & \\
+\hline
+& & \\
+{\tt UMFPACK\_Pt\_L}  &  (3) & $\m{P}\tr\m{Lx}=\m{b}$ \\
+{\tt UMFPACK\_L}      &  (4) & $\m{Lx}=\m{b}$ \\
+{\tt UMFPACK\_Lt\_P}  &  (5) & $\m{L}\he\m{Px}=\m{b}$ \\
+{\tt UMFPACK\_Lat\_P} &  (6) & $\m{L}\tr\m{Px}=\m{b}$ \\
+{\tt UMFPACK\_Lt}     &  (7) & $\m{L}\he\m{x}=\m{b}$ \\
+{\tt UMFPACK\_Lat}    &  (8) & $\m{L}\tr\m{x}=\m{b}$ \\
+& & \\
+\hline
+& & \\
+{\tt UMFPACK\_U\_Qt}  &  (9) & $\m{UQ}\tr\m{x}=\m{b}$ \\
+{\tt UMFPACK\_U}      & (10) & $\m{Ux}=\m{b}$ \\
+{\tt UMFPACK\_Q\_Ut}  & (11) & $\m{QU}\he\m{x}=\m{b}$ \\
+{\tt UMFPACK\_Q\_Uat} & (12) & $\m{QU}\tr\m{x}=\m{b}$ \\
+{\tt UMFPACK\_Ut}     & (13) & $\m{U}\he\m{x}=\m{b}$ \\
+{\tt UMFPACK\_Uat}    & (14) & $\m{U}\tr\m{x}=\m{b}$ \\
+& & \\
+\hline
+\end{tabular}
+}
+\end{center}
+\end{table}
+%-------------------------------------------------------------------------------
+\subsection{Primary routines: real/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+#include "umfpack.h"
+int status, sys, n, m, nz, Ap [n+1], Ai [nz] ;
+double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO], Ax [nz], X [n], B [n] ;
+void *Symbolic, *Numeric ;
+status = umfpack_di_symbolic (m, n, Ap, Ai, Ax, &Symbolic, Control, Info) ;
+status = umfpack_di_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;
+status = umfpack_di_solve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info) ;
+umfpack_di_free_symbolic (&Symbolic) ;
+umfpack_di_free_numeric (&Numeric) ;
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\subsection{Alternative routines: real/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+int Qinit [n], Wi [n] ;
+double W [5*n] ;
+umfpack_di_defaults (Control) ;
+status = umfpack_di_qsymbolic (m, n, Ap, Ai, Ax, Qinit, &Symbolic, Control, Info) ;
+status = umfpack_di_wsolve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info, Wi, W) ;
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\subsection{Matrix manipulation routines: real/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+int Ti [nz], Tj [nz], P [m], Q [n], Rp [m+1], Ri [nz], Map [nz] ;
+double Tx [nz], Rx [nz], Y [m], Z [m] ;
+status = umfpack_di_col_to_triplet (n, Ap, Tj) ;
+status = umfpack_di_triplet_to_col (m, n, nz, Ti, Tj, Tx, Ap, Ai, Ax, Map) ;
+status = umfpack_di_transpose (m, n, Ap, Ai, Ax, P, Q, Rp, Ri, Rx) ;
+status = umfpack_di_scale (Y, Z, Numeric) ;
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\subsection{Getting the contents of opaque objects: real/{\tt int}}
+%-------------------------------------------------------------------------------
+The {\tt filename} string should be large enough to hold the name of a file.
+{\footnotesize
+\begin{verbatim}
+int lnz, unz, Lp [m+1], Lj [lnz], Up [n+1], Ui [unz], do_recip ;
+double Lx [lnz], Ux [unz], D [min (m,n)], Rs [m], Mx [1], Ex [1] ;
+int nfr, nchains, P1 [m], Q1 [n], Front_npivcol [n+1], Front_parent [n+1], Front_1strow [n+1],
+Front_leftmostdesc [n+1], Chain_start [n+1], Chain_maxrows [n+1], Chain_maxcols [n+1] ;
+char filename [100] ;
+status = umfpack_di_get_lunz (&lnz, &unz, &m, &n, &nz_udiag, Numeric) ;
+status = umfpack_di_get_numeric (Lp, Lj, Lx, Up, Ui, Ux, P, Q, D,
+&do_recip, Rs, Numeric) ;
+status = umfpack_di_get_symbolic (&m, &n, &n1, &nz, &nfr, &nchains, P1, Q1,
+Front_npivcol, Front_parent, Front_1strow, Front_leftmostdesc,
+Chain_start, Chain_maxrows, Chain_maxcols, Symbolic) ;
+status = umfpack_di_load_numeric (&Numeric, filename) ;
+status = umfpack_di_save_numeric (Numeric, filename) ;
+status = umfpack_di_load_symbolic (&Symbolic, filename) ;
+status = umfpack_di_save_symbolic (Symbolic, filename) ;
+status = umfapck_di_get_determinant (Mx, Ex, Numeric, Info) ;
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\subsection{Reporting routines: real/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+umfpack_di_report_status (Control, status) ;
+umfpack_di_report_control (Control) ;
+umfpack_di_report_info (Control, Info) ;
+status = umfpack_di_report_matrix (m, n, Ap, Ai, Ax, 1, Control) ;
+status = umfpack_di_report_matrix (m, n, Rp, Ri, Rx, 0, Control) ;
+status = umfpack_di_report_numeric (Numeric, Control) ;
+status = umfpack_di_report_perm (m, P, Control) ;
+status = umfpack_di_report_perm (n, Q, Control) ;
+status = umfpack_di_report_symbolic (Symbolic, Control) ;
+status = umfpack_di_report_triplet (m, n, nz, Ti, Tj, Tx, Control) ;
+status = umfpack_di_report_vector (n, X, Control) ;
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\subsection{Primary routines: complex/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+double Az [nz], Xx [n], Xz [n], Bx [n], Bz [n] ;
+status = umfpack_zi_symbolic (m, n, Ap, Ai, Ax, Az, &Symbolic, Control, Info) ;
+status = umfpack_zi_numeric (Ap, Ai, Ax, Az, Symbolic, &Numeric, Control, Info) ;
+status = umfpack_zi_solve (sys, Ap, Ai, Ax, Az, Xx, Xz, Bx, Bz, Numeric, Control, Info) ;
+umfpack_zi_free_symbolic (&Symbolic) ;
+umfpack_zi_free_numeric (&Numeric) ;
+\end{verbatim}
+}
+The arrays {\tt Ax}, {\tt Bx}, and {\tt Xx} double in size if
+any imaginary argument ({\tt Az}, {\tt Xz}, or {\tt Bz}) is {\tt NULL}.
+%-------------------------------------------------------------------------------
+\subsection{Alternative routines: complex/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+double Wz [10*n] ;
+umfpack_zi_defaults (Control) ;
+status = umfpack_zi_qsymbolic (m, n, Ap, Ai, Ax, Az, Qinit, &Symbolic, Control, Info) ;
+status = umfpack_zi_wsolve (sys, Ap, Ai, Ax, Az, Xx, Xz, Bx, Bz, Numeric, Control, Info, Wi, Wz) ;
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\subsection{Matrix manipulation routines: complex/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+double Tz [nz], Rz [nz], Yx [m], Yz [m], Zx [m], Zz [m] ;
+status = umfpack_zi_col_to_triplet (n, Ap, Tj) ;
+status = umfpack_zi_triplet_to_col (m, n, nz, Ti, Tj, Tx, Tz, Ap, Ai, Ax, Az, Map) ;
+status = umfpack_zi_transpose (m, n, Ap, Ai, Ax, Az, P, Q, Rp, Ri, Rx, Rz, 1) ;
+status = umfpack_zi_transpose (m, n, Ap, Ai, Ax, Az, P, Q, Rp, Ri, Rx, Rz, 0) ;
+status = umfpack_zi_scale (Yx, Yz, Zx, Zz, Numeric) ;
+\end{verbatim}
+}
+The arrays {\tt Tx}, {\tt Rx}, {\tt Yx}, and {\tt Zx} double in size if
+any imaginary argument ({\tt Tz}, {\tt Rz}, {\tt Yz}, or {\tt Zz}) is {\tt NULL}.
+%-------------------------------------------------------------------------------
+\subsection{Getting the contents of opaque objects: complex/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+double Lz [lnz], Uz [unz], Dx [min (m,n)], Dz [min (m,n)], Mz [1] ;
+status = umfpack_zi_get_lunz (&lnz, &unz, &m, &n, &nz_udiag, Numeric) ;
+status = umfpack_zi_get_numeric (Lp, Lj, Lx, Lz, Up, Ui, Ux, Uz, P, Q, Dx, Dz,
+&do_recip, Rs, Numeric) ;
+status = umfpack_zi_get_symbolic (&m, &n, &n1, &nz, &nfr, &nchains, P1, Q1,
+Front_npivcol, Front_parent, Front_1strow, Front_leftmostdesc,
+Chain_start, Chain_maxrows, Chain_maxcols, Symbolic) ;
+status = umfpack_zi_load_numeric (&Numeric, filename) ;
+status = umfpack_zi_save_numeric (Numeric, filename) ;
+status = umfpack_zi_load_symbolic (&Symbolic, filename) ;
+status = umfpack_zi_save_symbolic (Symbolic, filename) ;
+status = umfapck_zi_get_determinant (Mx, Mz, Ex, Numeric, Info) ;
+\end{verbatim}
+}
+The arrays {\tt Lx}, {\tt Ux}, {\tt Dx}, and {\tt Mx} double in size if
+any imaginary argument ({\tt Lz}, {\tt Uz}, {\tt Dz}, or {\tt Mz}) is {\tt NULL}.
+%-------------------------------------------------------------------------------
+\subsection{Reporting routines: complex/{\tt int}}
+%-------------------------------------------------------------------------------
+{\footnotesize
+\begin{verbatim}
+umfpack_zi_report_status (Control, status) ;
+umfpack_zi_report_control (Control) ;
+umfpack_zi_report_info (Control, Info) ;
+status = umfpack_zi_report_matrix (m, n, Ap, Ai, Ax, Az, 1, Control) ;
+status = umfpack_zi_report_matrix (m, n, Rp, Ri, Rx, Rz, 0, Control) ;
+status = umfpack_zi_report_numeric (Numeric, Control) ;
+status = umfpack_zi_report_perm (m, P, Control) ;
+status = umfpack_zi_report_perm (n, Q, Control) ;
+status = umfpack_zi_report_symbolic (Symbolic, Control) ;
+status = umfpack_zi_report_triplet (m, n, nz, Ti, Tj, Tx, Tz, Control) ;
+status = umfpack_zi_report_vector (n, Xx, Xz, Control) ;
+\end{verbatim}
+}
+The arrays {\tt Ax}, {\tt Rx}, {\tt Tx}, and {\tt Xx} double in size if
+any imaginary argument ({\tt Az}, {\tt Rz}, {\tt Tz}, or {\tt Xz}) is {\tt NULL}.
+%-------------------------------------------------------------------------------
+\subsection{Utility routines}
+%-------------------------------------------------------------------------------
+These routines are the same in all four versions of UMFPACK.
+{\footnotesize
+\begin{verbatim}
+double t, s [2] ;
+t = umfpack_timer ( ) ;
+umfpack_tic (s) ;
+umfpack_toc (s) ;
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\subsection{AMD ordering routines}
+%-------------------------------------------------------------------------------
+UMFPACK makes use of the AMD ordering package for its symmetric ordering
+strategy.  You may also use these four user-callable routines in your own C
+programs.  You need to include the {\tt amd.h} file only if you make direct
+calls to the AMD routines themselves.  The {\tt int} versions are summarized
+below; {\tt long} versions are also available.  Refer to the AMD User Guide
+for more information, or to the file {\tt amd.h} which documents these routines.
+{\footnotesize
+\begin{verbatim}
+#include "amd.h"
+double amd_control [AMD_CONTROL], amd_info [AMD_INFO] ;
+amd_defaults (amd_control) ;
+status = amd_order (n, Ap, Ai, P, amd_control, amd_info) ;
+amd_control (amd_control) ;
+amd_info (amd_info) ;
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\section{Using UMFPACK in a Fortran program}
+%-------------------------------------------------------------------------------
+UMFPACK includes a basic Fortran 77 interface to some of the C-callable
+UMFPACK routines.
+Since interfacing C and Fortran programs is not portable, this interface might
+not work with all C and Fortran compilers.  Refer to Section~\ref{Install} for
+more details.  The following Fortran routines are provided.
+The list includes the C-callable routines that the Fortran interface
+routine calls.  Refer to the corresponding C routines in Section~\ref{C} for
+more details on what the Fortran routine does.
+\begin{itemize}
+\item {\tt umf4def}: sets the default control parameters
+({\tt umfpack\_di\_defaults}).
+\item {\tt umf4sym}: pre-ordering and symbolic factorization
+({\tt umfpack\_di\_symbolic}).
+\item {\tt umf4num}: numeric factorization
+({\tt umfpack\_di\_numeric}).
+\item {\tt umf4solr}: solve a linear system with iterative refinement
+({\tt umfpack\_di\_solve}).
+\item {\tt umf4sol}: solve a linear system without iterative refinement
+({\tt umfpack\_di\_solve}).  Sets {\tt Control [UMFPACK\_IRSTEP]}
+to zero, and does not require the matrix $\m{A}$.
+\item {\tt umf4scal}: scales a vector using UMFPACK's scale factors
+({\tt umfpack\_di\_scale}).
+\item {\tt umf4fnum}: free the {\tt Numeric} object
+({\tt umfpack\_di\_free\_numeric}).
+\item {\tt umf4fsym}: free the {\tt Symbolic} object
+({\tt umfpack\_di\_free\_symbolic}).
+\item {\tt umf4pcon}: prints the control parameters
+({\tt umfpack\_di\_report\_control}).
+\item {\tt umf4pinf}: print statistics
+({\tt umfpack\_di\_report\_info}).
+\item {\tt umf4snum}: save the {\tt Numeric} object to a file
+({\tt umfpack\_di\_save\_numeric}).
+\item {\tt umf4ssym}: save the {\tt Symbolic} object to a file
+({\tt umfpack\_di\_save\_symbolic}).
+\item {\tt umf4lnum}: load the {\tt Numeric} object from a file
+({\tt umfpack\_di\_load\_numeric}).
+\item {\tt umf4lsym}: load the {\tt Symbolic} object from a file
+({\tt umfpack\_di\_load\_symbolic}).
+\end{itemize}
+The matrix $\m{A}$ is passed to UMFPACK in compressed column form, with 0-based
+indices.  In Fortran, for an {\tt m}-by-{\tt n} matrix $\m{A}$ with {\tt nz}
+entries, the row indices of the first column (column 1) are in
+{\tt Ai (Ap(1)+1} \ldots {\tt Ap(2))}, with values in
+{\tt Ax (Ap(1)+1} \ldots {\tt Ap(2))}.  The last column (column {\tt n}) is in
+{\tt Ai (Ap(n)+1} \ldots {\tt Ap(n+1))} and
+{\tt Ax (Ap(n)+1} \ldots {\tt Ap(n+1))}.
+The number of entries in the matrix is thus {\tt nz = Ap (n+1)}.
+The row indices in {\tt Ai} are in the range 0 to {\tt m}-1.  They must be
+sorted, with no duplicate entries allowed.  None of the UMFPACK routines
+modify the input matrix $\m{A}$.
+The following definitions apply for the Fortran routines:
+{\footnotesize
+\begin{verbatim}
+integer m, n, Ap (n+1), Ai (nz), symbolic, numeric, filenum, status
+double precision Ax (nz), control (20), info (90), x (n), b (n)
+\end{verbatim}
+}
+UMFPACK's status is returned in either a {\tt status} argument, or in
+{\tt info (1)}.
+It is zero if UMFPACK was successful, 1 if the matrix is singular (this is a
+warning, not an error), and negative if an error occurred.
+Section~\ref{error_codes} summarizes the possible values of {\tt status}
+and {\tt info (1)}.
+See Table~\ref{sys} for a list of the values of the {\tt sys} argument.
+See Table~\ref{control} for a list of the control parameters (the
+Fortran usage is the same as the MATLAB usage for this array).
+For the {\tt Numeric} and {\tt Symbolic} handles, it is probably safe to
+assume that a Fortran {\tt integer} is sufficient to store a C pointer.  If
+that does not work, try defining {\tt numeric} and {\tt symbolic} in your
+Fortran program as integer arrays of size 2.  You will need to define them
+as {\tt integer*8} if you compile UMFPACK in the 64-bit mode.
+To avoid passing strings between C and Fortran in the load/save routines,
+a file number is passed instead, and the C interface constructs a file name
+(if {\tt filenum} is 42, the {\tt Numeric} file name is {\tt n42.umf}, and
+the {\tt Symbolic} file name is {\tt s42.umf}).
+The following is a summary of the calling sequence of each Fortran
+interface routine.  An example of their use is in the {\tt Demo/umf4hb.f}
+file.  That routine also includes an example of how to convert a 1-based
+sparse matrix into 0-based form.  For more details on the arguments of each
+routine, refer to the arguments of the same name in the corresponding
+C-callable routine, in Sections~\ref{Primary}~through~\ref{Utility}.
+The only exception is the {\tt control} argument of {\tt umf4sol},
+which sets {\tt control (8)} to zero to disable iterative refinement.
+Note that the solve routines do not overwrite {\tt b} with the solution,
+but return their solution in a different array, {\tt x}.
+{\footnotesize
+\begin{verbatim}
+call umf4def (control)
+call umf4sym (m, n, Ap, Ai, Ax, symbolic, control, info)
+call umf4num (Ap, Ai, Ax, symbolic, numeric, control, info)
+call umf4solr (sys, Ap, Ai, Ax, x, b, numeric, control, info)
+call umf4sol (sys, x, b, numeric, control, info)
+call umf4scal (x, b, numeric, status)
+call umf4fnum (numeric)
+call umf4fsym (symbolic)
+call umf4pcon (control)
+call umf4pinf (control)
+call umf4snum (numeric, filenum, status)
+call umf4ssym (symbolic, filenum, status)
+call umf4lnum (numeric, filenum, status)
+call umf4lsym (symbolic, filenum, status)
+\end{verbatim}
+}
+Access to the complex routines in UMFPACK is provided by the interface
+routines in {\tt umf4\_f77zwrapper.c}.  The following is a synopsis
+of each routine.  All the arguments are the same as the real versions,
+except {\tt Az}, {\tt xz}, and {\tt bz} are the imaginary parts of
+the matrix, solution, and right-hand-side, respectively.  The
+{\tt Ax}, {\tt x}, and {\tt b} are the real parts.
+{\footnotesize
+\begin{verbatim}
+call umf4zdef (control)
+call umf4zsym (m, n, Ap, Ai, Ax, Az, symbolic, control, info)
+call umf4znum (Ap, Ai, Ax, Az, symbolic, numeric, control, info)
+call umf4zsolr (sys, Ap, Ai, Ax, Az, x, xz, b, bz, numeric, control, info)
+call umf4zsol (sys, x, xz, b, bz, numeric, control, info)
+call umf4zscal (x, xz, b, bz, numeric, status)
+call umf4zfnum (numeric)
+call umf4zfsym (symbolic)
+call umf4zpcon (control)
+call umf4zpinf (control)
+call umf4zsnum (numeric, filenum, status)
+call umf4zssym (symbolic, filenum, status)
+call umf4zlnum (numeric, filenum, status)
+call umf4zlsym (symbolic, filenum, status)
+\end{verbatim}
+}
+The Fortran interface does not support the packed complex case.
+%-------------------------------------------------------------------------------
+\section{Installation}
+\label{Install}
+%-------------------------------------------------------------------------------
+%-------------------------------------------------------------------------------
+\subsection{Installing the C library}
+%-------------------------------------------------------------------------------
+The following discussion assumes you have the {\tt make} program, either in
+Unix, or in Windows with Cygwin\footnote{www.cygwin.com}.
+You can skip this section and go to next one if all you want to use is
+the UMFPACK and AMD mexFunctions in MATLAB.
+You will need to install both UMFPACK v4.4 and AMD v1.1 (or AMD v1.0) to use UMFPACK.
+The {\tt UMFPACK} and {\tt AMD} subdirectories must be placed side-by-side
+within the same directory.  AMD is a stand-alone package that
+is required by UMFPACK.  UMFPACK can be compiled without the
+BLAS \cite{DaydeDuff99,ACM679a,ATLAS,GotoVandeGeijn02},
+but your performance will be much less than what it should be.
+System-dependent configurations are in the {\tt AMD/Make}
+and {\tt UMFPACK/Make} directories (the \newline
+{\tt UMFPACK/Make} directory is actually
+just a symbolic link to {\tt AMD/Make}\footnote{Windows might not extract
+the symbolic link {\tt UMFPACK/Make} correctly.  If it doesn't, simply
+create the {\tt UMFPACK/Make} folder by copying it from {\tt AMD/Make}.}).
+You can edit the {\tt Make.include}
+files in either of those directories to customize the compilation.  The default
+settings will work on most systems, except that UMFPACK will be compiled so
+that it does not use the BLAS.  Sample configuration files are provided
+for Linux, Sun Solaris, SGI IRIX, IBM AIX, and the DEC/Compaq Alpha.
+To compile and install both packages,
+go to the {\tt UMFPACK} directory and type {\tt make}.  This will compile the
+libraries ({\tt AMD/Lib/libamd.a} and {\tt UMFPACK/Lib/libumfpack.a}).
+A demo of the AMD ordering routine will be compiled and tested in
+the {\tt AMD/Demo} directory, and five demo programs will then be
+compiled and tested in the {\tt UMFPACK/Demo} directory.
+The outputs of these demo programs will then be compared with output
+files in the distribution.  Expect to see a few differences, such as
+residual norms, compile-time control settings, and perhaps memory usage
+differences.  The AMD and UMFPACK mexFunctions for
+use in MATLAB will also be compiled.  If you do not have MATLAB 6.0 or
+later, type {\tt make lib} instead.
+If you have the GNU version of {\tt make}, the {\tt Source/GNUmakefile} and
+{\tt MATLAB/GNUmakefile} files are used.  These are much more concise than
+what the ``old'' version of {\tt make} can handle.  If you do not have
+GNU {\tt make}, the {\tt Source/Makefile} and {\tt MATLAB/Makefile} files
+are used instead.  Each UMFPACK source file is compiled into four
+versions ({\tt double} / complex, and {\tt int} / {\tt long}).  A proper
+old-style {\tt Makefile} is cumbersome in this case, so these two
+{\tt Makefile}'s have been constructed by brute force.  They ignore
+dependencies, and simply compile everything.  I highly recommend using GNU
+{\tt make} if you wish to modify UMFPACK.
+If you compile UMFPACK and AMD and then later change the {\tt Make.include}
+file or your system-specific configuration file such as {\tt Make.linux},
+then you should type {\tt make purge} and then {\tt make} to recompile.
+Here are the various parameters that you can control in your
+{\tt Make.include} file:
+\begin{itemize}
+\item {\tt CC = } your C compiler, such as {\tt cc}.
+\item {\tt RANLIB = } your system's {\tt ranlib} program, if needed.
+\item {\tt CFLAGS = } optimization flags, such as {\tt -O}.
+Add {\tt -DLP64} if you are compiling in 64-bit mode
+(32 bit {\tt int}'s, 64 bit {\tt long}'s, and 64 bit pointers).
+\item {\tt CONFIG = } configuration settings for the BLAS, memory allocation
+routines, and timing routines.
+\item {\tt LIB = } your libraries, such as {\tt -lm} or {\tt -lblas}.
+\item {\tt RM =} the command to delete a file.
+\item {\tt MV =} the command to rename a file.
+\item {\tt MEX =} the command to compile a MATLAB mexFunction.
+If you are using MATLAB 5, you need to add {\tt -DNBLAS} and
+{\tt -DNUTIL} to this command.  An example is provided in
+the {\tt Make/Make.include} file.
+\item {\tt F77 =} the command to compile a Fortran program (optional).
+\item {\tt F77FLAGS =} the Fortran compiler flags (optional).
+\item {\tt F77LIB =} the Fortran libraries (optional).
+\end{itemize}
+The {\tt CONFIG} string can include combinations of the following;
+most deal with how the BLAS are called:
+\begin{itemize}
+\item {\tt -DNBLAS} if you do not have any BLAS at all.
+\item {\tt -DCBLAS} if you have the C-BLAS \cite{ATLAS}.
+\item {\tt -DNSUNPERF} if you are on Solaris but do not have the Sun
+Performance Library (for the BLAS).
+\item {\tt -DNSCSL} if you on SGI IRIX but do not have the SCSL BLAS library.
+\item {\tt -DLONGBLAS} if your BLAS can take {\tt long} integer input
+arguments.  If not defined, then the {\tt umfpack\_*l\_*} versions of
+UMFPACK that use {\tt long} integers do not call the BLAS.
+This flag is set internally when using the Sun Performance BLAS
+or SGI's SCSL BLAS (both have 64-bit versions of the BLAS).
+\item Options for controlling how C calls the Fortran BLAS:
+{\tt -DBLAS\_BY\_VALUE}, {\tt -DBLAS\_NO\_UNDERSCORE},
+and {\tt -DBLAS\_CHAR\_ARG}.  These are set automatically for Windows,
+Sun Solaris, SGI Irix, Red Hat Linux, Compaq Alpha, and
+AIX (the IBM RS 6000).  They are ignored if you are using
+the C-BLAS interface to the BLAS.
+\item {\tt -DGETRUSAGE} if you have the {\tt getrusage} function.
+\item {\tt -DNUTIL} if you wish to compile the MATLAB-callable
+UMFPACK mexFunction with the {\tt mxMalloc}, {\tt mxRealloc}
+and {\tt mxFree} routines, instead of the undocumented (but
+superior) {\tt utMalloc}, {\tt utRealloc}, and {\tt utFree}
+routines.  The default is to use the {\tt ut*} routines on
+Unix, and the {\tt mx*} routines on Windows.
+\item {\tt -DNPOSIX} if you do not have the POSIX-compliant
+{\tt sysconf} and {\tt times} routines used by
+{\tt umfpack\_tic} and {\tt umfpack\_toc}.
+\item {\tt -DNRECIPROCAL} controls a trade-off between speed and accuracy.
+If defined (or if the pivot value itself is less than $10^{-12}$),
+then the pivot column is divided by the pivot value during numeric
+factorization.  Otherwise, it is multiplied by the reciprocal of the
+pivot, which is faster but can be less accurate.  The default is
+to multiply by the reciprocal unless the pivot value is small.
+This option also modifies how the rows of the matrix $\m{A}$ are
+scaled.  If {\tt -DNRECIPROCAL} is defined (or if any scale factor is
+less than $10^{-12}$), entries in the rows of $\m{A}$ are divided
+by the scale factors.  Otherwise, they are multiplied by the reciprocal.
+When compiling the complex routines with the GNU {\tt gcc} compiler, the
+pivot column is always divided by the pivot entry, because of a
+numerical accuracy issue encountered with {\tt gcc} version 3.2 with a
+few complex matrices on a Pentium 4M (running Linux).  You can still
+use {\tt -DNRECIPROCAL} to control how the scale factors
+for the rows of $\m{A}$ are applied.
+\item {\tt -DNO\_DIVIDE\_BY\_ZERO} controls how UMFPACK treats zeros
+on the diagonal of $\m{U}$, for a singular matrix $\m{A}$.
+If defined, then no division by
+zero is performed (a zero entry on the diagonal of $\m{U}$ is
+treated as if it were equal to one).  By default,
+UMFPACK will divide by zero.
+\item {\tt -DNO\_TIMER} controls whether or not timing routines
+are to be called.  If defined, no timers are used.
+Timers are included by default.
+\end{itemize}
+If a Fortran BLAS package is used you may see compiler warnings.  The BLAS
+routines
+{\tt dgemm}, {\tt dgemv}, {\tt dger}, {\tt dtrsm}, {\tt dtrsv}, {\tt dscal}
+and their corresponding complex versions are used.
+Header files are not provided for the Fortran
+BLAS.  You may safely ignore all of these warnings.
+I highly recommend the recent BLAS by Goto and van de Geijn
+\cite{GotoVandeGeijn02}.  Using this BLAS increased the performance
+of UMFPACK by up to 50\% on a Dell Latitude C840 laptop (2GHz Pentium 4M,
+512K L2 cache, 1GB main memory).  The peak performance of
+{\tt umfpack\_di\_numeric} with Goto and van de Geijn's BLAS is 1.6 Gflops
+on this computer.  In MATLAB, the peak performance of UMFPACK on
+a dense matrix (stored in sparse format) is 900 Mflops, as compared to
+1 Gflop for {\tt x = A}$\backslash${\tt b}
+when {\tt A} is stored as a regular full matrix.
+When you compile your program that uses the C-callable UMFPACK library,
+you need to link your program with both libraries
+({\tt UMFPACK/Lib/libumfpack.a} and {\tt AMD/Lib/libamd.a})
+and you need to tell your compiler to look in the
+directories {\tt UMFPACK/Include} and {\tt AMD/Include} for include
+files.  See {\tt UMFPACK/Demo/Makefile} for an example.
+You do not need to directly include any AMD include files in your
+program, unless you directly call AMD routines.  You only need the
+\begin{verbatim}
+#include "umfpack.h"
+\end{verbatim}
+statement, as described in Section~\ref{Synopsis}.
+If you would like to compile both 32-bit and 64-bit versions of the libraries,
+you will need to do it in two steps.  Modify your {\tt Make/Make.<arch>}
+file, and select the 32-bit option.  Type {\tt make} in the {\tt UMFPACK}
+directory, which creates the {\tt UMFPACK/Lib/libumfpack.a} and
+{\tt AMD/Lib/libamd.a} libraries.  Rename those two files.  Edit your
+{\tt Make/Make.<arch>} and select the 64-bit option.   Type {\tt make purge},
+and then {\tt make}, and you will create the 64-bit libraries.
+You can use the same {\tt umfpack.h} include file for both 32-bit and
+64-bit versions.  Simply link your program with the appropriate 32-bit
+or 64-bit compiled version of the UMFPACK and AMD libraries.
+Type {\tt make hb} in the {\tt UMFPACK/Demo/HB} directory
+to compile and run a C program that reads in and factorizes
+Harwell/Boeing matrices.  Note that this uses a stand-alone Fortran
+program to read in the Fortran-formatted Harwell/Boeing matrices and
+write them to a file which can be read by a C program.
+The {\tt umf\_multicompile.c} file has been added to assist in the
+compilation of UMFPACK in Microsoft Visual Studio, for Windows.
+%-------------------------------------------------------------------------------
+\subsection{Installing the MATLAB interface}
+%-------------------------------------------------------------------------------
+If all you want to do is use the UMFPACK mexFunction in MATLAB, you can skip
+the use of the {\tt make} command described above.  Simply type
+{\tt umfpack\_make} in MATLAB while in the {\tt UMFPACK/MATLAB} directory.
+You can also type {\tt amd\_make} in the {\tt AMD/MATLAB} directory
+to compile the stand-alone AMD mexFunction (this is not required to
+compile the UMFPACK mexFunction).  This works on any computer with MATLAB,
+including Windows.
+You will be prompted to select several configuration options, including
+whether or not to use the BLAS.
+MATLAB 5.3 (or earlier) does not include the BLAS, so you either have to
+compile UMFPACK without the BLAS (UMFPACK will be slow), or modify your
+{\tt <matlab>/bin/mexopts.sh} by adding your BLAS library
+to the {\tt CLIBS} string,
+where {\tt <matlab>} is the directory in which MATLAB is installed.
+If you are using Windows and the {\tt lcc} compiler bundled with
+MATLAB 6.1, then you may need to copy the
+{\tt UMFPACK}$\backslash${\tt MATLAB}$\backslash${\tt lcc\_lib}$\backslash${\tt libmwlapack.lib}
+file into the
+{\tt <matlab>}$\backslash${\tt extern}$\backslash${\tt lib}$\backslash${\tt win32}$\backslash${\tt lcc}$\backslash$
+directory.
+Next, type {\tt mex -setup}
+at the MATLAB prompt, and ask MATLAB to select the {\tt lcc} compiler.
+MATLAB 6.1 has built-in BLAS, but in that version of MATLAB the BLAS
+cannot be accessed by a mexFunction compiled by {\tt lcc} without first copying
+this file to the location listed above.
+If you have MATLAB 6.5 or later, you can probably skip this step.
+%-------------------------------------------------------------------------------
+\subsection{Installing the Fortran interface}
+%-------------------------------------------------------------------------------
+Once the 32-bit C-callable UMFPACK library is compiled, you can also compile
+the Fortran interface, by typing {\tt make fortran}.  This will create
+the {\tt umf4hb} program, test it, and compare the output with the
+file {\tt umf4hb.out} in the distribution.
+If you compiled UMFPACK in 64-bit mode, you need to use {\tt make fortran64}
+instead, which compiles the {\tt umf4hb64} program and compares its output
+with the file {\tt umf4hb64.out}.
+Refer to the comments in the {\tt Demo/umf4\_f77wrapper.c} file
+for more details.
+This interface is {\bf highly} non-portable, since it depends
+on how C and Fortran are interfaced.
+Because of this issue, the interface is included in the {\tt Demo} directory,
+and not as a primary part of the UMFPACK library.  The interface routines are
+not included in the compiled {\tt UMFPACK/Lib/libumfpack.a} library, but left
+as stand-alone compiled files ({\tt umf4\_f77wrapper.o} and
+{\tt umf4\_f77wrapper64.o} in the {\tt Demo} directory).
+You may need to modify the interface routines in the file
+{\tt umf4\_f77wrapper.c} if you are using compilers for which this interface
+has not been tested.
+%-------------------------------------------------------------------------------
+\subsection{Known Issues}
+%-------------------------------------------------------------------------------
+The Microsoft C or C++ compilers on a Pentium badly break the IEEE 754 standard,
+and do not treat NaN's properly.  According to IEEE 754, the expression
+{\tt (x != x)} is supposed to be true if and only if {\tt x} is NaN.  For
+non-compliant compilers in Windows that expression is always false, and another
+test must be used: {\tt (x < x)} is true if and only if {\tt x}
+is NaN.  For compliant compilers, {\tt (x < x)} is always false, for any
+value of {\tt x} (including NaN).
+To cover both cases, UMFPACK when running under Microsoft Windows
+defines the following macro, which is true if and only if {\tt x} is NaN,
+regardless of whether your compiler is compliant or not:
+\begin{verbatim}
+#define SCALAR_IS_NAN(x) (((x) != (x)) || ((x) < (x)))
+\end{verbatim}
+If your compiler breaks this test, then UMFPACK will fail catastrophically
+if it encounters a NaN.  You will not just see NaN's in your output; UMFPACK
+will probably crash with a segmentation fault.  In that case, you might try to
+see if the common (but non-ANSI C) routine {\tt isnan} is available, and modify
+the macro {\tt SCALAR\_IS\_NAN} in {\tt umf\_version.h} accordingly.  The
+simpler (and IEEE 754-compliant) test {\tt (x != x)} is always true with Linux
+on a PC, and on every Unix compiler I have tested.
+Some compilers will complain about the Fortran BLAS being defined implicitly.
+C prototypes for the BLAS are not used, except the C-BLAS.  Some compilers
+will complain about unrecognized {\tt \#pragma}'s.  You may safely ignore
+all of these warnings.
+%-------------------------------------------------------------------------------
+\section{Future work}
+\label{Future}
+%-------------------------------------------------------------------------------
+Here are a few features that are not in UMFPACK Version 4.4, 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.
+Future versions may return more statistics in the {\tt Info} array, and
+they may use more entries in the {\tt Control} array.
+These two arrays will probably become larger, since there are very few
+unused entries.  If they change in size, the constants
+{\tt UMFPACK\_CONTROL} and {\tt UMFPACK\_INFO} defined in {\tt umfpack.h}
+will be changed to reflect their new size.  Your C program should use
+these constants when declaring the size of these two arrays.  Do not
+define them as {\tt Control [20]} and {\tt Info [90]}.
+\item Forward/back solvers for the conventional row or column-form data
+structure for $\m{L}$ and $\m{U}$ (the output of
+{\tt umfpack\_*\_di\_get\_numeric}).  This would enable a separate
+solver that could be used to write a MATLAB mexFunction
+{\tt x = lu\_refine (A, b, L, U, P, Q, R)} that gives MATLAB access
+to the iterative refinement algorithm with sparse backward error
+analysis.  It would also be easier to handle sparse right-hand-sides
+in this data structure, and end up with good asymptotic run-time
+in this case
+(particularly for $\m{Lx}=\m{b}$; see \cite{GilbertPeierls88}).
+\item Complex absolute value computations could be
+based on FDLIBM (see \newline
+http://www.netlib.org/fdlibm),
+using the {\tt hypot(x,y)} routine.
+\item When using iterative refinement, the residual $\m{Ax}-\m{b}$ could be
+returned by {\tt umfpack\_solve}.
+\item The solve routines could handle multiple right-hand sides, and sparse
+right-hand sides.  See {\tt umfpack\_solve} for the MATLAB version
+of this feature.
+\item An option to redirect the error and diagnostic output.
+\item Permutation to block-triangular-form \cite{Duff78a} for the C-callable
+interface.  There are two routines in the ACM Collected
+Algorithms (529 and 575) \cite{Duff81b,Duff78b}
+that could be translated from Fortran
+to C and included in UMFPACK.  This would result in better performance
+for matrices from circuit simulation and
+chemical process engineering.  See {\tt umfpack\_btf.m} for the MATLAB
+version of this feature.  An upcoming package (KLU) will include
+this feature.
+\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.  You can do this at compile-time by
+modifying the definitions of {\tt ALLOCATE}, {\tt FREE}, and
+{\tt REALLOCATE} in the file {\tt umf\_internal.h}.  Be sure to document
+your changes carefully when you change UMFPACK source code.
+\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 A method that takes time proportional to the number of nonzeros in
+$\m{A}$ to compute the symbolic factorization \cite{GilbertNgPeyton94}.
+This would improve the performance of the symmetric and 2-by-2 strategies,
+and the unsymmetric strategy when dense rows are present.
+The current method takes
+time proportional to the number of nonzeros in the upper bound of $\m{U}$.
+The method used in UMFPACK exploits super-columns, however, so this
+bound is rarely reached.
+\item Other basic sparse matrix operations, such as sparse matrix
+multiplication, could be included.
+\item A more complete Fortran interface.
+\item A C++ interface.
+\item A parallel version using MPI.  This would require a large amount
+of effort.
+\end{enumerate}
+%-------------------------------------------------------------------------------
+\newpage
+\section{The primary UMFPACK routines}
+\label{Primary}
+%-------------------------------------------------------------------------------
+The include files are the same for all four 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}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_symbolic.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_numeric}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_numeric.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_solve}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_solve.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_free\_symbolic}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_free_symbolic.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_free\_numeric}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_free_numeric.h via sed
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\newpage
+\section{Alternative routines}
+\label{Alternative}
+%-------------------------------------------------------------------------------
+\subsection{umfpack\_*\_defaults}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_defaults.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_qsymbolic}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_qsymbolic.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_wsolve}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_wsolve.h via sed
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\newpage
+\section{Matrix manipulation routines}
+\label{Manipulate}
+%-------------------------------------------------------------------------------
+\subsection{umfpack\_*\_col\_to\_triplet}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_col_to_triplet.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_triplet\_to\_col}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_triplet_to_col.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_transpose}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_transpose.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_scale}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_scale.h via sed
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\newpage
+\section{Getting the contents of opaque objects}
+\label{Get}
+%-------------------------------------------------------------------------------
+\subsection{umfpack\_*\_get\_lunz}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_get_lunz.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_get\_numeric}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_get_numeric.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_get\_symbolic}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_get_symbolic.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_save\_numeric}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_save_numeric.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_load\_numeric}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_load_numeric.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_save\_symbolic}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_save_symbolic.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_load\_symbolic}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_load_symbolic.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_get\_determinant}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_get_determinant.h via sed
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\newpage
+\section{Reporting routines}
+\label{Report}
+%-------------------------------------------------------------------------------
+\subsection{umfpack\_*\_report\_status}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_status.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_report\_control}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_control.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_report\_info}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_info.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_report\_matrix}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_matrix.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_report\_numeric}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_numeric.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_report\_perm}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_perm.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_report\_symbolic}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_symbolic.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_report\_triplet}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_triplet.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_*\_report\_vector}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_report_vector.h via sed
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\newpage
+\section{Utility routines}
+\label{Utility}
+%-------------------------------------------------------------------------------
+\subsection{umfpack\_timer}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_timer.h via sed
+\end{verbatim}
+}
+\newpage
+\subsection{umfpack\_tic and umfpack\_toc}
+{\footnotesize
+\begin{verbatim}
+INCLUDE umfpack_tictoc.h via sed
+\end{verbatim}
+}
+%-------------------------------------------------------------------------------
+\newpage
+% References
+%-------------------------------------------------------------------------------
+\bibliographystyle{plain}
+\bibliography{UserGuide}
+\end{document}

Mercurial > octave-nkf

comparison liboctave/UMFPACK/UMFPACK/Doc/UserGuide.stex @ 5164:57077d0ddc8e