LFCQH
Computes the RH R factorization of a complex Hermitian positive definite matrix in band Hermitian storage mode and estimate its L1 condition number.
Required Arguments
A — Complex NCODA + 1 by N array containing the N by N positive definite band Hermitian matrix to be factored in band Hermitian storage mode. (Input)
NCODA — Number of upper or lower codiagonals of A. (Input)
FACT — Complex NCODA + 1 by N array containing the RH R factorization of the matrix A. (Output)
If A is not needed, A and FACT can share the same storage locations.
RCOND — Scalar containing an estimate of the reciprocal of the L1 condition number of A. (Output)
Optional Arguments
N — Order of the matrix. (Input)
Default: N = size (A,2).
LDA — Leading dimension of A exactly as specified in the dimension statement of the calling program. (Input)
Default: LDA = size (A,1).
LDFACT — Leading dimension of FACT exactly as specified in the dimension statement of the calling program. (Input)
Default: LDFACT = size (FACT,1).
FORTRAN 90 Interface
Generic: CALL LFCQH (A, NCODA, FACT, RCOND [, …])
Specific: The specific interface names are S_LFCQH and D_LFCQH.
FORTRAN 77 Interface
Single: CALL LFCQH (N, A, LDA, NCODA, FACT, LDFACT, RCOND)
Double: The double precision name is DLFCQH.
Description
Routine LFCQH computes an RH R Cholesky factorization and estimates the condition number of a complex Hermitian positive definite band coefficient matrix. R is an upper triangular band matrix.
The L1 condition number of the matrix A is defined to be κ(A) = ∥A ∥1∥A-1∥1. Since it is expensive to compute ∥A-1∥1, the condition number is only estimated. The estimation algorithm is the same as used by LINPACK and is described by Cline et al. (1979).
If the estimated condition number is greater than 1/ɛ (where ɛ is machine precision), a warning error is issued. This indicates that very small changes in A can cause very large changes in the solution x. Iterative refinement can sometimes find the solution to such a system.
LFCQH fails if any submatrix of R is not positive definite or if R has a zero diagonal element. These errors occur only if A either is very close to a singular matrix or is a matrix which is not positive definite.
The
RH R factors are returned in a form that is compatible with routines
LFIQH,
LFSQH and
LFDQH. To solve systems of equations with multiple right-hand-side vectors, use
LFCQH followed by either
LFIQH or
LFSQH called once for each right-hand side. The routine
LFDQH can be called to compute the determinant of the coefficient matrix after
LFCQH has performed the factorization.
LFCQH is based on the LINPACK routine CPBCO; see Dongarra et al. (1979).
Comments
1. Workspace may be explicitly provided, if desired, by use of L2CQH/DL2CQH. The reference is:
CALL L2CQH (N, A, LDA, NCODA, FACT, LDFACT, RCOND, WK)
The additional argument is:
WK — Complex work vector of length N.
2. Informational errors
Type | Code | Description |
|---|
3 | 1 | The input matrix is algorithmically singular. |
3 | 4 | The input matrix is not Hermitian. It has a diagonal entry with a small imaginary part. |
4 | 2 | The input matrix is not positive definite. |
4 | 4 | The input matrix is not Hermitian. It has a diagonal entry with an imaginary part |
Example
The inverse of a 5
× 5 band Hermitian matrix with one codiagonal is computed.
LFCQH is called to factor the matrix and to check for nonpositive definiteness or ill-conditioning.
LFIQH is called to determine the columns of the inverse.
USE LFCQH_INT
USE LFIQH_INT
USE UMACH_INT
USE WRCRN_INT
! Declare variables
INTEGER N, NCODA, LDA, LDFACT, NOUT
PARAMETER (N=5, NCODA=1, LDA=NCODA+1, LDFACT=LDA)
REAL RCOND
COMPLEX A(LDA,N), AINV(N,N), FACT(LDFACT,N), RES(N), RJ(N)
!
! Set values for A in band Hermitian form
!
! A = ( 0.0+0.0i -1.0+1.0i 1.0+2.0i 0.0+4.0i 1.0+1.0i )
! ( 2.0+0.0i 4.0+0.0i 10.0+0.0i 6.0+0.0i 9.0+0.0i )
!
DATA A/(0.0,0.0), (2.0,0.0), (-1.0,1.0), (4.0, 0.0), (1.0,2.0), &
(10.0,0.0), (0.0,4.0), (6.0,0.0), (1.0,1.0), (9.0,0.0)/
! Factor the matrix A
CALL LFCQH (A, NCODA, FACT, RCOND)
! Set up the columns of the identity
! matrix one at a time in RJ
RJ = (0.0E0,0.0E0)
DO 10 J=1, N
RJ(J) = (1.0E0,0.0E0)
! RJ is the J-th column of the identity
! matrix so the following LFIQH
! reference places the J-th column of
! the inverse of A in the J-th column
! of AINV
CALL LFIQH (A, NCODA, FACT, RJ, AINV(:,J), RES)
RJ(J) = (0.0E0,0.0E0)
10 CONTINUE
! Print the results
CALL UMACH (2, NOUT)
WRITE (NOUT,99999) RCOND, 1.0E0/RCOND
CALL WRCRN (’AINV’, AINV)
!
99999 FORMAT (’ RCOND = ’,F5.3,/,’ L1 Condition number = ’,F6.3)
END
Output
RCOND = 0.067
L1 Condition number = 14.961
AINV
1 2 3 4
1 ( 0.7166, 0.0000) ( 0.2166,-0.2166) (-0.0899,-0.0300) (-0.0207, 0.0622)
2 ( 0.2166, 0.2166) ( 0.4332, 0.0000) (-0.0599,-0.1198) (-0.0829, 0.0415)
3 (-0.0899, 0.0300) (-0.0599, 0.1198) ( 0.1797, 0.0000) ( 0.0000,-0.1244)
4 (-0.0207,-0.0622) (-0.0829,-0.0415) ( 0.0000, 0.1244) ( 0.2592, 0.0000)
5 ( 0.0092, 0.0046) ( 0.0138,-0.0046) (-0.0138,-0.0138) (-0.0288, 0.0288)
5
1 ( 0.0092,-0.0046)
2 ( 0.0138, 0.0046)
3 (-0.0138, 0.0138)
4 (-0.0288,-0.0288)
5 ( 0.1175, 0.0000)
