1*> \brief \b CTRT03 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8* Definition: 9* =========== 10* 11* SUBROUTINE CTRT03( UPLO, TRANS, DIAG, N, NRHS, A, LDA, SCALE, 12* CNORM, TSCAL, X, LDX, B, LDB, WORK, RESID ) 13* 14* .. Scalar Arguments .. 15* CHARACTER DIAG, TRANS, UPLO 16* INTEGER LDA, LDB, LDX, N, NRHS 17* REAL RESID, SCALE, TSCAL 18* .. 19* .. Array Arguments .. 20* REAL CNORM( * ) 21* COMPLEX A( LDA, * ), B( LDB, * ), WORK( * ), 22* $ X( LDX, * ) 23* .. 24* 25* 26*> \par Purpose: 27* ============= 28*> 29*> \verbatim 30*> 31*> CTRT03 computes the residual for the solution to a scaled triangular 32*> system of equations A*x = s*b, A**T *x = s*b, or A**H *x = s*b. 33*> Here A is a triangular matrix, A**T denotes the transpose of A, A**H 34*> denotes the conjugate transpose of A, s is a scalar, and x and b are 35*> N by NRHS matrices. The test ratio is the maximum over the number of 36*> right hand sides of 37*> norm(s*b - op(A)*x) / ( norm(op(A)) * norm(x) * EPS ), 38*> where op(A) denotes A, A**T, or A**H, and EPS is the machine epsilon. 39*> \endverbatim 40* 41* Arguments: 42* ========== 43* 44*> \param[in] UPLO 45*> \verbatim 46*> UPLO is CHARACTER*1 47*> Specifies whether the matrix A is upper or lower triangular. 48*> = 'U': Upper triangular 49*> = 'L': Lower triangular 50*> \endverbatim 51*> 52*> \param[in] TRANS 53*> \verbatim 54*> TRANS is CHARACTER*1 55*> Specifies the operation applied to A. 56*> = 'N': A *x = s*b (No transpose) 57*> = 'T': A**T *x = s*b (Transpose) 58*> = 'C': A**H *x = s*b (Conjugate transpose) 59*> \endverbatim 60*> 61*> \param[in] DIAG 62*> \verbatim 63*> DIAG is CHARACTER*1 64*> Specifies whether or not the matrix A is unit triangular. 65*> = 'N': Non-unit triangular 66*> = 'U': Unit triangular 67*> \endverbatim 68*> 69*> \param[in] N 70*> \verbatim 71*> N is INTEGER 72*> The order of the matrix A. N >= 0. 73*> \endverbatim 74*> 75*> \param[in] NRHS 76*> \verbatim 77*> NRHS is INTEGER 78*> The number of right hand sides, i.e., the number of columns 79*> of the matrices X and B. NRHS >= 0. 80*> \endverbatim 81*> 82*> \param[in] A 83*> \verbatim 84*> A is COMPLEX array, dimension (LDA,N) 85*> The triangular matrix A. If UPLO = 'U', the leading n by n 86*> upper triangular part of the array A contains the upper 87*> triangular matrix, and the strictly lower triangular part of 88*> A is not referenced. If UPLO = 'L', the leading n by n lower 89*> triangular part of the array A contains the lower triangular 90*> matrix, and the strictly upper triangular part of A is not 91*> referenced. If DIAG = 'U', the diagonal elements of A are 92*> also not referenced and are assumed to be 1. 93*> \endverbatim 94*> 95*> \param[in] LDA 96*> \verbatim 97*> LDA is INTEGER 98*> The leading dimension of the array A. LDA >= max(1,N). 99*> \endverbatim 100*> 101*> \param[in] SCALE 102*> \verbatim 103*> SCALE is REAL 104*> The scaling factor s used in solving the triangular system. 105*> \endverbatim 106*> 107*> \param[in] CNORM 108*> \verbatim 109*> CNORM is REAL array, dimension (N) 110*> The 1-norms of the columns of A, not counting the diagonal. 111*> \endverbatim 112*> 113*> \param[in] TSCAL 114*> \verbatim 115*> TSCAL is REAL 116*> The scaling factor used in computing the 1-norms in CNORM. 117*> CNORM actually contains the column norms of TSCAL*A. 118*> \endverbatim 119*> 120*> \param[in] X 121*> \verbatim 122*> X is COMPLEX array, dimension (LDX,NRHS) 123*> The computed solution vectors for the system of linear 124*> equations. 125*> \endverbatim 126*> 127*> \param[in] LDX 128*> \verbatim 129*> LDX is INTEGER 130*> The leading dimension of the array X. LDX >= max(1,N). 131*> \endverbatim 132*> 133*> \param[in] B 134*> \verbatim 135*> B is COMPLEX array, dimension (LDB,NRHS) 136*> The right hand side vectors for the system of linear 137*> equations. 138*> \endverbatim 139*> 140*> \param[in] LDB 141*> \verbatim 142*> LDB is INTEGER 143*> The leading dimension of the array B. LDB >= max(1,N). 144*> \endverbatim 145*> 146*> \param[out] WORK 147*> \verbatim 148*> WORK is COMPLEX array, dimension (N) 149*> \endverbatim 150*> 151*> \param[out] RESID 152*> \verbatim 153*> RESID is REAL 154*> The maximum over the number of right hand sides of 155*> norm(op(A)*x - s*b) / ( norm(op(A)) * norm(x) * EPS ). 156*> \endverbatim 157* 158* Authors: 159* ======== 160* 161*> \author Univ. of Tennessee 162*> \author Univ. of California Berkeley 163*> \author Univ. of Colorado Denver 164*> \author NAG Ltd. 165* 166*> \ingroup complex_lin 167* 168* ===================================================================== 169 SUBROUTINE CTRT03( UPLO, TRANS, DIAG, N, NRHS, A, LDA, SCALE, 170 $ CNORM, TSCAL, X, LDX, B, LDB, WORK, RESID ) 171* 172* -- LAPACK test routine -- 173* -- LAPACK is a software package provided by Univ. of Tennessee, -- 174* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 175* 176* .. Scalar Arguments .. 177 CHARACTER DIAG, TRANS, UPLO 178 INTEGER LDA, LDB, LDX, N, NRHS 179 REAL RESID, SCALE, TSCAL 180* .. 181* .. Array Arguments .. 182 REAL CNORM( * ) 183 COMPLEX A( LDA, * ), B( LDB, * ), WORK( * ), 184 $ X( LDX, * ) 185* .. 186* 187* ===================================================================== 188* 189* .. Parameters .. 190 REAL ONE, ZERO 191 PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) 192* .. 193* .. Local Scalars .. 194 INTEGER IX, J 195 REAL EPS, ERR, SMLNUM, TNORM, XNORM, XSCAL 196* .. 197* .. External Functions .. 198 LOGICAL LSAME 199 INTEGER ICAMAX 200 REAL SLAMCH 201 EXTERNAL LSAME, ICAMAX, SLAMCH 202* .. 203* .. External Subroutines .. 204 EXTERNAL CAXPY, CCOPY, CSSCAL, CTRMV 205* .. 206* .. Intrinsic Functions .. 207 INTRINSIC ABS, CMPLX, MAX, REAL 208* .. 209* .. Executable Statements .. 210* 211* Quick exit if N = 0 212* 213 IF( N.LE.0 .OR. NRHS.LE.0 ) THEN 214 RESID = ZERO 215 RETURN 216 END IF 217 EPS = SLAMCH( 'Epsilon' ) 218 SMLNUM = SLAMCH( 'Safe minimum' ) 219* 220* Compute the norm of the triangular matrix A using the column 221* norms already computed by CLATRS. 222* 223 TNORM = ZERO 224 IF( LSAME( DIAG, 'N' ) ) THEN 225 DO 10 J = 1, N 226 TNORM = MAX( TNORM, TSCAL*ABS( A( J, J ) )+CNORM( J ) ) 227 10 CONTINUE 228 ELSE 229 DO 20 J = 1, N 230 TNORM = MAX( TNORM, TSCAL+CNORM( J ) ) 231 20 CONTINUE 232 END IF 233* 234* Compute the maximum over the number of right hand sides of 235* norm(op(A)*x - s*b) / ( norm(op(A)) * norm(x) * EPS ). 236* 237 RESID = ZERO 238 DO 30 J = 1, NRHS 239 CALL CCOPY( N, X( 1, J ), 1, WORK, 1 ) 240 IX = ICAMAX( N, WORK, 1 ) 241 XNORM = MAX( ONE, ABS( X( IX, J ) ) ) 242 XSCAL = ( ONE / XNORM ) / REAL( N ) 243 CALL CSSCAL( N, XSCAL, WORK, 1 ) 244 CALL CTRMV( UPLO, TRANS, DIAG, N, A, LDA, WORK, 1 ) 245 CALL CAXPY( N, CMPLX( -SCALE*XSCAL ), B( 1, J ), 1, WORK, 1 ) 246 IX = ICAMAX( N, WORK, 1 ) 247 ERR = TSCAL*ABS( WORK( IX ) ) 248 IX = ICAMAX( N, X( 1, J ), 1 ) 249 XNORM = ABS( X( IX, J ) ) 250 IF( ERR*SMLNUM.LE.XNORM ) THEN 251 IF( XNORM.GT.ZERO ) 252 $ ERR = ERR / XNORM 253 ELSE 254 IF( ERR.GT.ZERO ) 255 $ ERR = ONE / EPS 256 END IF 257 IF( ERR*SMLNUM.LE.TNORM ) THEN 258 IF( TNORM.GT.ZERO ) 259 $ ERR = ERR / TNORM 260 ELSE 261 IF( ERR.GT.ZERO ) 262 $ ERR = ONE / EPS 263 END IF 264 RESID = MAX( RESID, ERR ) 265 30 CONTINUE 266* 267 RETURN 268* 269* End of CTRT03 270* 271 END 272