LFTZG
Computes the LU factorization of a complex general sparse matrix.
Required Arguments
A — Complex vector of length NZ containing the nonzero coefficients of the linear system. (Input)
IROW — Vector of length NZ containing the row numbers of the corresponding elements in A. (Input)
JCOL — Vector of length NZ containing the column numbers of the corresponding elements in A. (Input)
NFAC — On input, the dimension of vector FACT. (Input/Output)
On output, the number of nonzero coefficients in the triangular matrix L and U.
NL — The number of nonzero coefficients in the triangular matrix L excluding the diagonal elements. (Output)
FACT — Complex vector of length NFAC containing the nonzero elements of L (excluding the diagonals) in the first NL locations and the nonzero elements of U in NL + 1 to NFAC locations. (Output)
IRFAC — Vector of length NFAC containing the row numbers of the corresponding elements in FACT. (Output)
JCFAC — Vector of length NFAC containing the column numbers of the corresponding elements in FACT. (Output)
IPVT — Vector of length N containing the row pivoting information for the LU factorization. (Output)
JPVT — Vector of length N containing the column pivoting information for the LU factorization. (Output)
Optional Arguments
N — Number of equations. (Input)
Default: N = size (IPVT,1).
NZ — The number of nonzero coefficients in the linear system. (Input)
Default: NZ = size (A,1).
IPARAM — Parameter vector of length 6. (Input/Output)
Set IPARAM(1) to zero for default values of IPARAM and RPARAM. See Comment 3.
Default: IPARAM = 0.
RPARAM — Parameter vector of length 5. (Input/Output)
See Comment 3.
FORTRAN 90 Interface
Generic: CALL LFTZG (A, IROW, JCOL, NFAC, NL, FACT, IRFAC, JCFAC, IPVT, JPVT [, …])
Specific: The specific interface names are S_LFTZG and D_LFTZG.
FORTRAN 77 Interface
Single: CALL LFTZG (N, NZ, A, IROW, JCOL, IPARAM, RPARAM, NFAC, NL, FACT, IRFAC, JCFAC, IPVT, JPVT)
Double: The double precision name is DLFTZG.
Description
Consider the linear equation
Ax = b
where A is a complex n × n sparse matrix. The sparse coordinate format for the matrix A requires one complex and two integer vectors. The complex array a contains all the nonzeros in A. Let the number of nonzeros be nz. The two integer arrays irow and jcol, each of length nz, contain the row and column indices for these entries in A. That is
Airow(i), icol(i) = a(i), i = 1, …, nz
with all other entries in A zero.
The routine LFTZG performs an LU factorization of the coefficient matrix A. It uses by default a symmetric Markowitz strategy (Crowe et al. 1990) to choose pivots that most likely would reduce fill-ins while maintaining numerical stability. Different strategies are also provided as options for row oriented or column oriented problems. The algorithm can be expressed as
P AQ = LU
where P and Q are the row and column permutation matrices determined by the Markowitz strategy (Duff et al. 1986), and L and U are lower and upper triangular matrices, respectively.
Finally, the solution
x is obtained using
LFSZG by the following calculations:
1) Lz = Pb
2) Uy = z
3) x = Qy
Comments
1. Workspace may be explicitly provided, if desired, by use of L2TZG/DL2TZG. The reference is:
CALL L2TZG (N, NZ, A, IROW, JCOL, IPARAM, RPARAM, NFAC, NL, FACT, IRFAC, JCFAC, IPVT, JPVT, WK, LWK, IWK, LIWK)
The additional arguments are as follows:
WK — Complex work vector of length LWK.
LWK — The length of WK, LWK should be at least MAXNZ.
IWK — Integer work vector of length LIWK.
LIWK — The length of IWK, LIWK should be at least 15N + 4 * MAXNZ.
The workspace limit is determined by MAXNZ, where
MAXNZ = MIN0(LWK, INT(0.25(LIWK-15N)))
2. Informational errors
Type | Code | Description |
---|
3 | 1 | The coefficient matrix is numerically singular. |
3 | 2 | The growth factor is too large to continue. |
3. If the default parameters are desired for LFTZG, then set IPARAM(1) to zero and call the routine LFTZG. Otherwise, if any nondefault parameters are desired for IPARAM or RPARAM. then the following steps should be taken before calling LFTZG:
CALL L4LZG (IPARAM, RPARAM)
Set nondefault values for desired IPARAM, RPARAM elements.
Note that the call to L4LZG will set IPARAM and RPARAM to their default values so only nondefault values need to be set above. The arguments are as follows: |
IPARAM — Integer vector of length 6.
IPARAM(1) = Initialization flag.
IPARAM(2) = The pivoting strategy.
IPARAM(2) | Action |
1 | Markowitz row search |
2 | Markowitz column search |
3 | Symmetric Markowitz search |
IPARAM(3) = The number of rows which have least numbers of nonzero elements that will be searched for a pivotal element.
Default: 3.
IPARAM(4) = The maximal number of nonzero elements in A at any stage of the Gaussian elimination. (Output)
IPARAM(5) = The workspace limit.
IPARAM(5) | Action |
0 | Default limit, see Comment 1. |
integer | This integer value replaces the default workspace limit. When L2TZG is called, the values of LWK and LIWK are used instead of IPARAM(5). |
Default: 0.
IPARAM(6) = Not used in LFTZG.
RPARAM — Real vector of length 5.
RPARAM(1) = The upper limit on the growth factor. The computation stops when the growth factor exceeds the limit.
Default: 10.
RPARAM(2) = The stability factor. The absolute value of the pivotal element must be bigger than the largest element in absolute value in its row divided by RPARAM(2).
Default: 10.0.
RPARAM(3) = Drop-tolerance. Any element in the lower triangular factor L will be removed if its absolute value becomes smaller than the drop-tolerance at any stage of the Gaussian elimination.
Default: 0.0.
RPARAM(4) = The growth factor. It is calculated as the largest element in absolute value in A at any stage of the Gaussian elimination divided by the largest element in absolute value in the original A matrix. (Output)
Large value of the growth factor indicates that an appreciable error in the computed solution is possible.
RPARAM(5) = The value of the smallest pivotal element in absolute value. (Output)
If double precision is required, then DL4LZG is called and RPARAM is declared double precision.
Example
As an example, the following 6 × 6 matrix is factorized, and the outcome is printed:
The sparse coordinate form for A is given by:
USE LFTZG_INT
USE WRCRN_INT
USE WRIRN_INT
INTEGER N, NFAC, NZ
PARAMETER (N=6, NZ=15)
! SPECIFICATIONS FOR LOCAL VARIABLES
INTEGER IPVT(N), IRFAC(45), IROW(NZ), JCFAC(45),&
JCOL(NZ), JPVT(N), NL
COMPLEX A(NZ), FAC(45)
!
DATA A/(3.0,7.0), (3.0,2.0), (-3.0,0.0), (-1.0,3.0), (4.0,2.0),&
(10.0,7.0), (-5.0,4.0), (1.0,6.0), (-1.0,12.0), (-5.0,0.0),&
(12.0,2.0), (-2.0,8.0), (-2.0,-4.0), (-1.0,2.0), (-7.0,7.0)/
DATA IROW/6, 2, 2, 4, 3, 1, 5, 4, 6, 5, 5, 6, 4, 2, 5/
DATA JCOL/6, 2, 3, 5, 3, 1, 1, 4, 1, 4, 5, 2, 1, 4, 6/
DATA NFAC/45/
! Use default options
CALL LFTZG (A, IROW, JCOL, NFAC, NL, FACT, IRFAC, JCFAC, IPVT, JPVT)
!
CALL WRCRN (’fact’,FACT, 1, NFAC, 1)
CALL WRIRN (’ irfac ’,IRFAC, 1, NFAC, 1)
CALL WRIRN (’ jcfac ’,JCFAC, 1, NFAC, 1)
CALL WRIRN (’ p ’,IPVT, 1, N, 1)
CALL WRIRN (’ q ’,JPVT, 1, N, 1)
!
END
Output
fact
1 ( 0.50, 0.85)
2 ( 0.15, -0.41)
3 ( -0.60, 0.30)
4 ( 2.23, -1.97)
5 ( -0.15, 0.50)
6 ( -0.04, 0.26)
7 ( -0.32, -0.17)
8 ( -0.92, 7.46)
9 ( -6.71, -6.42)
10 ( 12.00, 2.00)
11 ( -1.00, 2.00)
12 ( -3.32, 0.21)
13 ( 3.00, 7.00)
14 ( -2.00, 8.00)
15 ( 10.00, 7.00)
16 ( 4.00, 2.00)
irfac
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
3 4 4 5 5 6 6 6 5 5 4 4 3 3 2 1
jcfac
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16
2 3 1 4 2 5 2 6 6 5 6 4 4 3 2 1
p
1 2 3 4 5 6
3 1 6 2 5 4
q
1 2 3 4 5 6
3 1 2 6 5 4