NAG CL Interfacef11mec (direct_​real_​gen_​lu)

Settings help

CL Name Style:

1Purpose

f11mec computes the $LU$ factorization of a real sparse matrix in compressed column (Harwell–Boeing), column-permuted format.

2Specification

 #include
 void f11mec (Integer n, const Integer irowix[], const double a[], Integer iprm[], double thresh, Integer nzlmx, Integer *nzlumx, Integer nzumx, Integer il[], double lval[], Integer iu[], double uval[], Integer *nnzl, Integer *nnzu, double *flop, NagError *fail)
The function may be called by the names: f11mec, nag_sparse_direct_real_gen_lu or nag_superlu_lu_factorize.

3Description

Given a real sparse matrix $A$, f11mec computes an $LU$ factorization of $A$ with partial pivoting, ${P}_{r}A{P}_{c}=LU$, where ${P}_{r}$ is a row permutation matrix (computed by f11mec), ${P}_{c}$ is a (supplied) column permutation matrix, $L$ is unit lower triangular and $U$ is upper triangular. The column permutation matrix, ${P}_{c}$, must be computed by a prior call to f11mdc. The matrix $A$ must be presented in the column permuted, compressed column (Harwell–Boeing) format.
The $LU$ factorization is output in the form of four one-dimensional arrays: integer arrays il and iu and real-valued arrays lval and uval. These describe the sparsity pattern and numerical values in the $L$ and $U$ matrices. The minimum required dimensions of these arrays cannot be given as a simple function of the size arguments (order and number of nonzero values) of the matrix $A$. This is due to unpredictable fill-in created by partial pivoting. f11mec will, on return, indicate which dimensions of these arrays were not adequate for the computation or (in the case of one of them) give a firm bound. You should then allocate more storage and try again.

4References

Demmel J W, Eisenstat S C, Gilbert J R, Li X S and Li J W H (1999) A supernodal approach to sparse partial pivoting SIAM J. Matrix Anal. Appl. 20 720–755
Demmel J W, Gilbert J R and Li X S (1999) An asynchronous parallel supernodal algorithm for sparse gaussian elimination SIAM J. Matrix Anal. Appl. 20 915–952

5Arguments

1: $\mathbf{n}$Integer Input
On entry: $n$, the order of the matrix $A$.
Constraint: ${\mathbf{n}}\ge 0$.
2: $\mathbf{irowix}\left[\mathit{dim}\right]$const Integer Input
Note: the dimension, dim, of the array irowix must be at least $\mathit{nnz}$, the number of nonzeros of the sparse matrix $A$.
On entry: the row index array of sparse matrix $A$. See Section 2.1.3 in the F11 Chapter Introduction.
3: $\mathbf{a}\left[\mathit{dim}\right]$const double Input
Note: the dimension, dim, of the array a must be at least $\mathit{nnz}$, the number of nonzeros of the sparse matrix $A$.
On entry: the array of nonzero values in the sparse matrix $A$.
4: $\mathbf{iprm}\left[7×{\mathbf{n}}\right]$Integer Input/Output
On entry: contains the column permutation which defines the permutation ${P}_{c}$ and associated data structures as computed by function f11mdc.
On exit: part of the array is modified to record the row permutation ${P}_{r}$ determined by pivoting.
5: $\mathbf{thresh}$double Input
On entry: the diagonal pivoting threshold, $t$. At step $j$ of the Gaussian elimination, if $|{A}_{jj}|\ge t\left(\underset{i\ge j}{\mathrm{max}}\phantom{\rule{0.25em}{0ex}}|{A}_{ij}|\right)$, use ${A}_{jj}$ as a pivot, otherwise use $\underset{i\ge j}{\mathrm{max}}\phantom{\rule{0.25em}{0ex}}|{A}_{ij}|$. A value of $t=1$ corresponds to partial pivoting, a value of $t=0$ corresponds to always choosing the pivot on the diagonal (unless it is zero).
Suggested value: ${\mathbf{thresh}}=1.0$. Smaller values may result in a faster factorization, but the benefits are likely to be small in most cases. It might be possible to use ${\mathbf{thresh}}=0.0$ if you are confident about the stability of the factorization, for example, if $A$ is diagonally dominant.
Constraint: $0.0\le {\mathbf{thresh}}\le 1.0$.
6: $\mathbf{nzlmx}$Integer Input
On entry: indicates the available size of array il. The dimension of il should be at least $7×{\mathbf{n}}+{\mathbf{nzlmx}}+4$. A good range for nzlmx that works for many problems is $\mathit{nnz}$ to $8×\mathit{nnz}$, where $\mathit{nnz}$ is the number of nonzeros in the sparse matrix $A$. If, on exit, ${\mathbf{fail}}\mathbf{.}\mathbf{code}=$ NE_NZLMX_TOO_SMALL, the given nzlmx was too small and you should attempt to provide more storage and call the function again.
Constraint: ${\mathbf{nzlmx}}\ge 1$.
7: $\mathbf{nzlumx}$Integer * Input/Output
On entry: indicates the available size of array lval. The dimension of lval should be at least nzlumx.
Constraint: ${\mathbf{nzlumx}}\ge 1$.
On exit: if ${\mathbf{fail}}\mathbf{.}\mathbf{code}=$ NE_NZLUMX_TOO_SMALL, the given nzlumx was too small and is reset to a value that will be sufficient. You should then provide the indicated storage and call the function again.
8: $\mathbf{nzumx}$Integer Input
On entry: indicates the available sizes of arrays iu and uval. The dimension of iu should be at least $2×{\mathbf{n}}+{\mathbf{nzumx}}+1$ and the dimension of uval should be at least nzumx. A good range for nzumx that works for many problems is $\mathit{nnz}$ to $8×\mathit{nnz}$, where $\mathit{nnz}$ is the number of nonzeros in the sparse matrix $A$. If, on exit, ${\mathbf{fail}}\mathbf{.}\mathbf{code}=$ NE_NZUMX_TOO_SMALL, the given nzumx was too small and you should attempt to provide more storage and call the function again.
Constraint: ${\mathbf{nzumx}}\ge 1$.
9: $\mathbf{il}\left[7×{\mathbf{n}}+{\mathbf{nzlmx}}+4\right]$Integer Output
On exit: encapsulates the sparsity pattern of matrix $L$.
10: $\mathbf{lval}\left[{\mathbf{nzlumx}}\right]$double Output
On exit: records the nonzero values of matrix $L$ and some of the nonzero values of matrix $U$.
11: $\mathbf{iu}\left[2×{\mathbf{n}}+{\mathbf{nzumx}}+1\right]$Integer Output
On exit: encapsulates the sparsity pattern of matrix $U$.
12: $\mathbf{uval}\left[{\mathbf{nzumx}}\right]$double Output
On exit: records some of the nonzero values of matrix $U$.
13: $\mathbf{nnzl}$Integer * Output
On exit: the number of nonzero values in the matrix $L$.
14: $\mathbf{nnzu}$Integer * Output
On exit: the number of nonzero values in the matrix $U$.
15: $\mathbf{flop}$double * Output
On exit: the number of floating-point operations performed.
16: $\mathbf{fail}$NagError * Input/Output
The NAG error argument (see Section 7 in the Introduction to the NAG Library CL Interface).

6Error Indicators and Warnings

NE_ALLOC_FAIL
Dynamic memory allocation failed.
See Section 3.1.2 in the Introduction to the NAG Library CL Interface for further information.
On entry, argument $⟨\mathit{\text{value}}⟩$ had an illegal value.
NE_INT
On entry, ${\mathbf{n}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{n}}\ge 0$.
On entry, ${\mathbf{nzlmx}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{nzlmx}}\ge 1$.
On entry, ${\mathbf{nzlumx}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{nzlumx}}\ge 1$.
On entry, ${\mathbf{nzumx}}=⟨\mathit{\text{value}}⟩$.
Constraint: ${\mathbf{nzumx}}\ge 1$.
NE_INTERNAL_ERROR
An internal error has occurred in this function. Check the function call and any array sizes. If the call is correct then please contact NAG for assistance.
See Section 7.5 in the Introduction to the NAG Library CL Interface for further information.
NE_NO_LICENCE
Your licence key may have expired or may not have been installed correctly.
See Section 8 in the Introduction to the NAG Library CL Interface for further information.
NE_NZLMX_TOO_SMALL
Insufficient nzlmx.
NE_NZLUMX_TOO_SMALL
Insufficient nzlumx.
NE_NZUMX_TOO_SMALL
Insufficient nzumx.
NE_REAL
On entry, ${\mathbf{thresh}}=⟨\mathit{\text{value}}⟩$.
Constraint: $0.0\le {\mathbf{thresh}}\le 1.0$.
NE_SINGULAR_MATRIX
The matrix is singular – no factorization possible.

7Accuracy

The computed factors $L$ and $U$ are the exact factors of a perturbed matrix $A+E$, where
 $|E| ≤ c (n) ε |L| |U| ,$
$c\left(n\right)$ is a modest linear function of $n$, and $\epsilon$ is the machine precision, when partial pivoting is used. If no partial pivoting is used, the factorization accuracy can be considerably worse. A call to f11mmc after f11mec can help determine the quality of the factorization.

8Parallelism and Performance

f11mec is threaded by NAG for parallel execution in multithreaded implementations of the NAG Library.
f11mec makes calls to BLAS and/or LAPACK routines, which may be threaded within the vendor library used by this implementation. Consult the documentation for the vendor library for further information.
Please consult the X06 Chapter Introduction for information on how to control and interrogate the OpenMP environment used within this function. Please also consult the Users' Note for your implementation for any additional implementation-specific information.

The total number of floating-point operations depends on the sparsity pattern of the matrix $A$.
A call to f11mec may be followed by calls to the functions:
• f11mfc to solve $AX=B$ or ${A}^{\mathrm{T}}X=B$;
• f11mgc to estimate the condition number of $A$;
• f11mmc to estimate the reciprocal pivot growth of the $LU$ factorization.

10Example

This example computes the $LU$ factorization of the matrix $A$, where
 $A=( 2.00 1.00 0 0 0 0 0 1.00 -1.00 0 4.00 0 1.00 0 1.00 0 0 0 1.00 2.00 0 -2.00 0 0 3.00 ) .$

10.1Program Text

Program Text (f11mece.c)

10.2Program Data

Program Data (f11mece.d)

10.3Program Results

Program Results (f11mece.r)