![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b CLA_GBRCOND_X 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download CLA_GBRCOND_X + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cla_gbrcond_x.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cla_gbrcond_x.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cla_gbrcond_x.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * REAL FUNCTION CLA_GBRCOND_X( TRANS, N, KL, KU, AB, LDAB, AFB, 00022 * LDAFB, IPIV, X, INFO, WORK, RWORK ) 00023 * 00024 * .. Scalar Arguments .. 00025 * CHARACTER TRANS 00026 * INTEGER N, KL, KU, KD, KE, LDAB, LDAFB, INFO 00027 * .. 00028 * .. Array Arguments .. 00029 * INTEGER IPIV( * ) 00030 * COMPLEX AB( LDAB, * ), AFB( LDAFB, * ), WORK( * ), 00031 * $ X( * ) 00032 * REAL RWORK( * ) 00033 * .. 00034 * 00035 * 00036 *> \par Purpose: 00037 * ============= 00038 *> 00039 *> \verbatim 00040 *> 00041 *> CLA_GBRCOND_X Computes the infinity norm condition number of 00042 *> op(A) * diag(X) where X is a COMPLEX vector. 00043 *> \endverbatim 00044 * 00045 * Arguments: 00046 * ========== 00047 * 00048 *> \param[in] TRANS 00049 *> \verbatim 00050 *> TRANS is CHARACTER*1 00051 *> Specifies the form of the system of equations: 00052 *> = 'N': A * X = B (No transpose) 00053 *> = 'T': A**T * X = B (Transpose) 00054 *> = 'C': A**H * X = B (Conjugate Transpose = Transpose) 00055 *> \endverbatim 00056 *> 00057 *> \param[in] N 00058 *> \verbatim 00059 *> N is INTEGER 00060 *> The number of linear equations, i.e., the order of the 00061 *> matrix A. N >= 0. 00062 *> \endverbatim 00063 *> 00064 *> \param[in] KL 00065 *> \verbatim 00066 *> KL is INTEGER 00067 *> The number of subdiagonals within the band of A. KL >= 0. 00068 *> \endverbatim 00069 *> 00070 *> \param[in] KU 00071 *> \verbatim 00072 *> KU is INTEGER 00073 *> The number of superdiagonals within the band of A. KU >= 0. 00074 *> \endverbatim 00075 *> 00076 *> \param[in] AB 00077 *> \verbatim 00078 *> AB is COMPLEX array, dimension (LDAB,N) 00079 *> On entry, the matrix A in band storage, in rows 1 to KL+KU+1. 00080 *> The j-th column of A is stored in the j-th column of the 00081 *> array AB as follows: 00082 *> AB(KU+1+i-j,j) = A(i,j) for max(1,j-KU)<=i<=min(N,j+kl) 00083 *> \endverbatim 00084 *> 00085 *> \param[in] LDAB 00086 *> \verbatim 00087 *> LDAB is INTEGER 00088 *> The leading dimension of the array AB. LDAB >= KL+KU+1. 00089 *> \endverbatim 00090 *> 00091 *> \param[in] AFB 00092 *> \verbatim 00093 *> AFB is COMPLEX array, dimension (LDAFB,N) 00094 *> Details of the LU factorization of the band matrix A, as 00095 *> computed by CGBTRF. U is stored as an upper triangular 00096 *> band matrix with KL+KU superdiagonals in rows 1 to KL+KU+1, 00097 *> and the multipliers used during the factorization are stored 00098 *> in rows KL+KU+2 to 2*KL+KU+1. 00099 *> \endverbatim 00100 *> 00101 *> \param[in] LDAFB 00102 *> \verbatim 00103 *> LDAFB is INTEGER 00104 *> The leading dimension of the array AFB. LDAFB >= 2*KL+KU+1. 00105 *> \endverbatim 00106 *> 00107 *> \param[in] IPIV 00108 *> \verbatim 00109 *> IPIV is INTEGER array, dimension (N) 00110 *> The pivot indices from the factorization A = P*L*U 00111 *> as computed by CGBTRF; row i of the matrix was interchanged 00112 *> with row IPIV(i). 00113 *> \endverbatim 00114 *> 00115 *> \param[in] X 00116 *> \verbatim 00117 *> X is COMPLEX array, dimension (N) 00118 *> The vector X in the formula op(A) * diag(X). 00119 *> \endverbatim 00120 *> 00121 *> \param[out] INFO 00122 *> \verbatim 00123 *> INFO is INTEGER 00124 *> = 0: Successful exit. 00125 *> i > 0: The ith argument is invalid. 00126 *> \endverbatim 00127 *> 00128 *> \param[in] WORK 00129 *> \verbatim 00130 *> WORK is COMPLEX array, dimension (2*N). 00131 *> Workspace. 00132 *> \endverbatim 00133 *> 00134 *> \param[in] RWORK 00135 *> \verbatim 00136 *> RWORK is REAL array, dimension (N). 00137 *> Workspace. 00138 *> \endverbatim 00139 * 00140 * Authors: 00141 * ======== 00142 * 00143 *> \author Univ. of Tennessee 00144 *> \author Univ. of California Berkeley 00145 *> \author Univ. of Colorado Denver 00146 *> \author NAG Ltd. 00147 * 00148 *> \date November 2011 00149 * 00150 *> \ingroup complexGBcomputational 00151 * 00152 * ===================================================================== 00153 REAL FUNCTION CLA_GBRCOND_X( TRANS, N, KL, KU, AB, LDAB, AFB, 00154 $ LDAFB, IPIV, X, INFO, WORK, RWORK ) 00155 * 00156 * -- LAPACK computational routine (version 3.4.0) -- 00157 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00158 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00159 * November 2011 00160 * 00161 * .. Scalar Arguments .. 00162 CHARACTER TRANS 00163 INTEGER N, KL, KU, KD, KE, LDAB, LDAFB, INFO 00164 * .. 00165 * .. Array Arguments .. 00166 INTEGER IPIV( * ) 00167 COMPLEX AB( LDAB, * ), AFB( LDAFB, * ), WORK( * ), 00168 $ X( * ) 00169 REAL RWORK( * ) 00170 * .. 00171 * 00172 * ===================================================================== 00173 * 00174 * .. Local Scalars .. 00175 LOGICAL NOTRANS 00176 INTEGER KASE, I, J 00177 REAL AINVNM, ANORM, TMP 00178 COMPLEX ZDUM 00179 * .. 00180 * .. Local Arrays .. 00181 INTEGER ISAVE( 3 ) 00182 * .. 00183 * .. External Functions .. 00184 LOGICAL LSAME 00185 EXTERNAL LSAME 00186 * .. 00187 * .. External Subroutines .. 00188 EXTERNAL CLACN2, CGBTRS, XERBLA 00189 * .. 00190 * .. Intrinsic Functions .. 00191 INTRINSIC ABS, MAX 00192 * .. 00193 * .. Statement Functions .. 00194 REAL CABS1 00195 * .. 00196 * .. Statement Function Definitions .. 00197 CABS1( ZDUM ) = ABS( REAL( ZDUM ) ) + ABS( AIMAG( ZDUM ) ) 00198 * .. 00199 * .. Executable Statements .. 00200 * 00201 CLA_GBRCOND_X = 0.0E+0 00202 * 00203 INFO = 0 00204 NOTRANS = LSAME( TRANS, 'N' ) 00205 IF ( .NOT. NOTRANS .AND. .NOT. LSAME(TRANS, 'T') .AND. .NOT. 00206 $ LSAME( TRANS, 'C' ) ) THEN 00207 INFO = -1 00208 ELSE IF( N.LT.0 ) THEN 00209 INFO = -2 00210 ELSE IF( KL.LT.0 .OR. KL.GT.N-1 ) THEN 00211 INFO = -3 00212 ELSE IF( KU.LT.0 .OR. KU.GT.N-1 ) THEN 00213 INFO = -4 00214 ELSE IF( LDAB.LT.KL+KU+1 ) THEN 00215 INFO = -6 00216 ELSE IF( LDAFB.LT.2*KL+KU+1 ) THEN 00217 INFO = -8 00218 END IF 00219 IF( INFO.NE.0 ) THEN 00220 CALL XERBLA( 'CLA_GBRCOND_X', -INFO ) 00221 RETURN 00222 END IF 00223 * 00224 * Compute norm of op(A)*op2(C). 00225 * 00226 KD = KU + 1 00227 KE = KL + 1 00228 ANORM = 0.0 00229 IF ( NOTRANS ) THEN 00230 DO I = 1, N 00231 TMP = 0.0E+0 00232 DO J = MAX( I-KL, 1 ), MIN( I+KU, N ) 00233 TMP = TMP + CABS1( AB( KD+I-J, J) * X( J ) ) 00234 END DO 00235 RWORK( I ) = TMP 00236 ANORM = MAX( ANORM, TMP ) 00237 END DO 00238 ELSE 00239 DO I = 1, N 00240 TMP = 0.0E+0 00241 DO J = MAX( I-KL, 1 ), MIN( I+KU, N ) 00242 TMP = TMP + CABS1( AB( KE-I+J, I ) * X( J ) ) 00243 END DO 00244 RWORK( I ) = TMP 00245 ANORM = MAX( ANORM, TMP ) 00246 END DO 00247 END IF 00248 * 00249 * Quick return if possible. 00250 * 00251 IF( N.EQ.0 ) THEN 00252 CLA_GBRCOND_X = 1.0E+0 00253 RETURN 00254 ELSE IF( ANORM .EQ. 0.0E+0 ) THEN 00255 RETURN 00256 END IF 00257 * 00258 * Estimate the norm of inv(op(A)). 00259 * 00260 AINVNM = 0.0E+0 00261 * 00262 KASE = 0 00263 10 CONTINUE 00264 CALL CLACN2( N, WORK( N+1 ), WORK, AINVNM, KASE, ISAVE ) 00265 IF( KASE.NE.0 ) THEN 00266 IF( KASE.EQ.2 ) THEN 00267 * 00268 * Multiply by R. 00269 * 00270 DO I = 1, N 00271 WORK( I ) = WORK( I ) * RWORK( I ) 00272 END DO 00273 * 00274 IF ( NOTRANS ) THEN 00275 CALL CGBTRS( 'No transpose', N, KL, KU, 1, AFB, LDAFB, 00276 $ IPIV, WORK, N, INFO ) 00277 ELSE 00278 CALL CGBTRS( 'Conjugate transpose', N, KL, KU, 1, AFB, 00279 $ LDAFB, IPIV, WORK, N, INFO ) 00280 ENDIF 00281 * 00282 * Multiply by inv(X). 00283 * 00284 DO I = 1, N 00285 WORK( I ) = WORK( I ) / X( I ) 00286 END DO 00287 ELSE 00288 * 00289 * Multiply by inv(X**H). 00290 * 00291 DO I = 1, N 00292 WORK( I ) = WORK( I ) / X( I ) 00293 END DO 00294 * 00295 IF ( NOTRANS ) THEN 00296 CALL CGBTRS( 'Conjugate transpose', N, KL, KU, 1, AFB, 00297 $ LDAFB, IPIV, WORK, N, INFO ) 00298 ELSE 00299 CALL CGBTRS( 'No transpose', N, KL, KU, 1, AFB, LDAFB, 00300 $ IPIV, WORK, N, INFO ) 00301 END IF 00302 * 00303 * Multiply by R. 00304 * 00305 DO I = 1, N 00306 WORK( I ) = WORK( I ) * RWORK( I ) 00307 END DO 00308 END IF 00309 GO TO 10 00310 END IF 00311 * 00312 * Compute the estimate of the reciprocal condition number. 00313 * 00314 IF( AINVNM .NE. 0.0E+0 ) 00315 $ CLA_GBRCOND_X = 1.0E+0 / AINVNM 00316 * 00317 RETURN 00318 * 00319 END