File: UserGuide.tex

package info (click to toggle)
umfpack 3.2-4
links: PTS
area: main
in suites: woody
size: 3,012 kB
ctags: 836
sloc: ansic: 15,486; makefile: 337; sh: 133
file content (5611 lines) | stat: -rw-r--r-- 216,077 bytes
%-------------------------------------------------------------------------------
% The UserGuide.m4tex file.  Processed into UserGuide.tex via m4.
%-------------------------------------------------------------------------------

% \documentstyle[12pt,fullpage]{article}
\documentclass[12pt]{article}
\usepackage{times}
\newcommand{\m}[1]{{\bf{#1}}}       % for matrices and vectors
\newcommand{\tr}{^{\sf T}}          % transpose
\usepackage[colorlinks,pdftex]{hyperref}
\begin{document}

\author{Timothy A. Davis \\
Dept. of Computer and Information Science and Engineering \\
Univ. of Florida, Gainesville, FL}
\title{UMFPACK Version 3.2 User Guide}
\date{January 1, 2002}
\maketitle

%-------------------------------------------------------------------------------
\begin{abstract}
    UMFPACK is a set of routines for solving unsymmetric sparse linear
    systems, \newline $\m{Ax}=\m{b}$, using the Unsymmetric MultiFrontal method
    and direct sparse LU factorization.  It is written in ANSI/ISO C, with a
    MATLAB (Version 6.0 or 6.1) interface.  UMFPACK relies on the Level-3 Basic
    Linear Algebra Subprograms (dense matrix multiply) for its performance.
\end{abstract}
%-------------------------------------------------------------------------------

Copyright\copyright 2002 by Timothy A. Davis, University of Florida,
davis@cise.ufl.edu.  All Rights Reserved.

{\bf UMFPACK License:}

    Your use or distribution of UMFPACK or any derivative code implies that
    you agree to this License.

    THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY
    EXPRESSED OR IMPLIED.  ANY USE IS AT YOUR OWN RISK.

    Permission is hereby granted to use or copy this program, provided
    that the Copyright, this License, and the Availability of the original
    version is retained on all copies.  User documentation of any code that
    uses this code or any derivative code must cite the Copyright, this
    License, the Availability note, and "Used by permission."  If this
    code or any derivative code is accessible from within MATLAB, then
    typing "help umfpack" must cite the Copyright, and "type umfpack"
    must also cite this License and the Availability note.  Permission to
    modify the code and to distribute modified code is granted, provided
    the Copyright, this License, and the Availability note are retained,
    and a notice that the code was modified is included.  This software
    was developed with support from the National Science Foundation, and
    is provided to you free of charge.

{\bf Acknowledgments:}

    This work was supported by the National Science Foundation, under
    grants DMS-9504974 and DMS-9803599.

%-------------------------------------------------------------------------------
\newpage
%-------------------------------------------------------------------------------

\tableofcontents

%-------------------------------------------------------------------------------
\newpage
\section{Overview}
%-------------------------------------------------------------------------------

UMFPACK Version 3.2 is a set of routines for solving systems of linear
equations, $\m{Ax}=\m{b}$, when $\m{A}$ is sparse and unsymmetric.  It is based
on the Unsymmetric MultiFrontal method \cite{DavisDuff97,DavisDuff99},
which factorizes $\m{PAQ}$ into the product $\m{LU}$, where $\m{L}$ and $\m{U}$
are lower and upper triangular, respectively, and $\m{P}$ are $\m{Q}$ are
permutation matrices.  Both $\m{P}$ and $\m{Q}$ are chosen to reduce fill-in
(new nonzeros in $\m{L}$ and $\m{U}$ that are not present in $\m{A}$).  The
permutation $\m{P}$ has the dual role of reducing fill-in and maintaining
numerical accuracy (via relaxed partial pivoting and row interchanges).

UMFPACK first finds a column pre-ordering that reduces fill-in, without regard
to numerical values, with a modified version of COLAMD
\cite{DavisGilbertLarimoreNg00_algo,DavisGilbertLarimoreNg00,Larimore98}.
The method finds a symmetric permutation $\m{Q}$ of the matrix $\m{A}\tr\m{A}$
(without forming $\m{A}\tr\m{A}$ explicitly).  This is a good choice for
$\m{Q}$, since the Cholesky factors of $\m{(AQ)\tr(AQ)}$ are an upper bound (in
terms of nonzero pattern) of the factor $\m{U}$ for the unsymmetric LU
factorization ($\m{PAQ}=\m{LU}$) regardless of the choice of $\m{P}$
\cite{GeorgeNg85,GeorgeNg87,GilbertNg93}.

Next, the method breaks the factorization of the matrix $\m{A}$ down into a
sequence of dense rectangular frontal matrices.  The frontal matrices are
related to each other by a supernodal column elimination tree, in which each
node in the tree represents one frontal matrix.  This analysis phase also
determines upper bounds on the memory usage, the floating-point operation count,
and the number of nonzeros in the LU factors.

UMFPACK factorizes each {\em chain} of frontal matrices in a single working
array, similar to how the unifrontal method \cite{dusc:96} factorizes the whole
matrix.  A chain of frontal matrices is a sequence of fronts where the parent
of front $i$ is $i$+1 in the supernodal column elimination tree.  UMFPACK is an
outer-product based, right-looking method.  At the $k$-th step of Gaussian
elimination, it represents the updated submatrix $\m{A}_k$ as an implicit
summation of a set of dense submatrices (referred to as {\em elements},
borrowing a phrase from finite-element methods) that arise when the frontal
matrices are factorized and their pivot rows and columns eliminated.

Each frontal matrix represents the elimination of one or more columns;
each column of $\m{A}$ will be eliminated in a specific frontal matrix,
and which frontal matrix will be used for each column is determined by
the pre-analysis phase.  The pre-analysis phase also determines the worst-case
size of each frontal matrix so that they can hold any candidate pivot column
and any candidate pivot row.  From the perspective of the analysis phase, any
candidate pivot column in the frontal matrix is identical (in terms of nonzero
pattern), and so is any row.  However, the numerical factorization phase has
more information than the analysis phase.  It uses this information to reorder
the columns within each frontal matrix to reduce fill-in.  Similarly, since
the number of nonzeros in each row and column are maintained (more precisely,
COLMMD-style approximate degrees \cite{GilbertMolerSchreiber}), a pivot row can
be selected based on sparsity-preserving criteria (low degree) as well as
numerical considerations (relaxed threshold partial pivoting).  This information
about row and column degrees is not available to left-looking methods such as
{\tt SuperLU} \cite{SuperLU99} or MATLAB's {\tt LU}
\cite{GilbertMolerSchreiber,GilbertPeierls88}.

More details of the method, including experimental results, are
described in \cite{Davis02,Davis02_algo}, available at
www.cise.ufl.edu/tech-reports.

%-------------------------------------------------------------------------------
\section{Availability}
%-------------------------------------------------------------------------------

UMFPACK Version 3.2 is available at www.cise.ufl.edu/research/sparse,
and has been submitted as a collected algorithm of the ACM
\cite{Davis02,Davis02_algo}.
It makes use of a modified version of COLAMD V2.0 by Timothy A.~Davis, Stefan
Larimore, John Gilbert, and Esmond Ng.  The original COLAMD V2.0 is available in
MATLAB V6.0 (or later), and at www.cise.ufl.edu/research/sparse.
These codes are also available in Netlib \cite{netlib} at
www.netlib.org.
Prior versions of UMFPACK, co-authored with Iain Duff, are available at
www.cise.ufl.edu/research/sparse and as MA38 (functionally
equivalent to Version 2.2.1) in the Harwell Subroutine Library.

%-------------------------------------------------------------------------------
\section{Using UMFPACK in MATLAB}
%-------------------------------------------------------------------------------

The easiest way to use UMFPACK is within MATLAB.  This discussion assumes that
you have MATLAB Version 6.0 or later (which includes the BLAS, and the
{\tt colamd} ordering routine).  To compile the UMFPACK mexFunction, you can
either type {\tt make umfpack} in the Unix system shell, or type
{\tt umfpack\_make} in MATLAB (which should work on any system, including
Windows).  See Section~\ref{Install} for more details on how to install UMFPACK.
Once installed, the UMFPACK mexFunction can analyze, factor, and solve linear
systems.  Table~\ref{matlab} summarizes some of the more common uses
of UMFPACK within MATLAB.

\begin{table}
\caption{Using UMFPACK's MATLAB interface}
\label{matlab}
\vspace{0.1in}
{\footnotesize
\begin{tabular}{l|l|l}
\hline
Function & Using UMFPACK & MATLAB 6.0 equivalent \\
\hline
 & & \\
\begin{minipage}[t]{1.5in}
Solve $\m{Ax}=\m{b}$.
\end{minipage}
&
\begin{minipage}[t]{2.2in}
\begin{verbatim}
x = umfpack (A,'\',b) ;
\end{verbatim}
\end{minipage}
&
\begin{minipage}[t]{2.2in}
\begin{verbatim}
x = A \ b ;
\end{verbatim}
\end{minipage}
 \\
 & & \\
\hline
 & & \\
\begin{minipage}[t]{1.5in}
Solve $\m{Ax}=\m{b}$ using a different column pre-ordering.
\end{minipage}
&
\begin{minipage}[t]{2.2in}
\begin{verbatim}

S = spones (A) ;
Q = symamd (S+S') ;
x = umfpack (A,Q,'\',b) ;
\end{verbatim}
\end{minipage}
&
\begin{minipage}[t]{2.2in}
\begin{verbatim}
spparms ('autommd',0) ; 
S = spones (A) ;
Q = symamd (S+S') ;
x = A (:,Q) \ b ;
x (Q) = x ;
spparms ('autommd',1) ;
\end{verbatim}
\end{minipage}
 \\
 & & \\
\hline
 & & \\
\begin{minipage}[t]{1.5in}
Solve $\m{A}\tr\m{x}\tr = \m{b}\tr$.
\end{minipage}
&
\begin{minipage}[t]{2.2in}
\begin{verbatim}
x = umfpack (b,'/',A) ;
\end{verbatim}
\end{minipage}
&
\begin{minipage}[t]{2.2in}
\begin{verbatim}
x = b / A ;
\end{verbatim}
\end{minipage}
 \\
 & & \\
\hline
 & & \\
\begin{minipage}[t]{1.5in}
Factorize $\m{A}$, then solve $\m{Ax}=\m{b}$.
\end{minipage}
&
\begin{minipage}[t]{2.2in}
\begin{verbatim}

[L,U,P,Q] = umfpack (A) ;
x = U \ (L \ (b (P))) ;
x (Q) = x ;
\end{verbatim}
\end{minipage}
&
\begin{minipage}[t]{2.2in}
\begin{verbatim}
Q = colamd (A) ;
[L,U,P] = lu (A (:,Q)) ;
x = U \ (L \ (P*b)) ;
x (Q) = x ;
\end{verbatim}
\end{minipage}
 \\
 & & \\
\hline
\end{tabular}
}
\end{table}

UMFPACK requires {\tt A} to be square, sparse,
nonsingular and not complex, and it requires {\tt b} to be a dense column vector
(and not complex).  This is more restrictive than what you can do with MATLAB's
backslash or {\tt LU}.  Future releases of UMFPACK may remove some of these
restrictions.

MATLAB's {\tt [L,U,P] = lu (A)} returns a lower triangular {\tt L}, an upper
triangular {\tt U}, and a permutation matrix {\tt P} such that {\tt P*A} is
equal to {\tt L*U}.  UMFPACK behaves differently; it returns {\tt P} and {\tt Q}
such that {\tt A (P,Q)} is equal to {\tt L*U}, where {\tt P} and {\tt Q} are
permutation vectors.  If you prefer permutation matrices, use the following
MATLAB code:

{\footnotesize
\begin{verbatim}
    [L,U,P,Q] = umfpack (A) ;
    n = size (A,1) ;
    I = speye (n) ;
    P = I (P,:) ;
    Q = I (:,Q) ;
\end{verbatim}
}

Now {\tt P*A*Q} is equal to {\tt L*U}.  Note that
        {\tt x = umfpack (A, '}$\backslash${\tt ', b)}
requires that {\tt b} be a dense column vector.
If you wish to use the LU factors from UMFPACK to solve a
linear system, $\m{Ax}=\m{b}$ where $\m{b}$ is a either a dense or sparse
matrix with more than one column, do this:

{\footnotesize
\begin{verbatim}
    [L,U,P,Q] = umfpack (A) ;
    x = U \ (L \ (b (P,:))) ;
    x (Q,:) = x ;
\end{verbatim}
}

The above two examples do not make use of the iterative refinement
that is built into
        {\tt x = umfpack (A, '}$\backslash${\tt ', b)}
however.

Since the numeric factorization refines its column pre-ordering, the {\tt Q}
in {\tt [L,U,P,Q] = umfpack (A)} and {\tt [Q,F,C] = umfpack (A, 'symbolic')}
will in general be different.

There are more options; you can provide your own column pre-ordering (in which
case UMFPACK does not call COLAMD), you can modify other control settings
(similar to the {\tt spparms} in MATLAB), and you get various statistics on
the analysis, factorization, and solution of the linear system.  Type
{\tt help umfpack\_details} and {\tt help umfpack\_report} in MATLAB for more
information.  Two demo m-files are provided.   Just type {\tt umfpack\_simple}
and {\tt umfpack\_demo} to run them.  They roughly correspond to the C programs
{\tt umfpack\_simple.c} and {\tt umfpack\_demo.c}.  You may want to type
{\tt more on} before running the {\tt umfpack\_simple} demo since it generates
lots of output.

A simple M-file ({\tt umfpack\_btf}) is provided that first permutes the matrix
to upper block triangular form, using MATLAB's {\tt dmperm} routine.  It
solves $\m{Ax}=\m{b}$; the LU factors are not returned.  Its usage is simple:
        {\tt x = umfpack\_btf (A, b)}.
Type {\tt help umfpack\_btf} for more options.

One issue you may encounter is how UMFPACK allocates its memory when being used
in a mexFunction.  One part of its working space is of variable size.   The
symbolic analysis phase determines an upper bound on the size of this memory,
but not all of this memory will typically be used in the numerical
factorization.  UMFPACK tries to allocate a decent amount of working space
(70\% of the upper bound, by default), with some elbow room so
that it can run more efficiently.  If this fails, it reduces its request and
uses less memory.  However, {\tt mxMalloc} aborts the {\tt umfpack} mexFunction
if it fails, so this strategy doesn't work in MATLAB.  The strategy works fine
when {\tt malloc} is used instead.  If you run out of memory in MATLAB, try
reducing {\tt Control(7)} to be less than 0.70, and try again.
Alternatively, set {\tt Control(7)} to 1.0 or 1.05 to avoid all reallocations
of memory.
Type {\tt help umfpack\_details} and {\tt umfpack\_report} for more
information, and refer to the {\tt Control [UMFPACK\_ALLOC\_INIT]} parameter
described in {\tt umfpack\_numeric} in Section~\ref{Primary}, below.

There is a solution to this problem, but it relies on undocumented internal
routines.  See the {\tt -DMATHWORKS} option in {\tt umf\_config.h} in
Section~\ref{Config} for details.

Memory allocation on a PC is notoriously bad, so I recommend setting
{\tt Control(7)} to a non-default value of 1.0 or even 1.05.
This will avoid most reallocations of memory.

%-------------------------------------------------------------------------------
\section{Using UMFPACK in a C program}
%-------------------------------------------------------------------------------

The C-callable UMFPACK library consists of 24 user-callable routines and one
include file.  Twenty-three of the routines come in dual versions, with
different sizes of integers.  All user-callable routine names begin with
{\tt umfpack\_} or {\tt umfpack\_l\_}; other routine names beginning with
{\tt umf\_} or {\tt umfl\_} are internal to the package, and should not be
called by the user.  The include file {\tt umfpack.h}, listed in
Section~\ref{Include}, must be included in any C program that uses UMFPACK.

%-------------------------------------------------------------------------------
\subsection{The size of an integer}
%-------------------------------------------------------------------------------

There are two versions of each user-callable routine (except for one routine).
The routine names starting with just {\tt umfpack\_} use {\tt int} integer
arguments; those starting with {\tt umfpack\_l\_} use {\tt long} integer
arguments.  If you compile UMFPACK in the standard ILP32 mode (32-bit
{\tt int}'s, {\tt long}'s, an pointers) then the versions are essentially
identical.  You will be able to solve problems using up to 4GB of memory.
If you compile UMFPACK in the standard LP64 mode, the size of an {\tt int}
remains 32-bits, but the size of a {\tt long} and a pointer both get promoted
to 64-bits.  In the LP64 mode, the {\tt umfpack\_l\_*} routines can solve huge
problems (not limited to 4GB), limited of course by the amount of available
memory.  The only drawback to the 64-bit mode is that few BLAS libraries
support 64-bit integers.  This limits the performance you will obtain.

Both versions are discussed below.  Use only one version for any one problem;
do not attempt to use one version to analyze the matrix and another version to
factorize the matrix, for example.

%-------------------------------------------------------------------------------
\subsection{Primary routines, and a simple example}
%-------------------------------------------------------------------------------

Five primary UMFPACK routines are required to solve $\m{Ax}=\m{b}$.  They are
fully described in Section~\ref{Primary}:

\begin{itemize}
\item {\tt umfpack\_symbolic}, {\tt umfpack\_l\_symbolic}:

    Pre-orders the columns of $\m{A}$ to reduce fill-in, based on its sparsity
    pattern only, finds the supernodal column elimination tree, and post-orders
    the tree.  Returns an opaque {\tt Symbolic} object as a {\tt void *}
    pointer.  The object contains the symbolic analysis and is needed for the
    numerical factorization.  This routine requires only $O(|\m{A}|)$ space,
    where $|\m{A}|$ is the number of nonzero entries in the matrix.  It computes
    upper bounds on the nonzeros in $\m{L}$ and $\m{U}$, the floating-point
    operations required, and the memory usage of {\tt umfpack\_numeric}.  The
    {\tt Symbolic} object is small; it contains just the column pre-ordering,
    the supernodal column elimination tree, and information about each frontal
    matrix, and is no larger than about $6n$ integers (where $\m{A}$ is
    $n$-by-$n$).  The matrix must be structurally non-singular (more precisely,
    each row and column must have at least one entry).

\item {\tt umfpack\_numeric}, {\tt umfpack\_l\_numeric}:

    Numerically factorizes a sparse matrix into $\m{PAQ}=\m{LU}$.  Requires the
    symbolic ordering and analysis computed by {\tt umfpack\_symbolic} or
    {\tt umfpack\_qsymbolic}.  Returns an opaque {\tt Numeric} object as a
    {\tt void *} pointer.  The object contains the numerical factorization and
    is used by {\tt umfpack\_solve}.  You can factorize a new matrix with a
    different values (but identical pattern) as the matrix analyzed by
    {\tt umfpack\_symbolic} or {\tt umfpack\_qsymbolic} by re-using the
    {\tt Symbolic} object (this feature is available when using UMFPACK in a
    C program, but not in MATLAB).  The matrix must be non-singular.

\item {\tt umfpack\_solve}, {\tt umfpack\_l\_solve}:

    Solves a sparse linear system ($\m{Ax}=\m{b}$, $\m{A}\tr\m{x}=\m{b}$, or
    systems involving just $\m{L}$ or $\m{U}$), using the numeric factorization
    computed by {\tt umfpack\_numeric}.  Iterative refinement with sparse
    backward error \cite{ardd:89} is used by default.

\item {\tt umfpack\_free\_symbolic}, {\tt umfpack\_l\_free\_symbolic}:

    Frees the {\tt Symbolic} object created by {\tt umfpack\_symbolic} or
    {\tt umfpack\_qsymbolic}.

\item {\tt umfpack\_free\_numeric}, {\tt umfpack\_l\_free\_numeric}:

    Frees the {\tt Numeric} object created by {\tt umfpack\_numeric}.

\end{itemize}

Be careful not to free a {\tt Symbolic} object with
{\tt umfpack\_free\_numeric}.  Nor should you attempt to free a {\tt Numeric}
object with {\tt umfpack\_free\_symbolic}.  Failure to free these objects will
lead to memory leaks.

The matrix $\m{A}$ is represented in compressed column form, which is
identical to the sparse matrix representation used by MATLAB.  It consists
of three arrays, where the matrix is {\tt n}-by{\tt n}, with {\tt nz} entries.
For the {\tt int} version of UMFPACK:

{\footnotesize
\begin{verbatim}
     int Ap [n+1] ;
     int Ai [nz] ;
     double Ax [nz] ;
\end{verbatim}
}

For the {\tt long} version of UMFPACK:

{\footnotesize
\begin{verbatim}
     long Ap [n+1] ;
     long Ai [nz] ;
     double Ax [nz] ;
\end{verbatim}
}

All nonzeros are entries, but an entry may be numerically zero.  The row indices
of entries in column {\tt j} are stored in
    {\tt Ai[Ap[j]} ... {\tt Ap[j+1]-1]}.
The corresponding numerical values are stored in
    {\tt Ax[Ap[j]} ... {\tt Ap[j+1]-1]}.

No duplicate row indices may be present, and the row indices in any given
column must be sorted in ascending order.  The first entry {\tt Ap [0]} must be
zero.  The total number of entries in the matrix is thus {\tt nz = Ap [n]}.
Except for the fact that extra zero entries can be included, there is thus a
unique compressed column representation of any given matrix $\m{A}$.

Here is a simple main program, {\tt umfpack\_simple.c}, that illustrates the
basic usage of UMFPACK.

{\footnotesize
\begin{verbatim}
#include <stdio.h>
#include "umfpack.h"

int    n = 5 ;
int    Ap [ ] = {0, 2, 5, 9, 10, 12} ;
int    Ai [ ] = { 0,  1,  0,   2,  4,  1,  2,  3,   4,  2,  1,  4} ;
double Ax [ ] = {2., 3., 3., -1., 4., 4., -3., 1., 2., 2., 6., 1.} ;
double b [ ] = {8., 45., -3., 3., 19.} ;
double x [5] ;

int main (int argc, char **argv)
{
    double *Control = (double *) NULL, *Info = (double *) NULL ;
    int i ;
    void *Symbolic, *Numeric ;
    (void) umfpack_symbolic (n, Ap, Ai, &Symbolic, Control, Info) ;
    (void) umfpack_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;
    umfpack_free_symbolic (&Symbolic) ;
    (void) umfpack_solve ("Ax=b", Ap, Ai, Ax, x, b, Numeric, Control, Info) ;
    umfpack_free_numeric (&Numeric) ;
    for (i = 0 ; i < n ; i++) printf ("x [%d] = %g\n", i, x [i]) ;
    return (0) ;
}
\end{verbatim}
}

It solves the same linear system as the {\tt umfpack\_simple.m} MATLAB m-file.
The {\tt Ap}, {\tt Ai}, and {\tt Ax} arrays represent the matrix
\[
\m{A} = \left[
\begin{array}{rrrrr}
 2 &  3 &  0 &  0 &  0 \\
 3 &  0 &  4 &  0 &  6 \\
 0 & -1 & -3 &  2 &  0 \\
 0 &  0 &  1 &  0 &  0 \\
 0 &  4 &  2 &  0 &  1 \\
\end{array}
\right].
\]
and the solution is $\m{x} = [1 \, 2 \, 3 \, 4 \, 5]\tr$.  The program uses
default control settings and does not return any statistics about the ordering,
factorization, or solution ({\tt Control} and {\tt Info} are both
{\tt (double *) NULL}).

%-------------------------------------------------------------------------------
\subsection{Alternative routines}
%-------------------------------------------------------------------------------

Three alternative routines are provided that modify UMFPACK's default
behavior.  They are fully described in Section~\ref{Alternative}:

\begin{itemize}
\item {\tt umfpack\_defaults}, {\tt umfpack\_l\_defaults}:

    Sets the default control parameters in the {\tt Control} array.  These can
    then be modified as desired before passing the array to the other UMFPACK
    routines.  Control parameters are fully described in Section~\ref{defaults}.
    One particular parameter deserves special notice.
    UMFPACK uses relaxed partial pivoting, where a candidate pivot entry is
    numerically acceptable if its magnitude is greater than or equal to a
    tolerance parameter times the magnitude of the largest entry in the same
    column.  The parameter {\tt Info [UMFPACK\_PIVOT\_TOLERANCE]} has a default
    value of 0.1.  This may be too small for some matrices, particularly for
    ill-conditioned or poorly scaled ones.  With the default pivot tolerance
    and default iterative refinement,
        {\tt x = umfpack (A, '}$\backslash${\tt ', b)}
    is just as accurate as
        {\tt x = A}$\backslash${\tt b}
    in MATLAB for nearly all matrices.

\item {\tt umfpack\_qsymbolic}, {\tt umfpack\_l\_qsymbolic}:

    An alternative to {\tt umfpack\_symbolic}.  Allows the user to specify his
    or her own column pre-ordering, rather than using the default COLAMD
    pre-ordering.

\item {\tt umfpack\_wsolve}, {\tt umfpack\_l\_wsolve}:

    An alternative to {\tt umfpack\_solve} which does not dynamically allocate
    any memory.  Requires the user to pass several additional size-$n$ work
    arrays.

\end{itemize}

%-------------------------------------------------------------------------------
\subsection{Matrix manipulation routines}
%-------------------------------------------------------------------------------

The compressed column data structure is compact, and simplifies the UMFPACK
routines that operate on the sparse matrix $\m{A}$.  However, it can be
inconvenient for the user to generate.  Section~\ref{Manipulate} presents the
details of routines for manipulating sparse matrices in {\em triplet} form,
compressed column form, and compressed row form (the transpose of the
compressed column form).  The triplet form of a matrix consists of three arrays.
For the {\tt int} version of UMFPACK:

{\footnotesize
\begin{verbatim}
     int Ti [nz] ;
     int Tj [nz] ;
     double Tx [nz] ;
\end{verbatim}
}

For the {\tt long} version:

{\footnotesize
\begin{verbatim}
     long Ti [nz] ;
     long Tj [nz] ;
     double Tx [nz] ;
\end{verbatim}
}

The {\tt k}-th triplet is $(i,j,a_{ij})$, where $i =$ {\tt Ti [k]},
$j =$ {\tt Tj [k]}, and $a_{ij} =$ {\tt Tx [k]}.  The triplets can be in any
order in the {\tt Ti}, {\tt Tj}, and {\tt Tx} arrays, and duplicate entries may
exist.  Any duplicate entries are summed when the triplet form is converted to
compressed column form.  This is a convenient way to create a matrix arising in
finite-element methods, for example.

Three routines are provided for manipulating sparse matrices:

\begin{itemize}
\item {\tt umfpack\_triplet\_to\_col}, {\tt umfpack\_l\_triplet\_to\_col}:

    Converts a triplet form of a matrix to compressed column form (ready for
    input to {\tt umfpack\_symbolic}, {\tt umfpack\_qsymbolic}, and
    {\tt umfpack\_numeric}).  Identical to {\tt A = spconvert (i,j,x)} in
    MATLAB, except that zero entries are not removed, so that the pattern of
    entries in the compressed column form of $\m{A}$ are fully under user
    control.  This is important if you want to factorize a new matrix with the
    {\tt Symbolic} object from a prior matrix with the same pattern as the new
    one.  MATLAB never stores explicitly zero entries.

\item {\tt umfpack\_col\_to\_triplet}, {\tt umfpack\_l\_col\_to\_triplet}:

    The opposite of {\tt umfpack\_triplet\_to\_col}.  Identical to
    {\tt [i,j,x] = find (A)} in MATLAB, except that numerically zero entries
    may be included.

\item {\tt umfpack\_transpose}, {\tt umfpack\_l\_transpose}:

    Transposes and optionally permutes a column form matrix \cite{Gustavson78}.
    Identical to {\tt B = A (P,Q)'} in MATLAB, except for the presence of
    numerically zero entries.

\end{itemize}

It is quite easy to add matrices in triplet form, transpose them, and permute
them.  See the discussion in {\tt umfpack\_triplet\_to\_col} in
Section~\ref{Manipulate} for more details.  All of the matrix manipulation
routines can correctly operate on singular matrices.

%-------------------------------------------------------------------------------
\subsection{Getting the contents of opaque objects}
%-------------------------------------------------------------------------------

There are cases where the user would like to do more with the LU factorization
of a matrix than solve a linear system.  The opaque {\tt Symbolic} and
{\tt Numeric} objects are just that - opaque.  In addition, the LU factors
are stored in the {\tt Numeric} object in a compact way that does not store
all of the row and column indices \cite{GeorgeLiu}.  These objects may not be
dereferenced by the user, and even if they were, it would be difficult for a
typical user to understand how the LU factors are stored.  Thus, three routines
are provided for copying their contents into user-provided arrays using simpler
data structures.  They are fully described in Section~\ref{Get}:

\begin{itemize}
\item {\tt umfpack\_get\_lunz}, {\tt umfpack\_l\_get\_lunz}:

    Returns the number of nonzeros in $\m{L}$ and $\m{U}$.

\item {\tt umfpack\_get\_numeric}, {\tt umfpack\_l\_get\_numeric}:

    Copies $\m{L}$, $\m{U}$, $\m{P}$, and $\m{Q}$ from the {\tt Numeric} object
    into arrays provided by the user.  The matrix $\m{L}$ is returned in
    compressed row form (with the column indices in each row sorted in ascending
    order).  The matrix $\m{U}$ is returned in compressed column form (also with
    sorted columns).  There are no explicit zero entries in $\m{L}$ and $\m{U}$,
    but such entries may exist in the {\tt Numeric} object.  The permutations
    $\m{P}$ and $\m{Q}$ are represented as permutation vectors, where
    {\tt P [k] = i} means that row {\tt i} of the original matrix is the
    {\tt k}-th pivot row (or the {\tt k}-th row of $\m{PAQ}$), and where
    {\tt Q [k] = j} means that column {\tt j} of the original matrix is the
    {\tt k}-th pivot column.  This is identical to how MATLAB uses permutation
    vectors.

\item {\tt umfpack\_get\_symbolic}, {\tt umfpack\_l\_get\_symbolic}:

    Copies the contents of the {\tt Symbolic} object (the initial column
    preordering, and supernodal column elimination tree, and information
    about each frontal matrix) into arrays provided by the user.

\end{itemize}

UMFPACK itself does not make use of the output of the {\tt umfpack\_get\_*}
routines; they are provided solely for returning the contents of the opaque
{\tt Symbolic} and {\tt Numeric} objects to the user.

%-------------------------------------------------------------------------------
\subsection{Reporting routines}
%-------------------------------------------------------------------------------

None of the UMFPACK routines discussed so far prints anything, even when an
error occurs.  UMFPACK provides you with nine routines for printing the input
and output arguments (including the {\tt Control} settings and {\tt Info}
statistics) of UMFPACK routines discussed above.  They are fully described in
Section~\ref{Report}:

\begin{itemize}
\item {\tt umfpack\_report\_status}, {\tt umfpack\_l\_report\_status}:

    Prints the status (return value) of other {\tt umfpack\_*} routines.

\item {\tt umfpack\_report\_info}, {\tt umfpack\_l\_report\_info}:

    Prints the statistics returned in the {\tt Info} array by
    {\tt umfpack\_*symbolic}, {\tt umfpack\_numeric}, and {\tt umfpack\_*solve}.

\item {\tt umfpack\_report\_control}, {\tt umfpack\_l\_report\_control}:

    Prints the {\tt Control} settings.

\item {\tt umfpack\_report\_matrix}, {\tt umfpack\_l\_report\_matrix}:

    Verifies and prints a compressed column-form or compressed row-form sparse
    matrix.

\item {\tt umfpack\_report\_triplet}, {\tt umfpack\_l\_report\_triplet}:

    Verifies and prints a matrix in triplet form.

\item {\tt umfpack\_report\_symbolic}, {\tt umfpack\_l\_report\_symbolic}:

    Verifies and prints a {\tt Symbolic} object.

\item {\tt umfpack\_report\_numeric}, {\tt umfpack\_l\_report\_numeric}:

    Verifies and prints a {\tt Numeric} object.

\item {\tt umfpack\_report\_perm}, {\tt umfpack\_l\_report\_perm}:

    Verifies and prints a permutation vector.

\item {\tt umfpack\_report\_vector}, {\tt umfpack\_l\_report\_vector}:

    Verifies and prints a real vector.

\end{itemize}

The {\tt umfpack\_report\_*} routines behave slightly differently when compiled
into the C-callable UMFPACK library than when used in the MATLAB mexFunction.
MATLAB stores its sparse matrices using the same compressed column data
structure discussed above, where row and column indices are in the range {\tt 0}
to {\tt n-1}, but it prints them as if they are in the range {\tt 1} to {\tt n}.
The UMFPACK mexFunction behaves the same way.

You can control how much the {\tt umfpack\_report\_*} routines print by
modifying the {\tt Control [UMFPACK\_PRL]} parameter.  Its default value is
{\tt UMFPACK\_DEFAULT\_PRL} which is equal to 1.  Here is a summary of how
the routines use this print level parameter:

\begin{itemize}
\item {\tt umfpack\_report\_status}, {\tt umfpack\_l\_report\_status}:

    No output if the print level is 0 or less, even when an error occurs.
    If 1, then error messages are printed, and nothing is printed if
    the status is {\tt UMFPACK\_OK}.  If 2 or more, then the status is always
    printed.  If 4 or more, then the UMFPACK Copyright is printed.
    If 6 or more, then the UMFPACK License is printed.  See also the first page
    of this User Guide for the Copyright and License.

\item {\tt umfpack\_report\_control}, {\tt umfpack\_l\_report\_control}:

    No output if the print level is 1 or less.  If 2 or more, all of
    {\tt Control} is printed.

\item {\tt umfpack\_report\_info}, {\tt umfpack\_l\_report\_info}:

    No output if the print level is 1 or less.  If 2 or more, all of
    {\tt Info} is printed.

\item all other {\tt umfpack\_report\_*} routines:

    If the print level is 2 or less, then these routines return silently without
    checking their inputs.  If 3 or more, the inputs are fully verified and a
    short status summary is printed.  If 4, then the first few entries of the
    input arguments are printed.  If 5, then all of the input arguments are
    printed.

\end{itemize}


%-------------------------------------------------------------------------------
\subsection{Utility routines}
%-------------------------------------------------------------------------------

UMFPACK includes a routine that returns the time used by the process,
{\tt umfpack\_timer}.  The routine uses either {\tt getrusage} (which is
preferred), or the ANSI C {\tt clock} routine if that is not available.
It is fully described in Section~\ref{Utility}.  It is the only routine
that is identical in both {\tt int} and {\tt long} versions (there is no
{\tt umfpack\_l\_timer} routine).

%-------------------------------------------------------------------------------
\subsection{Control parameters}
%-------------------------------------------------------------------------------

UMFPACK uses an optional {\tt double} array, {\tt Control}, to modify
its control parameters.  These may be modified by the user
(see {\tt umfpack\_defaults} and {\tt umfpack\_l\_defaults}).  Each
user-callable routine includes a complete description of how each control
setting modifies its behavior.  Table~\ref{control} summarizes the entire
contents of the {\tt Control} array.  Future versions may make use of
additional entries in the {\tt Control} array.
Note that ANSI C uses 0-based indexing, while MATLAB user's 1-based
indexing.  Thus, {\tt Control(1)} in MATLAB is the same as
{\tt Control[0]} or {\tt Control[UMFPACK\_PRL]} in ANSI C.

\begin{table}
\caption{UMFPACK Control parameters}
\label{control}
{\footnotesize
\begin{tabular}{llll}
\hline
MATLAB & ANSI C & default & description \\
\hline
\multicolumn{4}{l}{Used by reporting routines:} \\
{\tt Control(1)} & {\tt Control[UMFPACK\_PRL]} & 1 & printing level \\
\hline
% /* used in UMFPACK_*symbolic only: */
\multicolumn{4}{l}{Used by {\tt umfpack\_*symbolic:}} \\
{\tt Control(2)} & {\tt Control[UMFPACK\_DENSE\_ROW]} & 0.2 & dense row threshold \\
{\tt Control(3)} & {\tt Control[UMFPACK\_DENSE\_COL]} & 0.2 & dense column threshold \\
\hline
% /* used in UMFPACK_numeric only: */
\multicolumn{4}{l}{Used by {\tt umfpack\_*numeric:}} \\
{\tt Control(4)} & {\tt Control[UMFPACK\_PIVOT\_TOLERANCE]} & 0.1 & partial pivoting tolerance \\
{\tt Control(5)} & {\tt Control[UMFPACK\_BLOCK\_SIZE]} & 24 & BLAS block size \\
{\tt Control(6)} & {\tt Control[UMFPACK\_RELAXED\_AMALGAMATION]} & 0.25 & amalgamation \\
{\tt Control(7)} & {\tt Control[UMFPACK\_ALLOC\_INIT]} & 0.7 & initial memory allocation  \\
{\tt Control(13)} & {\tt Control[UMFPACK\_PIVOT\_OPTION]} & 0 & symmetric pivot preference \\
{\tt Control(14)} & {\tt Control[UMFPACK\_RELAXED2\_AMALGAMATION]} & 0.1 & amalgamation \\
{\tt Control(15)} & {\tt Control[UMFPACK\_RELAXED3\_AMALGAMATION]} & 0.125 & amalgamation \\
\hline
% /* used in UMFPACK_*solve only: */
\multicolumn{4}{l}{Used by {\tt umfpack\_*solve:}} \\
{\tt Control(8)} & {\tt Control[UMFPACK\_IRSTEP]} & 2 & max iter. refinement steps \\
\hline
%
% /* compile-time settings - Control [8..11] cannot be changed at run time: */
\multicolumn{4}{l}{Can only be changed at compile time:} \\
{\tt Control(9)} & {\tt Control[UMFPACK\_COMPILED\_WITH\_BLAS]} & - & true if BLAS is used \\
{\tt Control(10)} & {\tt Control[UMFPACK\_COMPILED\_FOR\_MATLAB]} & - & true for mexFunction \\
{\tt Control(11)} & {\tt Control[UMFPACK\_COMPILED\_WITH\_GETRUSAGE]} & - & true if {\tt getrusage} used \\
{\tt Control(12)} & {\tt Control[UMFPACK\_COMPILED\_IN\_DEBUG\_MODE]} & - & true if debug mode enabled \\
%
%/* Control [15...19] unused */
\hline
\end{tabular}
}
\end{table}

%-------------------------------------------------------------------------------
\subsection{Larger examples}
%-------------------------------------------------------------------------------

A full example of all user-callable UMFPACK routines (the {\tt int} routines)
is available in the C main program, {\tt umfpack\_demo.c} listed in
Section~\ref{Demo}.  A nearly identical program that uses the {\tt long} integer
version of UMFPACK is in {\tt umfpack\_l\_demo.c}.  Another example is
the UMFPACK mexFunction, {\tt umfpackmex.c}.  The mexFunction accesses only the
user-callable C interface to UMFPACK.  The only features that it does not use
are the support for the triplet form (MATLAB's sparse arrays are already in the
compressed column form) and the ability to reuse the {\tt Symbolic} object to
numerically factorize a matrix whose pattern is the same as a prior matrix
analyzed by {\tt umfpack\_symbolic} or {\tt umfpack\_qsymbolic}.  The latter is
an important feature, but the mexFunction does not return its opaque
{\tt Symbolic} and {\tt Numeric} objects to MATLAB.  Instead, it gets the
contents of these objects after extracting them via the {\tt umfpack\_get\_*}
routines, and returns them as MATLAB sparse matrices.

%-------------------------------------------------------------------------------
\section{Synopsis of all C-callable routines ({\tt int} version)}
%-------------------------------------------------------------------------------

Each subsection, below, summarizes the input variables, output variables, return
values, and calling sequences of the routines in one category.  Variables with
the same name as those already listed in a prior category have the same size
and type.

%-------------------------------------------------------------------------------
\subsection{Primary routines}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
#include "umfpack.h"
int status, n, nz, Ap [n+1], Ai [nz] ;
double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO], Ax [nz], X [n], B [n] ;
void *Symbolic, *Numeric ;
char *sys ;

status = umfpack_symbolic (n, Ap, Ai, &Symbolic, Control, Info) ;
status = umfpack_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;
status = umfpack_solve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info) ;
umfpack_free_symbolic (&Symbolic) ;
umfpack_free_numeric (&Numeric) ;
\end{verbatim}
}

%-------------------------------------------------------------------------------
\subsection{Alternative routines}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
int Qinit [n], Wi [n] ;
double W [n], Y [n], Z [n], S [n] ;

umfpack_defaults (Control) ;
status = umfpack_qsymbolic (n, Ap, Ai, Qinit, &Symbolic, Control, Info) ;
status = umfpack_wsolve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info,
    Wi, W, Y, Z, S) ;
\end{verbatim}
}

%-------------------------------------------------------------------------------
\subsection{Matrix manipulation routines}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
int Ti [nz], Tj [nz], Bp [n+1], Bi [max(n,nz)], P [n], Q [n], Cp [n+1], Ci [nz] ;
double Tx [nz], Cx [nz], Bx [nz] ;

status = umfpack_col_to_triplet (n, Ap, Tj) ;
status = umfpack_triplet_to_col (n, nz, Ti, Tj, Tx, Bp, Bi, Bx) ;
status = umfpack_transpose (n, Ap, Ai, Ax, P, Q, Cp, Ci, Cx) ;
\end{verbatim}
}

%-------------------------------------------------------------------------------
\subsection{Getting the contents of opaque objects}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
int lnz, unz, Lp [n+1], Li [lnz], Up [n+1], Ui [unz] ;
double Lx [lnz], Ux [unz] ;
int nfr, nchains, nsparse_col, Qtree [n], Front_npivots [n], Front_parent [n],
    Chain_start [n], Chain_maxrows [n], Chain_maxcols [n] ;

status = umfpack_get_lunz (&lnz, &unz, &n, Numeric) ;
status = umfpack_get_numeric (Lp, Li, Lx, Up, Ui, Ux, P, Q, Numeric) ;
status = umfpack_get_symbolic (&n, &nz, &nfr, &nchains, &nsparse_col,
    Qtree, Front_npivots, Front_parent, Chain_start, Chain_maxrows,
    Chain_maxcols, Symbolic) ;
\end{verbatim}
}

Note: the {\tt nsparse\_col} argument is no longer relevant.
It is always equal to {\tt n} in this version.

%-------------------------------------------------------------------------------
\subsection{Reporting routines}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
char *name, *form ;

umfpack_report_status (Control, status) ;
umfpack_report_control (Control) ;
umfpack_report_info (Control, Info) ;
status = umfpack_report_matrix (name, n, Ap, Ai, Ax, form, Control) ;
status = umfpack_report_numeric (name, Numeric, Control) ;
status = umfpack_report_perm (name, n, P, Control) ;
status = umfpack_report_symbolic (name, Symbolic, Control) ;
status = umfpack_report_triplet (name, n, nz, Ti, Tj, Tx, Control) ;
status = umfpack_report_vector (name, n, X, Control) ;
\end{verbatim}
}


%-------------------------------------------------------------------------------
\section{Synopsis of all C-callable routines ({\tt long} version)}
%-------------------------------------------------------------------------------

Each subsection, below, summarizes the input variables, output variables, return
values, and calling sequences of the routines in one category.  Variables with
the same name as those already listed in a prior category have the same size
and type.  Note that the include file, {\tt umfpack.h}, is the same for both
{\tt int} and {\tt long} versions of UMFPACK.

%-------------------------------------------------------------------------------
\subsection{Primary routines}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
#include "umfpack.h"
long status, n, nz, Ap [n+1], Ai [nz] ;
double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO], Ax [nz], X [n], B [n] ;
void *Symbolic, *Numeric ;
char *sys ;

status = umfpack_l_symbolic (n, Ap, Ai, &Symbolic, Control, Info) ;
status = umfpack_l_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;
status = umfpack_l_solve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info) ;
umfpack_l_free_symbolic (&Symbolic) ;
umfpack_l_free_numeric (&Numeric) ;
\end{verbatim}
}

%-------------------------------------------------------------------------------
\subsection{Alternative routines}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
long Qinit [n], Wi [n] ;
double W [n], Y [n], Z [n], S [n] ;

umfpack_l_defaults (Control) ;
status = umfpack_l_qsymbolic (n, Ap, Ai, Qinit, &Symbolic, Control, Info) ;
status = umfpack_l_wsolve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info,
    Wi, W, Y, Z, S) ;
\end{verbatim}
}

%-------------------------------------------------------------------------------
\subsection{Matrix manipulation routines}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
long Ti [nz], Tj [nz], Bp [n+1], Bi [max(n,nz)], P [n], Q [n], Cp [n+1], Ci [nz] ;
double Tx [nz], Cx [nz], Bx [nz] ;

status = umfpack_l_col_to_triplet (n, Ap, Tj) ;
status = umfpack_l_triplet_to_col (n, nz, Ti, Tj, Tx, Bp, Bi, Bx) ;
status = umfpack_l_transpose (n, Ap, Ai, Ax, P, Q, Cp, Ci, Cx) ;
\end{verbatim}
}

%-------------------------------------------------------------------------------
\subsection{Getting the contents of opaque objects}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
long lnz, unz, Lp [n+1], Li [lnz], Up [n+1], Ui [unz] ;
double Lx [lnz], Ux [unz] ;
long nfr, nchains, nsparse_col, Qtree [n], Front_npivots [n], Front_parent [n],
    Chain_start [n], Chain_maxrows [n], Chain_maxcols [n] ;

status = umfpack_l_get_lunz (&lnz, &unz, &n, Numeric) ;
status = umfpack_l_get_numeric (Lp, Li, Lx, Up, Ui, Ux, P, Q, Numeric) ;
status = umfpack_l_get_symbolic (&n, &nz, &nfr, &nchains, &nsparse_col,
    Qtree, Front_npivots, Front_parent, Chain_start, Chain_maxrows,
    Chain_maxcols, Symbolic) ;
\end{verbatim}
}

Note: the {\tt nsparse\_col} argument is no longer relevant.
It is always equal to {\tt n} in this version.

%-------------------------------------------------------------------------------
\subsection{Reporting routines}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
char *name, *form ;

umfpack_l_report_status (Control, status) ;
umfpack_l_report_control (Control) ;
umfpack_l_report_info (Control, Info) ;
status = umfpack_l_report_matrix (name, n, Ap, Ai, Ax, form, Control) ;
status = umfpack_l_report_numeric (name, Numeric, Control) ;
status = umfpack_l_report_perm (name, n, P, Control) ;
status = umfpack_l_report_symbolic (name, Symbolic, Control) ;
status = umfpack_l_report_triplet (name, n, nz, Ti, Tj, Tx, Control) ;
status = umfpack_l_report_vector (name, n, X, Control) ;
\end{verbatim}
}

%-------------------------------------------------------------------------------
\section{Synopsis of utility routines}
%-------------------------------------------------------------------------------

This routine is the same in both {\tt int} and {\tt long} versions of UMFPACK.

{\footnotesize
\begin{verbatim}
double t ;

t = umfpack_timer ( ) ;

\end{verbatim}
}

%-------------------------------------------------------------------------------
\section{Installation}
\label{Install}
%-------------------------------------------------------------------------------

UMFPACK comes with a {\tt Makefile} for compiling the C-callable {\tt umfpack.a}
library and the {\tt umfpack} mexFunction on Unix.  System-dependent
configurations are controlled by the {\tt Makefile}, and defined in
{\tt umf\_config.h} listed in Section~\ref{Config}.  You should not have to
modify {\tt umf\_config.h}.

To compile {\tt umfpack.a} on most Unix systems, all you need to do is to type
{\tt make}.  This will use the generic configuration, in {\tt Make.generic}.
The three demo programs will be executed, and the output of
{\tt umfpack\_demo.c} and {\tt umfpack\_l\_demo.c} will be compared with
{\tt umfpack\_demo.out} and {\tt umfpack\_l\_demo.out}.  These two demo programs
are identical, except that {\tt umfpack\_demo.c} uses the {\tt int} version,
while {\tt umfpack\_l\_demo.out} uses the {\tt long} version of UMFPACK.
Expect to see a few differences, such as residual norms, compile-time control
settings, and perhaps memory usage differences.  (The Compaq Alpha uses the
LP64 model by default, so if you're using that computer compare your output
with the 64-bit Solaris output in {\tt umfpack\_demo.out64} and
{\tt umfpack\_l\_demo.out64}).  The BLAS \cite{DaydeDuff99,ACM679a,ATLAS}
will not be used, so the performance of UMFPACK will not be as high as possible.
For better performance, edit the {\tt Makefile} and un-comment the
{\tt include Make.*} statement that is specific to your computer.  For example,

{\footnotesize
\begin{verbatim}
    # include Make.generic
    # include Make.linux
    # include Make.sgi
    include Make.solaris
    # include Make.alpha
    # include Make.rs6000
\end{verbatim}
}

will include the Solaris-specific configurations, which uses the Sun Performance
Library BLAS ({\tt sunperf}), and compiler optimizations that are different than
the generic settings.  If you change the {\tt Makefile} or your system-specific
{\tt Make.*} file, be sure to type {\tt make purge} before recompiling.  Here
are the various parameters that you can control in your {\tt Make.*} file;
more details are in {\tt umf\_config.h} listed in Section~\ref{Config}:

\begin{itemize}
\item {\tt CC = } your C compiler, usually, {\tt cc}.  If you don't modify
    this string at all in your {\tt Make.*}, then the {\tt make} program will
    use your default C compiler (if {\tt make} is installed properly).
\item {\tt RANLIB = } your system's {\tt ranlib} program, if needed.
\item {\tt CFLAGS = } optimization flags, such as {\tt -O}.
\item {\tt CONFIG = } configuration settings.
\item {\tt LIB = } your libraries, such as {\tt -lm} or {\tt -lblas}.
\end{itemize}

The {\tt CONFIG} string can include combinations of the following:
\begin{itemize}
\item {\tt -DNBLAS} if you do not have any BLAS at all.  By default,
    {\tt umf\_config.h} assumes you have some version of the BLAS.  The BLAS
    are de-selected in {\tt Make.generic} with the statement
    {\tt CONFIG = -DNBLAS}.
\item {\tt -DNCBLAS} if you do not have the C-BLAS \cite{ATLAS}.
    The interface to the C-BLAS is identical on any system (Unix or Windows).
    By default, {\tt umf\_config.h} assumes you have the C-BLAS,
    except for Solaris (which has {\tt sunperf}) and MATLAB,
    which has its own BLAS for compiling the MATLAB mexFunction on any system.
\item {\tt -DNSUNPERF} if you are on Solaris but do not have {\tt sunperf}.
\item {\tt -DLONGBLAS} if your BLAS can take {\tt long} integer input
    arguments.  If not defined, then the {\tt umfpack\_l\_*} version of
    UMFPACK that uses {\tt long} integers does not call the BLAS.
\item {\tt -DGETRUSAGE} if you have the {\tt getrusage} function.  This should
    exist on any UNIX system.
\item Options for controlling how C calls the Fortran BLAS:
    {\tt -DBLAS\_BY\_VALUE}, \newline {\tt -DBLAS\_NO\_UNDERSCORE},
    and {\tt -DBLAS\_CHAR\_ARG}.  These are set automatically for Sun Solaris,
    SGI Irix, Red Hat Linux, Compaq Alpha, and AIX (the IBM RS 6000).
\end{itemize}

To compile the {\tt umfpack} mexFunction on Unix, type {\tt make umfpack}.  The
MATLAB {\tt mex} command will select the appropriate compiler and compiler flags
for your system, and the BLAS internal to MATLAB will be used.  The
{\tt mexopts.sh} file in your UMFPACK directory has been modified from the
MATLAB default; the unmodified version is in {\tt mexopts.sh.orig} for
comparison.

If you're running Windows, and all you want to do is use UMFPACK in MATLAB,
then just type {\tt umfpack\_make} in MATLAB.  MATLAB Version 6.0 or higher is
required.  You won't be able to use the BLAS when compiling with the {\tt lcc}
compiler provided with MATLAB Version 6.0); you will get an error
stating that {\tt \_dgemm} is undefined.  There is no work-around for this
problem.  Either use a different C compiler, or don't use the BLAS.

%-------------------------------------------------------------------------------
\section{Future work}
\label{Future}
%-------------------------------------------------------------------------------

Here are a few features that are not in UMFPACK Version 3.2, in no particular
order.  They may appear in a future release of UMFPACK.  If you are interested,
let me know and I could consider including them:

\begin{enumerate}

\item Future versions may have different default {\tt Control} parameters.

\item a condition number estimator.  You can write your own in MATLAB by
    making a copy of the built-in MATLAB {\tt condest.m} routine and replacing
    {\tt LU} with {\tt umfpack}.  Be sure to do so only if your MATLAB license
    allows you to, and do not distribute the derivative MATLAB code without
    direct permission from The Mathworks, Inc.

\item an estimate of the 1-norm of $\m{PAQ}-\m{LU}$. \newline
    See {\tt ftp://ftp.mathworks.com/pub/contrib/v4/linalg/normest1.m}
    for a similar algorithm that computes the 1-norm estimate of
    $\sigma \m{I}+\m{AA}\tr-\m{LL}\tr$.  It can easily be modified to compute
    the 1-norm estimate of $\m{PAQ}-\m{LU}$.  See also \cite{DavisHager99}.

\item a complex version.

\item when using iterative refinement, the residual $\m{Ax}-\m{b}$ could be
    returned by {\tt umfpack\_solve} ({\tt umfpack\_wsolve} already does so,
    but this is not documented).

\item the solve routines could handle multiple right-hand sides, and sparse
    right-hand sides.

\item an option to redirect the error and diagnostic output to something
    other than standard output.

\item permutation to block-triangular-form \cite{Duff78a} for the C-callable
        interface.

\item the symbolic and numeric factorization could handle singular matrices,
    just like MATLAB's {\tt LU}.

\item the ability to use user-provided {\tt malloc}, {\tt free}, and
    {\tt realloc} memory allocation routines.  Note that UMFPACK makes very
    few calls to these routines.

\item the ability to use user-provided work arrays, so that {\tt malloc},
    {\tt free}, and {\tt realloc} realloc are not called.  The
    {\tt umfpack\_wsolve} routine is one example.

\item future versions may return more statistics in the {\tt Info} array, and
    they may use more entries in the {\tt Control} array.

\item use a method that takes time proportional to the number of nonzeros in
    $\m{A}$ to analyze $\m{A}$ when {\tt Qinit} is provided (or when
    {\tt Qinit} is not provided and {\tt umf\_colamd} ignores "dense" rows)
    \cite{GilbertNgPeyton94}.  The current method in {\tt umf\_analyze.c} takes
    time proportional to the number of nonzeros in the upper bound of $\m{U}$.

\item an option of extracting the diagonal of $\m{U}$ (or other subsets of
    $\m{L}$ and $\m{U}$) from the {\tt Numeric} object without having to extract
    the entire LU factorization.

\item a Fortran interface (this would probably require modifying UMFPACK to use
    user-provided work arrays).

\item a C++ interface.

\item a parallel version using MPI.

\end{enumerate}


%-------------------------------------------------------------------------------
\newpage
\section{The primary UMFPACK routines}
\label{Primary}
%-------------------------------------------------------------------------------

The include files are the same for both {\tt int} and {\tt long} versions of
UMFPACK.  The generic integer type is {\tt Int}, which is an {\tt int} or
{\tt long}, depending on which version of UMFPACK you are using.  

\subsection{umfpack\_symbolic and umfpack\_l\_symbolic}

{\footnotesize
\begin{verbatim}
int umfpack_symbolic
(
    int n,
    const int Ap [ ],
    const int Ai [ ],
    void **Symbolic,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO]
) ;

long umfpack_l_symbolic
(
    long n,
    const long Ap [ ],
    const long Ai [ ],
    void **Symbolic,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO]
) ;

int Syntax:

    #include "umfpack.h"
    void *Symbolic ;
    int n, *Ap, *Ai, status ;
    double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO] ;
    status = umfpack_symbolic (n, Ap, Ai, &Symbolic, Control, Info) ;

long Syntax:

    #include "umfpack.h"
    void *Symbolic ;
    long n, *Ap, *Ai, status ;
    double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO] ;
    status = umfpack_l_symbolic (n, Ap, Ai, &Symbolic, Control, Info) ;

Purpose:

    Given nonzero pattern of a sparse matrix A in column-oriented form,
    umfpack_symbolic performs a column pre-ordering to reduce fill-in
    (using UMF_colamd, modified from colamd V2.0 for UMFPACK), and a symbolic
    factorization.  This is required before the matrix can be numerically
    factorized with umfpack_numeric.  If you wish to bypass the UMF_colamd
    pre-ordering, use umfpack_qsymbolic instead.

Returns:

    The status code is returned.  See Info [UMFPACK_STATUS], below.

Arguments:

    Int n ;             Input argument, not modified.

        A is an n-by-n matrix.  Restriction: n > 0.

    Int Ap [n+1] ;      Input argument, not modified.

        Ap is an integer array of size n+1.  On input, it holds the "pointers"
        for the column form of the sparse matrix A.  Column j of the matrix A
        is held in Ai [(Ap [j]) ... (Ap [j+1]-1)].  The first entry, Ap [0],
        must be zero, and Ap [j] < Ap [j+1] must hold for all j in the range
        0 to n-1.  The value nz = Ap [n] is thus the total number of entries
        in the pattern of the matrix A.  nz must be greater than zero.

    Int Ai [nz] ;       Input argument, not modified, of size nz = Ap [n].

        The nonzero pattern (row indices) for column j is stored in
        Ai [(Ap [j]) ... (Ap [j+1]-1)].  The row indices in a given column j
        must be in ascending order, and no duplicate row indices may be present.
        Row indices must be in the range 0 to n-1 (the matrix is 0-based).
        See umfpack_triplet_to_col for how to sort the columns of a matrix
        and sum up the duplicate entries.  See umfpack_report_matrix for how to
        print the matrix A.

    void **Symbolic ;   Output argument.

        **Symbolic is the address of a (void *) pointer variable in the user's
        calling routine (see Syntax, above).  On input, the contents of this
        variable are not defined.  On output, this variable holds a (void *)
        pointer to the Symbolic object (if successful), or (void *) NULL if
        a failure occurred.

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_DENSE_ROW]:  rows with more than
            max (16, Control [UMFPACK_DENSE_ROW] * 16 * sqrt (n))
            entries (after "dense" columns are removed) are ignored in the
            column pre-ordering, UMF_colamd.  Default: 0.2.

        Control [UMFPACK_DENSE_COL]:  columns with more than
            max (16, Control [UMFPACK_DENSE_COL] * 16 * sqrt (n))
            entries are placed placed last in the column pre-ordering by
            UMF_colamd.  Default: 0.2.

        Control [UMFPACK_BLOCK_SIZE]:  the block size to use for Level-3 BLAS
            in the subsequent numerical factorization (umfpack_numeric).
            A value less than 1 is treated as 1.  Default: 24.  Modifying this
            parameter affects when updates are applied to the working frontal
            matrix, and can indirectly affect fill-in and operation count.
            As long as the block size is large enough (8 or so), this parameter
            has modest effect on performance.  In Version 3.0, this parameter
            was an input to umfpack_numeric, and had a default value of 16.
            On a Sun UltraSparc, a block size of 24 is better for larger
            matrices (16 is better for smaller ones, but not by much).  In the
            current version, it is required in the symbolic analysis phase, and
            is thus an input to this phase.

    double Info [UMFPACK_INFO] ;        Output argument, not defined on input.

        Contains statistics about the symbolic analysis.  If a (double *) NULL
        pointer is passed, then no statistics are returned in Info (this is not
        an error condition).  The entire Info array is cleared (all entries set
        to -1) and then the following statistics are computed:

        Info [UMFPACK_STATUS]: status code.  This is also the return value,
            whether or not Info is present.

            UMFPACK_OK

                Each column of the input matrix contained row indices
                in increasing order, with no duplicates.  Only in this case
                does umfpack_symbolic compute a valid symbolic factorization.
                For the other cases below, no Symbolic object is created
                (*Symbolic is (void *) NULL).

            UMFPACK_ERROR_jumbled_matrix

                Columns of input matrix were jumbled (unsorted columns or
                duplicate entries).

            UMFPACK_ERROR_n_nonpositive

                n is less than or equal to zero.

            UMFPACK_ERROR_singular_matrix

                Matrix is singular.

            UMFPACK_ERROR_nz_negative

                Number of entries in the matrix is negative.

            UMFPACK_ERROR_Ap0_nonzero

                Ap [0] is nonzero.

            UMFPACK_ERROR_col_length_negative

                A column has a negative number of entries.

            UMFPACK_ERROR_row_index_out_of_bounds

                A row index is out of bounds.

            UMFPACK_ERROR_out_of_memory

                Insufficient memory to perform the symbolic analysis.

            UMFPACK_ERROR_argument_missing

                One or more required arguments (Ap and/or Ai) is missing.

            UMFPACK_ERROR_internal_error

                Something very serious went wrong.  This is a bug.
                Please contact the author (davis@cise.ufl.edu).

        Info [UMFPACK_N]:  the value of the input argument n.

        Info [UMFPACK_NZ]:  the number of entries in the input matrix (Ap [n]).

        Info [UMFPACK_SIZE_OF_UNIT]:  the number of bytes in a Unit,
            for memory usage statistics below.

        Info [UMFPACK_SIZE_OF_INT]:  the number of bytes in an int.

        Info [UMFPACK_SIZE_OF_LONG]:  the number of bytes in a long.

        Info [UMFPACK_SIZE_OF_POINTER]:  the number of bytes in a void *
            pointer.

        Info [UMFPACK_SIZE_OF_ENTRY]:  the number of bytes in a numerical entry.

        Info [UMFPACK_NDENSE_ROW]:  number of "dense" rows in A.  These rows are
            ignored when the column pre-ordering is computed in UMF_colamd.
            If > 0, then the matrix had to be re-analyzed by UMF_analyze, which
            does not ignore these rows.  Note that all rows are stored in the
            same data structure, regardless of whether they are "sparse",
            "dense", or "empty".

        Info [UMFPACK_NEMPTY_ROW]:  number of "empty" rows in A.  These are
            rows whose entries are all in "dense" columns.  Any given row
            is classified as either "dense" or "empty" or "sparse".

        Info [UMFPACK_NDENSE_COL]:  number of "dense" columns in A.
            These columns are ordered last in the factorization.
            Any given column is classified as either "dense" or "empty" or
            "sparse".  All columns are stored in the same data structure,
            however (Version 3.0 stored dense columns in a separate dense
            array, but this is no longer true for Version 3.1 and following).

        Info [UMFPACK_NEMPTY_COL]:  number of "empty" columns in A.  These are
            columns whose entries are all in "dense" rows.  These columns are
            ordered last in the factorization, along with "dense" columns.

        Info [UMFPACK_SYMBOLIC_DEFRAG]:  number of garbage collections
            performed in UMF_colamd, the column pre-ordering routine, and in
            UMF_analyze, which is called if UMF_colamd isn't, or if UMF_colamd
            ignores one or more "dense" rows.

        Info [UMFPACK_SYMBOLIC_PEAK_MEMORY]:  the amount of memory (in Units)
            required for umfpack_symbolic to complete.  This is roughly
            2.2*nz + (20 to 25)*n integers, depending on the matrix.  This
            count includes the size of the Symbolic object itself, which is
            reported in Info [UMFPACK_SYMBOLIC_SIZE].

        Info [UMFPACK_SYMBOLIC_SIZE]: the final size of the Symbolic object (in
            Units).  This is fairly small, roughly (1 to 6)*n integers,
            depending on the matrix.

        Info [UMFPACK_VARIABLE_INIT_ESTIMATE]: the Numeric object contains two
            components.  The first is fixed in size (O (n) information, plus
            the "dense" part of the LU factors).   The second part holds the
            sparse LU factors and the contribution blocks from factorized
            frontal matrices.  This part changes in size during factorization.
            Info [UMFPACK_VARIABLE_INIT_ESTIMATE] is the exact size (in Units)
            required for this second variable-sized part in order for the
            numerical factorization to start.

        Info [UMFPACK_VARIABLE_PEAK_ESTIMATE]: the estimated peak size (in
            Units) of the variable-sized part of the Numeric object.  This is
            usually an upper bound, but that is not guaranteed.

        Info [UMFPACK_VARIABLE_FINAL_ESTIMATE]: the estimated final size (in
            Units) of the variable-sized part of the Numeric object.  This is
            usually an upper bound, but that is not guaranteed.  It holds just
            the sparse LU factors.

        Info [UMFPACK_NUMERIC_SIZE_ESTIMATE]:  an estimate of the final size (in
            Units) of the entire Numeric object (both fixed-size and variable-
            sized parts), which holds the LU factorization (including the L, U,
            P and Q matrices).

        Info [UMFPACK_PEAK_MEMORY_ESTIMATE]:  an estimate of the total amount of
            memory (in Units) required by umfpack_symbolic and umfpack_numeric
            to perform both the symbolic and numeric factorization.  This is the
            larger of the amount of memory needed in umfpack_numeric itself, and
            the amount of memory needed in umfpack_symbolic
            (Info [UMFPACK_SYMBOLIC_PEAK_MEMORY]).  The count includes the size
            of both the Symbolic and Numeric objects themselves.

        Info [UMFPACK_FLOPS_ESTIMATE]:  an estimate of the total floating-point
            operations required to factorize the matrix.  This is a "true"
            theoretical estimate of the number of flops that would be performed
            by a flop-parsimonious sparse LU algorithm.  It assumes that no
            extra flops are performed except for what is strictly required to
            compute the LU factorization.  It ignores, for example, the flops
            performed by umfpack_numeric to add contribution blocks of frontal
            matrices together.  If L and U are the upper bound on the pattern
            of the factors, then this flop count estimate can be represented in
            Matlab as:

                Lnz = full (sum (spones (L))) - 1 ;     % nz in each col of L
                Unz = full (sum (spones (U')))' - 1 ;   % nz in each row of U
                flops = 2*Lnz*Unz + sum (Lnz) ;

            The flop counts include add, subtract, multiply, and divide.  They
            exclude max, absolute value computations, and comparisons.

            The actual "true flop" count found by umfpack_numeric will be less
            than this estimate.

        Info [UMFPACK_LNZ_ESTIMATE]:  an estimate of the number of nonzeros in
            L, including the diagonal.  Since L is unit-diagonal, the diagonal
            of L is not stored.  This estimate is a strict upper bound on the
            actual nonzeros in L to be computed by umfpack_numeric.

        Info [UMFPACK_UNZ_ESTIMATE]:  an estimate of the number of nonzeros in
            U, including the diagonal.  This estimate is a strict upper bound on
            the actual nonzeros in U to be computed by umfpack_numeric.

        Info [UMFPACK_SYMBOLIC_TIME]:  The time taken by umfpack_symbolic, in
            seconds.  In the ANSI C version, this may be invalid if the time
            taken is more than about 36 minutes, because of wrap-around in
            the ANSI C clock ( ) function.  Compile UMFPACK with -DGETRUSAGE
            if you have the more accurate getrusage ( ) function.

        At the start of umfpack_symbolic, all of Info is set of -1, and then
        after that only the above listed Info [...] entries are accessed.
        Future versions might modify different parts of Info.
\end{verbatim}
}

\newpage
\subsection{umfpack\_numeric and umfpack\_l\_numeric}

{\footnotesize
\begin{verbatim}
int umfpack_numeric
(
    const int Ap [ ],
    const int Ai [ ],
    const double Ax [ ],
    void *Symbolic,
    void **Numeric,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO]
) ;

long umfpack_l_numeric
(
    const long Ap [ ],
    const long Ai [ ],
    const double Ax [ ],
    void *Symbolic,
    void **Numeric,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO]
) ;

int Syntax:

    #include "umfpack.h"
    void *Symbolic, *Numeric ;
    int *Ap, *Ai, status ;
    double *Ax, Control [UMFPACK_CONTROL], Info [UMFPACK_INFO] ;
    status = umfpack_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;

long Syntax:

    #include "umfpack.h"
    void *Symbolic, *Numeric ;
    long *Ap, *Ai, status ;
    double *Ax, Control [UMFPACK_CONTROL], Info [UMFPACK_INFO] ;
    status = umfpack_l_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;

Purpose:

    Given a sparse matrix A in column-oriented form, and a symbolic analysis
    computed by umfpack_symbolic, the umfpack_numeric routine performs the
    numerical factorization, PAQ=LU, where P and Q are permutation matrices
    (represented as permutation vectors), L is unit-lower triangular, and U
    is upper triangular.  This is required before the system Ax=b (or other
    related linear systems) can be solved.  umfpack_numeric can be called
    multiple times for each call to umfpack_symbolic, to factorize a sequence
    of matrices with identical nonzero pattern.  Simply compute the Symbolic
    object once, with umfpack_*symbolic, and reuse it for subsequent matrices.
    umfpack_numeric safely detects if the pattern changes, and sets an
    appropriate error code.

Returns:

    The status code is returned.  See Info [UMFPACK_STATUS], below.

Arguments:

    Int Ap [n+1] ;      Input argument, not modified.

        This must be identical to the Ap array passed to umfpack_symbolic.
        The value of n is what was passed to umfpack_symbolic (this is held in
        the Symbolic object).

    Int Ai [nz] ;       Input argument, not modified, of size nz = Ap [n].

        This must be identical to the Ai array passed to umfpack_symbolic.

        Not all changes to Ai and Ap are detected; if the matrix has the same
        number of nonzeros and can be factorized in the existing frontal
        matrices as defined in the Symbolic object, then umfpack_numeric will
        not complain, and will successfully factorize the matrix and return a
        valid Numeric object.

    double Ax [nz] ;    Input argument, not modified, of size nz = Ap [n].

        The numerical values of the sparse matrix A.  The nonzero pattern (row
        indices) for column j is stored in Ai [(Ap [j]) ... (Ap [j+1]-1)], and
        the corresponding numerical values are stored in
        Ax [(Ap [j]) ... (Ap [j+1]-1)].

    void *Symbolic ;    Input argument, not modified.

        The Symbolic object, which holds the symbolic factorization computed by
        umfpack_symbolic.  The Symbolic object is not modified by
        umfpack_numeric.

    void **Numeric ;    Output argument.

        **Numeric is the addres of a (void *) pointer variable in the user's
        calling routine (see Syntax, above).  On input, the contents of this
        variable are not defined.  On output, this variable holds a (void *)
        pointer to the Numeric object (if successful), or (void *) NULL if
        a failure occurred.

    double Control [UMFPACK_CONTROL] ;   Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PIVOT_TOLERANCE]:  relative pivot tolerance for
            threshold partial pivoting with row interchanges.  In any given
            column, an entry is numerically acceptable if it is greater than or
            equal to Control [UMFPACK_PIVOT_TOLERANCE] times the largest
            absolute value in the column.  A value of 1.0 gives true partial
            pivoting.  A value of zero is treated as 1.0.  Default: 0.1.
            Smaller values tend to lead to sparser LU factors, but the solution
            to the linear system can become inaccurate.  Larger values can lead
            to a more accurate solution (but not always), and usually an
            increase in the total work.

        Control [UMFPACK_RELAXED_AMALGAMATION]:  This controls the creation of
            "elements" (small dense submatrices) that are formed when a frontal
            matrix is factorized.  A new element is created if the current one,
            plus the new pivot, contains "too many" explicitly zero numerical
            entries.  The two elements are merged if the number of extra zero
            entries is < Control [UMFPACK_RELAXED_AMALGAMATION] times the
            size of the merged element.  A lower setting
            decreases fill-in, but run-time and memory usage can increase.
            A larger setting increases fill-in (because the extra zeros are
            treated as normal entries during pivot selection), but this can
            lead to an increase in run-time but (paradoxically) a decrease in
            memory usage (one merged elements can take less space than two
            separate elements).  Except for the initial column ordering,
            this parameter has the most impact on the run-time, fill-in,
            operation count, and memory usage.
            Default: 0.25, which is fine for nearly all matrices.
            (For nearly all matrices, different values of this parameter can
            decrease the run-time by at most 5%, but can also dramatically
            increase the run time for some matrices).

        Control [UMFPACK_RELAXED2_AMALGAMATION]:  This, along with the block
            size (Control [UMFPACK_BLOCK_SIZE]), controls how often the
            pending updates are applied when the next pivot entry resides in
            the current frontal matrix.  If the number of zero entries in the
            LU part of the current frontal matrix would exceed this parameter
            times the size of the LU part, then the pending updates are applied
            before the next pivot is included in the frontal matrix.
            Default: 0.20 (that is, more than 10% zero entries causes the
            pending updates to be applied).  This input parameter is new
            since Version 3.1.

        Control [UMFPACK_RELAXED3_AMALGAMATION]:  This, along with the block
            size (Control [UMFPACK_BLOCK_SIZE]), controls how often the
            pending updates are applied when the next pivot entry does NOT reside
            in the current frontal matrix.  If the number of zero entries in the
            LU part of the current frontal matrix would exceed this parameter
            times the size of the LU part, then the pending updates are applied
            before the next pivot is included in the frontal matrix.
            Default: 0.10 (that is, more than 10% zero entries causes the
            pending updates to be applied).  This input parameter is new
            since Version 3.1.

        Control [UMFPACK_ALLOC_INIT]: When umfpack_numeric starts, it allocates
            memory for the Numeric object.  Part of this is of fixed size
            (approximately n double's + 12*n integers).
            The remainder is of variable size, which grows to hold the LU
            factors and the frontal matrices created during factorization.
            A estimate of the upper bound is computed by umfpack_symbolic, and
            returned by umfpack_*symbolic in
            Info [UMFPACK_VARIABLE_PEAK_ESTIMATE].  umfpack_numeric initially
            allocates space for the variable-sized part equal to this estimate
            times Control [UMFPACK_ALLOC_INIT].  Typically, umfpack_numeric
            needs only about half the estimated memory space, so a setting of
            0.5 or 0.6 often provides enough memory for umfpack_numeric to
            factorize the matrix with no subsequent increases in the size of
            this block.  A value less than zero is treated as zero (in which
            case, just the bare minimum amount of memory needed to start the
            factorization is initially allocated).  The bare initial memory
            required is returned by umfpack_*symbolic in
            Info [UMFPACK_VARIABLE_INIT_ESTIMATE] (which in fact not an
            estimate, but exact).  If the variable-size part of the Numeric
            object is found to be too small sometime after numerical
            factorization has started, the memory is increased in size by a
            factor of 1.2.   If this fails, the request is reduced by a factor
            of 0.95 until it succeeds, or until it determines that no increase
            in size is possible.  Garbage collection then occurs.  These two
            factors (1.2 and 0.95) are fixed control parameters defined in
            umf_internal.h and cannot be changed at run-time.  You would need
            to edit umf_internal.h to modify them.  If you do this, be sure that
            the two factors are greater than 1 and less than 1, respectively.

            The strategy of attempting to malloc a working space, and re-trying
            with a smaller space, may not work under Matlab, since mxMalloc
            aborts the mexFunction if it fails.  I may try to address this is
            issue in a future release - in the mean time, decrease
            Control [UMFPACK_ALLOC_INIT] if you run out of memory in Matlab.

            Default initial allocation size: 0.7.  Thus, with the default
            control settings, the upper-bound is reached after two reallocations
            (0.7 * 1.2 * 1.2 = 1.008).

            Changing this parameter has no affect on fill-in or operation count.
            It has a small impact on run-time (the extra time required to do
            the garbage collection and memory reallocation).

        Control [UMFPACK_PIVOT_OPTION]:  If this is nonzero, then entries on
            the diagonal of A are given preference over off-diagonal entries.
            This can improve the fill-in on matrices with symmmetric nonzero
            pattern.  Default: 0 (do not give preference to the diagonal of A).
            This parameter was added for UMFPACK Version 3.1.

    double Info [UMFPACK_INFO] ;        Output argument.

        Contains statistics about the numeric factorization.  If a
        (double *) NULL pointer is passed, then no statistics are returned in
        Info (this is not an error condition).  The following statistics are
        computed in umfpack_numeric:

        Info [UMFPACK_STATUS]: status code.  This is also the return value,
            whether or not Info is present.

            UMFPACK_OK

                Numeric factorization was successful.  Only in this case
                does umfpack_numeric compute a valid numeric factorization.
                For the other cases below, no Numeric object is created
                (*Numeric is (void *) NULL).

            UMFPACK_ERROR_out_of_memory

                Insufficient memory to complete the numeric factorization.

            UMFPACK_ERROR_argument_missing

                One or more required arguments (Ap, Ai, and/or Ax) are missing.

            UMFPACK_ERROR_singular_matrix

                The input matrix is singular.

            UMFPACK_ERROR_invalid_Symbolic_object

                Symbolic object provided as input is invalid.

            UMFPACK_ERROR_different_pattern

                The pattern (Ap and/or Ai) has changed since the call to
                umfpack_*symbolic which produced the Symbolic object.

        Info [UMFPACK_N]:  the value of n stored in the Symbolic object.

        Info [UMFPACK_NZ]:  the number of entries in the input matrix.
            This value is obtained from the Symbolic object.

        Info [UMFPACK_SIZE_OF_UNIT]:  the number of bytes in a Unit, for memory
            usage statistics below.

        Info [UMFPACK_VARIABLE_INIT]: the initial size (in Units) of the
            variable-sized part of the Numeric object.  If this differs from
            Info [UMFPACK_VARIABLE_INIT_ESTIMATE], then the pattern (Ap and/or
            Ai) has changed since the last call to umfpack_*symbolic, which is
            an error condition.

        Info [UMFPACK_VARIABLE_PEAK]: the peak size (in Units) of the
            variable-sized part of the Numeric object.  This size is the amount
            of space actually used inside the block of memory, not the space
            allocated via UMF_malloc.  You can reduce UMFPACK's memory
            requirements by setting Control [UMFPACK_ALLOC_INIT] to the ratio
            Info [UMFPACK_VARIABLE_PEAK] / Info[UMFPACK_VARIABLE_PEAK_ESTIMATE].
            This will ensure that no memory reallocations occur (you may want to
            add 0.001 to make sure that integer roundoff does not lead to a
            memory size that is 1 Unit too small; otherwise, garbage collection
            and reallocation will occur).

        Info [UMFPACK_VARIABLE_FINAL]: the final size (in Units) of the
            variable-sized part of the Numeric object.  It holds just the
            sparse LU factors.

        Info [UMFPACK_NUMERIC_SIZE]:  the actual final size (in Units) of the
            entire Numeric object, including the final size of the variable
            part of the object.  Info [UMFPACK_NUMERIC_SIZE_ESTIMATE],
            an estimate, was computed by umfpack_symbolic.  The estimate is
            normally an upper bound on the actual final size, but this is not
            guaranteed.

        Info [UMFPACK_PEAK_MEMORY]:  the actual peak memory usage (in Units) of
            both umfpack_symbolic and umfpack_numeric.  An estimate,
            Info [UMFPACK_PEAK_MEMORY_ESTIMATE], was computed by
            umfpack_symbolic.  The estimate is normally an upper bound on the
            actual peak usage, but this is not guaranteed.  With testing on
            hundreds of matrix arising in real applications, I have never
            observed a matrix where this estimate or the Numeric size estimate
            was less than the actual result, but this is theoretically possible.
            Please send me one if you find such a matrix.

        Info [UMFPACK_FLOPS]:  the actual count of the (useful) floating-point
            operations performed.  An estimate, Info [UMFPACK_FLOPS_ESTIMATE],
            was computed by umfpack_symbolic.  The estimate is guaranteed to be
            an upper bound on this flop count.  The flop count excludes
            "useless" flops on zero values, flops performed during the pivot
            search (for tentative updates and assembly of candidate columns),
            and flops performed to add frontal matrices together.  It does
            include the flops performed to factorize the "dense" and "empty"
            columns.

        Info [UMFPACK_LNZ]: the actual nonzero entries in final factor L,
            including the diagonal.  This excludes any zero entries in L,
            although some of these are stored in the Numeric object.  It does
            include entries in "dense" or "empty" columns.  The
            Info [UMFPACK_LU_ENTRIES] statistic does account for all
            explicitly stored zeros, however.  Info [UMFPACK_LNZ_ESTIMATE],
            an estimate, was computed by umfpack_symbolic.  The estimate is
            guaranteed to be an upper bound on Info [UMFPACK_LNZ].

        Info [UMFPACK_UNZ]: the actual nonzero entries in final factor U,
            including the diagonal.  This excludes any zero entries in U,
            although some of these are stored in the Numeric object.  It does
            include entries in "dense" or "empty" columns.  The
            Info [UMFPACK_LU_ENTRIES] statistic does account for all
            explicitly stored zeros, however.  Info [UMFPACK_UNZ_ESTIMATE],
            an estimate, was computed by umfpack_symbolic.  The estimate is
            guaranteed to be an upper bound on Info [UMFPACK_UNZ].

        Info [UMFPACK_NUMERIC_DEFRAG]:  The number of garbage collections
            performed during umfpack_numeric, to compact contents of the
            variable-sized workspace used by umfpack_numeric.  No estimate was
            computed by umfpack_symbolic.  In the current version of UMFPACK,
            garbage collection is performed and then the memory is reallocated,
            so this statistic is the same as Info [UMFPACK_NUMERIC_REALLOC],
            below.  It may differ in future releases.

        Info [UMFPACK_NUMERIC_REALLOC]:  The number of times that the Numeric
            object was increased in size from its initial size.  A rough upper
            bound on the peak size of the Numeric object was computed by
            umfpack_symbolic, so reallocations should be rare.  However, if
            umfpack_numeric is unable to allocate that much storage, it reduces
            its request until either the allocation succeeds, or until it gets
            too small to do anything with.  If the memory that it finally got
            was small, but usable, then the reallocation count could be high.
            No estimate of this count was computed by umfpack_symbolic.

        Info [UMFPACK_NUMERIC_COSTLY_REALLOC]:  The number of times that the
            system realloc ( ) library routine had to move the workspace.
            Realloc can sometimes increase the size of a block of memory
            without moving it, which is much faster.  This statistic will
            always be <= Info [UMFPACK_NUMERIC_REALLOC].  If your memory space
            is fragmented, then the number of "costly" realloc's will be equal
            to Info [UMFPACK_NUMERIC_REALLOC].

        Info [UMFPACK_COMPRESSED_PATTERN]:  The number of integers used to
            represent the pattern of "sparse" part L and U.  The "sparse" part
            of L and U excludes entries on the diagonal, which is stored
            separately.  It excludes entries in the "dense" and "empty"
            columns.  Those are stored together in a single dense array of
            size n by (Info [UMFPACK_NDENSE_COL] + Info [UMFPACK_NEMPTY_COL]),
            and no integers are required to represent their pattern.

        Info [UMFPACK_LU_ENTRIES]:  The total number of numerical values that
            are stored for the LU factors, including the dense array for "dense"
            and "empty" columns.  Some of the values may be explicitly zero.

        Info [UMFPACK_NUMERIC_TIME]:  The time taken by umfpack_numeric, in
            seconds.  In the ANSI C version, this may be invalid if the time
            taken is more than about 36 minutes, because of wrap-around in
            the ANSI C clock ( ) function.  Compile UMFPACK with -DGETRUSAGE
            if you have the more accurate getrusage ( ) function.

        Only the above listed Info [...] entries are accessed.  The remaining
        entries of Info are not accessed or modified by umfpack_numeric.
        Future versions might modify different parts of Info.
\end{verbatim}
}

\newpage
\subsection{umfpack\_solve and umfpack\_l\_solve}

{\footnotesize
\begin{verbatim}
int umfpack_solve
(
    const char sys [ ],
    const int Ap [ ],
    const int Ai [ ],
    const double Ax [ ],
    double X [ ],
    const double B [ ],
    void *Numeric,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO]
) ;

long umfpack_l_solve
(
    const char sys [ ],
    const long Ap [ ],
    const long Ai [ ],
    const double Ax [ ],
    double X [ ],
    const double B [ ],
    void *Numeric,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO]
) ;

int Syntax:

    #include "umfpack.h"
    void *Numeric ;
    int status, *Ap, *Ai ;
    char *sys ;
    double *B, *X, *Ax, Info [UMFPACK_INFO], Control [UMFPACK_CONTROL] ;
    status = umfpack_solve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info) ;

long Syntax:

    #include "umfpack.h"
    void *Numeric ;
    long status, *Ap, *Ai ;
    char *sys ;
    double *B, *X, *Ax, Info [UMFPACK_INFO], Control [UMFPACK_CONTROL] ;
    status = umfpack_l_solve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info);

Purpose:

    Given LU factors computed by umfpack_numeric (PAQ=LU) and the
    right-hand-side, B, solve a linear system for the solution X.  Iterative
    refinement is optionally performed.  This routine dynamically allocates
    workspace of size O(n).

Returns:

    The status code is returned.  See Info [UMFPACK_STATUS], below.

Arguments:

    char sys [ ] ;              Input argument, not modified.

        A string that defines which system to solve.  Iterative refinement can
        be optionally performed when the sys argument is:

            "Ax=b"
            "A'x=b"

        For these values of the sys argument, iterative refinement is not
        performed (Control [UMFPACK_IRSTEP], Ap, Ai, and Ax are ignored):

            "P'Lx=b"
            "L'Px=b"
            "UQ'x=b"
            "QU'x=b"
            "Lx=b"
            "L'x=b"
            "Ux=b"
            "U'x=b"

    Int Ap [n+1] ;              Input argument, not modified.
    Int Ai [nz] ;               Input argument, not modified.
    double Ax [nz] ;            Input argument, not modified.

        If iterative refinement is requested (Control [UMFPACK_IRSTEP] >= 1 and
        Ax=b or A'x=b is being solved), then these arrays must be identical to
        the same ones passed to umfpack_numeric.  The umfpack_solve routine
        does not check the contents of these three arguments, so the results
        are undefined if Ap, Ai, and/or Ax are modified between the calls the
        umfpack_numeric and umfpack_solve.  These three arrays do not need to
        be present (NULL pointers can be passed) if Control [UMFPACK_IRSTEP] is
        zero, or if a system other than Ax=b or A'x=b is being solved.

    double X [n] ;              Output argument.

        The solution to the linear system.

    double B [n] ;              Input argument, not modified.

        The right-hand side vector, b, stored as a conventional array of size n.
        This routine does not solve for multiple right-hand-sides, nor does it
        allow b to be stored in a sparse-column form.

    void *Numeric ;             Input argument, not modified.

        Numeric must point to a valid Numeric object, computed by
        umfpack_numeric.

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_IRSTEP]:  The maximum number of iterative refinement
            steps to attempt.  A value less than zero is treated as zero.  If
            less than 1, or if Ax=b or A'x=b is not being solved, then the Ap,
            Ai, and Ax arguments are not accessed.  Default: 2.

    double Info [UMFPACK_INFO] ;        Output argument.

        Contains statistics about the solution factorization.  If a
        (double *) NULL pointer is passed, then no statistics are returned in
        Info (this is not an error condition).  The following statistics are
        computed in umfpack_solve:

        Info [UMFPACK_STATUS]: status code.  This is also the return value,
            whether or not Info is present.

            UMFPACK_OK

                The linear system was successfully solved.

            UMFPACK_ERROR_out_of_memory

                Insufficient memory to solve the linear system.

            UMFPACK_ERROR_argument_missing

                One or more required arguments are missing.  The B, X, and
                sys arguments are always required.  Info and Control are
                not required.  Ap, Ai, and Ax are required if Ax=b or A'x=b
                is to be solve and the (default) iterative refinement is
                requested.

            UMFPACK_ERROR_invalid_Numeric_object

                The Numeric object is not valid.

        Info [UMFPACK_N]:  the value of n stored in the Numeric object.

        Info [UMFPACK_NZ]:  the number of entries in the input matrix, Ap [n],
            if iterative refinement is requested (sys is "Ax=b" or "A'x=b"
            and Control [UMFPACK_IRSTEP] >= 1).

        Info [UMFPACK_IR_TAKEN]:  The number of iterative refinement steps
            effectively taken.  The number of steps attempted may be one more
            than this; the refinement algorithm backtracks if the last
            refinement step worsens the solution.  This is set to -1 if
            iterative refinement was not requested.

        Info [UMFPACK_IR_ATTEMPTED]:   The number of iterative refinement steps
            attempted.  The number of times a linear system was solved is one
            more than this (once for the initial Ax=b, and once for each Ay=r
            solved for each iterative refinement step attempted).  This
            statistic is set to -1 if iterative refinement was not requested.

        Info [UMFPACK_OMEGA1]:  sparse backward error estimate, omega1, if
            iterative refinement was performed, or -1 if iterative refinement
            not performed.

        Info [UMFPACK_OMEGA2]:  sparse backward error estimate, omega2, if
            iterative refinement was performed, or -1 if iterative refinement
            not performed.

        Info [UMFPACK_SOLVE_FLOPS]:  the number of floating point operations
            performed to solve the linear system.  This includes the work
            taken for all iterative refinement steps, including the backtrack
            (if any).

        Info [UMFPACK_SOLVE_TIME]:  The time taken by umfpack_solve, in
            seconds.  In the ANSI C version, this may be invalid if the time
            taken is more than about 36 minutes, because of wrap-around in
            the ANSI C clock ( ) function.  Compile UMFPACK with -DGETRUSAGE
            if you have the more accurate getrusage ( ) function.

        Only the above listed Info [...] entries are accessed.  The remaining
        entries of Info are not accessed or modified by umfpack_solve.
        Future versions might modify different parts of Info.
\end{verbatim}
}

\newpage
\subsection{umfpack\_free\_symbolic and umfpack\_l\_free\_symbolic}

{\footnotesize
\begin{verbatim}
void umfpack_free_symbolic
(
    void **Symbolic
) ;

void umfpack_l_free_symbolic
(
    void **Symbolic
) ;

int Syntax:

    #include "umfpack.h"
    void *Symbolic ;
    umfpack_free_symbolic (&Symbolic) ;

long Syntax:

    #include "umfpack.h"
    void *Symbolic ;
    umfpack_l_free_symbolic (&Symbolic) ;

Purpose:

    Deallocates the Symbolic object and sets the Symbolic handle to NULL.
    This routine is the only valid way of destroying the Symbolic object;
    any other action (such as using "free (Symbolic) ;" or not freeing Symbolic
    at all) will lead to memory leaks.

Arguments:

    void **Symbolic ;           Input argument, deallocated and Symbolic is
                                set to (void *) NULL on output.

        Symbolic must point to a valid Symbolic object, computed by
        umfpack_symbolic.  No action is taken if Symbolic is a (void *) NULL
        pointer.
\end{verbatim}
}

\newpage
\subsection{umfpack\_free\_numeric and umfpack\_l\_free\_numeric}

{\footnotesize
\begin{verbatim}
void umfpack_free_numeric
(
    void **Numeric
) ;

void umfpack_l_free_numeric
(
    void **Numeric
) ;

int Syntax:

    #include "umfpack.h"
    void *Numeric ;
    umfpack_free_numeric (&Numeric) ;

long Syntax:

    #include "umfpack.h"
    void *Numeric ;
    umfpack_l_free_numeric (&Numeric) ;

Purpose:

    Deallocates the Numeric object and sets the Numeric handle to NULL.
    This routine is the only valid way of destroying the Numeric object;
    any other action (such as using "free (Numeric) ;" or not freeing Numeric
    at all) will lead to memory leaks.

Arguments:

    void **Numeric ;            Input argument, deallocated and Numeric is
                                set to (void *) NULL on output.

        Numeric must point to a valid Numeric object, computed by
        umfpack_numeric.  No action is taken if Numeric is a (void *) NULL
        pointer.
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
\section{Alternatives routines}
\label{Alternative}
%-------------------------------------------------------------------------------

\subsection{umfpack\_defaults and umfpack\_l\_defaults}
\label{defaults}

{\footnotesize
\begin{verbatim}
void umfpack_defaults
(
    double Control [UMFPACK_CONTROL]
) ;

void umfpack_l_defaults
(
    double Control [UMFPACK_CONTROL]
) ;

int Syntax:

    #include "umfpack.h"
    double Control [UMFPACK_CONTROL] ;
    umfpack_defaults (Control) ;

long Syntax:

    #include "umfpack.h"
    double Control [UMFPACK_CONTROL] ;
    umfpack_l_defaults (Control) ;

Purpose:

    Sets the default control parameter settings.

Arguments:

    double Control [UMFPACK_CONTROL] ;  Output argument.

        Control is set to the default control parameter settings.  You can
        then modify individual settings by changing specific entries in the
        Control array.  If Control is a (double *) NULL pointer, then
        umfpack_defaults returns silently (no error is generated, since
        passing a NULL pointer for Control to any UMFPACK routine is valid).
\end{verbatim}
}

\newpage
\subsection{umfpack\_qsymbolic and umfpack\_l\_qsymbolic}

{\footnotesize
\begin{verbatim}
int umfpack_qsymbolic
(
    int n,
    const int Ap [ ],
    const int Ai [ ],
    const int Qinit [ ],
    void **Symbolic,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO]
) ;

long umfpack_l_qsymbolic
(
    long n,
    const long Ap [ ],
    const long Ai [ ],
    const long Qinit [ ],
    void **Symbolic,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO]
) ;

int Syntax:

    #include "umfpack.h"
    void *Symbolic ;
    int n, *Ap, *Ai, *Qinit, status ;
    double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO] ;
    status = umfpack_qsymbolic (n, Ap, Ai, Qinit, &Symbolic, Control, Info) ;

long Syntax:

    #include "umfpack.h"
    void *Symbolic ;
    long n, *Ap, *Ai, *Qinit, status ;
    double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO] ;
    status = umfpack_l_qsymbolic (n, Ap, Ai, Qinit, &Symbolic, Control, Info) ;

Purpose:

    Given the nonzero pattern of a sparse matrix A in column-oriented form, and
    a sparsity preserving column preordering Qinit, umfpack_qsymbolic performs
    the symbolic factorization of A*Qinit (or A (:,Qinit) in Matlab notation).
    It also computes the column elimination tree post-ordering.  This is
    identical to umfpack_symbolic, except that colamd is not called and the
    user input column order Qinit is used instead.  Note that in general, the
    Qinit passed to umfpack_qsymbolic will differ from the final Q found in
    umfpack_numeric, because of the column etree postordering done in
    umfpack_qsymbolic and sparsity-preserving modifications made within each
    frontal matrix during umfpack_numeric.

    *** WARNING ***  A poor choice of Qinit can easily cause umfpack_numeric to
    use a huge amount of memory and do a lot of work.  The "default" symbolic
    analysis method is umfpack_symbolic, not this routine.  If you use this
    routine, the performance of UMFPACK is your responsibility;  UMFPACK will
    not try to second-guess a poor choice of Qinit.  If you are unsure about
    the quality of your Qinit, then call both umfpack_symbolic and
    umfpack_qsymbolic, and pick the one with lower estimates of work and
    memory usage (Info [UMFPACK_FLOPS_ESTIMATE] and
    Info [UMFPACK_PEAK_MEMORY_ESTIMATE]).  Don't forget to call
    umfpack_free_symbolic to free the Symbolic object that you don't need.

Returns:

    The value of Info [UMFPACK_STATUS]; see below.

Arguments:

    All arguments are the same as umfpack_symbolic, except for the following:

    Int Qinit [n] ;             Input argument, not modified.

        The user's fill-reducing initial column preordering.  This must be a
        permutation of 0..n-1.  If Qinit [k] = j, then column j is the kth
        column of the matrix A (:,Qinit) to be factorized.  If Qinit is an
        (Int *) NULL pointer, then colamd is called instead.  In fact,

        Symbolic = umfpack_symbolic (n, Ap, Ai, Control, Info) ;

        is identical to

        Symbolic = umfpack_qsymbolic (n, Ap, Ai, (Int *) NULL, Control, Info) ;

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        Identical to umfpack_symbolic if Qinit is (Int *) NULL.  Otherwise,
        if Qinit is present, it is identical to umfpack_symbolic except for the
        following:

        Control [UMFPACK_DENSE_ROW]:  ignored.

        Control [UMFPACK_DENSE_COL]:  Let j be the leftmost column in
            A (:,Qinit) with more entries than the value determined by the
            dense column control parameter (see umfpack_symbolic), or j=n if
            there is no such column.  Columns j through n-1 are all treated as
            "dense", and factorized in a (n-j)-by-n dense array.  When
            determining Qinit, be sure the "dense" columns of A (:,Qinit) are
            as far to the right as possible.

    double Info [UMFPACK_INFO] ;        Output argument, not defined on input.

        Identical to umfpack_symbolic if Qinit is (Int *) NULL.  Otherwise,
        if Qinit is present, it is identical to umfpack_symbolic except for the
        following:

        Info [UMFPACK_NDENSE_ROW]:  zero
        Info [UMFPACK_NEMPTY_ROW]:  zero
        Info [UMFPACK_NDENSE_COL]:  n-j, where j is defined above.
        Info [UMFPACK_NEMPTY_COL]:  zero
\end{verbatim}
}

\newpage
\subsection{umfpack\_wsolve and umfpack\_l\_wsolve}

{\footnotesize
\begin{verbatim}
int umfpack_wsolve
(
    const char sys [ ],
    const int Ap [ ],
    const int Ai [ ],
    const double Ax [ ],
    double X [ ],
    const double B [ ],
    void *Numeric,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO],
    int Wi [ ],
    double W [ ],
    double Y [ ],
    double Z [ ],
    double S [ ]
) ;

long umfpack_l_wsolve
(
    const char sys [ ],
    const long Ap [ ],
    const long Ai [ ],
    const double Ax [ ],
    double X [ ],
    const double B [ ],
    void *Numeric,
    const double Control [UMFPACK_CONTROL],
    double Info [UMFPACK_INFO],
    long Wi [ ],
    double W [ ],
    double Y [ ],
    double Z [ ],
    double S [ ]
) ;

int Syntax:

    #include "umfpack.h"
    void *Numeric ;
    int status, *Ap, *Ai ;
    char *sys ;
    double *B, *X, *Ax, Info [UMFPACK_INFO], Control [UMFPACK_CONTROL] ;
    int *Wi ;
    double *W, *Y, *Z, *S ;
    status = umfpack_wsolve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info,
        Wi, W, Y, Z, S) ;

long Syntax:

    #include "umfpack.h"
    void *Numeric ;
    long status, *Ap, *Ai ;
    char *sys ;
    double *B, *X, *Ax, Info [UMFPACK_INFO], Control [UMFPACK_CONTROL] ;
    long *Wi ;
    double *W, *Y, *Z, *S ;
    status = umfpack_l_wsolve (sys, Ap, Ai, Ax, X, B, Numeric, Control, Info,
        Wi, W, Y, Z, S) ;

Purpose:

    Given LU factors computed by umfpack_numeric (PAQ=LU) and the
    right-hand-side, B, solve a linear system for the solution X.  Iterative
    refinement is optionally performed.  This routine is identical to
    umfpack_solve, except that it does not dynamically allocate any workspace.
    When you have many linear systems to solve, this routine is slightly faster
    than umfpack_solve, since the workspace (Wi, W, Y, Z, and S) needs to be
    allocated only once, prior to calling umfpack_wsolve.

Returns:

    The status code is returned.  See Info [UMFPACK_STATUS], below.

Arguments:

    char sys [ ] ;              Input argument, not modified.
    Int Ap [n+1] ;              Input argument, not modified.
    Int Ai [nz] ;               Input argument, not modified.
    double Ax [nz] ;            Input argument, not modified.
    double X [n] ;              Output argument.
    double B [n] ;              Input argument, not modified.
    void *Numeric ;             Input argument, not modified.
    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.
    double Info [UMFPACK_INFO] ;        Output argument.

        The above arguments are identical to umfpack_solve, except that the
        error code UMFPACK_ERROR_out_of_memory will not be returned in
        Info [UMFPACK_STATUS], since umfpack_wsolve does not allocate any
        memory.

    Int Wi [2*n] ;              Workspace.
    double W [n] ;              Workspace.
    double Y [n] ;              Workspace, only needed for iterative refinement.
    double Z [n] ;              Workspace, only needed for iterative refinement.
    double S [n] ;              Workspace, only needed for iterative refinement.

        The Wi, W, Y, Z, and S arguments are workspace used by umfpack_wsolve.
        Their contents are undefined on output.
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
\section{Matrix manipulation routines}
\label{Manipulate}
%-------------------------------------------------------------------------------

\subsection{umfpack\_col\_to\_triplet and umfpack\_l\_col\_to\_triplet}

{\footnotesize
\begin{verbatim}
int umfpack_col_to_triplet
(
    int n,
    const int Ap [ ],
    int Tj [ ]
) ;

long umfpack_l_col_to_triplet
(
    long n,
    const long Ap [ ],
    long Tj [ ]
) ;

int Syntax:

    #include "umfpack.h"
    int n, *Tj, *Ap, status ;
    status = umfpack_col_to_triplet (n, Ap, Tj) ;

long Syntax:

    #include "umfpack.h"
    long n, *Tj, *Ap, status ;
    status = umfpack_l_col_to_triplet (n, Ap, Tj) ;

Purpose:

    Converts a column-oriented matrix to a triplet form.  Only the column
    pointers, Ap, are required, and only the column indices of the triplet form
    are constructed.   This routine is the opposite of umfpack_triplet_to_col.
    The matrix may be singular.

Returns:

    UMFPACK_OK if successful
    UMFPACK_ERROR_argument_missing if Ap or Tj is missing
    UMFPACK_ERROR_n_nonpositive if n <= 0
    UMFPACK_ERROR_Ap0_nonzero if Ap [0] != 0
    UMFPACK_ERROR_nz_negative if Ap [n] < 0
    UMFPACK_ERROR_col_length_negative if Ap [j] > Ap [j+1] for any j in the
        range 0 to n-1.
    Empty rows, unsorted columns, and duplicate entries do not cause an error
    (these would only be evident by examining Ai).  Empty columns are OK.

Arguments:

    Int n ;             Input argument, not modified.

        A is an n-by-n matrix.  Restriction: n > 0.

    Int Ap [n+1] ;      Input argument, not modified.

        The column pointers of the column-oriented form of the matrix.  See
        umfpack_symbolic for a description.  The number of entries in
        the matrix is nz = Ap [n].  Restrictions on Ap are the same as those
        for umfpack_transpose.  Ap [0] must be zero, nz must be >= 0, and
        Ap [j] <= Ap [j+1] and Ap [j] <= Ap [n] must be true for all j in the
        range 0 to n-1.  Empty columns are OK (that is, Ap [j] may equal
        Ap [j+1] for any j in the range 0 to n-1).

    Int Tj [nz] ;       Output argument.

        Tj is an integer array of size nz on input, where nz = Ap [n].
        Suppose the column-form of the matrix is held in Ap, Ai, and Ax
        (see umfpack_symbolic for a description).  Then on output, the triplet
        form of the same matrix is held in Ai (row indices), Tj (column
        indices), and Ax (numerical values).  Note, however, that this routine
        does not require Ai and Ax in order to do the conversion.
\end{verbatim}
}

\newpage
\subsection{umfpack\_triplet\_to\_col and umfpack\_l\_triplet\_to\_col}

{\footnotesize
\begin{verbatim}
int umfpack_triplet_to_col
(
    int n,
    int nz,
    const int Ti [ ],
    const int Tj [ ],
    const double Tx [ ],
    int Bp [ ],
    int Bi [ ],
    double Bx [ ]
) ;

long umfpack_l_triplet_to_col
(
    long n,
    long nz,
    const long Ti [ ],
    const long Tj [ ],
    const double Tx [ ],
    long Bp [ ],
    long Bi [ ],
    double Bx [ ]
) ;

int Syntax:

    #include "umfpack.h"
    int n, nz, *Ti, *Tj, *Bp, *Bi, status ;
    double *Tx, *Bx ;
    status = umfpack_triplet_to_col (n, nz, Ti, Tj, Tx, Bp, Bi, Bx) ;

long Syntax:

    #include "umfpack.h"
    long n, nz, *Ti, *Tj, *Bp, *Bi, status ;
    double *Tx, *Bx ;
    status = umfpack_l_triplet_to_col (n, nz, Ti, Tj, Tx, Bp, Bi, Bx) ;

Purpose:

    Converts a sparse matrix from "triplet" form to compressed-column form.

    The triplet form of a matrix is a very simple data structure for basic
    sparse matrix operations.  For example, suppose you wish to factorize a
    matrix A coming from a finite element method, in which A is a sum of
    dense submatrices, A = E1 + E2 + E3 + ... .  The entries in each element
    matrix Ei can be concatenated together in the three triplet arrays, and
    any overlap between the elements will be correctly summed by
    umfpack_triplet_to_col.

    Transposing a matrix in triplet form is simple; just interchange the
    use of Ti and Tj.

    Permuting a matrix in triplet form is also simple.  If you want the matrix
    PAQ, or A (P,Q) in Matlab notation, where P [k] = i means that row i of
    A is the kth row of PAQ and Q [k] = j means that column j of A is the kth
    column of PAQ, then do the following.  First, create inverse permutations
    Pinv and Qinv such that Pinv [i] = k if P [k] = i and Qinv [j] = k if
    Q [k] = j.  Next, for the mth triplet (Ti [m], Tj [m], Tx [m]), replace
    Ti [m] with Pinv [Ti [m]] and replace Tj [m] with Qinv [Tj [m]].

    If you have a column-form matrix with duplicate entries or unsorted
    columns, you can sort it and sum up the duplicates by first converting it
    to triplet form with umfpack_col_to_triplet, and then coverting it back
    with umfpack_triplet_to_col.

    You can print the input triplet form with umfpack_report_triplet, and
    the output matrix with umfpack_report_matrix.

    The matrix may be singular (nz can be zero, and empty rows and/or columns
    may exist).

Returns:

    UMFPACK_OK if successful.
    UMFPACK_ERROR_argument_missing if Bp, Bi, Ti, and/or Tj are missing.
    UMFPACK_ERROR_n_nonpositive if n <= 0.
    UMFPACK_ERROR_nz_negative if nz < 0.
    UMFPACK_ERROR_invalid_triplet if for any k, Ti [k] and/or Tj [k] are not in
        the range 0 to n-1.
    UMFPACK_ERROR_out_of_memory if unable to allocate sufficient workspace.

Arguments:

    Int n ;             Input argument, not modified.

        A is an n-by-n matrix.  Restriction: n > 0.  All row and column indices
        in the triplet form must be in the range 0 to n-1.

    Int nz ;            Input argument, not modified.

        The number of entries in the triplet form of the matrix.  Restriction:
        nz >= 0.

    Int Ti [nz] ;       Input argument, not modified.
    Int Tj [nz] ;       Input argument, not modified.
    double Tx [nz] ;    Input argument, not modified.

        Ti, Tj, and Tx hold the "triplet" form of a sparse matrix.  The kth
        nonzero entry is in row i = Ti [k], column j = Tj [k], and has a
        numerical value of a_ij = Tx [k].  The row and column indices i and j
        must be in the range 0 to n-1.  Duplicate entries may be present; they
        are summed in the output matrix.  This is not an error condition.
        The "triplets" may be in any order.  Tx is optional; if Tx and/or Bx
        are not present (a (double *) NULL pointer), then Bx is not computed.

    Int Bp [n+1] ;      Output argument, not modified.

        Bp is an integer array of size n+1 on input.
        On output, Bp holds the "pointers" for the column form of the sparse
        matrix A.  Column j of the matrix A is held in
        Bi [(Bp [j]) ... (Bp [j+1]-1)].  The first entry, Bp [0], is zero, and
        Bp [j] <= Bp [j+1] holds for all j in the range 0 to n-1.  The value
        nz2 = Bp [n] is thus the total number of entries in the pattern of the
        matrix A.  Equivalently, the number of duplicate triplets is
        nz - Bp [n].

    Int Bi [max(n,nz2)] ;       Output argument, not modified.

        Bi is an integer array of size max (n,nz2) on input, where nz2 <= nz.
        Bi is also used as workspace during the conversion, and for this use
        the size of Bi must also be at least n.

        The nonzero pattern (row indices) for column j is stored in
        Bi [(Bp [j]) ... (Bp [j+1]-1)].  The row indices in a given column j
        are in ascending order, and no duplicate row indices are present.
        Row indices are in the range 0 to n-1 (the matrix is 0-based).

    double Bx [nz2] ;   Output argument, not modified, of size nz2 = Bp [n].

        Bx is a double array of size nz2 on input, where nz2 <= nz.  Bx is
        optional; if Tx and/or Bx are not present (a (double *) NULL pointer),
        then Bx is not computed.  If present, Bx holds the numerical values of
        the sparse matrix A.  The nonzero pattern (row indices) for column j is
        stored in Bi [(Bp [j]) ... (Bp [j+1]-1)], and the corresponding
        numerical values are stored in Bx [(Bp [j]) ... (Bp [j+1]-1)].
\end{verbatim}
}

\newpage
\subsection{umfpack\_transpose and umfpack\_l\_transpose}

{\footnotesize
\begin{verbatim}
int umfpack_transpose
(
    int n,
    const int Ap [ ],
    const int Ai [ ],
    const double Ax [ ],
    const int P [ ],
    const int Q [ ],
    int Cp [ ],
    int Ci [ ],
    double Cx [ ]
) ;

long umfpack_l_transpose
(
    long n,
    const long Ap [ ],
    const long Ai [ ],
    const double Ax [ ],
    const long P [ ],
    const long Q [ ],
    long Cp [ ],
    long Ci [ ],
    double Cx [ ]
) ;

int Syntax:

    #include "umfpack.h"
    int n, status, *Ap, *Ai, *P, *Q, *Cp, *Ci ;
    double *Ax, *Cx ;
    status = umfpack_transpose (n, Ap, Ai, Ax, P, Q, Cp, Ci, Cx) ;

long Syntax:

    #include "umfpack.h"
    long n, status, *Ap, *Ai, *P, *Q, *Cp, *Ci ;
    double *Ax, *Cx ;
    status = umfpack_l_transpose (n, Ap, Ai, Ax, P, Q, Cp, Ci, Cx) ;

Purpose:

    Transposes and optionally permutes a sparse matrix in row or column-form,
    C = (PAQ)'.  In Matlab notation, C = (A (P,Q))'.  Alternatively, this
    routine can be viewed as converting A (P,Q) from column-form to row-form,
    or visa versa.  Empty rows and columns may exist.
    The matrix A may be singular.

    umfpack_transpose is useful if you want to factorize A' instead of A.
    Factorizing A' instead of A can be much better, particularly if AA' is much
    sparser than A'A.  You can still solve Ax=b if you factorize A', by solving
    with the sys argument "A'x=b" in umfpack_*solve.

    The input A and output C can be printed with umfpack_report matrix, and the
    permutation vectors can be printed with umfpack_report_perm.

Returns:

    UMFPACK_OK if successful.
    UMFPACK_ERROR_out_of_memory if umfpack_transpose fails to allocate a
        size-n workspace.
    UMFPACK_ERROR_argument_missing if Ai, Ap, Ci, and/or Cp are missing.
    UMFPACK_ERROR_n_nonpositive if n <= 0.
    UMFPACK_ERROR_invalid_permutation if P and/or Q are invalid.
    UMFPACK_ERROR_nz_negative if Ap [n] < 0.
    UMFPACK_ERROR_Ap0_nonzero if Ap [0] != 0.
    UMFPACK_ERROR_col_length_negative if Ap [j] > Ap [j+1] for any j in the
        range 0 to n-1.
    UMFPACK_ERROR_row_index_out_of_bounds if any row index i is < 0 or >= n.
    UMFPACK_ERROR_jumbled_matrix if the row indices in any column are not in
        ascending order.

Arguments:

    Int n ;             Input argument, not modified.

        A is an n-by-n matrix.  Restriction: n > 0.

    Int Ap [n+1] ;      Input argument, not modified.

        The column pointers of the column-oriented form of the matrix A.  See
        umfpack_symbolic for a description.  The number of entries in
        the matrix is nz = Ap [n].  Ap [0] must be zero, Ap [n] must be > 0, and
        Ap [j] <= Ap [j+1] and Ap [j] <= Ap [n] must be true for all j in the
        range 0 to n-1.  Empty columns are OK (that is, Ap [j] may equal Ap
        [j+1] for any j in the range 0 to n-1).

    Int Ai [nz] ;       Input argument, not modified, of size nz = Ap [n].

        The nonzero pattern (row indices) for column j is stored in
        Ai [(Ap [j]) ... (Ap [j+1]-1)].  The row indices in a given column j
        must be in ascending order, and no duplicate row indices may be present.
        Row indices must be in the range 0 to n-1 (the matrix is 0-based).

    double Ax [nz] ;    Input argument, not modified, of size nz = Ap [n].

        If present, these are the numerical values of the sparse matrix A.
        The nonzero pattern (row indices) for column j is stored in
        Ai [(Ap [j]) ... (Ap [j+1]-1)], and the corresponding numerical values
        are stored in Ax [(Ap [j]) ... (Ap [j+1]-1)].  If Ax and/or Cx are not
        present, then the output Cx [...] is not computed, and only the pattern
        is transposed.  This is not an error condition.

    Int P [n] ;         Input argument, not modified.

        The permutation vector P is defined as P [k] = i, where the original
        row i of A is the kth row of PAQ.  If you want to use the identity
        permutation for P, simply pass (Int *) NULL for P.  This is not an error
        condition.

    Int Q [n] ;         Input argument, not modified.

        The permutation vector Q is defined as Q [k] = j, where the original
        column j of A is the kth column of PAQ.  If you want to use the identity
        permutation for Q, simply pass (Int *) NULL for Q.  This is not an error
        condition.

    Int Cp [n+1] ;      Output argument.

        The column pointers of the matrix C = (A (P,Q))', in the same form
        as the column pointers Ap for the matrix A.

    Int Ci [nz] ;       Output argument.

        The row indices of the matrix C = (A (P,Q))', in the same form
        as the row indices Ai for the matrix A.

    double Cx [nz] ;    Output argument.

        If present, these are the numerical values of the sparse matrix C,
        in the same form as the values Ax of the matrix A.  If Ax and/or Cx
        are not present, then the output Cx [...] is not computed, and only
        the pattern is transposed.  This is not an error condition.
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
\section{Getting the contents of opaque objects}
\label{Get}
%-------------------------------------------------------------------------------

\subsection{umfpack\_get\_lunz and umfpack\_l\_get\_lunz}

{\footnotesize
\begin{verbatim}
int umfpack_get_lunz
(
    int *lnz,
    int *unz,
    int *n,
    void *Numeric
) ;

long umfpack_l_get_lunz
(
    long *lnz,
    long *unz,
    long *n,
    void *Numeric
) ;

int Syntax:

    #include "umfpack.h"
    void *Numeric ;
    int status, lnz, unz, n ;
    status = umfpack_get_lunz (&lnz, &unz, &n, Numeric) ;

long Syntax:

    #include "umfpack.h"
    void *Numeric ;
    long status, lnz, unz, n ;
    status = umfpack_l_get_lunz (&lnz, &unz, &n, Numeric) ;

Purpose:

    Determines the size and number of nonzeros in the LU factors held by the
    Numeric object.  These are also the sizes of the output arrays required
    by umfpack_get_numeric.

    This routine may seem redundant, since n is a value originally passed to
    umfpack_symbolic by the user, and lnz and unz are available from the Info
    array.  However, the Info array is not always returned (its use is
    optional).  This routine is also useful in the context of many sparse linear
    systems, with many Numeric handles.  The user could store an array of
    Numeric objects in an array of (void *) pointers, for example.  The Info
    array is re-initialized each time an UMFPACK routine is called, and thus the
    lnz and unz information could be lost.  Lnz and unz can differ from
    different calls to umfpack_numeric with different numerical values (Ax),
    even when using the same Symbolic object.   This routine allows the LU
    factors to be extracted from the Numeric object (with umfpack_get_numeric)
    without the use of the corresponding Info array.

Returns:

    UMFPACK_OK if successful.
    UMFPACK_ERROR_invalid_Numeric_object if Numeric is not a valid object.
    UMFPACK_ERROR_argument_missing if lnz, unz, or n are (Int *) NULL

Arguments:

    Int *lnz ;          Output argument.

        The number of nonzeros in L, including the diagonal (which is all
        one's).  This value is the required size of the Li and Lx arrays as
        computed by umfpack_get_numeric.  The value of lnz is identical to
        Info [UMFPACK_LNZ], if that value was returned by umfpack_numeric.

    Int *unz ;          Output argument.

        The number of nonzeros in U, including the diagonal.  This value is the
        required size of the Ui and Ux arrays as computed by
        umfpack_get_numeric.  The value of unz is identical to
        Info [UMFPACK_UNZ], if that value was returned by umfpack_numeric.

    Int *n ;            Output argument.

        The order of the L and U matrices.  The size of Lp and Up, as required
        by umfpack_get_numeric, is n+1, and the size of P and Q are n.  The
        value of n is identical to that passed to umfpack_symbolic.

    void *Numeric ;     Input argument, not modified.

        Numeric must point to a valid Numeric object, computed by
        umfpack_numeric.
\end{verbatim}
}

\newpage
\subsection{umfpack\_get\_numeric and umfpack\_l\_get\_numeric}

{\footnotesize
\begin{verbatim}
int umfpack_get_numeric
(
    int Lp [ ],
    int Li [ ],
    double Lx [ ],
    int Up [ ],
    int Ui [ ],
    double Ux [ ],
    int P [ ],
    int Q [ ],
    void *Numeric
) ;

long umfpack_l_get_numeric
(
    long Lp [ ],
    long Li [ ],
    double Lx [ ],
    long Up [ ],
    long Ui [ ],
    double Ux [ ],
    long P [ ],
    long Q [ ],
    void *Numeric
) ;

int Syntax:

    #include "umfpack.h"
    void *Numeric ;
    int *Lp, *Li, *Up, *Ui, *P, *Q, status ;
    double *Lx, *Ux ;
    status = umfpack_get_numeric (Lp, Li, Lx, Up, Ui, Ux, P, Q, Numeric) ;

long Syntax:

    #include "umfpack.h"
    void *Numeric ;
    long *Lp, *Li, *Up, *Ui, *P, *Q, status ;
    double *Lx, *Ux ;
    status = umfpack_l_get_numeric (Lp, Li, Lx, Up, Ui, Ux, P, Q, Numeric) ;

Purpose:

    This routine copies the LU factors and permutation vectors from the Numeric
    object into user-accessible arrays.  This routine is not needed to solve a
    linear system.  Note that the output arrays Lp, Li, Lx, Up, Ui, Ux, P, and
    Q are not allocated by umfpack_get_numeric; they must exist on input.

Returns:

    Returns UMFPACK_OK if successful.  Returns UMFPACK_ERROR_out_of_memory
    if insufficient memory is available for the 2*n integer workspace that
    UMFPACK_get_numeric allocates to construct L and/or U.  Returns
    UMFPACK_ERROR_invalid_Numeric_object if the Numeric object provided as
    input is invalid.

Arguments:

    Int Lp [n+1] ;      Output argument.
    Int Li [lnz] ;      Output argument.
    double Lx [lnz] ;   Output argument.

        The matrix L is returned in compressed-row form.  The column indices
        of row i and corresponding numerical values are in:

            Li [Lp [i] ... Lp [i+1]-1]
            Lx [Lp [i] ... Lp [i+1]-1]

        respectively.  Each row is stored in sorted order, from low column
        indices to higher.  The last entry in each row is the diagonal, which
        is numerically equal to one.  The sizes of Lp, Li, and Lx are returned
        by umfpack_get_lunz.  If any one of the Lp, Li, or Lx arrays are not
        present, the L matrix is not returned.  This is not an error condition.
        Thus, if you do not want the L matrix returned, simply pass
        (Int *) NULL for Lp and Li, and (double *) NULL for Lx.
        The L matrix can be printed if n, Lp, Li, and Lx are passed to
        umfpack_report_matrix (using the "row" form).

    Int Up [n+1] ;      Output argument.
    Int Ui [unz] ;      Output argument.
    double Ux [unz] ;   Output argument.

        The matrix U is returned in compressed-column form.  The row indices
        of column j and corresponding numerical values are in

            Ui [Up [j] ... Up [j+1]-1]
            Ux [Up [j] ... Up [j+1]-1]

        respectively.  Each column is stored in sorted order, from low row
        indices to higher.  The last entry in each column is the diagonal.  The
        sizes of Up, Ui, and Ux are returned by umfpack_get_lunz.  If any one of
        the Up, Ui, or Ux arrays are not present, the U matrix is not returned.
        This is not an error condition.  Thus, if you do not want the U matrix
        returned, simply pass (Int *) NULL for Up and Ui, and (double *) NULL
        for Ux.  The U matrix can be printed if n, Up, Ui, and Ux are passed to
        umfpack_report_matrix (using the "column" form).

    Int P [n] ;         Output argument.

        The permutation vector P is defined as P [k] = i, where the original
        row i of A is the kth pivot row in PAQ.  If you do not want the P vector
        to be returned, simply pass (Int *) NULL for P.  This is not an error
        condition.  You can print P and Q with umfpack_report_perm.

    Int Q [n] ;         Output argument.

        The permutation vector Q is defined as Q [k] = j, where the original
        column j of A is the kth pivot column in PAQ.  If you not want the Q
        vector to be returned, simply pass (Int *) NULL for Q.  This is not
        an error condition.  Note that Q is not necessarily identical to
        Qtree, the column preordering held in the Symbolic object.  Refer to
        the description of Qtree and Front_npivots in umfpack_get_symbolic for
        details.

    void *Numeric ;     Input argument, not modified.

        Numeric must point to a valid Numeric object, computed by
        umfpack_numeric.
\end{verbatim}
}

\newpage
\subsection{umfpack\_get\_symbolic and umfpack\_l\_get\_symbolic}

{\footnotesize
\begin{verbatim}
int umfpack_get_symbolic
(
    int *n,
    int *nz,
    int *nfr,
    int *nchains,
    int *nsparse_col,
    int Qtree [ ],
    int Front_npivots [ ],
    int Front_parent [ ],
    int Chain_start [ ],
    int Chain_maxrows [ ],
    int Chain_maxcols [ ],
    void *Symbolic
) ;

long umfpack_l_get_symbolic
(
    long *n,
    long *nz,
    long *nfr,
    long *nchains,
    long *nsparse_col,
    long Qtree [ ],
    long Front_npivots [ ],
    long Front_parent [ ],
    long Chain_start [ ],
    long Chain_maxrows [ ],
    long Chain_maxcols [ ],
    void *Symbolic
) ;

int Syntax:

    #include "umfpack.h"
    int status, n, nz, nfr, nchains, nsparse_col, *Qtree,
        *Front_npivots, *Front_parent, *Chain_start, *Chain_maxrows,
        *Chain_maxcols ;
    void *Symbolic ;
    status = umfpack_get_symbolic (&n, &nz, &nfr, &nchains, &nsparse_col,
        Qtree, Front_npivots, Front_parent, Chain_start, Chain_maxrows,
        Chain_maxcols, Symbolic) ;

long Syntax:

    #include "umfpack.h"
    long status, n, nz, nfr, nchains, nsparse_col, *Qtree,
        *Front_npivots, *Front_parent, *Chain_start, *Chain_maxrows,
        *Chain_maxcols ;
    void *Symbolic ;
    status = umfpack_l_get_symbolic (&n, &nz, &nfr, &nchains, &nsparse_col,
        Qtree, Front_npivots, Front_parent, Chain_start, Chain_maxrows,
        Chain_maxcols, Symbolic) ;

Purpose:

    Copies the contents of the Symbolic object into simple integer arrays
    accessible to the user.  This routine is not needed to factorize
    and/or solve a sparse linear system using UMFPACK.  Note that the output
    arrays Qtree, Front_npivots, Front_parent, Chain_start, Chain_maxrows,
    and Chain_maxcols are not allocated by umfpack_get_symbolic;  they must
    exist on input.

    The Symbolic object is small.  It size, in integers, is
    (3*nchains + n + 2*nfr + 20), which is no greater than 6*n+20.  The object
    holds the initial column permutation, the supernodal column elimination
    tree, and information about each frontal matrix.  You can print it with
    umfpack_report_symbolic.

Returns:

    Returns UMFPACK_OK if successful, UMFPACK_ERROR_invalid_Symbolic_object
    if Symbolic is an invalid object.

Arguments:

    Note that if any of the output arguments are (Int *) NULL pointers, then
    that argument is not returned.  This is not an error condition.  Thus,
    if you do not want a particular component of the Symbolic object to be
    returned to you, simply pass a (Int *) NULL pointer for that particular
    output argument.

    Int *n ;            Output argument.

        The dimension of the matrix A analyzed by the call to umfpack_symbolic
        that generated the Symbolic object.

    Int *nz ;           Output argument.

        The number of nonzeros in A.

    Int *nfr ;  Output argument.

        The number of frontal matrices that will be used by umfpack_numeric
        to factorize the matrix A.  One or more pivots are contained in each
        frontal matrix, and the total number of pivots in the frontal matrices
        is n (see the description of nsparse_col, below).
        Thus, nfr is in the range 1 to n.

    Int *nchains ;      Output argument.

        The frontal matrices are related to one another by the supernodal
        column elimination tree.  Each node in this tree is one frontal matrix.
        The tree is partitioned into a set of disjoint paths, and a frontal
        matrix chain is one path in this tree.  Each chain is factorized using
        a unifrontal technique, with a single working array that holds each
        frontal matrix in the chain, one at a time.  nchains is in the range
        1 to nfr.

    Int *nsparse_col ;  Output argument.

        This is equal to n.  It differed from n in Version 3.0.

    Int Qtree [n] ;     Output argument.

        The initial column permutation.  If Qtree [k] = j, then this means that
        column j is the kth pivot column in the preordered matrix.
        Qtree is not necessarily the same as the final
        column permutation Q, computed by umfpack_numeric.  The numeric
        factorization may reorder the pivot columns within each frontal matrix
        to reduce fill-in.

    Int Front_npivots [nfr] ;   Output argument.

        This array should be of size at least n, in order to guarantee that it
        will be large enough to hold the output.  Only the first nfr entries
        are used, however.  The kth frontal matrix holds Front_npivots [k] pivot
        columns.  Thus, the first frontal matrix, front 0, is used to factorize
        the first Front_npivots [0] columns; these correspond to the original
        columns Qtree [0] through Qtree [Front_npivots [0]-1].  The next frontal
        matrix is used to factorize the next Front_npivots [1] columns, which
        are thus the original columns Qtree [Front_npivots [0]] through
        Qtree [Front_npivots [0] + Front_npivots [1] - 1], and so on.
        The sum of Front_npivots [0..nfr-1] is equal to n.

        Any modifications that umfpack_numeric makes to the initial column
        permutation are constrained to within each frontal matrix.  Thus,
        for the first frontal matrix, Q [0] through Q [Front_npivots [0]-1] is
        some permutation of the columns Qtree [0] through
        Qtree [Front_npivots [0]-1].  For second frontal matrix,
        Q [Front_npivots [0]] through Q [Front_npivots [0] + Front_npivots[1]-1]
        is some permutation of the same portion of Qtree, and so on.  All pivot
        columns are numerically factorized within the frontal matrix originally
        determined by the symbolic factorization; there is no delayed pivoting
        across frontal matrices.

    Int Front_parent [nfr] ;    Output argument.

        This array should be of size at least n, in order to guarantee that it
        will be large enough to hold the output.  Only the first nfr entries
        are used, however.  Front_parent [0..nfr-1] holds the supernodal column
        elimination tree.  Each node in the tree corresponds to a single frontal
        matrix.  The parent of node f is Front_parent [f].

    Int Chain_start [nchains+1] ;       Output argument.

        This array should be of size at least n+1, in order to guarantee that it
        will be large enough to hold the output.  Only the first nchains+1
        entries are used, however.  The kth frontal matrix chain consists of
        frontal matrices Chain_start [k] through Chain_start [k+1]-1.  Thus,
        Chain_start [0] is always 0, and Chain_start [nchains] is the total
        number of frontal matrices, nfr.  For two adjacent fronts f and f+1
        within a single chain, f+1 is always the parent of f (that is,
        Front_parent [f] = f+1).

    Int Chain_maxrows [nchains] ;       Output argument.
    Int Chain_maxcols [nchains] ;       Output argument.

        These arrays should be of size at least n, in order to guarantee that
        they will be large enough to hold the output.  Only the first nchains
        entries of Chain_maxrows and Chain_maxcols are used, however.  The kth
        frontal matrix chain requires a single working array of dimension
        Chain_maxrows [k] by Chain_maxcols [k], for the unifrontal technique
        that factorizes the frontal matrix chain.  Since the symbolic
        factorization only provides an upper bound on the size of each frontal
        matrix, not all of the working array is necessarily used during the
        numerical factorization.

        Note that the upper bound on the number of rows and columns of each
        frontal matrix is computed by umfpack_symbolic, but all that is
        required by umfpack_numeric is the maximum of these two sets of
        values for each frontal matrix chain.  Thus, the size of each
        individual frontal matrix is not preserved in the Symbolic object.

    void *Symbolic ;                    Input argument, not modified.

        The Symbolic object, which holds the symbolic factorization computed by
        umfpack_symbolic.  The Symbolic object is not modified by
        umfpack_get_symbolic.
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
\section{Reporting routines}
\label{Report}
%-------------------------------------------------------------------------------

\subsection{umfpack\_report\_status and umfpack\_l\_report\_status}

{\footnotesize
\begin{verbatim}
void umfpack_report_status
(
    const double Control [UMFPACK_CONTROL],
    int status
) ;

void umfpack_l_report_status
(
    const double Control [UMFPACK_CONTROL],
    long status
) ;

int Syntax:

    #include "umfpack.h"
    double Control [UMFPACK_CONTROL] ;
    int status ;
    umfpack_report_status (Control, status) ;

long Syntax:

    #include "umfpack.h"
    double Control [UMFPACK_CONTROL] ;
    long status ;
    umfpack_l_report_status (Control, status) ;

Purpose:

    Prints the status (return value) of other umfpack_* routines.

Arguments:

    double Control [UMFPACK_CONTROL] ;   Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            0 or less: no output, even when an error occurs
            1: error messages only
            2 or more: print status, whether or not an error occured
            4 or more: also print the UMFPACK Copyright
            6 or more: also print the UMFPACK License
            Default: 1

    Int status ;                        Input argument, not modified.

        The return value from another umfpack_* routine.
\end{verbatim}
}

\newpage
\subsection{umfpack\_report\_control and umfpack\_l\_report\_control}

{\footnotesize
\begin{verbatim}
void umfpack_report_control
(
    const double Control [UMFPACK_CONTROL]
) ;

void umfpack_l_report_control
(
    const double Control [UMFPACK_CONTROL]
) ;

int Syntax:

    #include "umfpack.h"
    double Control [UMFPACK_CONTROL] ;
    umfpack_report_control (Control) ;

long Syntax:

    #include "umfpack.h"
    double Control [UMFPACK_CONTROL] ;
    umfpack_l_report_control (Control) ;

Purpose:

    Prints the current control settings.  Note that with the default print
    level, nothing is printed.  Does nothing if Control is (double *) NULL.

Arguments:

    double Control [UMFPACK_CONTROL] ;   Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            1 or less: no output
            2 or more: print all of Control
            Default: 1
\end{verbatim}
}

\newpage
\subsection{umfpack\_report\_info and umfpack\_l\_report\_info}

{\footnotesize
\begin{verbatim}
void umfpack_report_info
(
    const double Control [UMFPACK_CONTROL],
    const double Info [UMFPACK_INFO]
) ;

void umfpack_l_report_info
(
    const double Control [UMFPACK_CONTROL],
    const double Info [UMFPACK_INFO]
) ;

int Syntax:

    #include "umfpack.h"
    double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO] ;
    umfpack_report_info (Control, Info) ;

long Syntax:

    #include "umfpack.h"
    double Control [UMFPACK_CONTROL], Info [UMFPACK_INFO] ;
    umfpack_l_report_info (Control, Info) ;

Purpose:

    Reports statistics from the umfpack_*symbolic, umfpack_numeric, and
    umfpack_*solve routines.

Arguments:

    double Control [UMFPACK_CONTROL] ;   Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            0 or less: no output, even when an error occurs
            1: error messages only
            2 or more: error messages, and print all of Info
            Default: 1

    double Info [UMFPACK_INFO] ;                Input argument, not modified.

        Info is an output argument of several UMFPACK routines.
        The contents of Info are printed on standard output.
\end{verbatim}
}

\newpage
\subsection{umfpack\_report\_matrix and umfpack\_l\_report\_matrix}

{\footnotesize
\begin{verbatim}
int umfpack_report_matrix
(
    const char name [ ],
    int n,
    const int Ap [ ],
    const int Ai [ ],
    const double Ax [ ],
    const char form [ ],
    const double Control [UMFPACK_CONTROL]
) ;

long umfpack_l_report_matrix
(
    const char name [ ],
    long n,
    const long Ap [ ],
    const long Ai [ ],
    const double Ax [ ],
    const char form [ ],
    const double Control [UMFPACK_CONTROL]
) ;

int Syntax:

    #include "umfpack.h"
    int n, *Ap, *Ai, status ;
    double *Ax, Control [UMFPACK_CONTROL] ;
    status = umfpack_report_matrix ("A", n, Ap, Ai, Ax, "column", Control) ;
or:
    status = umfpack_report_matrix ("A", n, Ap, Ai, Ax, "row", Control) ;

long Syntax:

    #include "umfpack.h"
    long n, *Ap, *Ai, status ;
    double *Ax, Control [UMFPACK_CONTROL] ;
    status = umfpack_l_report_matrix ("A", n, Ap, Ai, Ax, "column", Control) ;
or:
    status = umfpack_l_report_matrix ("A", n, Ap, Ai, Ax, "row", Control) ;

Purpose:

    Verifies and prints a row or column-oriented sparse matrix.

Returns:

    UMFPACK_OK if Control [UMFPACK_PRL] <= 2 (the input is not checked).

    Otherwise:

    UMFPACK_OK if the matrix is valid and non-singular.

    UMFPACK_ERROR_singular_matrix if the matrix is structurally singular but
        otherwise valid.  It has one or more rows or columns with no entries.
        This test is made without considering the numerical values, but by just
        looking at the pattern of the entries.  Thus, all structurally singular
        matrices are numerically singular, but not all numerically singular
        matrices are structurally singular.  The matrix may be operated on by
        the matrix manipulation routines (umfpack_transpose,
        umfpack_col_to_triplet) but it may not be analyzed by umfpack_*symbolic
        or factorized by umfpack_numeric).

    UMFPACK_ERROR_n_nonpositive if n <= 0.
    UMFPACK_ERROR_argument_missing if Ap and/or Ai are missing.
    UMFPACK_ERROR_nz_negative if Ap [n] < 0.
    UMFPACK_ERROR_Ap0_nonzero if Ap [0] is not zero.
    UMFPACK_ERROR_col_length_negative if Ap [j+1] < Ap [j] for any j in the
        range 0 to n-1.
    UMFPACK_ERROR_out_of_memory if out of memory.
    UMFPACK_ERROR_row_index_out_of_bounds if any row index in Ai is not in
        the range 0 to n-1.
    UMFPACK_ERROR_jumbled_matrix if the row indices in any column are not in
        ascending order, or contain duplicates.

Arguments:

    char name [ ] ;     Input argument, not modified.

        The name of the matrix.  This is optional; no name is printed if
        a (char *) NULL pointer is passed.

    Int n ;             Input argument, not modified.

        A is an n-by-n matrix.  Restriction: n > 0.

    Int Ap [n+1] ;      Input argument, not modified.

        Ap is an integer array of size n+1.  If the form argument is "column",
        then on input, it holds the "pointers" for the column form of the
        sparse matrix A.  The row indices of column j of the matrix A are held
        in Ai [(Ap [j]) ... (Ap [j+1]-1)].  If form is "row", then Ap holds the
        row pointers.  The column indices of row j of the matrix are held
        in Ai [(Ap [j]) ... (Ap [j+1]-1)].

        The first entry, Ap [0], must be zero, and Ap [j] <= Ap [j+1] must hold
        for all j in the range 0 to n-1.  The value nz = Ap [n] is thus the
        total number of entries in the pattern of the matrix A.
        Restriction:  Ap [0] == 0 and Ap [n] > 0.

    Int Ai [nz] ;       Input argument, not modified, of size nz = Ap [n].

        If form is "column", then the nonzero pattern (row indices) for column
        j is stored in Ai [(Ap [j]) ... (Ap [j+1]-1)].  Row indices must be in
        the range 0 to n-1 (the matrix is 0-based).

        If form is "row", then the nonzero pattern (column indices) for row
        j is stored in Ai [(Ap [j]) ... (Ap [j+1]-1)]. Column indices must be in
        the range 0 to n-1 (the matrix is 0-based).

    double Ax [nz] ;    Input argument, not modified, of size nz = Ap [n].

        The numerical values of the sparse matrix A.

        If form is "row", then the nonzero pattern (row indices) for column j
        is stored in Ai [(Ap [j]) ... (Ap [j+1]-1)], and the corresponding
        numerical values are stored in Ax [(Ap [j]) ... (Ap [j+1]-1)].

        If form is "column", then the nonzero pattern (column indices) for row j
        is stored in Ai [(Ap [j]) ... (Ap [j+1]-1)], and the corresponding
        numerical values are stored in Ax [(Ap [j]) ... (Ap [j+1]-1)].

        No numerical values are printed if Ax is a (double *) NULL pointer.

    char *form ;        Input argument, not modified.

        The matrix is in row-oriented form if form is "row".  Otherwise,
        the matrix is in column-oriented form.  The form argument may be
        (char *) NULL, in which case the matrix is in column-oriented form.

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            2 or less: no output.  returns silently without checking anything.
            3: fully check input, and print a short summary of its status
            4: as 3, but print first few entries of the input
            5: as 3, but print all of the input
            Default: 1
\end{verbatim}
}

\newpage
\subsection{umfpack\_report\_numeric and umfpack\_l\_report\_numeric}

{\footnotesize
\begin{verbatim}
int umfpack_report_numeric
(
    const char name [ ],
    void *Numeric,
    const double Control [UMFPACK_CONTROL]
) ;

long umfpack_l_report_numeric
(
    const char name [ ],
    void *Numeric,
    const double Control [UMFPACK_CONTROL]
) ;

int Syntax:

    #include "umfpack.h"
    void *Numeric ;
    double Control [UMFPACK_CONTROL] ;
    int status ;
    status = umfpack_report_numeric ("Numeric", Numeric, Control) ;

long Syntax:

    #include "umfpack.h"
    void *Numeric ;
    double Control [UMFPACK_CONTROL] ;
    long status ;
    status = umfpack_l_report_numeric ("Numeric", Numeric, Control) ;

Purpose:

    Verifies and prints a Numeric object.  This routine checks the object more
    carefully than the computational routines.  Normally, this check is not
    required, since umfpack_numeric either returns (void *) NULL, or a valid
    Numeric object.  However, if you suspect that your own code has corrupted
    the Numeric object (by overruning memory bounds, for example), then this
    routine might be able to detect a corrupted Numeric object.  Since this
    is a complex object, not all such user-generated errors are guaranteed to
    be caught by this routine.

Returns:

    UMFPACK_OK if Control [UMFPACK_PRL] <= 2 (the input is not checked).

    Otherwise:

    UMFPACK_OK if the Numeric object is valid.
    UMFPACK_ERROR_invalid_Numeric_object if the Numeric object is invalid.
    UMFPACK_ERROR_out_of_memory if out of memory.

Arguments:

    char name [ ] ;     Input argument, not modified.

        The name of the Numeric object.  This is optional; no name is printed
        if a (char *) NULL pointer is passed.

    void *Numeric ;                     Input argument, not modified.

        The Numeric object, which holds the numeric factorization computed by
        umfpack_numeric.

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            2 or less: no output.  returns silently without checking anything.
            3: fully check input, and print a short summary of its status
            4: as 3, but print first few entries of the input
            5: as 3, but print all of the input
            Default: 1
\end{verbatim}
}

\newpage
\subsection{umfpack\_report\_perm and umfpack\_l\_report\_perm}

{\footnotesize
\begin{verbatim}
int umfpack_report_perm
(
    const char name [ ],
    int n,
    const int P [ ],
    const double Control [UMFPACK_CONTROL]
) ;

long umfpack_l_report_perm
(
    const char name [ ],
    long n,
    const long P [ ],
    const double Control [UMFPACK_CONTROL]
) ;

int Syntax:

    #include "umfpack.h"
    int n, *P, status ;
    double Control [UMFPACK_CONTROL] ;
    status = umfpack_report_perm ("P", n, P, Control) ;

long Syntax:

    #include "umfpack.h"
    long n, *P, status ;
    double Control [UMFPACK_CONTROL] ;
    status = umfpack_l_report_perm ("P", n, P, Control) ;

Purpose:

    Verifies and prints a permutation vector.

Returns:

    UMFPACK_OK if Control [UMFPACK_PRL] <= 2 (the input is not checked).

    Otherwise:

    UMFPACK_OK if the permutation vector is valid (this includes that case
        when P is (Int *) NULL, which is not an error condition.

    UMFPACK_ERROR_n_nonpositive if n <= 0.
    UMFPACK_ERROR_out_of_memory if out of memory.
    UMFPACK_ERROR_invalid_permutation if P is not a valid permutation vector.

Arguments:

    char name [ ] ;     Input argument, not modified.

        The name of the permutation vector.  This is optional; no name is
        printed if a (char *) NULL pointer is passed.

    Int n ;             Input argument, not modified.

        P is an integer vector of size n.  Restriction: n > 0.

    Int P [n] ;         Input argument, not modified.

        A permutation vector of size n.  If P is not present (a (Int *) NULL
        pointer, then P is assumed to be the identity permutation.  This is
        consistent with its use as an input argument to umfpack_qsymbolic.
        If P is present, the entries in P must range between 0 and n-1, and no
        duplicates may exists.

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            2 or less: no output.  returns silently without checking anything.
            3: fully check input, and print a short summary of its status
            4: as 3, but print first few entries of the input
            5: as 3, but print all of the input
            Default: 1
\end{verbatim}
}

\newpage
\subsection{umfpack\_report\_symbolic and umfpack\_l\_report\_symbolic}

{\footnotesize
\begin{verbatim}
int umfpack_report_symbolic
(
    const char name [ ],
    void *Symbolic,
    const double Control [UMFPACK_CONTROL]
) ;

long umfpack_l_report_symbolic
(
    const char name [ ],
    void *Symbolic,
    const double Control [UMFPACK_CONTROL]
) ;

int Syntax:

    #include "umfpack.h"
    void *Symbolic ;
    double Control [UMFPACK_CONTROL] ;
    int status ;
    status = umfpack_report_symbolic ("Symbolic", Symbolic, Control) ;

long Syntax:

    #include "umfpack.h"
    void *Symbolic ;
    double Control [UMFPACK_CONTROL] ;
    long status ;
    status = umfpack_l_report_symbolic ("Symbolic", Symbolic, Control) ;

Purpose:

    Verifies and prints a Symbolic object.  This routine checks the object more
    carefully than the computational routines.  Normally, this check is not
    required, since umfpack_*symbolic either returns (void *) NULL, or a valid
    Symbolic object.  However, if you suspect that your own code has corrupted
    the Symbolic object (by overruning memory bounds, for example), then this
    routine might be able to detect a corrupted Symbolic object.  Since this is
    a complex object, not all such user-generated errors are guaranteed to be
    caught by this routine.

Returns:

    UMFPACK_OK if Control [UMFPACK_PRL] is <= 2 (no inputs are checked).

    Otherwise:

    UMFPACK_OK if the Symbolic object is valid.
    UMFPACK_ERROR_invalid_Symbolic_object if the Symbolic object is invalid.
    UMFPACK_ERROR_out_of_memory if out of memory.

Arguments:

    char name [ ] ;     Input argument, not modified.

        The name of the Symbolic object.  This is optional; no name is printed
        if a (char *) NULL pointer is passed.

    void *Symbolic ;                    Input argument, not modified.

        The Symbolic object, which holds the symbolic factorization computed by
        umfpack_symbolic.

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            2 or less: no output.  returns silently without checking anything.
            3: fully check input, and print a short summary of its status
            4: as 3, but print first few entries of the input
            5: as 3, but print all of the input
            Default: 1
\end{verbatim}
}

\newpage
\subsection{umfpack\_report\_triplet and umfpack\_l\_report\_triplet}

{\footnotesize
\begin{verbatim}
int umfpack_report_triplet
(
    const char name [ ],
    int n,
    int nz,
    const int Ti [ ],
    const int Tj [ ],
    const double Tx [ ],
    const double Control [UMFPACK_CONTROL]
) ;

long umfpack_l_report_triplet
(
    const char name [ ],
    long n,
    long nz,
    const long Ti [ ],
    const long Tj [ ],
    const double Tx [ ],
    const double Control [UMFPACK_CONTROL]
) ;

int Syntax:

    #include "umfpack.h"
    int n, nz, *Ti, *Tj, status ;
    double *Tx, Control [UMFPACK_CONTROL] ;
    status = umfpack_report_triplet ("Triplet", n, nz, Ti, Tj, Tx, Control) ;

long Syntax:

    #include "umfpack.h"
    long n, nz, *Ti, *Tj, status ;
    double *Tx, Control [UMFPACK_CONTROL] ;
    status = umfpack_l_report_triplet ("Triplet", n, nz, Ti, Tj, Tx, Control) ;

Purpose:

    Verifies and prints a matrix in triplet form.

Returns:

    UMFPACK_OK if Control [UMFPACK_PRL] <= 2 (the input is not checked).

    Otherwise:

    UMFPACK_OK if the Triplet matrix is OK.
    UMFPACK_ERROR_argument_missing if Ti and/or Tj are missing.
    UMFPACK_ERROR_n_nonpositive if n <= 0.
    UMFPACK_ERROR_nz_negative if nz < 0.
    UMFPACK_ERROR_invalid_triplet if any row or column index in Ti and/or Tj
        is not in the range 0 to n-1.

Arguments:

    char name [ ] ;     Input argument, not modified.

        The name of the matrix.  This is optional; no name is printed if
        a (char *) NULL pointer is passed.

    Int n ;             Input argument, not modified.

        A is an n-by-n matrix.

    Int nz ;            Input argument, not modified.

        The number of entries in the triplet form of the matrix.

    Int Ti [nz] ;       Input argument, not modified.
    Int Tj [nz] ;       Input argument, not modified.
    double Tx [nz] ;    Input argument, not modified.

        Ti, Tj, and Tx hold the "triplet" form of a sparse matrix.  The kth
        nonzero entry is in row i = Ti [k], column j = Tj [k], and has a
        numerical value of a_ij = Tx [k].  The row and column indices i and j
        must be in the range 0 to n-1.  Duplicate entries
        may be present; they are summed in the output matrix.  This is not an
        error condition.  The "triplets" may be in any order.  Tx is optional;
        if Tx is not present (a (double *) NULL pointer), then the numerical
        values are not printed.

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            2 or less: no output.  returns silently without checking anything.
            3: fully check input, and print a short summary of its status
            4: as 3, but print first few entries of the input
            5: as 3, but print all of the input
            Default: 1
\end{verbatim}
}

\newpage
\subsection{umfpack\_report\_vector and umfpack\_l\_report\_vector}

{\footnotesize
\begin{verbatim}
int umfpack_report_vector
(
    const char name [ ],
    int n,
    const double X [ ],
    const double Control [UMFPACK_CONTROL]
) ;

long umfpack_l_report_vector
(
    const char name [ ],
    long n,
    const double X [ ],
    const double Control [UMFPACK_CONTROL]
) ;

int Syntax:

    #include "umfpack.h"
    int n, status ;
    double *X, Control [UMFPACK_CONTROL] ;
    status = umfpack_report_vector ("X", n, X, Control) ;

long Syntax:

    #include "umfpack.h"
    long n, status ;
    double *X, Control [UMFPACK_CONTROL] ;
    status = umfpack_l_report_vector ("X", n, X, Control) ;

Purpose:

    Verifies and prints a real vector.

Returns:

    UMFPACK_OK if Control [UMFPACK_PRL] <= 2 (the input is not checked).

    Otherwise:
    
    UMFPACK_OK if the vector is valid.
    UMFPACK_ERROR_argument_missing if X is missing.
    UMFPACK_ERROR_n_nonpositive if n <= 0.

Arguments:

    char name [ ] ;     Input argument, not modified.

        The name of the vector.  This is optional; no name is
        printed if a (char *) NULL pointer is passed.

    Int n ;             Input argument, not modified.

        X is a real vector of size n.  Restriction: n > 0.

    double X [n] ;              Input argument, not modified.

        A real vector of size n.  X must not be (double *) NULL.

    double Control [UMFPACK_CONTROL] ;  Input argument, not modified.

        If a (double *) NULL pointer is passed, then the default control
        settings are used.  Otherwise, the settings are determined from the
        Control array.  See umfpack_defaults on how to fill the Control
        array with the default settings.  The following Control parameters
        are used:

        Control [UMFPACK_PRL]:  printing level.

            2 or less: no output.  returns silently without checking anything.
            3: fully check input, and print a short summary of its status
            4: as 3, but print first few entries of the input
            5: as 3, but print all of the input
            Default: 1
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
\section{Utility routines}
\label{Utility}
%-------------------------------------------------------------------------------

\subsection{umfpack\_timer}

{\footnotesize
\begin{verbatim}
double umfpack_timer ( void ) ;

Syntax (for both int and long versions):

    #include "umfpack.h"
    double t ;
    t = umfpack_timer ( ) ;

Purpose:

    Returns the CPU time used by the process.  Includes both "user" and "system"
    time (the latter is time spent by the system on behalf of the process, and
    is thus charged to the process).

    This routine uses the Unix getrusage ( ) routine, if available.  It is not
    subject to overflow.  If getrusage ( ) is not available, the portable ANSI
    C clock ( ) routine is used instead.  Unfortunately, clock ( ) overflows
    if the CPU time exceeds 2147 seconds (about 36 minutes) when
    sizeof (clock_t) is 4 bytes.  If you have getrusage ( ), be sure to compile
    UMFPACK with the -DGETRUSAGE flag set; see umf_config.h and the User Guide
    for details.

Arguments:

    None.
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
\section{{\tt umfpack.h} include file}
\label{Include}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
/*
    This is the umfpack.h include file, and should be included in all user code
    that uses UMFPACK.  Do not include any of the umf_* header files in user
    code.  All routines in UMFPACK starting with "umfpack_" are user-callable
    (the 24 routines listed below).  All other routines are prefixed "umf_",
    and are not user-callable.
*/

#ifndef UMFPACK_H
#define UMFPACK_H

/* -------------------------------------------------------------------------- */
/* size of Info and Control arrays */
/* -------------------------------------------------------------------------- */

#define UMFPACK_INFO 90
#define UMFPACK_CONTROL 20

/* -------------------------------------------------------------------------- */
/* User-callable routines */
/* -------------------------------------------------------------------------- */

/* Primary routines: */
#include "umfpack_symbolic.h"
#include "umfpack_numeric.h"
#include "umfpack_solve.h"
#include "umfpack_free_symbolic.h"
#include "umfpack_free_numeric.h"

/* Alternative routines: */
#include "umfpack_defaults.h"
#include "umfpack_qsymbolic.h"
#include "umfpack_wsolve.h"

/* Matrix manipulation routines: */
#include "umfpack_triplet_to_col.h"
#include "umfpack_col_to_triplet.h"
#include "umfpack_transpose.h"

/* Getting the contents of the Symbolic and Numeric opaque objects: */
#include "umfpack_get_lunz.h"
#include "umfpack_get_numeric.h"
#include "umfpack_get_symbolic.h"

/* Reporting routines (the above 14 routines print nothing): */
#include "umfpack_report_status.h"
#include "umfpack_report_info.h"
#include "umfpack_report_control.h"
#include "umfpack_report_matrix.h"
#include "umfpack_report_triplet.h"
#include "umfpack_report_symbolic.h"
#include "umfpack_report_numeric.h"
#include "umfpack_report_perm.h"
#include "umfpack_report_vector.h"

/* Utility routines: */
#include "umfpack_timer.h"


/* -------------------------------------------------------------------------- */
/* Version, copyright, and license */
/* -------------------------------------------------------------------------- */

#define UMFPACK_VERSION "UMFPACK V3.2"

#define UMFPACK_COPYRIGHT \
"UMFPACK:  Copyright (c) 2002 by Timothy A. Davis, University of Florida,\n" \
"davis@cise.ufl.edu.  All Rights Reserved.\n"

#define UMFPACK_LICENSE \
"\nUMFPACK License:\n" \
"\n" \
"    Your use or distribution of UMFPACK or any derivative code implies that\n"\
"    you agree to this License.\n" \
"\n" \
"    THIS MATERIAL IS PROVIDED AS IS, WITH ABSOLUTELY NO WARRANTY\n" \
"    EXPRESSED OR IMPLIED.  ANY USE IS AT YOUR OWN RISK.\n" \
"\n" \
"    Permission is hereby granted to use or copy this program, provided\n" \
"    that the Copyright, this License, and the Availability of the original\n" \
"    version is retained on all copies.  User documentation of any code that\n"\
"    uses this code or any derivative code must cite the Copyright, this\n" \
"    License, the Availability note, and \"Used by permission.\"  If this\n" \
"    code or any derivative code is accessible from within MATLAB, then\n" \
"    typing \"help umfpack\" must cite the Copyright, and \"type umfpack\"\n" \
"    must also cite this License and the Availability note.  Permission to\n" \
"    modify the code and to distribute modified code is granted, provided\n" \
"    the Copyright, this License, and the Availability note are retained,\n" \
"    and a notice that the code was modified is included.  This software\n" \
"    was developed with support from the National Science Foundation, and\n" \
"    is provided to you free of charge.\n" \
"\n" \
"Availability:  http://www.cise.ufl.edu/research/sparse\n" \
"\n"


/* -------------------------------------------------------------------------- */
/* contents of Info */
/* -------------------------------------------------------------------------- */

/* Note that umfpack_report.m must coincide with these definitions. */

/* returned by all routines that use Info: */
#define UMFPACK_STATUS 0
#define UMFPACK_N 1
#define UMFPACK_NZ 2

/* computed in UMFPACK_*symbolic and UMFPACK_numeric: */
#define UMFPACK_SIZE_OF_UNIT 3

/* computed in UMFPACK_*symbolic: */
#define UMFPACK_SIZE_OF_INT 4
#define UMFPACK_SIZE_OF_LONG 5
#define UMFPACK_SIZE_OF_POINTER 6
#define UMFPACK_SIZE_OF_ENTRY 7
#define UMFPACK_NDENSE_ROW 8
#define UMFPACK_NEMPTY_ROW 9
#define UMFPACK_NDENSE_COL 10
#define UMFPACK_NEMPTY_COL 11
#define UMFPACK_SYMBOLIC_DEFRAG 12
#define UMFPACK_SYMBOLIC_PEAK_MEMORY 13
#define UMFPACK_SYMBOLIC_SIZE 14
#define UMFPACK_SYMBOLIC_TIME 15

/* Info [16..19] unused */

/* estimates computed in UMFPACK_*symbolic: */
#define UMFPACK_NUMERIC_SIZE_ESTIMATE 20
#define UMFPACK_PEAK_MEMORY_ESTIMATE 21
#define UMFPACK_FLOPS_ESTIMATE 22
#define UMFPACK_LNZ_ESTIMATE 23
#define UMFPACK_UNZ_ESTIMATE 24
#define UMFPACK_VARIABLE_INIT_ESTIMATE 25
#define UMFPACK_VARIABLE_PEAK_ESTIMATE 26
#define UMFPACK_VARIABLE_FINAL_ESTIMATE 27
#define UMFPACK_MAX_FRONT_SIZE_ESTIMATE 28

/* Info [29..39] unused */

/* exact values, (estimates shown above) computed in UMFPACK_numeric: */
#define UMFPACK_NUMERIC_SIZE 40
#define UMFPACK_PEAK_MEMORY 41
#define UMFPACK_FLOPS 42
#define UMFPACK_LNZ 43
#define UMFPACK_UNZ 44
#define UMFPACK_VARIABLE_INIT 45
#define UMFPACK_VARIABLE_PEAK 46
#define UMFPACK_VARIABLE_FINAL 47
#define UMFPACK_MAX_FRONT_SIZE 48

/* Info [49..59] unused */

/* computed in UMFPACK_numeric: */
#define UMFPACK_NUMERIC_DEFRAG 60
#define UMFPACK_NUMERIC_REALLOC 61
#define UMFPACK_NUMERIC_COSTLY_REALLOC 62
#define UMFPACK_COMPRESSED_PATTERN 63
#define UMFPACK_LU_ENTRIES 64
#define UMFPACK_NUMERIC_TIME 65

/* Info [66..79] unused */

/* computed in UMFPACK_solve: */
#define UMFPACK_IR_TAKEN 80
#define UMFPACK_IR_ATTEMPTED 81
#define UMFPACK_OMEGA1 82
#define UMFPACK_OMEGA2 83
#define UMFPACK_SOLVE_FLOPS 84
#define UMFPACK_SOLVE_TIME 85

/* Info [86..89] unused */

/* Unused parts of Info may be used in future versions of UMFPACK. */


/* -------------------------------------------------------------------------- */
/* contents of Control */
/* -------------------------------------------------------------------------- */

/* used in all UMFPACK_report_* routines: */
#define UMFPACK_PRL 0

/* used in UMFPACK_*symbolic only: */
#define UMFPACK_DENSE_ROW 1
#define UMFPACK_DENSE_COL 2

/* used in UMFPACK_numeric only: */
#define UMFPACK_PIVOT_TOLERANCE 3
#define UMFPACK_BLOCK_SIZE 4
#define UMFPACK_RELAXED_AMALGAMATION 5
#define UMFPACK_ALLOC_INIT 6
#define UMFPACK_PIVOT_OPTION 12
#define UMFPACK_RELAXED2_AMALGAMATION 13
#define UMFPACK_RELAXED3_AMALGAMATION 14

/* used in UMFPACK_*solve only: */
#define UMFPACK_IRSTEP 7

/* compile-time settings - Control [8..11] cannot be changed at run time: */
#define UMFPACK_COMPILED_WITH_BLAS 8
#define UMFPACK_COMPILED_FOR_MATLAB 9
#define UMFPACK_COMPILED_WITH_GETRUSAGE 10
#define UMFPACK_COMPILED_IN_DEBUG_MODE 11

/* Control [15...19] unused */

/* Unused parts of Control may be used in future versions of UMFPACK. */


/* -------------------------------------------------------------------------- */
/* default values of Control [0..7,12..13]: */
/* -------------------------------------------------------------------------- */

/* Note that the default block sized changed for Version 3.1 and following. */

#define UMFPACK_DEFAULT_PRL 1
#define UMFPACK_DEFAULT_DENSE_ROW 0.2
#define UMFPACK_DEFAULT_DENSE_COL 0.2
#define UMFPACK_DEFAULT_PIVOT_TOLERANCE 0.1
#define UMFPACK_DEFAULT_BLOCK_SIZE 24
#define UMFPACK_DEFAULT_RELAXED_AMALGAMATION 0.25
#define UMFPACK_DEFAULT_RELAXED2_AMALGAMATION 0.1
#define UMFPACK_DEFAULT_RELAXED3_AMALGAMATION 0.125
#define UMFPACK_DEFAULT_ALLOC_INIT 0.7
#define UMFPACK_DEFAULT_IRSTEP 2
#define UMFPACK_DEFAULT_PIVOT_OPTION 0

/* default values of Control [0..7,12..13] may change in future versions */
/* of UMFPACK. */

/* -------------------------------------------------------------------------- */
/* status codes */
/* -------------------------------------------------------------------------- */

#define UMFPACK_OK (0)
#define UMFPACK_ERROR_out_of_memory (-1)
#define UMFPACK_ERROR_singular_matrix (-2)
#define UMFPACK_ERROR_invalid_Numeric_object (-3)
#define UMFPACK_ERROR_invalid_Symbolic_object (-4)
#define UMFPACK_ERROR_argument_missing (-5)
#define UMFPACK_ERROR_n_nonpositive (-6)
#define UMFPACK_ERROR_nz_negative (-7)
#define UMFPACK_ERROR_jumbled_matrix (-8)
#define UMFPACK_ERROR_Ap0_nonzero (-9)
#define UMFPACK_ERROR_row_index_out_of_bounds (-10)
#define UMFPACK_ERROR_different_pattern (-11)
#define UMFPACK_ERROR_col_length_negative (-12)
#define UMFPACK_ERROR_invalid_system (-13)
#define UMFPACK_ERROR_invalid_triplet (-14)
#define UMFPACK_ERROR_invalid_permutation (-15)
#define UMFPACK_ERROR_problem_too_large (-16)
#define UMFPACK_ERROR_internal_error (-911)

#endif /* UMFPACK_H */
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
\section{Demo C main program, {\tt umfpack\_demo.c}}
\label{Demo}
%-------------------------------------------------------------------------------

The {\tt umfpack\_l\_demo.c} is identical except that all the routine names
are changed, and {\tt long}'s are used instead of {\tt int}'s.

{\footnotesize
\begin{verbatim}
/*
  A demo of UMFPACK Version 3.2:  See umfpack_demo.m for a (roughly)
  equivalent Matlab version.  The only difference is that while the Matlab
  umfpack mexFunction provides separate access to umfpack_symbolic, via
  umfpack (A, 'symbolic'), it does not use its output for a subsequent
  numerical factorization.  Thus, you will find that the output of this
  program and the Matlab diary are slightly different.  The Matlab output also
  uses 1-based matrix row and column indices, not 0-based (the internal
  represention is the same).

  First, factor and solve a 5-by-5 system, Ax=b, using default parameters,

      [ 2  3  0  0  0 ]      [  8 ]                  [ 1 ]
      [ 3  0  4  0  6 ]      [ 45 ]                  [ 2 ]
  A = [ 0 -1 -3  2  0 ], b = [ -3 ]. Solution is x = [ 3 ].
      [ 0  0  1  0  0 ]      [  3 ]                  [ 4 ]
      [ 0  4  2  0  1 ]      [ 19 ]                  [ 5 ]

  Then solve A'x=b, with solution:
        x = [  1.8158  1.4561 1.5000 -24.8509 10.2632 ]'
  using the factors of A.   Modify one entry (A (1,4) = 0, where the row and
  column indices range from 0 to 4, obtaining the system:

      [ 2  3  0  0  0 ]      [  8 ]                  [ 11.0    ]
      [ 3  0  4  0  0 ]      [ 45 ]                  [ -4.6667 ]
  A = [ 0 -1 -3  2  0 ], b = [ -3 ]. Solution is x = [  3.0    ].
      [ 0  0  1  0  0 ]      [  3 ]                  [  0.6667 ]
      [ 0  4  2  0  1 ]      [ 19 ]                  [ 31.6667 ]

  The pattern of A has not changed (it has explicitly zero entry), so a
  reanalysis with umfpack_symbolic does not need to be done (the Matlab
  umfpack_demo.m will need to redo it, because the Matlab caller is not
  provided with the Symbolic object).  Refactorize (with umfpack_numeric),
  and solve Ax=b.  Note that the pivot ordering has changed.  Next, change all
  of the entries in A, but not the pattern.  The system becomes

      [ 2 13  0  0  0 ]      [  8 ]                  [  8.5012 ]
      [ 2  0 23  0 39 ]      [ 45 ]                  [ -0.6925 ]
  A = [ 0  7 15 30  0 ], b = [ -3 ]. Solution is x = [  0.1667 ].
      [ 0  0 18  0  0 ]      [  3 ]                  [ -0.0218 ]
      [ 0 10 18  0 37 ]      [ 19 ]                  [  0.6196 ]

  Finally, compute B = A', and do the symbolic and numeric factorization of B.
  Factorizing A' can sometimes be better than factorizing A itself (less work
  and memory usage).  Solve B'x=b twice; the solution is the same as the
  solution to Ax=b for the above A.
*/

/* -------------------------------------------------------------------------- */
/* definitions */
/* -------------------------------------------------------------------------- */

#include <stdio.h>
#include <stdlib.h>
#include "umfpack.h"

#define ABS(x) ((x) >= 0 ? (x) : -(x))
#define MAX(a,b) (((a) > (b)) ? (a) : (b))
#ifndef TRUE
#define TRUE (1)
#endif
#ifndef FALSE
#define FALSE (0)
#endif

/* -------------------------------------------------------------------------- */
/* triplet form of the matrix.  The triplets can be in any order. */
/* -------------------------------------------------------------------------- */

static int    n = 5 ;
static int    nz = 12 ;
static int    Arow [ ] = { 0,  4,  1,  1,   2,   2,  0,  1,  2,  3,  4,  4} ;
static int    Acol [ ] = { 0,  4,  0,  2,   1,   2,  1,  4,  3,  2,  1,  2} ;
static double Aval [ ] = {2., 1., 3., 4., -1., -3., 3., 6., 2., 1., 4., 2.} ;
static double b [ ] = {8., 45., -3., 3., 19.} ;
static double x [5] ;
static double r [5] ;

/* -------------------------------------------------------------------------- */
/* error: print a message and exit */
/* -------------------------------------------------------------------------- */

static void error
(
    char *message
)
{
    printf ("\n\n====== error: %s =====\n\n", message) ;
    exit (1) ;
}


/* -------------------------------------------------------------------------- */
/* resid: compute the residual, r = Ax-b or r = A'x=b and return maxnorm (r) */
/* -------------------------------------------------------------------------- */

static double resid
(
    int n,
    int Ap [ ],
    int Ai [ ],
    double Ax [ ],
    double x [ ],
    double r [ ],
    int transpose
)
{
    int i, j, p ;
    double norm ;

    for (i = 0 ; i < n ; i++)
    {
        r [i] = -b [i] ;
    }
    if (transpose)
    {
        for (j = 0 ; j < n ; j++)
        {
            for (p = Ap [j] ; p < Ap [j+1] ; p++)
            {
                i = Ai [p] ;
                r [j] += Ax [p] * x [i] ;
            }
        }
    }
    else
    {
        for (j = 0 ; j < n ; j++)
        {
            for (p = Ap [j] ; p < Ap [j+1] ; p++)
            {
                i = Ai [p] ;
                r [i] += Ax [p] * x [j] ;
            }
        }
    }
    norm = 0. ;
    for (i = 0 ; i < n ; i++)
    {
        norm = MAX (norm, ABS (r [i])) ;
    }
    return (norm) ;
}


/* -------------------------------------------------------------------------- */
/* main program */
/* -------------------------------------------------------------------------- */

/*ARGSUSED0*/   /* argc and argv are unused */
int main
(
    int argc,
    char **argv
)
{
    double Info [UMFPACK_INFO], Control [UMFPACK_CONTROL], *Ax, *Bx, *Lx, *Ux,
        *W, *Y, *Z, *S, t ;
    int *Ap, *Ai, *Bp, *Bi, row, col, p, lnz, unz, nn, *Lp, *Li, *Ui, *Up,
        *P, *Q, *Lj, i, j, k, anz, nfr, nchains, nsparse_col, *Qtree, fnpiv,
        status, *Front_npivots, *Front_parent, *Chain_start, *Wi,
        *Chain_maxrows, *Chain_maxcols ;
    void *Symbolic, *Numeric ;

    /* ---------------------------------------------------------------------- */
    /* initializations */
    /* ---------------------------------------------------------------------- */

    t = umfpack_timer ( ) ;

    printf ("\n%s demo:\n", UMFPACK_VERSION) ;

    /* get the default control parameters */
    umfpack_defaults (Control) ;

    /* change the default print level for this demo */
    /* (otherwise, nothing will print) */
    Control [UMFPACK_PRL] = 6 ;

    /* print the license agreement */
    umfpack_report_status (Control, UMFPACK_OK) ;
    Control [UMFPACK_PRL] = 5 ;

    /* print the control parameters */
    umfpack_report_control (Control) ;

    /* ---------------------------------------------------------------------- */
    /* print A and b, and convert A to column-form */
    /* ---------------------------------------------------------------------- */

    /* print the right-hand-side */
    (void) umfpack_report_vector ("b", n, b, Control) ;

    /* print the triplet form of the matrix */
    (void) umfpack_report_triplet ("A", n, nz, Arow, Acol, Aval, Control) ;

    /* convert to column form */
    Ap = (int *) malloc ((n+1) * sizeof (int)) ;
    Ai = (int *) malloc (nz * sizeof (int)) ;
    Ax = (double *) malloc (nz * sizeof (double)) ;
    if (!Ap || !Ai || !Ax)
    {
        error ("out of memory") ;
    }
    status = umfpack_triplet_to_col (n, nz, Arow, Acol, Aval, Ap, Ai, Ax) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_status (Control, status) ;
        error ("umfpack_triplet_to_col failed") ;
    }

    /* print the column-form of A */
    (void) umfpack_report_matrix ("A", n, Ap, Ai, Ax, "column", Control) ;

    /* ---------------------------------------------------------------------- */
    /* symbolic factorization */
    /* ---------------------------------------------------------------------- */

    status = umfpack_symbolic (n, Ap, Ai, &Symbolic, Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_info (Control, Info) ;
        umfpack_report_status (Control, status) ;
        error ("umfpack_symbolic failed") ;
    }

    /* print the symbolic factorization */
    (void) umfpack_report_symbolic ("Symbolic factorization of A",
        Symbolic, Control) ;

    /* ---------------------------------------------------------------------- */
    /* numeric factorization */
    /* ---------------------------------------------------------------------- */

    status = umfpack_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_info (Control, Info) ;
        umfpack_report_status (Control, status) ;
        error ("umfpack_numeric failed") ;
    }

    /* print the numeric factorization */
    (void) umfpack_report_numeric ("Numeric factorization of A",
        Numeric, Control) ;

    /* ---------------------------------------------------------------------- */
    /* solve Ax=b */
    /* ---------------------------------------------------------------------- */

    status = umfpack_solve ("Ax=b", Ap, Ai, Ax, x, b, Numeric, Control, Info) ;
    umfpack_report_info (Control, Info) ;
    umfpack_report_status (Control, status) ;
    if (status != UMFPACK_OK)
    {
        error ("umfpack_solve failed") ;
    }
    (void) umfpack_report_vector ("x (solution of Ax=b)", n, x, Control) ;
    printf ("maxnorm of residual: %g\n\n", resid (n, Ap, Ai, Ax, x, r, FALSE)) ;

    /* ---------------------------------------------------------------------- */
    /* solve A'x=b */
    /* ---------------------------------------------------------------------- */

    status = umfpack_solve ("A'x=b", Ap, Ai, Ax, x, b, Numeric, Control, Info) ;
    umfpack_report_info (Control, Info) ;
    if (status != UMFPACK_OK)
    {
        error ("umfpack_solve failed") ;
    }
    (void) umfpack_report_vector ("x (solution of A'x=b)", n, x, Control) ;
    printf ("maxnorm of residual: %g\n\n", resid (n, Ap, Ai, Ax, x, r, TRUE)) ;

    /* ---------------------------------------------------------------------- */
    /* modify one numerical value in the column-form of A */
    /* ---------------------------------------------------------------------- */

    /* change A (1,4), look for row index 1 in column 4. */
    row = 1 ;
    col = 4 ;
    for (p = Ap [col] ; p < Ap [col+1] ; p++)
    {
        if (row == Ai [p])
        {
            printf ("\nchanging A (%d,%d) from %g", row, col, Ax [p]) ;
            Ax [p] = 0.0 ;
            printf (" to %g\n", Ax [p]) ;
            break ;
        }
    }
    (void) umfpack_report_matrix ("modified A", n, Ap, Ai, Ax, "column",
        Control) ;

    /* ---------------------------------------------------------------------- */
    /* redo the numeric factorization */
    /* ---------------------------------------------------------------------- */

    /* The pattern (Ap and Ai) hasn't changed, so the symbolic factorization */
    /* doesn't have to be redone, no matter how much we change Ax. */

    /* We don't need the Numeric object any more, so free it. */
    umfpack_free_numeric (&Numeric) ;

    /* Note that a memory leak would have occured if the old Numeric */
    /* had not been free'd with umfpack_free_numeric above. */
    status = umfpack_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_info (Control, Info) ;
        umfpack_report_status (Control, status) ;
        error ("umfpack_numeric failed") ;
    }
    (void) umfpack_report_numeric ("Numeric factorization of modified A",
        Numeric, Control) ;

    /* ---------------------------------------------------------------------- */
    /* solve Ax=b, with the modified A */
    /* ---------------------------------------------------------------------- */

    status = umfpack_solve ("Ax=b", Ap, Ai, Ax, x, b, Numeric, Control, Info) ;
    umfpack_report_info (Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_status (Control, status) ;
        error ("umfpack_solve failed") ;
    }
    (void) umfpack_report_vector ("x (with modified A)", n, x, Control) ;
    printf ("maxnorm of residual: %g\n\n", resid (n, Ap, Ai, Ax, x, r, FALSE)) ;

    /* ---------------------------------------------------------------------- */
    /* modify all of the numerical values of A, but not the pattern */
    /* ---------------------------------------------------------------------- */

    for (col = 0 ; col < n ; col++)
    {
        for (p = Ap [col] ; p < Ap [col+1] ; p++)
        {
            row = Ai [p] ;
            printf ("changing A (%d,%d) from %g", row, col, Ax [p]) ;
            Ax [p] = Ax [p] + col*10 - row ;
            printf (" to %g\n", Ax [p]) ;
        }
    }
    (void) umfpack_report_matrix ("completely modified A (same pattern)",
        n, Ap, Ai, Ax, "column", Control) ;

    /* ---------------------------------------------------------------------- */
    /* redo the numeric factorization */
    /* ---------------------------------------------------------------------- */

    umfpack_free_numeric (&Numeric) ;
    status = umfpack_numeric (Ap, Ai, Ax, Symbolic, &Numeric, Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_info (Control, Info) ;
        umfpack_report_status (Control, status) ;
        error ("umfpack_numeric failed") ;
    }
    (void) umfpack_report_numeric (
    "Numeric factorization of completely modified A", Numeric, Control) ;

    /* ---------------------------------------------------------------------- */
    /* solve Ax=b, with the modified A */
    /* ---------------------------------------------------------------------- */

    status = umfpack_solve ("Ax=b", Ap, Ai, Ax, x, b, Numeric, Control, Info) ;
    umfpack_report_info (Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_status (Control, status) ;
        error ("umfpack_solve failed") ;
    }
    (void) umfpack_report_vector ("x (with completely modified A)",
        n, x, Control) ;
    printf ("maxnorm of residual: %g\n\n", resid (n, Ap, Ai, Ax, x, r, FALSE)) ;

    /* ---------------------------------------------------------------------- */
    /* free the symbolic and numeric factorization */
    /* ---------------------------------------------------------------------- */

    umfpack_free_symbolic (&Symbolic) ;
    umfpack_free_numeric (&Numeric) ;

    /* ---------------------------------------------------------------------- */
    /* B = transpose of A */
    /* ---------------------------------------------------------------------- */

    Bp = (int *) malloc ((n+1) * sizeof (int)) ;
    Bi = (int *) malloc (nz * sizeof (int)) ;
    Bx = (double *) malloc (nz * sizeof (double)) ;
    if (!Bp || !Bi || !Bx)
    {
        error ("out of memory") ;
    }
    status = umfpack_transpose (n, Ap, Ai, Ax, (int *) NULL, (int *) NULL,
        Bp, Bi, Bx) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_status (Control, status) ;
        error ("umfpack_transpose failed") ;
    }
    (void) umfpack_report_matrix ("B (transpose of A)",
        n, Bp, Bi, Bx, "column", Control) ;

    /* ---------------------------------------------------------------------- */
    /* symbolic factorization of B */
    /* ---------------------------------------------------------------------- */

    status = umfpack_symbolic (n, Bp, Bi, &Symbolic, Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_info (Control, Info) ;
        umfpack_report_status (Control, status) ;
        error ("umfpack_symbolic failed") ;
    }
    (void) umfpack_report_symbolic ("Symbolic factorization of B",
        Symbolic, Control) ;

    /* ---------------------------------------------------------------------- */
    /* copy the contents of Symbolic into user arrays print them */
    /* ---------------------------------------------------------------------- */

    printf ("\nGet the contents of the Symbolic object for B:\n") ;
    printf ("(compare with umfpack_report_symbolic output, above)\n") ;
    Qtree = (int *) malloc (n * sizeof (int)) ;
    Front_npivots = (int *) malloc (n * sizeof (int)) ;
    Front_parent = (int *) malloc (n * sizeof (int)) ;
    Chain_start = (int *) malloc ((n+1) * sizeof (int)) ;
    Chain_maxrows = (int *) malloc (n * sizeof (int)) ;
    Chain_maxcols = (int *) malloc (n * sizeof (int)) ;
    if (!Qtree || !Front_npivots || !Front_parent || !Chain_start ||
        !Chain_maxrows || !Chain_maxcols)
    {
        error ("out of memory") ;
    }

    status = umfpack_get_symbolic (&nn, &anz, &nfr, &nchains, &nsparse_col,
        Qtree, Front_npivots, Front_parent, Chain_start,
        Chain_maxrows, Chain_maxcols, Symbolic) ;

    printf ("From the Symbolic object, B is of dimension n = %d\n", nn) ;
    printf ("   with nz = %d, number of fronts = %d,\n", nz, nfr) ;
    printf ("   number of frontal matrix chains = %d\n", nchains) ;

    printf ("\nPivot columns in each front, and parent of each front:\n") ;
    k = 0 ;
    for (i = 0 ; i < nfr ; i++)
    {
        fnpiv = Front_npivots [i] ;
        printf ("    Front %d: parent front: %d number of pivots: %d\n",
                i, Front_parent [i], fnpiv) ;
        for (j = 0 ; j < fnpiv ; j++)
        {
            col = Qtree [k] ;
            printf (
                "        %d-th pivot column is column %d in original matrix\n",
                k, col) ;
            k++ ;
        }
    }

    printf ("\nNote that the column ordering, above, will be refined\n") ;
    printf ("in the numeric factorization below.  The assignment of pivot\n") ;
    printf ("columns to frontal matrices will always remain unchanged.\n") ;

    printf ("\nTotal number of pivot columns in frontal matrices: %d\n", k) ;

    printf ("\nFrontal matrix chains:\n") ;
    for (j = 0 ; j < nchains ; j++)
    {
        printf ("   Frontal matrices %d to %d are factorized in a single\n",
            Chain_start [j], Chain_start [j+1] - 1) ;
        printf ("        working array of size %d-by-%d\n",
            Chain_maxrows [j], Chain_maxcols [j]) ;
    }

    /* ---------------------------------------------------------------------- */
    /* numeric factorization of B */
    /* ---------------------------------------------------------------------- */

    status = umfpack_numeric (Bp, Bi, Bx, Symbolic, &Numeric, Control, Info) ;
    if (status != UMFPACK_OK)
    {
        error ("umfpack_numeric failed") ;
    }
    (void) umfpack_report_numeric ("Numeric factorization of B",
        Numeric, Control) ;

    /* ---------------------------------------------------------------------- */
    /* extract the LU factors of B and print them */
    /* ---------------------------------------------------------------------- */

    if (umfpack_get_lunz (&lnz, &unz, &nn, Numeric) != UMFPACK_OK)
    {
        error ("umfpack_get_lunz failed") ;
    }
    Lp = (int *) malloc ((n+1) * sizeof (int)) ;
    Li = (int *) malloc (lnz * sizeof (int)) ;
    Lx = (double *) malloc (lnz * sizeof (double)) ;
    Up = (int *) malloc ((n+1) * sizeof (int)) ;
    Ui = (int *) malloc (unz * sizeof (int)) ;
    Ux = (double *) malloc (unz * sizeof (double)) ;
    P = (int *) malloc (n * sizeof (int)) ;
    Q = (int *) malloc (n * sizeof (int)) ;
    if (!Lp || !Li || !Lx || !Up || !Ui || !Ux || !P || !Q)
    {
        error ("out of memory") ;
    }
    status = umfpack_get_numeric (Lp, Li, Lx, Up, Ui, Ux, P, Q, Numeric) ;
    if (status != UMFPACK_OK)
    {
        error ("umfpack_get_numeric failed") ;
    }
    (void) umfpack_report_matrix ("L (lower triangular factor of B)",
        n, Lp, Li, Lx, "row", Control) ;
    (void) umfpack_report_matrix ("U (upper triangular factor of B)",
        n, Up, Ui, Ux, "column", Control) ;
    (void) umfpack_report_perm ("P", n, P, Control) ;
    (void) umfpack_report_perm ("Q", n, Q, Control) ;

    /* ---------------------------------------------------------------------- */
    /* convert L to triplet form and print it */
    /* ---------------------------------------------------------------------- */

    printf ("\nConverting L to triplet form, and printing it:\n") ;
    Lj = (int *) malloc (lnz * sizeof (int)) ;
    if (!Lj)
    {
        error ("out of memory") ;
    }
    if (umfpack_col_to_triplet (n, Lp, Lj) != UMFPACK_OK)
    {
        error ("umfpack_col_to_triplet failed") ;
    }
    (void) umfpack_report_triplet ("L, in triplet form", n, lnz, Li, Lj, Lx,
        Control) ;

    /* ---------------------------------------------------------------------- */
    /* solve B'x=b */
    /* ---------------------------------------------------------------------- */

    status = umfpack_solve ("A'x=b", Bp, Bi, Bx, x, b, Numeric, Control, Info) ;
    umfpack_report_info (Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_status (Control, status) ;
        error ("umfpack_solve failed") ;
    }
    (void) umfpack_report_vector ("x (solution of B'x=b)", n, x, Control) ;
    printf ("maxnorm of residual: %g\n\n", resid (n, Bp, Bi, Bx, x, r, TRUE)) ;

    /* ---------------------------------------------------------------------- */
    /* solve B'x=b again, using umfpack_wsolve instead */
    /* ---------------------------------------------------------------------- */

    printf ("\nSolving B'x=b again, using umfpack_wsolve instead:\n") ;
    Wi = (int *) malloc (n * sizeof (int)) ;
    W = (double *) malloc (n * sizeof (double)) ;
    Y = (double *) malloc (n * sizeof (double)) ;
    Z = (double *) malloc (n * sizeof (double)) ;
    S = (double *) malloc (n * sizeof (double)) ;
    if (!Wi || !W || !Y || !Z || !S)
    {
        error ("out of memory") ;
    }

    status = umfpack_wsolve ("A'x=b", Bp, Bi, Bx, x, b, Numeric, Control, Info,
        Wi, W, Y, Z, S) ;
    umfpack_report_info (Control, Info) ;
    if (status != UMFPACK_OK)
    {
        umfpack_report_status (Control, status) ;
        error ("umfpack_wsolve failed") ;
    }
    (void) umfpack_report_vector ("x (solution of B'x=b)", n, x, Control) ;
    printf ("maxnorm of residual: %g\n\n", resid (n, Bp, Bi, Bx, x, r, TRUE)) ;

    /* ---------------------------------------------------------------------- */
    /* free everything */
    /* ---------------------------------------------------------------------- */

    /* This is not strictly required since the process is exiting and the */
    /* system will reclaim the memory anyway.  It's useful, though, just as */
    /* a list of what is currently malloc'ed by this program.  Plus, it's */
    /* always a good habit to explicitly free whatever you malloc. */

    free (Ap) ;
    free (Ai) ;
    free (Ax) ;

    free (Bp) ;
    free (Bi) ;
    free (Bx) ;

    free (Qtree) ;
    free (Front_npivots) ;
    free (Front_parent) ;
    free (Chain_start) ;
    free (Chain_maxrows) ;
    free (Chain_maxcols) ;

    free (Lp) ;
    free (Li) ;
    free (Lx) ;

    free (Up) ;
    free (Ui) ;
    free (Ux) ;

    free (P) ;
    free (Q) ;

    free (Lj) ;

    free (Wi) ;
    free (W) ;
    free (Y) ;
    free (Z) ;
    free (S) ;

    umfpack_free_symbolic (&Symbolic) ;
    umfpack_free_numeric (&Numeric) ;

    /* ---------------------------------------------------------------------- */
    /* print the total time spent in this demo */
    /* ---------------------------------------------------------------------- */

    t = umfpack_timer ( ) - t ;
    printf ("\numfpack demo complete.  Total time: %5.2f (seconds)\n", t) ;
    return (0) ;
}
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
\section{Configuration file {\tt umf\_config.h}}
\label{Config}
%-------------------------------------------------------------------------------

{\footnotesize
\begin{verbatim}
/*
    This file controls the compile-time configuration of UMFPACK.  Modify the
    Makefile, the architecture-dependent Make.* file, and this file if
    necessary, to control these options.  The following compile-time flags are
    available:

        -DNBLAS

            BLAS mode.  If -DNBLAS is set, then no BLAS will be used.  Vanilla
            C code will be used instead.  This is portable, and easier to
            install, but you won't get the best performance.

            If -DNBLAS is not set, then externally-available BLAS routines
            (dgemm, dger, and dgemv or the equivalent C-BLAS routines) will be
            used.  This will give you the best performance, but perhaps at the
            expense of portability.

            The default is to use the BLAS, for both the C-callable umfpack.a
            library and the Matlab mexFunction.  If you have trouble installing
            UMFPACK, set -DNBLAS.

        -DNCBLAS

            If -DNCBLAS is set, then the C-BLAS will not be called.  This is the
            default when compiling the Matlab mexFunction, or when compiling
            umfpack.a on Sun Solaris or SGI IRIX.

            If -DNCBLAS is not set, then the C-BLAS interface to the BLAS is
            used.  If your vendor-supplied BLAS library does not have a C-BLAS
            interface, you can obtain the ATLAS BLAS, available at
            http://www.netlib.org/atlas.

            Using the C-BLAS is the default when compiling umfpack.a on all
            architectures except Sun Solaris (the Sun Performance Library is
            somewhat faster).  The ANSI C interface to the BLAS is fully
            portable.

            This flag is ignored if -DNBLAS is set.

        -DLONGBLAS

            If not defined, then the BLAS are not called in the long integer
            version of UMFPACK (the umfpack_l_* routines).  The most common
            definitions of the BLAS, unfortunately, use int arguments, and
            are thus not suitable for use in the LP64 model.  Only the Sun
            Performance Library, as far as I can tell, has a version of the
            BLAS that allows long integer (64-bit) input arguments.  This
            flag is set automatically in Sun Solaris if you are using the
            Sun Performance BLAS.  You can set it yourself, too, if your BLAS
            routines can take long integer input arguments.

        -DNSUNPERF

            Applies only to Sun Solaris.  If -DNSUNPERF is set, then the Sun
            Performance Library BLAS will not be used.

            The Sun Performance Library BLAS is used by default when compiling
            the C-callable umfpack.a library on Sun Solaris.

            This flag is ignored if -DNBLAS is set.

        -DNSCSL

            Applies only to SGI IRIX.  If -DSCSL is set, then the SGI SCSL
            Scientific Library BLAS will not be used.

            The SGI SCSL Scientific Library BLAS is used by default when
            compiling the C-callable umfpack.a library on SGI IRIX.

            This flag is ignored if -DNBLAS is set.

        -DGETRUSAGE

            If -DGETRUSAGE is set, then your system's getrusage routine will be
            used for getting the process CPU time.  Otherwise the ANSI C clock
            routine will be used.  The default is to use getrusage on Sun
            Solaris, SGI Irix, Linux, and AIX (IBM RS 6000) and to use clock on
            all other architectures.

    C-to-Fortran interface, for the Fortran BLAS (these are set automatically
        for the C-BLAS or Sun Performance BLAS):

        -DBLAS_BY_VALUE         if scalars are passed by value, not reference
        -DBLAS_NO_UNDERSCORE    if no underscore should be appended
        -DBLAS_CHAR_ARG         if BLAS options are single char's, not strings

    You should normally not set these flags yourself:

        -DMATLAB_MEX_FILE

            This flag is turned on when compiling the umfpack mexFunction for
            use in Matlab.  When compiling the Matlab mexFunction, the Matlab
            BLAS are used by default (this is set in the Makefile).

        -DMATHWORKS

            This flag is turned on when compiling umfpack as a built-in routine
            in Matlab.  It can also be used when compiling umfpack as a
            mexFunction.  Internal routines utMalloc, utFree, utRealloc, and
            utPrintf are used, and the "util.h" file is included.  This avoids
            the problem discussed in the User Guide regarding memory allocation
            in Matlab.  utMalloc returns NULL on failure, instead of terminating
            the mexFunction (which is what mxMalloc does).  However, the ut*
            routines are not documented by The MathWorks, Inc., so I cannot
            guarantee that you will always be able to use them.

        -DNDEBUG

            Debugging mode (if NDEBUG is not defined).  The default, of course,
            is no debugging.  Turning on debugging takes some work (see below).

*/



/* ========================================================================== */
/* === NDEBUG =============================================================== */
/* ========================================================================== */

/*
    UMFPACK will be exceedingly slow when running in debug mode.  The next three
    lines ensure that debugging is turned off.  If you want to compile UMFPACK
    in debugging mode, you must comment out the three lines below:
*/
#ifndef NDEBUG
#define NDEBUG
#endif

/*
    Next, you must either remove the -DNDEBUG option in the Makefile, or simply
    add the following line:
#undef NDEBUG
*/


/* ========================================================================== */
/* === Memory allocator ===================================================== */
/* ========================================================================== */

/* The Matlab mexFunction uses Matlab's memory manager, while the C-callable */
/* umfpack.a library uses the ANSI C malloc, free, and realloc routines. */

#ifdef MATLAB_MEX_FILE
#include "matrix.h"
#define ALLOCATE mxMalloc
#define FREE mxFree
#define REALLOCATE mxRealloc
#else
#ifdef MATHWORKS
#include "util.h"
/* Compiling UMFPACK as a built-in routine. */
/* Since UMFPACK carefully checks for out-of-memory after every allocation, */
/* we can use ut* routines here. */
#define ALLOCATE utMalloc
#define FREE utFree
#define REALLOCATE utRealloc
#else
#define ALLOCATE malloc
#define FREE free
#define REALLOCATE realloc
#endif
#endif


/* ========================================================================== */
/* === PRINTF macro ========================================================= */
/* ========================================================================== */

/* All output goes through the PRINTF macro.  Printing occurs only from the */
/* UMFPACK_report_* routines. */

#ifdef MATLAB_MEX_FILE
#include "mex.h"
#define PRINTF(params) { (void) mexPrintf params ; }
#else
#ifdef MATHWORKS
/* Already #included "util.h" above in Memory allocator section */
#define PRINTF(params) { (void) utPrintf params ; }
#else
#define PRINTF(params) { (void) printf params ; }
#endif
#endif


/* ========================================================================== */
/* === 0-based or 1-based printing ========================================== */
/* ========================================================================== */

#if defined (MATLAB_MEX_FILE) || defined (MATHWORKS)
/* In Matlab, matrices are 1-based to the user, but 0-based internally. */
/* One is added to all row and column indices when printing matrices */
/* in UMFPACK_report_*.  */
#define INDEX(i) ((i)+1)
#else
/* In ANSI C, matrices are 0-based and indices are reported as such. */
#define INDEX(i) (i)
#endif


/* ========================================================================== */
/* === Architecture ========================================================= */
/* ========================================================================== */

#if defined (__sun)
#define UMFPACK_ARCHITECTURE "Sun Solaris"
#endif

#if defined (__sgi)
#define UMFPACK_ARCHITECTURE "SGI Irix"
#endif

#if defined (__linux)
#define UMFPACK_ARCHITECTURE "Linux"
#endif

#if defined (_AIX)
#define UMFPACK_ARCHITECTURE "IBM AIX"
#endif

#if defined (__alpha)
#define UMFPACK_ARCHITECTURE "Compaq Alpha"
#endif

/* ========================================================================== */
/* === Timer ================================================================ */
/* ========================================================================== */

/*
    If you have the getrusage routine (all Unix systems I've test do), then use
    that.  Otherwise, use the ANSI C clock function.   Note that on many
    systems, the ANSI clock function wraps around after only 2147 seconds, or
    about 36 minutes.  BE CAREFUL:  if you compare the run time of UMFPACK with
    other sparse matrix packages, be sure to use the same timer.  See
    umfpack_timer.c for the timer used by UMFPACK.
*/

/* Sun Solaris, SGI Irix, Linux, Compaq Alpha, and IBM RS 6000 all have */
/* getrusage.  It's in BSD unix, so perhaps all unix systems have it. */
#if defined (__sun) || defined (__sgi) || defined (__linux) \
|| defined (__alpha) || defined (_AIX)
#define GETRUSAGE
#endif


/* ========================================================================== */
/* === BLAS ================================================================= */
/* ========================================================================== */

/*
    Determine if the BLAS exists for the long integer version.  It exists if
    LONGBLAS is defined in the Makefile, or if using the BLAS from the
    Sun Performance Library, or SGI's SCSL Scientific Library.
*/

#if !defined (MATLAB_MEX_FILE) && defined (__sun) && !defined (NSUNPERF)
#define BLAS_SUNPERF
#ifndef LONGBLAS
#define LONGBLAS
#endif
#endif

#if !defined (MATLAB_MEX_FILE) && defined (__sgi) && !defined (NSCSL)
#define BLAS_SCSL
#ifndef LONGBLAS
#define LONGBLAS
#endif
#endif

#if defined (DLONG) && !defined (LONGBLAS) && !defined (NBLAS)
#define NBLAS
#endif

/* -------------------------------------------------------------------------- */

#ifndef NBLAS

/*
    If the compile-time flag -DNBLAS is defined, then the BLAS are not used,
    portable vanilla C code is used instead, and the remainder of this file
    is ignored.

    Using the BLAS is much faster, but how C calls the Fortran BLAS is
    machine-dependent and thus can cause portability problems.  Thus, use
    -DNBLAS to ensure portability (at the expense of speed).

    Preferences:

        *** The best interface to use, regardless of the option you select
            below, is the standard C-BLAS interface.  Not all vendor-supplied
            BLAS libraries use this interface (as of April 2001).  The only
            problem with this interface is that it does not extend to the LP64
            model.

        1) most preferred:  use the optimized vendor-supplied library (such as
            the Sun Performance Library, or IBM's ESSL).  This is often the
            fastest, but might not be portable and might not always be
            available.  When compiling a Matlab mexFunction it might be
            difficult get the mex compiler script to recognize the vendor-
            supplied BLAS (I was not able get my mexFunction to use the
            Sun Performance Library BLAS, for example, because of linking
            difficulties).

        2) When compiling the UMFPACK mexFunction to use UMFPACK in Matlab, use
            the BLAS provided by The Mathworks, Inc.  This assumes you are using
            Matlab V6 or higher, since the BLAS are not incorporated in V5 or
            earlier versions.  On my Sun workstation, the Matlab BLAS gave
            slightly worse performance than the Sun Perf. BLAS.  The advantage
            of using the Matlab BLAS is that it's available on any computer that
            has Matlab V6 or higher.  I have not tried using Matlab BLAS outside
            of a mexFunction in a stand-alone C code, but Matlab (V6) allows for
            this.  This is well worth trying if you have Matlab and don't want
            to bother installing the ATLAS BLAS (option 3a, below).  The only
            glitch to this is that Matlab does not provide a portable interface
            to the BLAS (an underscore is required for some but not all
            architectures).  These variations are taken into account in the
            mexopts.sh file provided with UMFPACK.

        3) Use a portable high-performance BLAS library:

            (a) The ATLAS BLAS, available at http://www.netlib.org/atlas,
                by R. Clint Whaley, Antoine Petitet, and Jack Dongarra.
                This has a standard C interface, and thus the interface to it is
                fully portable.  Its performance rivals, and sometimes exceeds,
                the vendor-supplied BLAS on many computers.

            (b) The Fortran RISC BLAS by Michel Dayde', Iain Duff, Antoine
                Petitet, and Abderrahim Qrichi Aniba, available via anonymous
                ftp to ftp.enseeiht.fr in the pub/numerique/BLAS/RISC directory,
                See M. J. Dayde' and I. S. Duff, "The RISC BLAS:  A blocked
                implementation of level 3 BLAS for RISC processors, ACM Trans.
                Math. Software, vol. 25, no. 3., Sept. 1999.  This will give
                you good performance, but with the same C-to-Fortran portability
                problems as option (1).

        4) Use UMFPACK's built-in vanilla C code by setting -DNBLAS at compile
            time.  The key advantage is portability, which is guaranteed if you
            have an ANSI C compliant compiler.  You also don't need to download
            any other package - UMFPACK is stand-alone.  No Fortran is used
            anywhere in UMFPACK.  UMFPACK will be much slower than when using
            options (1) through (3), however.

        5) least preferred:  use the standard Fortran implementation of the
            BLAS, also available at Netlib (http://www.netlib.org/blas).  This
            will be no faster than option (4), and not portable because of
            C-to-Fortran calling conventions.  Don't bother trying option (5).

    The mechanics of how C calls the BLAS on various computers are as follows:

        * C-BLAS (from the ATLAS library, for example):

            The same interface is used on all computers.  This is the default
            (except on Sun Solaris, or when compiling the Matlab mexFunction).
            SGI Irix provides a C-BLAS interface to its vendor-supplied BLAS.

        * Defaults for calling the Fortran BLAS:
            add underscore, pass scalars by reference, use string arguments.

        * The Fortran BLAS on Sun Solaris (when compiling the Matlab mexFunction
            or when using the Fortran RISC BLAS), SGI, Linux:
            use defaults.

        * Sun Solaris (when using the C-callable Sun Performance library):
            no underscore, pass scalars by value, use character arguments.

        * The Fortran BLAS (ESSL Library) on the IBM RS 6000:
            no underscore, pass scalars by reference, use string arguments.

        * The Fortran BLAS on the HP PA:
            no underscore, pass scalars by reference, use string arguments.
            This has not been tested.  For the umfpack.a library, I recommend
            using the C-BLAS in the ATLAS library instead.  The Matlab
            mexFunction needs to have the -DBLAS_NO_UNDERSCORE compile-time
            flag set.  I've modified the mexopts.sh file to do this, but have
            not tested it.

        * The Fortran BLAS on Windows:
            no underscore, pass scalars by reference, use string arguments.
            This has not been tested.  For the umfpack.a library, I recommend
            using the C-BLAS in the ATLAS library instead.  The Mathworks-
            provided lcc compiler prepends an underscore to all C routine names.
            Thus, dgemm becomes _dgemm.  However, the Mathworks BLAS library has
            dgemm, not _dgemm.  I've contacted Mathworks and so far there is
            no work-around for this problem.  Use another compiler.  If you must
            use lcc then either do not use the BLAS in Matlab or use the C-BLAS.

*/



#ifndef NCBLAS

/* -------------------------------------------------------------------------- */
/* use the C-BLAS (any computer) */
/* -------------------------------------------------------------------------- */

/* This is the default, except for Solaris and IRIX umfpack.a, and for the */
/* mexFunction on any architecture. */

/* If you use the ATLAS C-BLAS, then be sure to set the -I flag to */
/* -I/path/ATLAS/include, where /path/ATLAS is the ATLAS installation */
/* directory.  Note that UMFPACK uses column-major storage for its dense */
/* matrices, but these are not visible to the user. */

#include "cblas.h"

#define BLAS_DGEMM_ROUTINE cblas_dgemm
#define BLAS_DGEMV_ROUTINE cblas_dgemv
#define BLAS_DGER_ROUTINE  cblas_dger

#define BLAS_NO_TRANSPOSE CblasNoTrans
#define BLAS_TRANSPOSE CblasTrans

/* This argument is present only for the C-BLAS: */
#define BLAS_COLUMN_MAJOR_ORDER CblasColMajor,

#define BLAS_SCALAR(n) n
#define BLAS_INT_SCALAR(n) n

#else

/* No such argument when not using the C-BLAS */
#define BLAS_COLUMN_MAJOR_ORDER


/* -------------------------------------------------------------------------- */
/* use Fortran (or other architecture-specific) BLAS */
/* -------------------------------------------------------------------------- */

/* Determine which architecture we're on and set options accordingly. */

#ifdef BLAS_SUNPERF
/* Sun Solaris sunperf library - the default for Solaris umfpack.a */
#include <sunperf.h>
#define BLAS_BY_VALUE
#define BLAS_NO_UNDERSCORE
#define BLAS_CHAR_ARG
#endif

#ifdef BLAS_SCSL
/* SGI SCSL library - the default for SGI umfpack.a */
#if defined (LP64)
#include <scsl_blas_i8.h>
#else
#include <scsl_blas.h>
#endif
#define BLAS_BY_VALUE
#define BLAS_NO_UNDERSCORE
#endif

/* The IBM RS 6000 does not add the underscore */
#if defined (_AIX)
#define BLAS_NO_UNDERSCORE
#endif

/*
    Add your own architecture-dependent settings here.
    For example, to call the Fortran BLAS on Windows, or HP PA:

    #if defined (__win32) || defined (__hppa)
    #define BLAS_NO_UNDERSCORE
    #endif
*/


/* -------------------------------------------------------------------------- */
/* BLAS names */
/* -------------------------------------------------------------------------- */

#if defined (LP64) && defined (BLAS_SUNPERF)

/* 64-bit sunperf BLAS, for Sun Solaris only */
#define BLAS_DGEMM_ROUTINE dgemm_64
#define BLAS_DGEMV_ROUTINE dgemv_64
#define BLAS_DGER_ROUTINE  dger_64

#else

/* naming convention (use underscore, or not) */
#ifdef BLAS_NO_UNDERSCORE
#define BLAS_DGEMM_ROUTINE dgemm
#define BLAS_DGEMV_ROUTINE dgemv
#define BLAS_DGER_ROUTINE  dger
#else
/* default:  add underscore */
#define BLAS_DGEMM_ROUTINE dgemm_
#define BLAS_DGEMV_ROUTINE dgemv_
#define BLAS_DGER_ROUTINE  dger_
#endif

#endif /* LP64 && BLAS_SUNPERF */

/* -------------------------------------------------------------------------- */
/* BLAS scalars */
/* -------------------------------------------------------------------------- */

/* pass scalars by value or by reference */
#ifdef BLAS_BY_VALUE
#define BLAS_SCALAR(n) n
#else
/* default:  pass scalars by reference */
#define BLAS_SCALAR(n) &(n)
#endif

#if defined (BLAS_SCSL) && defined (LP64)
#define BLAS_INT_SCALAR(n) ((long long) n)
#else
#define BLAS_INT_SCALAR(n) BLAS_SCALAR(n)
#endif

/* -------------------------------------------------------------------------- */
/* BLAS strings */
/* -------------------------------------------------------------------------- */

/* pass strings or single characters */
#ifdef BLAS_CHAR_ARG
#define BLAS_NO_TRANSPOSE 'N'
#define BLAS_TRANSPOSE 'T'
#else
/* default:  use string arguments */
#define BLAS_NO_TRANSPOSE "N"
#define BLAS_TRANSPOSE "T"
#endif


#endif /* NCBLAS */



/* -------------------------------------------------------------------------- */
/* BLAS macros, for all interfaces */
/* -------------------------------------------------------------------------- */


/* C = C - A*B, where A is m-by-k, B is k-by-n, and leading dimension is d */
#define BLAS_DGEMM(m,n,k,A,B,C,d) \
{ \
    double alpha = -1.0 ; \
    double beta = 1.0 ; \
    (void) BLAS_DGEMM_ROUTINE (BLAS_COLUMN_MAJOR_ORDER \
        BLAS_NO_TRANSPOSE, BLAS_NO_TRANSPOSE, \
        BLAS_INT_SCALAR (m), BLAS_INT_SCALAR (n), BLAS_INT_SCALAR (k), \
        BLAS_SCALAR (alpha), \
        A, BLAS_INT_SCALAR (d), B, BLAS_INT_SCALAR (d), BLAS_SCALAR (beta), \
        C, BLAS_INT_SCALAR (d)) ; \
}


/* A = A - x*y', where A is m-by-n with leading dimension d */
#define BLAS_DGER(m,n,x,y,A,d) \
{ \
    double alpha = -1.0 ; \
    Int incx = 1 ; \
    (void) BLAS_DGER_ROUTINE (BLAS_COLUMN_MAJOR_ORDER \
        BLAS_INT_SCALAR (m), BLAS_INT_SCALAR (n), BLAS_SCALAR (alpha), \
        x, BLAS_INT_SCALAR (incx), y, BLAS_INT_SCALAR (d), A, BLAS_INT_SCALAR (d)) ; \
}


/* y = y - A'*x, where A is m-by-n with leading dimension d, */
/* and x and y are row vectors with stride d */
#define BLAS_DGEMV_ROW(m,n,A,x,y,d) \
{ \
    double alpha = -1.0 ; \
    double beta = 1.0 ; \
    (void) BLAS_DGEMV_ROUTINE (BLAS_COLUMN_MAJOR_ORDER \
        BLAS_TRANSPOSE, \
        BLAS_INT_SCALAR (m), BLAS_INT_SCALAR (n), BLAS_SCALAR (alpha), \
        A, BLAS_INT_SCALAR (d), x, BLAS_INT_SCALAR (d), BLAS_SCALAR (beta), \
        y, BLAS_INT_SCALAR (d)) ; \
}


/* y = y - A*x, where A is m-by-n with leading dimension d, */
/* and x and y are column vectors with stride 1 */
#define BLAS_DGEMV_COL(m,n,A,x,y,d) \
{ \
    double alpha = -1.0 ; \
    double beta = 1.0 ; \
    Int incx = 1 ; \
    Int incy = 1 ; \
    (void) BLAS_DGEMV_ROUTINE (BLAS_COLUMN_MAJOR_ORDER \
        BLAS_NO_TRANSPOSE, \
        BLAS_INT_SCALAR (m), BLAS_INT_SCALAR (n), BLAS_SCALAR (alpha), \
        A, BLAS_INT_SCALAR (d), x, BLAS_INT_SCALAR (incx), BLAS_SCALAR (beta), \
        y, BLAS_INT_SCALAR (incy)) ; \
}


#endif  /* NBLAS */
\end{verbatim}
}

%-------------------------------------------------------------------------------
\newpage
% References
%-------------------------------------------------------------------------------

\bibliographystyle{plain}
\bibliography{UserGuide}

\end{document}