LNFZD
Computes the numerical Cholesky factorization of a sparse Hermitian matrix A.
Required Arguments
A — Complex vector of length NZ containing the nonzero coefficients of the lower triangle of the linear system. (Input)
IROW — Vector of length NZ containing the row numbers of the corresponding elements in the lower triangle of A. (Input)
JCOL — Vector of length NZ containing the column numbers of the corresponding elements in the lower triangle of A. (Input)
MAXSUB — Number of subscripts contained in array NZSUB as output from subroutine LSCXD/DLSCXD. (Input)
NZSUB — Vector of length MAXSUB containing the row subscripts for the nonzeros in the Cholesky factor in compressed format as output from subroutine LSCXD/DLSCXD. (Input)
INZSUB — Vector of length N + 1 containing pointers for NZSUB as output from subroutine LSCXD/DLSCXD. (Input)
The row subscripts for the nonzeros in column J are stored from location INZSUB(J) to INZSUB(J + 1) - 1.
MAXNZ — Length of RLNZ as output from subroutine LSCXD/DLSCXD. (Input)
ILNZ — Vector of length N + 1 containing pointers to the Cholesky factor as output from subroutine LSCXD/DLSCXD. (Input)
The row subscripts for the nonzeros in column J of the factor are stored from location ILNZ(J) to ILNZ(J + 1) − 1. (ILNZ , NZSUB, INZSUB) sets up the compressed data structure in column ordered form for the Cholesky factor.
IPER — Vector of length N containing the permutation as output from subroutine LSCXD/DLSCXD. (Input)
INVPER — Vector of length N containing the inverse permutation as output from subroutine LSCXD/DLSCXD. (Input)
ISPACE — The storage space needed for the stack of frontal matrices as output from subroutine LSCXD/DLSCXD. (Input)
DIAGNL — Complex vector of length N containing the diagonal of the factor. (Output)
RLNZ — Complex vector of length MAXNZ containing the strictly lower triangle nonzeros of the Cholesky factor. (Output)
RPARAM — Parameter vector containing factorization information. (Output)
RPARAM (1) = smallest diagonal element in absolute value.
RPARAM (2) = largest diagonal element in absolute value.
Optional Arguments
N — Number of equations. (Input)
Default: N = size (IPER,1).
NZ — The number of nonzero coefficients in the linear system. (Input)
Default: NZ = size (A,1).
IJOB — Integer parameter selecting factorization method. (Input)
IJOB = 1 yields factorization in sparse column format.
IJOB = 2 yields factorization using multifrontal method.
Default: IJOB = 1.
ITWKSP — The total workspace needed. (Input)
If the default is desired, set ITWKSP to zero. See Comment 1 for the default.
Default: ITWKSP = 0.
FORTRAN 90 Interface
Generic: CALL LNFZD (A, IROW, JCOL, MAXSUB, NZSUB, INZSUB, MAXNZ, ILNZ, IPER, INVPER, ISPACE, DIAGNL, RLNZ, RPARAM [, …])
Specific: The specific interface names are S_LNFZD and D_LNFZD.
FORTRAN 77 Interface
Single: CALL LNFZD (N, NZ, A, IROW, JCOL, IJOB, MAXSUB, NZSUB, INZSUB, MAXNZ, ILNZ, IPER, INVPER, ISPACE, ITWKSP, DIAGNL, RLNZ, RPARAM)
Double: The double precision name is DLNFZD.
Description
Consider the linear equation
where A is sparse, positive definite and Hermitian. The sparse coordinate format for the matrix A requires one complex and two integer vectors. The complex array a contains all the nonzeros in the lower triangle of A including the diagonal. 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
irow(i) ≥ jcol(i) i = 1, …, nz
with all other entries in the lower triangle of A zero.
The routine LNFZD produces the Cholesky factorization of P APT given the symbolic factorization of A which is computed by LSCXD. That is, this routine computes L which satisfies
P APT= LLH
The diagonal of L is stored in DIAGNL and the strictly lower triangular part of L is stored in compressed subscript form in R = RLNZ as follows. The nonzeros in the jth column of L are stored in locations R(i), …, R(i + k) where i = ILNZ(j) and k = ILNZ(j + 1)− ILNZ(j) − 1. The row subscripts are stored in the vector NZSUB from locations INZSUB(j) to INZSUB(j) + k.
The numerical computations can be carried out in one of two ways. The first method
(when IJOB = 2) performs the factorization using a multifrontal technique. This option requires more storage but in certain cases will be faster. The multifrontal method is based on the routines in Liu (1987). For detailed description of this method, see Liu (1990), also Duff and Reid (1983, 1984), Ashcraft (1987), Ashcraft et al. (1987), and Liu (1986, 1989). The second method (when IJOB = 1) is fully described in George and Liu (1981). This is just the standard factorization method based on the sparse compressed storage scheme.
Comments
1. Workspace may be explicitly provided by use of L2FZD/DL2FZD. The reference is:
CALL L2FZD (N, NZ, A, IROW, JCOL, IJOB, MAXSUB, NZSUB, INZSUB, MAXNZ, ILNZ, IPER, INVPER, ISPACE, DIAGNL, RLNZ, RPARAM, 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 N + 3NZ.
IWK — Integer work vector of length LIWK.
LIWK — The length of IWK, LIWK should be at least 2N.
Note that the parameter ITWKSP is not an argument to this routine.
2. Informational errors
Type |
Code |
Description |
4 |
1 |
The coefficient matrix is not positive definite. |
4 |
2 |
A column without nonzero elements has been found in the coefficient matrix. |
Example
As an example, consider the 3 × 3 linear system:
The number of nonzeros in the lower triangle of A is nz = 5. The sparse coordinate form for the lower triangle of A is given by:
or equivalently by
We first call LSCXD to produce the symbolic information needed to pass on to LNFZD. Then call LNFZD to factor this matrix. The results are displayed below.
USE LNFZD_INT
USE LSCXD_INT
USE WRCRN_INT
INTEGER N, NZ, NRLNZ
PARAMETER (N=3, NZ=5, NRLNZ=5)
!
INTEGER IJOB, ILNZ(N+1), INVPER(N), INZSUB(N+1), IPER(N),&
IROW(NZ), ISPACE, JCOL(NZ), MAXNZ, MAXSUB,&
NZSUB(3*NZ)
REAL RPARAM(2)
COMPLEX A(NZ), DIAGNL(N), RLNZ(NRLNZ)
!
DATA A/(2.0,0.0), (4.0,0.0), (10.0,0.0), (-1.0,-1.0), (1.0,-2.0)/
DATA IROW/1, 2, 3, 2, 3/
DATA JCOL/1, 2, 3, 1, 2/
! Select minimum degree ordering
! for multifrontal method
IJOB = 3
MAXSUB = 3*NZ
CALL LSCXD (IROW, JCOL, NZSUB, INZSUB, MAXNZ, ILNZ, INVPER, &
IJOB=IJOB, MAXSUB=MAXSUB)
! Check if NRLNZ is large enough
IF (NRLNZ .GE. MAXNZ) THEN
! Choose multifrontal method
IJOB = 2
CALL LNFZD (A, IROW, JCOL, MAXSUB, NZSUB, INZSUB, MAXNZ, &
ILNZ, IPER, INVPER, ISPACE, DIAGNL, RLNZ, RPARAM, &
IJOB=IJOB)
! Print results
CALL WRCRN (’ diagnl ’, DIAGNL, 1, N, 1)
CALL WRCRN (’ rlnz ’, RLNZ, 1, MAXNZ, 1)
END IF
!
END
diagnl
1 2 3
( 1.414, 0.000) ( 1.732, 0.000) ( 2.887, 0.000)
rlnz
1 2
(-0.707,-0.707) ( 0.577,-1.155)