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