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