 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 REAL FUNCTION CLA_GERCOND_C( TRANS, N, A, LDA, AF, LDAF, IPIV, C, 00002 $ CAPPLY, INFO, WORK, RWORK ) 00003 * 00004 * -- LAPACK routine (version 3.2.1) -- 00005 * -- Contributed by James Demmel, Deaglan Halligan, Yozo Hida and -- 00006 * -- Jason Riedy of Univ. of California Berkeley. -- 00007 * -- April 2009 -- 00008 * 00009 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00010 * -- Univ. of California Berkeley and NAG Ltd. -- 00011 * 00012 IMPLICIT NONE 00013 * .. 00014 * .. Scalar Aguments .. 00015 CHARACTER TRANS 00016 LOGICAL CAPPLY 00017 INTEGER N, LDA, LDAF, INFO 00018 * .. 00019 * .. Array Arguments .. 00020 INTEGER IPIV( * ) 00021 COMPLEX A( LDA, * ), AF( LDAF, * ), WORK( * ) 00022 REAL C( * ), RWORK( * ) 00023 * .. 00024 * 00025 * Purpose 00026 * ======= 00027 * 00028 * CLA_GERCOND_C computes the infinity norm condition number of 00029 * op(A) * inv(diag(C)) where C is a REAL vector. 00030 * 00031 * Arguments 00032 * ========= 00033 * 00034 * TRANS (input) CHARACTER*1 00035 * Specifies the form of the system of equations: 00036 * = 'N': A * X = B (No transpose) 00037 * = 'T': A**T * X = B (Transpose) 00038 * = 'C': A**H * X = B (Conjugate Transpose = Transpose) 00039 * 00040 * N (input) INTEGER 00041 * The number of linear equations, i.e., the order of the 00042 * matrix A. N >= 0. 00043 * 00044 * A (input) COMPLEX array, dimension (LDA,N) 00045 * On entry, the N-by-N matrix A 00046 * 00047 * LDA (input) INTEGER 00048 * The leading dimension of the array A. LDA >= max(1,N). 00049 * 00050 * AF (input) COMPLEX array, dimension (LDAF,N) 00051 * The factors L and U from the factorization 00052 * A = P*L*U as computed by CGETRF. 00053 * 00054 * LDAF (input) INTEGER 00055 * The leading dimension of the array AF. LDAF >= max(1,N). 00056 * 00057 * IPIV (input) INTEGER array, dimension (N) 00058 * The pivot indices from the factorization A = P*L*U 00059 * as computed by CGETRF; row i of the matrix was interchanged 00060 * with row IPIV(i). 00061 * 00062 * C (input) REAL array, dimension (N) 00063 * The vector C in the formula op(A) * inv(diag(C)). 00064 * 00065 * CAPPLY (input) LOGICAL 00066 * If .TRUE. then access the vector C in the formula above. 00067 * 00068 * INFO (output) INTEGER 00069 * = 0: Successful exit. 00070 * i > 0: The ith argument is invalid. 00071 * 00072 * WORK (input) COMPLEX array, dimension (2*N). 00073 * Workspace. 00074 * 00075 * RWORK (input) REAL array, dimension (N). 00076 * Workspace. 00077 * 00078 * ===================================================================== 00079 * 00080 * .. Local Scalars .. 00081 LOGICAL NOTRANS 00082 INTEGER KASE, I, J 00083 REAL AINVNM, ANORM, TMP 00084 COMPLEX ZDUM 00085 * .. 00086 * .. Local Arrays .. 00087 INTEGER ISAVE( 3 ) 00088 * .. 00089 * .. External Functions .. 00090 LOGICAL LSAME 00091 EXTERNAL LSAME 00092 * .. 00093 * .. External Subroutines .. 00094 EXTERNAL CLACN2, CGETRS, XERBLA 00095 * .. 00096 * .. Intrinsic Functions .. 00097 INTRINSIC ABS, MAX, REAL, AIMAG 00098 * .. 00099 * .. Statement Functions .. 00100 REAL CABS1 00101 * .. 00102 * .. Statement Function Definitions .. 00103 CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) ) 00104 * .. 00105 * .. Executable Statements .. 00106 CLA_GERCOND_C = 0.0E+0 00107 * 00108 INFO = 0 00109 NOTRANS = LSAME( TRANS, 'N' ) 00110 IF ( .NOT. NOTRANS .AND. .NOT. LSAME( TRANS, 'T' ) .AND. .NOT. 00111 $ LSAME( TRANS, 'C' ) ) THEN 00112 ELSE IF( N.LT.0 ) THEN 00113 INFO = -2 00114 END IF 00115 IF( INFO.NE.0 ) THEN 00116 CALL XERBLA( 'CLA_GERCOND_C', -INFO ) 00117 RETURN 00118 END IF 00119 * 00120 * Compute norm of op(A)*op2(C). 00121 * 00122 ANORM = 0.0E+0 00123 IF ( NOTRANS ) THEN 00124 DO I = 1, N 00125 TMP = 0.0E+0 00126 IF ( CAPPLY ) THEN 00127 DO J = 1, N 00128 TMP = TMP + CABS1( A( I, J ) ) / C( J ) 00129 END DO 00130 ELSE 00131 DO J = 1, N 00132 TMP = TMP + CABS1( A( I, J ) ) 00133 END DO 00134 END IF 00135 RWORK( I ) = TMP 00136 ANORM = MAX( ANORM, TMP ) 00137 END DO 00138 ELSE 00139 DO I = 1, N 00140 TMP = 0.0E+0 00141 IF ( CAPPLY ) THEN 00142 DO J = 1, N 00143 TMP = TMP + CABS1( A( J, I ) ) / C( J ) 00144 END DO 00145 ELSE 00146 DO J = 1, N 00147 TMP = TMP + CABS1( A( J, I ) ) 00148 END DO 00149 END IF 00150 RWORK( I ) = TMP 00151 ANORM = MAX( ANORM, TMP ) 00152 END DO 00153 END IF 00154 * 00155 * Quick return if possible. 00156 * 00157 IF( N.EQ.0 ) THEN 00158 CLA_GERCOND_C = 1.0E+0 00159 RETURN 00160 ELSE IF( ANORM .EQ. 0.0E+0 ) THEN 00161 RETURN 00162 END IF 00163 * 00164 * Estimate the norm of inv(op(A)). 00165 * 00166 AINVNM = 0.0E+0 00167 * 00168 KASE = 0 00169 10 CONTINUE 00170 CALL CLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE ) 00171 IF( KASE.NE.0 ) THEN 00172 IF( KASE.EQ.2 ) THEN 00173 * 00174 * Multiply by R. 00175 * 00176 DO I = 1, N 00177 WORK( I ) = WORK( I ) * RWORK( I ) 00178 END DO 00179 * 00180 IF (NOTRANS) THEN 00181 CALL CGETRS( 'No transpose', N, 1, AF, LDAF, IPIV, 00182 $ WORK, N, INFO ) 00183 ELSE 00184 CALL CGETRS( 'Conjugate transpose', N, 1, AF, LDAF, IPIV, 00185 $ WORK, N, INFO ) 00186 ENDIF 00187 * 00188 * Multiply by inv(C). 00189 * 00190 IF ( CAPPLY ) THEN 00191 DO I = 1, N 00192 WORK( I ) = WORK( I ) * C( I ) 00193 END DO 00194 END IF 00195 ELSE 00196 * 00197 * Multiply by inv(C**H). 00198 * 00199 IF ( CAPPLY ) THEN 00200 DO I = 1, N 00201 WORK( I ) = WORK( I ) * C( I ) 00202 END DO 00203 END IF 00204 * 00205 IF ( NOTRANS ) THEN 00206 CALL CGETRS( 'Conjugate transpose', N, 1, AF, LDAF, IPIV, 00207 $ WORK, N, INFO ) 00208 ELSE 00209 CALL CGETRS( 'No transpose', N, 1, AF, LDAF, IPIV, 00210 $ WORK, N, INFO ) 00211 END IF 00212 * 00213 * Multiply by R. 00214 * 00215 DO I = 1, N 00216 WORK( I ) = WORK( I ) * RWORK( I ) 00217 END DO 00218 END IF 00219 GO TO 10 00220 END IF 00221 * 00222 * Compute the estimate of the reciprocal condition number. 00223 * 00224 IF( AINVNM .NE. 0.0E+0 ) 00225 $ CLA_GERCOND_C = 1.0E+0 / AINVNM 00226 * 00227 RETURN 00228 * 00229 END
1.7.3 
