1*> \brief \b DCHKTR 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 DCHKTR( DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, 12* THRESH, TSTERR, NMAX, A, AINV, B, X, XACT, 13* WORK, RWORK, IWORK, NOUT ) 14* 15* .. Scalar Arguments .. 16* LOGICAL TSTERR 17* INTEGER NMAX, NN, NNB, NNS, NOUT 18* DOUBLE PRECISION THRESH 19* .. 20* .. Array Arguments .. 21* LOGICAL DOTYPE( * ) 22* INTEGER IWORK( * ), NBVAL( * ), NSVAL( * ), NVAL( * ) 23* DOUBLE PRECISION A( * ), AINV( * ), B( * ), RWORK( * ), 24* $ WORK( * ), X( * ), XACT( * ) 25* .. 26* 27* 28*> \par Purpose: 29* ============= 30*> 31*> \verbatim 32*> 33*> DCHKTR tests DTRTRI, -TRS, -RFS, and -CON, and DLATRS 34*> \endverbatim 35* 36* Arguments: 37* ========== 38* 39*> \param[in] DOTYPE 40*> \verbatim 41*> DOTYPE is LOGICAL array, dimension (NTYPES) 42*> The matrix types to be used for testing. Matrices of type j 43*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) = 44*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used. 45*> \endverbatim 46*> 47*> \param[in] NN 48*> \verbatim 49*> NN is INTEGER 50*> The number of values of N contained in the vector NVAL. 51*> \endverbatim 52*> 53*> \param[in] NVAL 54*> \verbatim 55*> NVAL is INTEGER array, dimension (NN) 56*> The values of the matrix column dimension N. 57*> \endverbatim 58*> 59*> \param[in] NNB 60*> \verbatim 61*> NNB is INTEGER 62*> The number of values of NB contained in the vector NBVAL. 63*> \endverbatim 64*> 65*> \param[in] NBVAL 66*> \verbatim 67*> NBVAL is INTEGER array, dimension (NNB) 68*> The values of the blocksize NB. 69*> \endverbatim 70*> 71*> \param[in] NNS 72*> \verbatim 73*> NNS is INTEGER 74*> The number of values of NRHS contained in the vector NSVAL. 75*> \endverbatim 76*> 77*> \param[in] NSVAL 78*> \verbatim 79*> NSVAL is INTEGER array, dimension (NNS) 80*> The values of the number of right hand sides NRHS. 81*> \endverbatim 82*> 83*> \param[in] THRESH 84*> \verbatim 85*> THRESH is DOUBLE PRECISION 86*> The threshold value for the test ratios. A result is 87*> included in the output file if RESULT >= THRESH. To have 88*> every test ratio printed, use THRESH = 0. 89*> \endverbatim 90*> 91*> \param[in] TSTERR 92*> \verbatim 93*> TSTERR is LOGICAL 94*> Flag that indicates whether error exits are to be tested. 95*> \endverbatim 96*> 97*> \param[in] NMAX 98*> \verbatim 99*> NMAX is INTEGER 100*> The leading dimension of the work arrays. 101*> NMAX >= the maximum value of N in NVAL. 102*> \endverbatim 103*> 104*> \param[out] A 105*> \verbatim 106*> A is DOUBLE PRECISION array, dimension (NMAX*NMAX) 107*> \endverbatim 108*> 109*> \param[out] AINV 110*> \verbatim 111*> AINV is DOUBLE PRECISION array, dimension (NMAX*NMAX) 112*> \endverbatim 113*> 114*> \param[out] B 115*> \verbatim 116*> B is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 117*> where NSMAX is the largest entry in NSVAL. 118*> \endverbatim 119*> 120*> \param[out] X 121*> \verbatim 122*> X is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 123*> \endverbatim 124*> 125*> \param[out] XACT 126*> \verbatim 127*> XACT is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 128*> \endverbatim 129*> 130*> \param[out] WORK 131*> \verbatim 132*> WORK is DOUBLE PRECISION array, dimension 133*> (NMAX*max(3,NSMAX)) 134*> \endverbatim 135*> 136*> \param[out] RWORK 137*> \verbatim 138*> RWORK is DOUBLE PRECISION array, dimension 139*> (max(NMAX,2*NSMAX)) 140*> \endverbatim 141*> 142*> \param[out] IWORK 143*> \verbatim 144*> IWORK is INTEGER array, dimension (NMAX) 145*> \endverbatim 146*> 147*> \param[in] NOUT 148*> \verbatim 149*> NOUT is INTEGER 150*> The unit number for output. 151*> \endverbatim 152* 153* Authors: 154* ======== 155* 156*> \author Univ. of Tennessee 157*> \author Univ. of California Berkeley 158*> \author Univ. of Colorado Denver 159*> \author NAG Ltd. 160* 161*> \date November 2011 162* 163*> \ingroup double_lin 164* 165* ===================================================================== 166 SUBROUTINE DCHKTR( DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, 167 $ THRESH, TSTERR, NMAX, A, AINV, B, X, XACT, 168 $ WORK, RWORK, IWORK, NOUT ) 169* 170* -- LAPACK test routine (version 3.4.0) -- 171* -- LAPACK is a software package provided by Univ. of Tennessee, -- 172* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 173* November 2011 174* 175* .. Scalar Arguments .. 176 LOGICAL TSTERR 177 INTEGER NMAX, NN, NNB, NNS, NOUT 178 DOUBLE PRECISION THRESH 179* .. 180* .. Array Arguments .. 181 LOGICAL DOTYPE( * ) 182 INTEGER IWORK( * ), NBVAL( * ), NSVAL( * ), NVAL( * ) 183 DOUBLE PRECISION A( * ), AINV( * ), B( * ), RWORK( * ), 184 $ WORK( * ), X( * ), XACT( * ) 185* .. 186* 187* ===================================================================== 188* 189* .. Parameters .. 190 INTEGER NTYPE1, NTYPES 191 PARAMETER ( NTYPE1 = 10, NTYPES = 18 ) 192 INTEGER NTESTS 193 PARAMETER ( NTESTS = 9 ) 194 INTEGER NTRAN 195 PARAMETER ( NTRAN = 3 ) 196 DOUBLE PRECISION ONE, ZERO 197 PARAMETER ( ONE = 1.0D0, ZERO = 0.0D0 ) 198* .. 199* .. Local Scalars .. 200 CHARACTER DIAG, NORM, TRANS, UPLO, XTYPE 201 CHARACTER*3 PATH 202 INTEGER I, IDIAG, IMAT, IN, INB, INFO, IRHS, ITRAN, 203 $ IUPLO, K, LDA, N, NB, NERRS, NFAIL, NRHS, NRUN 204 DOUBLE PRECISION AINVNM, ANORM, DUMMY, RCOND, RCONDC, RCONDI, 205 $ RCONDO, SCALE 206* .. 207* .. Local Arrays .. 208 CHARACTER TRANSS( NTRAN ), UPLOS( 2 ) 209 INTEGER ISEED( 4 ), ISEEDY( 4 ) 210 DOUBLE PRECISION RESULT( NTESTS ) 211* .. 212* .. External Functions .. 213 LOGICAL LSAME 214 DOUBLE PRECISION DLANTR 215 EXTERNAL LSAME, DLANTR 216* .. 217* .. External Subroutines .. 218 EXTERNAL ALAERH, ALAHD, ALASUM, DCOPY, DERRTR, DGET04, 219 $ DLACPY, DLARHS, DLATRS, DLATTR, DTRCON, DTRRFS, 220 $ DTRT01, DTRT02, DTRT03, DTRT05, DTRT06, DTRTRI, 221 $ DTRTRS, XLAENV 222* .. 223* .. Scalars in Common .. 224 LOGICAL LERR, OK 225 CHARACTER*32 SRNAMT 226 INTEGER INFOT, IOUNIT 227* .. 228* .. Common blocks .. 229 COMMON / INFOC / INFOT, IOUNIT, OK, LERR 230 COMMON / SRNAMC / SRNAMT 231* .. 232* .. Intrinsic Functions .. 233 INTRINSIC MAX 234* .. 235* .. Data statements .. 236 DATA ISEEDY / 1988, 1989, 1990, 1991 / 237 DATA UPLOS / 'U', 'L' / , TRANSS / 'N', 'T', 'C' / 238* .. 239* .. Executable Statements .. 240* 241* Initialize constants and the random number seed. 242* 243 PATH( 1: 1 ) = 'Double precision' 244 PATH( 2: 3 ) = 'TR' 245 NRUN = 0 246 NFAIL = 0 247 NERRS = 0 248 DO 10 I = 1, 4 249 ISEED( I ) = ISEEDY( I ) 250 10 CONTINUE 251* 252* Test the error exits 253* 254 IF( TSTERR ) 255 $ CALL DERRTR( PATH, NOUT ) 256 INFOT = 0 257 CALL XLAENV( 2, 2 ) 258* 259 DO 120 IN = 1, NN 260* 261* Do for each value of N in NVAL 262* 263 N = NVAL( IN ) 264 LDA = MAX( 1, N ) 265 XTYPE = 'N' 266* 267 DO 80 IMAT = 1, NTYPE1 268* 269* Do the tests only if DOTYPE( IMAT ) is true. 270* 271 IF( .NOT.DOTYPE( IMAT ) ) 272 $ GO TO 80 273* 274 DO 70 IUPLO = 1, 2 275* 276* Do first for UPLO = 'U', then for UPLO = 'L' 277* 278 UPLO = UPLOS( IUPLO ) 279* 280* Call DLATTR to generate a triangular test matrix. 281* 282 SRNAMT = 'DLATTR' 283 CALL DLATTR( IMAT, UPLO, 'No transpose', DIAG, ISEED, N, 284 $ A, LDA, X, WORK, INFO ) 285* 286* Set IDIAG = 1 for non-unit matrices, 2 for unit. 287* 288 IF( LSAME( DIAG, 'N' ) ) THEN 289 IDIAG = 1 290 ELSE 291 IDIAG = 2 292 END IF 293* 294 DO 60 INB = 1, NNB 295* 296* Do for each blocksize in NBVAL 297* 298 NB = NBVAL( INB ) 299 CALL XLAENV( 1, NB ) 300* 301*+ TEST 1 302* Form the inverse of A. 303* 304 CALL DLACPY( UPLO, N, N, A, LDA, AINV, LDA ) 305 SRNAMT = 'DTRTRI' 306 CALL DTRTRI( UPLO, DIAG, N, AINV, LDA, INFO ) 307* 308* Check error code from DTRTRI. 309* 310 IF( INFO.NE.0 ) 311 $ CALL ALAERH( PATH, 'DTRTRI', INFO, 0, UPLO // DIAG, 312 $ N, N, -1, -1, NB, IMAT, NFAIL, NERRS, 313 $ NOUT ) 314* 315* Compute the infinity-norm condition number of A. 316* 317 ANORM = DLANTR( 'I', UPLO, DIAG, N, N, A, LDA, RWORK ) 318 AINVNM = DLANTR( 'I', UPLO, DIAG, N, N, AINV, LDA, 319 $ RWORK ) 320 IF( ANORM.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN 321 RCONDI = ONE 322 ELSE 323 RCONDI = ( ONE / ANORM ) / AINVNM 324 END IF 325* 326* Compute the residual for the triangular matrix times 327* its inverse. Also compute the 1-norm condition number 328* of A. 329* 330 CALL DTRT01( UPLO, DIAG, N, A, LDA, AINV, LDA, RCONDO, 331 $ RWORK, RESULT( 1 ) ) 332* 333* Print the test ratio if it is .GE. THRESH. 334* 335 IF( RESULT( 1 ).GE.THRESH ) THEN 336 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 337 $ CALL ALAHD( NOUT, PATH ) 338 WRITE( NOUT, FMT = 9999 )UPLO, DIAG, N, NB, IMAT, 339 $ 1, RESULT( 1 ) 340 NFAIL = NFAIL + 1 341 END IF 342 NRUN = NRUN + 1 343* 344* Skip remaining tests if not the first block size. 345* 346 IF( INB.NE.1 ) 347 $ GO TO 60 348* 349 DO 40 IRHS = 1, NNS 350 NRHS = NSVAL( IRHS ) 351 XTYPE = 'N' 352* 353 DO 30 ITRAN = 1, NTRAN 354* 355* Do for op(A) = A, A**T, or A**H. 356* 357 TRANS = TRANSS( ITRAN ) 358 IF( ITRAN.EQ.1 ) THEN 359 NORM = 'O' 360 RCONDC = RCONDO 361 ELSE 362 NORM = 'I' 363 RCONDC = RCONDI 364 END IF 365* 366*+ TEST 2 367* Solve and compute residual for op(A)*x = b. 368* 369 SRNAMT = 'DLARHS' 370 CALL DLARHS( PATH, XTYPE, UPLO, TRANS, N, N, 0, 371 $ IDIAG, NRHS, A, LDA, XACT, LDA, B, 372 $ LDA, ISEED, INFO ) 373 XTYPE = 'C' 374 CALL DLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 375* 376 SRNAMT = 'DTRTRS' 377 CALL DTRTRS( UPLO, TRANS, DIAG, N, NRHS, A, LDA, 378 $ X, LDA, INFO ) 379* 380* Check error code from DTRTRS. 381* 382 IF( INFO.NE.0 ) 383 $ CALL ALAERH( PATH, 'DTRTRS', INFO, 0, 384 $ UPLO // TRANS // DIAG, N, N, -1, 385 $ -1, NRHS, IMAT, NFAIL, NERRS, 386 $ NOUT ) 387* 388* This line is needed on a Sun SPARCstation. 389* 390 IF( N.GT.0 ) 391 $ DUMMY = A( 1 ) 392* 393 CALL DTRT02( UPLO, TRANS, DIAG, N, NRHS, A, LDA, 394 $ X, LDA, B, LDA, WORK, RESULT( 2 ) ) 395* 396*+ TEST 3 397* Check solution from generated exact solution. 398* 399 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 400 $ RESULT( 3 ) ) 401* 402*+ TESTS 4, 5, and 6 403* Use iterative refinement to improve the solution 404* and compute error bounds. 405* 406 SRNAMT = 'DTRRFS' 407 CALL DTRRFS( UPLO, TRANS, DIAG, N, NRHS, A, LDA, 408 $ B, LDA, X, LDA, RWORK, 409 $ RWORK( NRHS+1 ), WORK, IWORK, 410 $ INFO ) 411* 412* Check error code from DTRRFS. 413* 414 IF( INFO.NE.0 ) 415 $ CALL ALAERH( PATH, 'DTRRFS', INFO, 0, 416 $ UPLO // TRANS // DIAG, N, N, -1, 417 $ -1, NRHS, IMAT, NFAIL, NERRS, 418 $ NOUT ) 419* 420 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 421 $ RESULT( 4 ) ) 422 CALL DTRT05( UPLO, TRANS, DIAG, N, NRHS, A, LDA, 423 $ B, LDA, X, LDA, XACT, LDA, RWORK, 424 $ RWORK( NRHS+1 ), RESULT( 5 ) ) 425* 426* Print information about the tests that did not 427* pass the threshold. 428* 429 DO 20 K = 2, 6 430 IF( RESULT( K ).GE.THRESH ) THEN 431 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 432 $ CALL ALAHD( NOUT, PATH ) 433 WRITE( NOUT, FMT = 9998 )UPLO, TRANS, 434 $ DIAG, N, NRHS, IMAT, K, RESULT( K ) 435 NFAIL = NFAIL + 1 436 END IF 437 20 CONTINUE 438 NRUN = NRUN + 5 439 30 CONTINUE 440 40 CONTINUE 441* 442*+ TEST 7 443* Get an estimate of RCOND = 1/CNDNUM. 444* 445 DO 50 ITRAN = 1, 2 446 IF( ITRAN.EQ.1 ) THEN 447 NORM = 'O' 448 RCONDC = RCONDO 449 ELSE 450 NORM = 'I' 451 RCONDC = RCONDI 452 END IF 453 SRNAMT = 'DTRCON' 454 CALL DTRCON( NORM, UPLO, DIAG, N, A, LDA, RCOND, 455 $ WORK, IWORK, INFO ) 456* 457* Check error code from DTRCON. 458* 459 IF( INFO.NE.0 ) 460 $ CALL ALAERH( PATH, 'DTRCON', INFO, 0, 461 $ NORM // UPLO // DIAG, N, N, -1, -1, 462 $ -1, IMAT, NFAIL, NERRS, NOUT ) 463* 464 CALL DTRT06( RCOND, RCONDC, UPLO, DIAG, N, A, LDA, 465 $ RWORK, RESULT( 7 ) ) 466* 467* Print the test ratio if it is .GE. THRESH. 468* 469 IF( RESULT( 7 ).GE.THRESH ) THEN 470 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 471 $ CALL ALAHD( NOUT, PATH ) 472 WRITE( NOUT, FMT = 9997 )NORM, UPLO, N, IMAT, 473 $ 7, RESULT( 7 ) 474 NFAIL = NFAIL + 1 475 END IF 476 NRUN = NRUN + 1 477 50 CONTINUE 478 60 CONTINUE 479 70 CONTINUE 480 80 CONTINUE 481* 482* Use pathological test matrices to test DLATRS. 483* 484 DO 110 IMAT = NTYPE1 + 1, NTYPES 485* 486* Do the tests only if DOTYPE( IMAT ) is true. 487* 488 IF( .NOT.DOTYPE( IMAT ) ) 489 $ GO TO 110 490* 491 DO 100 IUPLO = 1, 2 492* 493* Do first for UPLO = 'U', then for UPLO = 'L' 494* 495 UPLO = UPLOS( IUPLO ) 496 DO 90 ITRAN = 1, NTRAN 497* 498* Do for op(A) = A, A**T, and A**H. 499* 500 TRANS = TRANSS( ITRAN ) 501* 502* Call DLATTR to generate a triangular test matrix. 503* 504 SRNAMT = 'DLATTR' 505 CALL DLATTR( IMAT, UPLO, TRANS, DIAG, ISEED, N, A, 506 $ LDA, X, WORK, INFO ) 507* 508*+ TEST 8 509* Solve the system op(A)*x = b. 510* 511 SRNAMT = 'DLATRS' 512 CALL DCOPY( N, X, 1, B, 1 ) 513 CALL DLATRS( UPLO, TRANS, DIAG, 'N', N, A, LDA, B, 514 $ SCALE, RWORK, INFO ) 515* 516* Check error code from DLATRS. 517* 518 IF( INFO.NE.0 ) 519 $ CALL ALAERH( PATH, 'DLATRS', INFO, 0, 520 $ UPLO // TRANS // DIAG // 'N', N, N, 521 $ -1, -1, -1, IMAT, NFAIL, NERRS, NOUT ) 522* 523 CALL DTRT03( UPLO, TRANS, DIAG, N, 1, A, LDA, SCALE, 524 $ RWORK, ONE, B, LDA, X, LDA, WORK, 525 $ RESULT( 8 ) ) 526* 527*+ TEST 9 528* Solve op(A)*X = b again with NORMIN = 'Y'. 529* 530 CALL DCOPY( N, X, 1, B( N+1 ), 1 ) 531 CALL DLATRS( UPLO, TRANS, DIAG, 'Y', N, A, LDA, 532 $ B( N+1 ), SCALE, RWORK, INFO ) 533* 534* Check error code from DLATRS. 535* 536 IF( INFO.NE.0 ) 537 $ CALL ALAERH( PATH, 'DLATRS', INFO, 0, 538 $ UPLO // TRANS // DIAG // 'Y', N, N, 539 $ -1, -1, -1, IMAT, NFAIL, NERRS, NOUT ) 540* 541 CALL DTRT03( UPLO, TRANS, DIAG, N, 1, A, LDA, SCALE, 542 $ RWORK, ONE, B( N+1 ), LDA, X, LDA, WORK, 543 $ RESULT( 9 ) ) 544* 545* Print information about the tests that did not pass 546* the threshold. 547* 548 IF( RESULT( 8 ).GE.THRESH ) THEN 549 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 550 $ CALL ALAHD( NOUT, PATH ) 551 WRITE( NOUT, FMT = 9996 )'DLATRS', UPLO, TRANS, 552 $ DIAG, 'N', N, IMAT, 8, RESULT( 8 ) 553 NFAIL = NFAIL + 1 554 END IF 555 IF( RESULT( 9 ).GE.THRESH ) THEN 556 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 557 $ CALL ALAHD( NOUT, PATH ) 558 WRITE( NOUT, FMT = 9996 )'DLATRS', UPLO, TRANS, 559 $ DIAG, 'Y', N, IMAT, 9, RESULT( 9 ) 560 NFAIL = NFAIL + 1 561 END IF 562 NRUN = NRUN + 2 563 90 CONTINUE 564 100 CONTINUE 565 110 CONTINUE 566 120 CONTINUE 567* 568* Print a summary of the results. 569* 570 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS ) 571* 572 9999 FORMAT( ' UPLO=''', A1, ''', DIAG=''', A1, ''', N=', I5, ', NB=', 573 $ I4, ', type ', I2, ', test(', I2, ')= ', G12.5 ) 574 9998 FORMAT( ' UPLO=''', A1, ''', TRANS=''', A1, ''', DIAG=''', A1, 575 $ ''', N=', I5, ', NB=', I4, ', type ', I2, ', 576 $ test(', I2, ')= ', G12.5 ) 577 9997 FORMAT( ' NORM=''', A1, ''', UPLO =''', A1, ''', N=', I5, ',', 578 $ 11X, ' type ', I2, ', test(', I2, ')=', G12.5 ) 579 9996 FORMAT( 1X, A, '( ''', A1, ''', ''', A1, ''', ''', A1, ''', ''', 580 $ A1, ''',', I5, ', ... ), type ', I2, ', test(', I2, ')=', 581 $ G12.5 ) 582 RETURN 583* 584* End of DCHKTR 585* 586 END 587