|
LAPACK 3.3.0
|
|
LAPACK 3.3.0
|
00001 SUBROUTINE CCHKTB( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, 00002 $ NMAX, AB, AINV, B, X, XACT, WORK, RWORK, NOUT ) 00003 * 00004 * -- LAPACK test routine (version 3.1) -- 00005 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 LOGICAL TSTERR 00010 INTEGER NMAX, NN, NNS, NOUT 00011 REAL THRESH 00012 * .. 00013 * .. Array Arguments .. 00014 LOGICAL DOTYPE( * ) 00015 INTEGER NSVAL( * ), NVAL( * ) 00016 REAL RWORK( * ) 00017 COMPLEX AB( * ), AINV( * ), B( * ), WORK( * ), X( * ), 00018 $ XACT( * ) 00019 * .. 00020 * 00021 * Purpose 00022 * ======= 00023 * 00024 * CCHKTB tests CTBTRS, -RFS, and -CON, and CLATBS. 00025 * 00026 * Arguments 00027 * ========= 00028 * 00029 * DOTYPE (input) LOGICAL array, dimension (NTYPES) 00030 * The matrix types to be used for testing. Matrices of type j 00031 * (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) = 00032 * .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used. 00033 * 00034 * NN (input) INTEGER 00035 * The number of values of N contained in the vector NVAL. 00036 * 00037 * NVAL (input) INTEGER array, dimension (NN) 00038 * The values of the matrix column dimension N. 00039 * 00040 * NNS (input) INTEGER 00041 * The number of values of NRHS contained in the vector NSVAL. 00042 * 00043 * NSVAL (input) INTEGER array, dimension (NNS) 00044 * The values of the number of right hand sides NRHS. 00045 * 00046 * THRESH (input) REAL 00047 * The threshold value for the test ratios. A result is 00048 * included in the output file if RESULT >= THRESH. To have 00049 * every test ratio printed, use THRESH = 0. 00050 * 00051 * TSTERR (input) LOGICAL 00052 * Flag that indicates whether error exits are to be tested. 00053 * 00054 * NMAX (input) INTEGER 00055 * The leading dimension of the work arrays. 00056 * NMAX >= the maximum value of N in NVAL. 00057 * 00058 * AB (workspace) COMPLEX array, dimension (NMAX*NMAX) 00059 * 00060 * AINV (workspace) COMPLEX array, dimension (NMAX*NMAX) 00061 * 00062 * B (workspace) COMPLEX array, dimension (NMAX*NSMAX) 00063 * where NSMAX is the largest entry in NSVAL. 00064 * 00065 * X (workspace) COMPLEX array, dimension (NMAX*NSMAX) 00066 * 00067 * XACT (workspace) COMPLEX array, dimension (NMAX*NSMAX) 00068 * 00069 * WORK (workspace) COMPLEX array, dimension 00070 * (NMAX*max(3,NSMAX)) 00071 * 00072 * RWORK (workspace) REAL array, dimension 00073 * (max(NMAX,2*NSMAX)) 00074 * 00075 * NOUT (input) INTEGER 00076 * The unit number for output. 00077 * 00078 * ===================================================================== 00079 * 00080 * .. Parameters .. 00081 INTEGER NTYPE1, NTYPES 00082 PARAMETER ( NTYPE1 = 9, NTYPES = 17 ) 00083 INTEGER NTESTS 00084 PARAMETER ( NTESTS = 8 ) 00085 INTEGER NTRAN 00086 PARAMETER ( NTRAN = 3 ) 00087 REAL ONE, ZERO 00088 PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) 00089 * .. 00090 * .. Local Scalars .. 00091 CHARACTER DIAG, NORM, TRANS, UPLO, XTYPE 00092 CHARACTER*3 PATH 00093 INTEGER I, IDIAG, IK, IMAT, IN, INFO, IRHS, ITRAN, 00094 $ IUPLO, J, K, KD, LDA, LDAB, N, NERRS, NFAIL, 00095 $ NIMAT, NIMAT2, NK, NRHS, NRUN 00096 REAL AINVNM, ANORM, RCOND, RCONDC, RCONDI, RCONDO, 00097 $ SCALE 00098 * .. 00099 * .. Local Arrays .. 00100 CHARACTER TRANSS( NTRAN ), UPLOS( 2 ) 00101 INTEGER ISEED( 4 ), ISEEDY( 4 ) 00102 REAL RESULT( NTESTS ) 00103 * .. 00104 * .. External Functions .. 00105 LOGICAL LSAME 00106 REAL CLANTB, CLANTR 00107 EXTERNAL LSAME, CLANTB, CLANTR 00108 * .. 00109 * .. External Subroutines .. 00110 EXTERNAL ALAERH, ALAHD, ALASUM, CCOPY, CERRTR, CGET04, 00111 $ CLACPY, CLARHS, CLASET, CLATBS, CLATTB, CTBCON, 00112 $ CTBRFS, CTBSV, CTBT02, CTBT03, CTBT05, CTBT06, 00113 $ CTBTRS 00114 * .. 00115 * .. Scalars in Common .. 00116 LOGICAL LERR, OK 00117 CHARACTER*32 SRNAMT 00118 INTEGER INFOT, IOUNIT 00119 * .. 00120 * .. Common blocks .. 00121 COMMON / INFOC / INFOT, IOUNIT, OK, LERR 00122 COMMON / SRNAMC / SRNAMT 00123 * .. 00124 * .. Intrinsic Functions .. 00125 INTRINSIC CMPLX, MAX, MIN 00126 * .. 00127 * .. Data statements .. 00128 DATA ISEEDY / 1988, 1989, 1990, 1991 / 00129 DATA UPLOS / 'U', 'L' / , TRANSS / 'N', 'T', 'C' / 00130 * .. 00131 * .. Executable Statements .. 00132 * 00133 * Initialize constants and the random number seed. 00134 * 00135 PATH( 1: 1 ) = 'Complex precision' 00136 PATH( 2: 3 ) = 'TB' 00137 NRUN = 0 00138 NFAIL = 0 00139 NERRS = 0 00140 DO 10 I = 1, 4 00141 ISEED( I ) = ISEEDY( I ) 00142 10 CONTINUE 00143 * 00144 * Test the error exits 00145 * 00146 IF( TSTERR ) 00147 $ CALL CERRTR( PATH, NOUT ) 00148 INFOT = 0 00149 * 00150 DO 140 IN = 1, NN 00151 * 00152 * Do for each value of N in NVAL 00153 * 00154 N = NVAL( IN ) 00155 LDA = MAX( 1, N ) 00156 XTYPE = 'N' 00157 NIMAT = NTYPE1 00158 NIMAT2 = NTYPES 00159 IF( N.LE.0 ) THEN 00160 NIMAT = 1 00161 NIMAT2 = NTYPE1 + 1 00162 END IF 00163 * 00164 NK = MIN( N+1, 4 ) 00165 DO 130 IK = 1, NK 00166 * 00167 * Do for KD = 0, N, (3N-1)/4, and (N+1)/4. This order makes 00168 * it easier to skip redundant values for small values of N. 00169 * 00170 IF( IK.EQ.1 ) THEN 00171 KD = 0 00172 ELSE IF( IK.EQ.2 ) THEN 00173 KD = MAX( N, 0 ) 00174 ELSE IF( IK.EQ.3 ) THEN 00175 KD = ( 3*N-1 ) / 4 00176 ELSE IF( IK.EQ.4 ) THEN 00177 KD = ( N+1 ) / 4 00178 END IF 00179 LDAB = KD + 1 00180 * 00181 DO 90 IMAT = 1, NIMAT 00182 * 00183 * Do the tests only if DOTYPE( IMAT ) is true. 00184 * 00185 IF( .NOT.DOTYPE( IMAT ) ) 00186 $ GO TO 90 00187 * 00188 DO 80 IUPLO = 1, 2 00189 * 00190 * Do first for UPLO = 'U', then for UPLO = 'L' 00191 * 00192 UPLO = UPLOS( IUPLO ) 00193 * 00194 * Call CLATTB to generate a triangular test matrix. 00195 * 00196 SRNAMT = 'CLATTB' 00197 CALL CLATTB( IMAT, UPLO, 'No transpose', DIAG, ISEED, 00198 $ N, KD, AB, LDAB, X, WORK, RWORK, INFO ) 00199 * 00200 * Set IDIAG = 1 for non-unit matrices, 2 for unit. 00201 * 00202 IF( LSAME( DIAG, 'N' ) ) THEN 00203 IDIAG = 1 00204 ELSE 00205 IDIAG = 2 00206 END IF 00207 * 00208 * Form the inverse of A so we can get a good estimate 00209 * of RCONDC = 1/(norm(A) * norm(inv(A))). 00210 * 00211 CALL CLASET( 'Full', N, N, CMPLX( ZERO ), 00212 $ CMPLX( ONE ), AINV, LDA ) 00213 IF( LSAME( UPLO, 'U' ) ) THEN 00214 DO 20 J = 1, N 00215 CALL CTBSV( UPLO, 'No transpose', DIAG, J, KD, 00216 $ AB, LDAB, AINV( ( J-1 )*LDA+1 ), 1 ) 00217 20 CONTINUE 00218 ELSE 00219 DO 30 J = 1, N 00220 CALL CTBSV( UPLO, 'No transpose', DIAG, N-J+1, 00221 $ KD, AB( ( J-1 )*LDAB+1 ), LDAB, 00222 $ AINV( ( J-1 )*LDA+J ), 1 ) 00223 30 CONTINUE 00224 END IF 00225 * 00226 * Compute the 1-norm condition number of A. 00227 * 00228 ANORM = CLANTB( '1', UPLO, DIAG, N, KD, AB, LDAB, 00229 $ RWORK ) 00230 AINVNM = CLANTR( '1', UPLO, DIAG, N, N, AINV, LDA, 00231 $ RWORK ) 00232 IF( ANORM.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN 00233 RCONDO = ONE 00234 ELSE 00235 RCONDO = ( ONE / ANORM ) / AINVNM 00236 END IF 00237 * 00238 * Compute the infinity-norm condition number of A. 00239 * 00240 ANORM = CLANTB( 'I', UPLO, DIAG, N, KD, AB, LDAB, 00241 $ RWORK ) 00242 AINVNM = CLANTR( 'I', UPLO, DIAG, N, N, AINV, LDA, 00243 $ RWORK ) 00244 IF( ANORM.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN 00245 RCONDI = ONE 00246 ELSE 00247 RCONDI = ( ONE / ANORM ) / AINVNM 00248 END IF 00249 * 00250 DO 60 IRHS = 1, NNS 00251 NRHS = NSVAL( IRHS ) 00252 XTYPE = 'N' 00253 * 00254 DO 50 ITRAN = 1, NTRAN 00255 * 00256 * Do for op(A) = A, A**T, or A**H. 00257 * 00258 TRANS = TRANSS( ITRAN ) 00259 IF( ITRAN.EQ.1 ) THEN 00260 NORM = 'O' 00261 RCONDC = RCONDO 00262 ELSE 00263 NORM = 'I' 00264 RCONDC = RCONDI 00265 END IF 00266 * 00267 *+ TEST 1 00268 * Solve and compute residual for op(A)*x = b. 00269 * 00270 SRNAMT = 'CLARHS' 00271 CALL CLARHS( PATH, XTYPE, UPLO, TRANS, N, N, KD, 00272 $ IDIAG, NRHS, AB, LDAB, XACT, LDA, 00273 $ B, LDA, ISEED, INFO ) 00274 XTYPE = 'C' 00275 CALL CLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 00276 * 00277 SRNAMT = 'CTBTRS' 00278 CALL CTBTRS( UPLO, TRANS, DIAG, N, KD, NRHS, AB, 00279 $ LDAB, X, LDA, INFO ) 00280 * 00281 * Check error code from CTBTRS. 00282 * 00283 IF( INFO.NE.0 ) 00284 $ CALL ALAERH( PATH, 'CTBTRS', INFO, 0, 00285 $ UPLO // TRANS // DIAG, N, N, KD, 00286 $ KD, NRHS, IMAT, NFAIL, NERRS, 00287 $ NOUT ) 00288 * 00289 CALL CTBT02( UPLO, TRANS, DIAG, N, KD, NRHS, AB, 00290 $ LDAB, X, LDA, B, LDA, WORK, RWORK, 00291 $ RESULT( 1 ) ) 00292 * 00293 *+ TEST 2 00294 * Check solution from generated exact solution. 00295 * 00296 CALL CGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00297 $ RESULT( 2 ) ) 00298 * 00299 *+ TESTS 3, 4, and 5 00300 * Use iterative refinement to improve the solution 00301 * and compute error bounds. 00302 * 00303 SRNAMT = 'CTBRFS' 00304 CALL CTBRFS( UPLO, TRANS, DIAG, N, KD, NRHS, AB, 00305 $ LDAB, B, LDA, X, LDA, RWORK, 00306 $ RWORK( NRHS+1 ), WORK, 00307 $ RWORK( 2*NRHS+1 ), INFO ) 00308 * 00309 * Check error code from CTBRFS. 00310 * 00311 IF( INFO.NE.0 ) 00312 $ CALL ALAERH( PATH, 'CTBRFS', INFO, 0, 00313 $ UPLO // TRANS // DIAG, N, N, KD, 00314 $ KD, NRHS, IMAT, NFAIL, NERRS, 00315 $ NOUT ) 00316 * 00317 CALL CGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00318 $ RESULT( 3 ) ) 00319 CALL CTBT05( UPLO, TRANS, DIAG, N, KD, NRHS, AB, 00320 $ LDAB, B, LDA, X, LDA, XACT, LDA, 00321 $ RWORK, RWORK( NRHS+1 ), 00322 $ RESULT( 4 ) ) 00323 * 00324 * Print information about the tests that did not 00325 * pass the threshold. 00326 * 00327 DO 40 K = 1, 5 00328 IF( RESULT( K ).GE.THRESH ) THEN 00329 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00330 $ CALL ALAHD( NOUT, PATH ) 00331 WRITE( NOUT, FMT = 9999 )UPLO, TRANS, 00332 $ DIAG, N, KD, NRHS, IMAT, K, RESULT( K ) 00333 NFAIL = NFAIL + 1 00334 END IF 00335 40 CONTINUE 00336 NRUN = NRUN + 5 00337 50 CONTINUE 00338 60 CONTINUE 00339 * 00340 *+ TEST 6 00341 * Get an estimate of RCOND = 1/CNDNUM. 00342 * 00343 DO 70 ITRAN = 1, 2 00344 IF( ITRAN.EQ.1 ) THEN 00345 NORM = 'O' 00346 RCONDC = RCONDO 00347 ELSE 00348 NORM = 'I' 00349 RCONDC = RCONDI 00350 END IF 00351 SRNAMT = 'CTBCON' 00352 CALL CTBCON( NORM, UPLO, DIAG, N, KD, AB, LDAB, 00353 $ RCOND, WORK, RWORK, INFO ) 00354 * 00355 * Check error code from CTBCON. 00356 * 00357 IF( INFO.NE.0 ) 00358 $ CALL ALAERH( PATH, 'CTBCON', INFO, 0, 00359 $ NORM // UPLO // DIAG, N, N, KD, KD, 00360 $ -1, IMAT, NFAIL, NERRS, NOUT ) 00361 * 00362 CALL CTBT06( RCOND, RCONDC, UPLO, DIAG, N, KD, AB, 00363 $ LDAB, RWORK, RESULT( 6 ) ) 00364 * 00365 * Print the test ratio if it is .GE. THRESH. 00366 * 00367 IF( RESULT( 6 ).GE.THRESH ) THEN 00368 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00369 $ CALL ALAHD( NOUT, PATH ) 00370 WRITE( NOUT, FMT = 9998 ) 'CTBCON', NORM, UPLO, 00371 $ DIAG, N, KD, IMAT, 6, RESULT( 6 ) 00372 NFAIL = NFAIL + 1 00373 END IF 00374 NRUN = NRUN + 1 00375 70 CONTINUE 00376 80 CONTINUE 00377 90 CONTINUE 00378 * 00379 * Use pathological test matrices to test CLATBS. 00380 * 00381 DO 120 IMAT = NTYPE1 + 1, NIMAT2 00382 * 00383 * Do the tests only if DOTYPE( IMAT ) is true. 00384 * 00385 IF( .NOT.DOTYPE( IMAT ) ) 00386 $ GO TO 120 00387 * 00388 DO 110 IUPLO = 1, 2 00389 * 00390 * Do first for UPLO = 'U', then for UPLO = 'L' 00391 * 00392 UPLO = UPLOS( IUPLO ) 00393 DO 100 ITRAN = 1, NTRAN 00394 * 00395 * Do for op(A) = A, A**T, and A**H. 00396 * 00397 TRANS = TRANSS( ITRAN ) 00398 * 00399 * Call CLATTB to generate a triangular test matrix. 00400 * 00401 SRNAMT = 'CLATTB' 00402 CALL CLATTB( IMAT, UPLO, TRANS, DIAG, ISEED, N, KD, 00403 $ AB, LDAB, X, WORK, RWORK, INFO ) 00404 * 00405 *+ TEST 7 00406 * Solve the system op(A)*x = b 00407 * 00408 SRNAMT = 'CLATBS' 00409 CALL CCOPY( N, X, 1, B, 1 ) 00410 CALL CLATBS( UPLO, TRANS, DIAG, 'N', N, KD, AB, 00411 $ LDAB, B, SCALE, RWORK, INFO ) 00412 * 00413 * Check error code from CLATBS. 00414 * 00415 IF( INFO.NE.0 ) 00416 $ CALL ALAERH( PATH, 'CLATBS', INFO, 0, 00417 $ UPLO // TRANS // DIAG // 'N', N, N, 00418 $ KD, KD, -1, IMAT, NFAIL, NERRS, 00419 $ NOUT ) 00420 * 00421 CALL CTBT03( UPLO, TRANS, DIAG, N, KD, 1, AB, LDAB, 00422 $ SCALE, RWORK, ONE, B, LDA, X, LDA, 00423 $ WORK, RESULT( 7 ) ) 00424 * 00425 *+ TEST 8 00426 * Solve op(A)*x = b again with NORMIN = 'Y'. 00427 * 00428 CALL CCOPY( N, X, 1, B, 1 ) 00429 CALL CLATBS( UPLO, TRANS, DIAG, 'Y', N, KD, AB, 00430 $ LDAB, B, SCALE, RWORK, INFO ) 00431 * 00432 * Check error code from CLATBS. 00433 * 00434 IF( INFO.NE.0 ) 00435 $ CALL ALAERH( PATH, 'CLATBS', INFO, 0, 00436 $ UPLO // TRANS // DIAG // 'Y', N, N, 00437 $ KD, KD, -1, IMAT, NFAIL, NERRS, 00438 $ NOUT ) 00439 * 00440 CALL CTBT03( UPLO, TRANS, DIAG, N, KD, 1, AB, LDAB, 00441 $ SCALE, RWORK, ONE, B, LDA, X, LDA, 00442 $ WORK, RESULT( 8 ) ) 00443 * 00444 * Print information about the tests that did not pass 00445 * the threshold. 00446 * 00447 IF( RESULT( 7 ).GE.THRESH ) THEN 00448 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00449 $ CALL ALAHD( NOUT, PATH ) 00450 WRITE( NOUT, FMT = 9997 )'CLATBS', UPLO, TRANS, 00451 $ DIAG, 'N', N, KD, IMAT, 7, RESULT( 7 ) 00452 NFAIL = NFAIL + 1 00453 END IF 00454 IF( RESULT( 8 ).GE.THRESH ) THEN 00455 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00456 $ CALL ALAHD( NOUT, PATH ) 00457 WRITE( NOUT, FMT = 9997 )'CLATBS', UPLO, TRANS, 00458 $ DIAG, 'Y', N, KD, IMAT, 8, RESULT( 8 ) 00459 NFAIL = NFAIL + 1 00460 END IF 00461 NRUN = NRUN + 2 00462 100 CONTINUE 00463 110 CONTINUE 00464 120 CONTINUE 00465 130 CONTINUE 00466 140 CONTINUE 00467 * 00468 * Print a summary of the results. 00469 * 00470 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS ) 00471 * 00472 9999 FORMAT( ' UPLO=''', A1, ''', TRANS=''', A1, '' 00473 ', $ DIAG=''', A1, ''', N=', I5, ', KD=', I5, ', NRHS=', I5, 00474 $ ', type ', I2, ', test(', I2, ')=', G12.5 ) 00475 9998 FORMAT( 1X, A, '( ''', A1, ''', ''', A1, ''', ''', A1, ''',', 00476 $ I5, ',', I5, ', ... ), type ', I2, ', test(', I2, ')=', 00477 $ G12.5 ) 00478 9997 FORMAT( 1X, A, '( ''', A1, ''', ''', A1, ''', ''', A1, ''', ''', 00479 $ A1, ''',', I5, ',', I5, ', ... ), type ', I2, ', test(', 00480 $ I1, ')=', G12.5 ) 00481 RETURN 00482 * 00483 * End of CCHKTB 00484 * 00485 END
1.7.3 
