CHFAC
Computes an upper triangular factorization of a real symmetric nonnegative definite matrix.
Required Arguments
A ‑ N by N symmetric nonnegative definite matrix for which an upper triangular factorization is desired. (Input)
Only elements in the upper triangle of A are referenced.
IRANK ‑ Rank of A. (Output)
N ‑ IRANK is the number of effective zero pivots.
R ‑ N by N upper triangular matrix containing the R matrix from a Cholesky decomposition RTR of A. (Output)
The elements of the appropriate rows of R are set to 0.0 if linear dependence of the columns of A is declared. (There are N ‑ IRANK rows of R whose elements are set to 0.0.) If A is not needed, then R and A can share the same storage locations.
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 in the calling program. (Input)
Default: LDA = size (A,1).
TOL ‑ Tolerance used in determining linear dependence. (Input)
TOL = 100 * AMACH(4) is a common choice. See documentation for routine AMACH.
Default: TOL = 1.19e –5 for single precision and 2.2d –14 for double precision.
LDR ‑ Leading dimension of R exactly as specified in the dimension statement in the calling program. (Input)
Default: LDR = size (R,1).
FORTRAN 90 Interface
Generic: CALL CHFAC (A, IRANK, R [, …])
Specific: The specific interface names are S_CHFAC and D_CHFAC.
FORTRAN 77 Interface
Single: CALL CHFAC (N, A, LDA, TOL, IRANK, R, LDR)
Double: The double precision name is DCHFAC.
Description
Routine CHFAC computes a Cholesky factorization RTR = A of an n × n symmetric nonnegative definite matrix A. The matrix R is taken to be an upper triangular matrix. The diagonal elements of R are taken to be nonnegative. If A is singular and has rank r, n ‑ r rows of R have all their elements taken to be zero.
The algorithm is based on the work of Healy (1968). The algorithm proceeds sequentially by columns. The i‑th column is declared to be linearly dependent on the first i ‑ 1 columns if
where ɛ (stored in TOL) is the input tolerance. When a linear dependence is declared, all elements in the i‑th row of R are set to zero.
Modifications due to Farebrother and Berry (1974) and Barrett and Healy (1978) for checking for matrices that are not nonnegative definite are also incorporated. Routine CHFAC declares A to not be nonnegative definite and issues an error message with an error code of 1 if either of the
following conditions is satisfied:
Healy’s (1968) algorithm and CHFAC permit the matrices A and R to occupy the same storage locations. Barrett and Healy (1978) in their remark neglect this fact. Routine CHFAC uses
for aii in the Condition 2 above to remedy this problem.
Comments
1. Informational error
Type | Code | Description |
|---|
3 | 1 | The input matrix is not nonnegative definite within the tolerance defined by TOL. |
2. Elements of row i of R are set to 0.0 if a linear dependence is declared. Linear dependence is declared if
Example
A Cholesky factorization of a 5 × 5 symmetric nonnegative definite matrix is computed. Maindonald (1984, pages 85‑86) discusses in detail the computations for this problem.
! SPECIFICATIONS FOR PARAMETERS
USE IMSL_LIBRARIES
IMPLICIT NONE
INTEGER LDA, LDR, N, J
PARAMETER (N=5, LDA=N, LDR=N)
!
INTEGER IRANK, NOUT
REAL A(LDA,N), R(LDR,N), TOL
!
DATA (A(1,J),J=1,N)/36.0, 12.0, 30.0, 6.0, 18.0/
DATA (A(2,J),J=1,N)/12.0, 20.0, 2.0, 10.0, 22.0/
DATA (A(3,J),J=1,N)/30.0, 2.0, 29.0, 1.0, 7.0/
DATA (A(4,J),J=1,N)/ 6.0, 10.0, 1.0, 14.0, 20.0/
DATA (A(5,J),J=1,N)/ 8.0, 22.0, 7.0, 20.0, 40.0/
!
CALL CHFAC (A, IRANK, R)
CALL UMACH (2, NOUT)
WRITE (NOUT,*) 'IRANK = ', IRANK
CALL WRRRN ('R', R)
END
Output
IRANK = 4
R
1 2 3 4 5
1 6.000 2.000 5.000 1.000 3.000
2 0.000 4.000 -2.000 2.000 4.000
3 0.000 0.000 0.000 0.000 0.000
4 0.000 0.000 0.000 3.000 3.000
5 0.000 0.000 0.000 0.000 2.449
