00001 SUBROUTINE CGTCON( NORM, N, DL, D, DU, DU2, IPIV, ANORM, RCOND, 00002 $ WORK, INFO ) 00003 * 00004 * -- LAPACK routine (version 3.2) -- 00005 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00006 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00007 * November 2006 00008 * 00009 * Modified to call CLACN2 in place of CLACON, 10 Feb 03, SJH. 00010 * 00011 * .. Scalar Arguments .. 00012 CHARACTER NORM 00013 INTEGER INFO, N 00014 REAL ANORM, RCOND 00015 * .. 00016 * .. Array Arguments .. 00017 INTEGER IPIV( * ) 00018 COMPLEX D( * ), DL( * ), DU( * ), DU2( * ), WORK( * ) 00019 * .. 00020 * 00021 * Purpose 00022 * ======= 00023 * 00024 * CGTCON estimates the reciprocal of the condition number of a complex 00025 * tridiagonal matrix A using the LU factorization as computed by 00026 * CGTTRF. 00027 * 00028 * An estimate is obtained for norm(inv(A)), and the reciprocal of the 00029 * condition number is computed as RCOND = 1 / (ANORM * norm(inv(A))). 00030 * 00031 * Arguments 00032 * ========= 00033 * 00034 * NORM (input) CHARACTER*1 00035 * Specifies whether the 1-norm condition number or the 00036 * infinity-norm condition number is required: 00037 * = '1' or 'O': 1-norm; 00038 * = 'I': Infinity-norm. 00039 * 00040 * N (input) INTEGER 00041 * The order of the matrix A. N >= 0. 00042 * 00043 * DL (input) COMPLEX array, dimension (N-1) 00044 * The (n-1) multipliers that define the matrix L from the 00045 * LU factorization of A as computed by CGTTRF. 00046 * 00047 * D (input) COMPLEX array, dimension (N) 00048 * The n diagonal elements of the upper triangular matrix U from 00049 * the LU factorization of A. 00050 * 00051 * DU (input) COMPLEX array, dimension (N-1) 00052 * The (n-1) elements of the first superdiagonal of U. 00053 * 00054 * DU2 (input) COMPLEX array, dimension (N-2) 00055 * The (n-2) elements of the second superdiagonal of U. 00056 * 00057 * IPIV (input) INTEGER array, dimension (N) 00058 * The pivot indices; for 1 <= i <= n, row i of the matrix was 00059 * interchanged with row IPIV(i). IPIV(i) will always be either 00060 * i or i+1; IPIV(i) = i indicates a row interchange was not 00061 * required. 00062 * 00063 * ANORM (input) REAL 00064 * If NORM = '1' or 'O', the 1-norm of the original matrix A. 00065 * If NORM = 'I', the infinity-norm of the original matrix A. 00066 * 00067 * RCOND (output) REAL 00068 * The reciprocal of the condition number of the matrix A, 00069 * computed as RCOND = 1/(ANORM * AINVNM), where AINVNM is an 00070 * estimate of the 1-norm of inv(A) computed in this routine. 00071 * 00072 * WORK (workspace) COMPLEX array, dimension (2*N) 00073 * 00074 * INFO (output) INTEGER 00075 * = 0: successful exit 00076 * < 0: if INFO = -i, the i-th argument had an illegal value 00077 * 00078 * ===================================================================== 00079 * 00080 * .. Parameters .. 00081 REAL ONE, ZERO 00082 PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) 00083 * .. 00084 * .. Local Scalars .. 00085 LOGICAL ONENRM 00086 INTEGER I, KASE, KASE1 00087 REAL AINVNM 00088 * .. 00089 * .. Local Arrays .. 00090 INTEGER ISAVE( 3 ) 00091 * .. 00092 * .. External Functions .. 00093 LOGICAL LSAME 00094 EXTERNAL LSAME 00095 * .. 00096 * .. External Subroutines .. 00097 EXTERNAL CGTTRS, CLACN2, XERBLA 00098 * .. 00099 * .. Intrinsic Functions .. 00100 INTRINSIC CMPLX 00101 * .. 00102 * .. Executable Statements .. 00103 * 00104 * Test the input arguments. 00105 * 00106 INFO = 0 00107 ONENRM = NORM.EQ.'1' .OR. LSAME( NORM, 'O' ) 00108 IF( .NOT.ONENRM .AND. .NOT.LSAME( NORM, 'I' ) ) THEN 00109 INFO = -1 00110 ELSE IF( N.LT.0 ) THEN 00111 INFO = -2 00112 ELSE IF( ANORM.LT.ZERO ) THEN 00113 INFO = -8 00114 END IF 00115 IF( INFO.NE.0 ) THEN 00116 CALL XERBLA( 'CGTCON', -INFO ) 00117 RETURN 00118 END IF 00119 * 00120 * Quick return if possible 00121 * 00122 RCOND = ZERO 00123 IF( N.EQ.0 ) THEN 00124 RCOND = ONE 00125 RETURN 00126 ELSE IF( ANORM.EQ.ZERO ) THEN 00127 RETURN 00128 END IF 00129 * 00130 * Check that D(1:N) is non-zero. 00131 * 00132 DO 10 I = 1, N 00133 IF( D( I ).EQ.CMPLX( ZERO ) ) 00134 $ RETURN 00135 10 CONTINUE 00136 * 00137 AINVNM = ZERO 00138 IF( ONENRM ) THEN 00139 KASE1 = 1 00140 ELSE 00141 KASE1 = 2 00142 END IF 00143 KASE = 0 00144 20 CONTINUE 00145 CALL CLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE ) 00146 IF( KASE.NE.0 ) THEN 00147 IF( KASE.EQ.KASE1 ) THEN 00148 * 00149 * Multiply by inv(U)*inv(L). 00150 * 00151 CALL CGTTRS( 'No transpose', N, 1, DL, D, DU, DU2, IPIV, 00152 $ WORK, N, INFO ) 00153 ELSE 00154 * 00155 * Multiply by inv(L')*inv(U'). 00156 * 00157 CALL CGTTRS( 'Conjugate transpose', N, 1, DL, D, DU, DU2, 00158 $ IPIV, WORK, N, INFO ) 00159 END IF 00160 GO TO 20 00161 END IF 00162 * 00163 * Compute the estimate of the reciprocal condition number. 00164 * 00165 IF( AINVNM.NE.ZERO ) 00166 $ RCOND = ( ONE / AINVNM ) / ANORM 00167 * 00168 RETURN 00169 * 00170 * End of CGTCON 00171 * 00172 END
1.7.2