1*> \brief \b ZCHKBB 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 ZCHKBB( NSIZES, MVAL, NVAL, NWDTHS, KK, NTYPES, DOTYPE, 12* NRHS, ISEED, THRESH, NOUNIT, A, LDA, AB, LDAB, 13* BD, BE, Q, LDQ, P, LDP, C, LDC, CC, WORK, 14* LWORK, RWORK, RESULT, INFO ) 15* 16* .. Scalar Arguments .. 17* INTEGER INFO, LDA, LDAB, LDC, LDP, LDQ, LWORK, NOUNIT, 18* $ NRHS, NSIZES, NTYPES, NWDTHS 19* DOUBLE PRECISION THRESH 20* .. 21* .. Array Arguments .. 22* LOGICAL DOTYPE( * ) 23* INTEGER ISEED( 4 ), KK( * ), MVAL( * ), NVAL( * ) 24* DOUBLE PRECISION BD( * ), BE( * ), RESULT( * ), RWORK( * ) 25* COMPLEX*16 A( LDA, * ), AB( LDAB, * ), C( LDC, * ), 26* $ CC( LDC, * ), P( LDP, * ), Q( LDQ, * ), 27* $ WORK( * ) 28* .. 29* 30* 31*> \par Purpose: 32* ============= 33*> 34*> \verbatim 35*> 36*> ZCHKBB tests the reduction of a general complex rectangular band 37*> matrix to real bidiagonal form. 38*> 39*> ZGBBRD factors a general band matrix A as Q B P* , where * means 40*> conjugate transpose, B is upper bidiagonal, and Q and P are unitary; 41*> ZGBBRD can also overwrite a given matrix C with Q* C . 42*> 43*> For each pair of matrix dimensions (M,N) and each selected matrix 44*> type, an M by N matrix A and an M by NRHS matrix C are generated. 45*> The problem dimensions are as follows 46*> A: M x N 47*> Q: M x M 48*> P: N x N 49*> B: min(M,N) x min(M,N) 50*> C: M x NRHS 51*> 52*> For each generated matrix, 4 tests are performed: 53*> 54*> (1) | A - Q B PT | / ( |A| max(M,N) ulp ), PT = P' 55*> 56*> (2) | I - Q' Q | / ( M ulp ) 57*> 58*> (3) | I - PT PT' | / ( N ulp ) 59*> 60*> (4) | Y - Q' C | / ( |Y| max(M,NRHS) ulp ), where Y = Q' C. 61*> 62*> The "types" are specified by a logical array DOTYPE( 1:NTYPES ); 63*> if DOTYPE(j) is .TRUE., then matrix type "j" will be generated. 64*> Currently, the list of possible types is: 65*> 66*> The possible matrix types are 67*> 68*> (1) The zero matrix. 69*> (2) The identity matrix. 70*> 71*> (3) A diagonal matrix with evenly spaced entries 72*> 1, ..., ULP and random signs. 73*> (ULP = (first number larger than 1) - 1 ) 74*> (4) A diagonal matrix with geometrically spaced entries 75*> 1, ..., ULP and random signs. 76*> (5) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP 77*> and random signs. 78*> 79*> (6) Same as (3), but multiplied by SQRT( overflow threshold ) 80*> (7) Same as (3), but multiplied by SQRT( underflow threshold ) 81*> 82*> (8) A matrix of the form U D V, where U and V are orthogonal and 83*> D has evenly spaced entries 1, ..., ULP with random signs 84*> on the diagonal. 85*> 86*> (9) A matrix of the form U D V, where U and V are orthogonal and 87*> D has geometrically spaced entries 1, ..., ULP with random 88*> signs on the diagonal. 89*> 90*> (10) A matrix of the form U D V, where U and V are orthogonal and 91*> D has "clustered" entries 1, ULP,..., ULP with random 92*> signs on the diagonal. 93*> 94*> (11) Same as (8), but multiplied by SQRT( overflow threshold ) 95*> (12) Same as (8), but multiplied by SQRT( underflow threshold ) 96*> 97*> (13) Rectangular matrix with random entries chosen from (-1,1). 98*> (14) Same as (13), but multiplied by SQRT( overflow threshold ) 99*> (15) Same as (13), but multiplied by SQRT( underflow threshold ) 100*> \endverbatim 101* 102* Arguments: 103* ========== 104* 105*> \param[in] NSIZES 106*> \verbatim 107*> NSIZES is INTEGER 108*> The number of values of M and N contained in the vectors 109*> MVAL and NVAL. The matrix sizes are used in pairs (M,N). 110*> If NSIZES is zero, ZCHKBB does nothing. NSIZES must be at 111*> least zero. 112*> \endverbatim 113*> 114*> \param[in] MVAL 115*> \verbatim 116*> MVAL is INTEGER array, dimension (NSIZES) 117*> The values of the matrix row dimension M. 118*> \endverbatim 119*> 120*> \param[in] NVAL 121*> \verbatim 122*> NVAL is INTEGER array, dimension (NSIZES) 123*> The values of the matrix column dimension N. 124*> \endverbatim 125*> 126*> \param[in] NWDTHS 127*> \verbatim 128*> NWDTHS is INTEGER 129*> The number of bandwidths to use. If it is zero, 130*> ZCHKBB does nothing. It must be at least zero. 131*> \endverbatim 132*> 133*> \param[in] KK 134*> \verbatim 135*> KK is INTEGER array, dimension (NWDTHS) 136*> An array containing the bandwidths to be used for the band 137*> matrices. The values must be at least zero. 138*> \endverbatim 139*> 140*> \param[in] NTYPES 141*> \verbatim 142*> NTYPES is INTEGER 143*> The number of elements in DOTYPE. If it is zero, ZCHKBB 144*> does nothing. It must be at least zero. If it is MAXTYP+1 145*> and NSIZES is 1, then an additional type, MAXTYP+1 is 146*> defined, which is to use whatever matrix is in A. This 147*> is only useful if DOTYPE(1:MAXTYP) is .FALSE. and 148*> DOTYPE(MAXTYP+1) is .TRUE. . 149*> \endverbatim 150*> 151*> \param[in] DOTYPE 152*> \verbatim 153*> DOTYPE is LOGICAL array, dimension (NTYPES) 154*> If DOTYPE(j) is .TRUE., then for each size in NN a 155*> matrix of that size and of type j will be generated. 156*> If NTYPES is smaller than the maximum number of types 157*> defined (PARAMETER MAXTYP), then types NTYPES+1 through 158*> MAXTYP will not be generated. If NTYPES is larger 159*> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES) 160*> will be ignored. 161*> \endverbatim 162*> 163*> \param[in] NRHS 164*> \verbatim 165*> NRHS is INTEGER 166*> The number of columns in the "right-hand side" matrix C. 167*> If NRHS = 0, then the operations on the right-hand side will 168*> not be tested. NRHS must be at least 0. 169*> \endverbatim 170*> 171*> \param[in,out] ISEED 172*> \verbatim 173*> ISEED is INTEGER array, dimension (4) 174*> On entry ISEED specifies the seed of the random number 175*> generator. The array elements should be between 0 and 4095; 176*> if not they will be reduced mod 4096. Also, ISEED(4) must 177*> be odd. The random number generator uses a linear 178*> congruential sequence limited to small integers, and so 179*> should produce machine independent random numbers. The 180*> values of ISEED are changed on exit, and can be used in the 181*> next call to ZCHKBB to continue the same random number 182*> sequence. 183*> \endverbatim 184*> 185*> \param[in] THRESH 186*> \verbatim 187*> THRESH is DOUBLE PRECISION 188*> A test will count as "failed" if the "error", computed as 189*> described above, exceeds THRESH. Note that the error 190*> is scaled to be O(1), so THRESH should be a reasonably 191*> small multiple of 1, e.g., 10 or 100. In particular, 192*> it should not depend on the precision (single vs. double) 193*> or the size of the matrix. It must be at least zero. 194*> \endverbatim 195*> 196*> \param[in] NOUNIT 197*> \verbatim 198*> NOUNIT is INTEGER 199*> The FORTRAN unit number for printing out error messages 200*> (e.g., if a routine returns IINFO not equal to 0.) 201*> \endverbatim 202*> 203*> \param[in,out] A 204*> \verbatim 205*> A is DOUBLE PRECISION array, dimension 206*> (LDA, max(NN)) 207*> Used to hold the matrix A. 208*> \endverbatim 209*> 210*> \param[in] LDA 211*> \verbatim 212*> LDA is INTEGER 213*> The leading dimension of A. It must be at least 1 214*> and at least max( NN ). 215*> \endverbatim 216*> 217*> \param[out] AB 218*> \verbatim 219*> AB is DOUBLE PRECISION array, dimension (LDAB, max(NN)) 220*> Used to hold A in band storage format. 221*> \endverbatim 222*> 223*> \param[in] LDAB 224*> \verbatim 225*> LDAB is INTEGER 226*> The leading dimension of AB. It must be at least 2 (not 1!) 227*> and at least max( KK )+1. 228*> \endverbatim 229*> 230*> \param[out] BD 231*> \verbatim 232*> BD is DOUBLE PRECISION array, dimension (max(NN)) 233*> Used to hold the diagonal of the bidiagonal matrix computed 234*> by ZGBBRD. 235*> \endverbatim 236*> 237*> \param[out] BE 238*> \verbatim 239*> BE is DOUBLE PRECISION array, dimension (max(NN)) 240*> Used to hold the off-diagonal of the bidiagonal matrix 241*> computed by ZGBBRD. 242*> \endverbatim 243*> 244*> \param[out] Q 245*> \verbatim 246*> Q is COMPLEX*16 array, dimension (LDQ, max(NN)) 247*> Used to hold the unitary matrix Q computed by ZGBBRD. 248*> \endverbatim 249*> 250*> \param[in] LDQ 251*> \verbatim 252*> LDQ is INTEGER 253*> The leading dimension of Q. It must be at least 1 254*> and at least max( NN ). 255*> \endverbatim 256*> 257*> \param[out] P 258*> \verbatim 259*> P is COMPLEX*16 array, dimension (LDP, max(NN)) 260*> Used to hold the unitary matrix P computed by ZGBBRD. 261*> \endverbatim 262*> 263*> \param[in] LDP 264*> \verbatim 265*> LDP is INTEGER 266*> The leading dimension of P. It must be at least 1 267*> and at least max( NN ). 268*> \endverbatim 269*> 270*> \param[out] C 271*> \verbatim 272*> C is COMPLEX*16 array, dimension (LDC, max(NN)) 273*> Used to hold the matrix C updated by ZGBBRD. 274*> \endverbatim 275*> 276*> \param[in] LDC 277*> \verbatim 278*> LDC is INTEGER 279*> The leading dimension of U. It must be at least 1 280*> and at least max( NN ). 281*> \endverbatim 282*> 283*> \param[out] CC 284*> \verbatim 285*> CC is COMPLEX*16 array, dimension (LDC, max(NN)) 286*> Used to hold a copy of the matrix C. 287*> \endverbatim 288*> 289*> \param[out] WORK 290*> \verbatim 291*> WORK is COMPLEX*16 array, dimension (LWORK) 292*> \endverbatim 293*> 294*> \param[in] LWORK 295*> \verbatim 296*> LWORK is INTEGER 297*> The number of entries in WORK. This must be at least 298*> max( LDA+1, max(NN)+1 )*max(NN). 299*> \endverbatim 300*> 301*> \param[out] RWORK 302*> \verbatim 303*> RWORK is DOUBLE PRECISION array, dimension (max(NN)) 304*> \endverbatim 305*> 306*> \param[out] RESULT 307*> \verbatim 308*> RESULT is DOUBLE PRECISION array, dimension (4) 309*> The values computed by the tests described above. 310*> The values are currently limited to 1/ulp, to avoid 311*> overflow. 312*> \endverbatim 313*> 314*> \param[out] INFO 315*> \verbatim 316*> INFO is INTEGER 317*> If 0, then everything ran OK. 318*> 319*>----------------------------------------------------------------------- 320*> 321*> Some Local Variables and Parameters: 322*> ---- ----- --------- --- ---------- 323*> ZERO, ONE Real 0 and 1. 324*> MAXTYP The number of types defined. 325*> NTEST The number of tests performed, or which can 326*> be performed so far, for the current matrix. 327*> NTESTT The total number of tests performed so far. 328*> NMAX Largest value in NN. 329*> NMATS The number of matrices generated so far. 330*> NERRS The number of tests which have exceeded THRESH 331*> so far. 332*> COND, IMODE Values to be passed to the matrix generators. 333*> ANORM Norm of A; passed to matrix generators. 334*> 335*> OVFL, UNFL Overflow and underflow thresholds. 336*> ULP, ULPINV Finest relative precision and its inverse. 337*> RTOVFL, RTUNFL Square roots of the previous 2 values. 338*> The following four arrays decode JTYPE: 339*> KTYPE(j) The general type (1-10) for type "j". 340*> KMODE(j) The MODE value to be passed to the matrix 341*> generator for type "j". 342*> KMAGN(j) The order of magnitude ( O(1), 343*> O(overflow^(1/2) ), O(underflow^(1/2) ) 344*> \endverbatim 345* 346* Authors: 347* ======== 348* 349*> \author Univ. of Tennessee 350*> \author Univ. of California Berkeley 351*> \author Univ. of Colorado Denver 352*> \author NAG Ltd. 353* 354*> \ingroup complex16_eig 355* 356* ===================================================================== 357 SUBROUTINE ZCHKBB( NSIZES, MVAL, NVAL, NWDTHS, KK, NTYPES, DOTYPE, 358 $ NRHS, ISEED, THRESH, NOUNIT, A, LDA, AB, LDAB, 359 $ BD, BE, Q, LDQ, P, LDP, C, LDC, CC, WORK, 360 $ LWORK, RWORK, RESULT, INFO ) 361* 362* -- LAPACK test routine (input) -- 363* -- LAPACK is a software package provided by Univ. of Tennessee, -- 364* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 365* 366* .. Scalar Arguments .. 367 INTEGER INFO, LDA, LDAB, LDC, LDP, LDQ, LWORK, NOUNIT, 368 $ NRHS, NSIZES, NTYPES, NWDTHS 369 DOUBLE PRECISION THRESH 370* .. 371* .. Array Arguments .. 372 LOGICAL DOTYPE( * ) 373 INTEGER ISEED( 4 ), KK( * ), MVAL( * ), NVAL( * ) 374 DOUBLE PRECISION BD( * ), BE( * ), RESULT( * ), RWORK( * ) 375 COMPLEX*16 A( LDA, * ), AB( LDAB, * ), C( LDC, * ), 376 $ CC( LDC, * ), P( LDP, * ), Q( LDQ, * ), 377 $ WORK( * ) 378* .. 379* 380* ===================================================================== 381* 382* .. Parameters .. 383 COMPLEX*16 CZERO, CONE 384 PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ), 385 $ CONE = ( 1.0D+0, 0.0D+0 ) ) 386 DOUBLE PRECISION ZERO, ONE 387 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 ) 388 INTEGER MAXTYP 389 PARAMETER ( MAXTYP = 15 ) 390* .. 391* .. Local Scalars .. 392 LOGICAL BADMM, BADNN, BADNNB 393 INTEGER I, IINFO, IMODE, ITYPE, J, JCOL, JR, JSIZE, 394 $ JTYPE, JWIDTH, K, KL, KMAX, KU, M, MMAX, MNMAX, 395 $ MNMIN, MTYPES, N, NERRS, NMATS, NMAX, NTEST, 396 $ NTESTT 397 DOUBLE PRECISION AMNINV, ANORM, COND, OVFL, RTOVFL, RTUNFL, ULP, 398 $ ULPINV, UNFL 399* .. 400* .. Local Arrays .. 401 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KMAGN( MAXTYP ), 402 $ KMODE( MAXTYP ), KTYPE( MAXTYP ) 403* .. 404* .. External Functions .. 405 DOUBLE PRECISION DLAMCH 406 EXTERNAL DLAMCH 407* .. 408* .. External Subroutines .. 409 EXTERNAL DLAHD2, DLASUM, XERBLA, ZBDT01, ZBDT02, ZGBBRD, 410 $ ZLACPY, ZLASET, ZLATMR, ZLATMS, ZUNT01 411* .. 412* .. Intrinsic Functions .. 413 INTRINSIC ABS, DBLE, MAX, MIN, SQRT 414* .. 415* .. Data statements .. 416 DATA KTYPE / 1, 2, 5*4, 5*6, 3*9 / 417 DATA KMAGN / 2*1, 3*1, 2, 3, 3*1, 2, 3, 1, 2, 3 / 418 DATA KMODE / 2*0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0, 419 $ 0, 0 / 420* .. 421* .. Executable Statements .. 422* 423* Check for errors 424* 425 NTESTT = 0 426 INFO = 0 427* 428* Important constants 429* 430 BADMM = .FALSE. 431 BADNN = .FALSE. 432 MMAX = 1 433 NMAX = 1 434 MNMAX = 1 435 DO 10 J = 1, NSIZES 436 MMAX = MAX( MMAX, MVAL( J ) ) 437 IF( MVAL( J ).LT.0 ) 438 $ BADMM = .TRUE. 439 NMAX = MAX( NMAX, NVAL( J ) ) 440 IF( NVAL( J ).LT.0 ) 441 $ BADNN = .TRUE. 442 MNMAX = MAX( MNMAX, MIN( MVAL( J ), NVAL( J ) ) ) 443 10 CONTINUE 444* 445 BADNNB = .FALSE. 446 KMAX = 0 447 DO 20 J = 1, NWDTHS 448 KMAX = MAX( KMAX, KK( J ) ) 449 IF( KK( J ).LT.0 ) 450 $ BADNNB = .TRUE. 451 20 CONTINUE 452* 453* Check for errors 454* 455 IF( NSIZES.LT.0 ) THEN 456 INFO = -1 457 ELSE IF( BADMM ) THEN 458 INFO = -2 459 ELSE IF( BADNN ) THEN 460 INFO = -3 461 ELSE IF( NWDTHS.LT.0 ) THEN 462 INFO = -4 463 ELSE IF( BADNNB ) THEN 464 INFO = -5 465 ELSE IF( NTYPES.LT.0 ) THEN 466 INFO = -6 467 ELSE IF( NRHS.LT.0 ) THEN 468 INFO = -8 469 ELSE IF( LDA.LT.NMAX ) THEN 470 INFO = -13 471 ELSE IF( LDAB.LT.2*KMAX+1 ) THEN 472 INFO = -15 473 ELSE IF( LDQ.LT.NMAX ) THEN 474 INFO = -19 475 ELSE IF( LDP.LT.NMAX ) THEN 476 INFO = -21 477 ELSE IF( LDC.LT.NMAX ) THEN 478 INFO = -23 479 ELSE IF( ( MAX( LDA, NMAX )+1 )*NMAX.GT.LWORK ) THEN 480 INFO = -26 481 END IF 482* 483 IF( INFO.NE.0 ) THEN 484 CALL XERBLA( 'ZCHKBB', -INFO ) 485 RETURN 486 END IF 487* 488* Quick return if possible 489* 490 IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 .OR. NWDTHS.EQ.0 ) 491 $ RETURN 492* 493* More Important constants 494* 495 UNFL = DLAMCH( 'Safe minimum' ) 496 OVFL = ONE / UNFL 497 ULP = DLAMCH( 'Epsilon' )*DLAMCH( 'Base' ) 498 ULPINV = ONE / ULP 499 RTUNFL = SQRT( UNFL ) 500 RTOVFL = SQRT( OVFL ) 501* 502* Loop over sizes, widths, types 503* 504 NERRS = 0 505 NMATS = 0 506* 507 DO 160 JSIZE = 1, NSIZES 508 M = MVAL( JSIZE ) 509 N = NVAL( JSIZE ) 510 MNMIN = MIN( M, N ) 511 AMNINV = ONE / DBLE( MAX( 1, M, N ) ) 512* 513 DO 150 JWIDTH = 1, NWDTHS 514 K = KK( JWIDTH ) 515 IF( K.GE.M .AND. K.GE.N ) 516 $ GO TO 150 517 KL = MAX( 0, MIN( M-1, K ) ) 518 KU = MAX( 0, MIN( N-1, K ) ) 519* 520 IF( NSIZES.NE.1 ) THEN 521 MTYPES = MIN( MAXTYP, NTYPES ) 522 ELSE 523 MTYPES = MIN( MAXTYP+1, NTYPES ) 524 END IF 525* 526 DO 140 JTYPE = 1, MTYPES 527 IF( .NOT.DOTYPE( JTYPE ) ) 528 $ GO TO 140 529 NMATS = NMATS + 1 530 NTEST = 0 531* 532 DO 30 J = 1, 4 533 IOLDSD( J ) = ISEED( J ) 534 30 CONTINUE 535* 536* Compute "A". 537* 538* Control parameters: 539* 540* KMAGN KMODE KTYPE 541* =1 O(1) clustered 1 zero 542* =2 large clustered 2 identity 543* =3 small exponential (none) 544* =4 arithmetic diagonal, (w/ singular values) 545* =5 random log (none) 546* =6 random nonhermitian, w/ singular values 547* =7 (none) 548* =8 (none) 549* =9 random nonhermitian 550* 551 IF( MTYPES.GT.MAXTYP ) 552 $ GO TO 90 553* 554 ITYPE = KTYPE( JTYPE ) 555 IMODE = KMODE( JTYPE ) 556* 557* Compute norm 558* 559 GO TO ( 40, 50, 60 )KMAGN( JTYPE ) 560* 561 40 CONTINUE 562 ANORM = ONE 563 GO TO 70 564* 565 50 CONTINUE 566 ANORM = ( RTOVFL*ULP )*AMNINV 567 GO TO 70 568* 569 60 CONTINUE 570 ANORM = RTUNFL*MAX( M, N )*ULPINV 571 GO TO 70 572* 573 70 CONTINUE 574* 575 CALL ZLASET( 'Full', LDA, N, CZERO, CZERO, A, LDA ) 576 CALL ZLASET( 'Full', LDAB, N, CZERO, CZERO, AB, LDAB ) 577 IINFO = 0 578 COND = ULPINV 579* 580* Special Matrices -- Identity & Jordan block 581* 582* Zero 583* 584 IF( ITYPE.EQ.1 ) THEN 585 IINFO = 0 586* 587 ELSE IF( ITYPE.EQ.2 ) THEN 588* 589* Identity 590* 591 DO 80 JCOL = 1, N 592 A( JCOL, JCOL ) = ANORM 593 80 CONTINUE 594* 595 ELSE IF( ITYPE.EQ.4 ) THEN 596* 597* Diagonal Matrix, singular values specified 598* 599 CALL ZLATMS( M, N, 'S', ISEED, 'N', RWORK, IMODE, 600 $ COND, ANORM, 0, 0, 'N', A, LDA, WORK, 601 $ IINFO ) 602* 603 ELSE IF( ITYPE.EQ.6 ) THEN 604* 605* Nonhermitian, singular values specified 606* 607 CALL ZLATMS( M, N, 'S', ISEED, 'N', RWORK, IMODE, 608 $ COND, ANORM, KL, KU, 'N', A, LDA, WORK, 609 $ IINFO ) 610* 611 ELSE IF( ITYPE.EQ.9 ) THEN 612* 613* Nonhermitian, random entries 614* 615 CALL ZLATMR( M, N, 'S', ISEED, 'N', WORK, 6, ONE, 616 $ CONE, 'T', 'N', WORK( N+1 ), 1, ONE, 617 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, KL, 618 $ KU, ZERO, ANORM, 'N', A, LDA, IDUMMA, 619 $ IINFO ) 620* 621 ELSE 622* 623 IINFO = 1 624 END IF 625* 626* Generate Right-Hand Side 627* 628 CALL ZLATMR( M, NRHS, 'S', ISEED, 'N', WORK, 6, ONE, 629 $ CONE, 'T', 'N', WORK( M+1 ), 1, ONE, 630 $ WORK( 2*M+1 ), 1, ONE, 'N', IDUMMA, M, NRHS, 631 $ ZERO, ONE, 'NO', C, LDC, IDUMMA, IINFO ) 632* 633 IF( IINFO.NE.0 ) THEN 634 WRITE( NOUNIT, FMT = 9999 )'Generator', IINFO, N, 635 $ JTYPE, IOLDSD 636 INFO = ABS( IINFO ) 637 RETURN 638 END IF 639* 640 90 CONTINUE 641* 642* Copy A to band storage. 643* 644 DO 110 J = 1, N 645 DO 100 I = MAX( 1, J-KU ), MIN( M, J+KL ) 646 AB( KU+1+I-J, J ) = A( I, J ) 647 100 CONTINUE 648 110 CONTINUE 649* 650* Copy C 651* 652 CALL ZLACPY( 'Full', M, NRHS, C, LDC, CC, LDC ) 653* 654* Call ZGBBRD to compute B, Q and P, and to update C. 655* 656 CALL ZGBBRD( 'B', M, N, NRHS, KL, KU, AB, LDAB, BD, BE, 657 $ Q, LDQ, P, LDP, CC, LDC, WORK, RWORK, 658 $ IINFO ) 659* 660 IF( IINFO.NE.0 ) THEN 661 WRITE( NOUNIT, FMT = 9999 )'ZGBBRD', IINFO, N, JTYPE, 662 $ IOLDSD 663 INFO = ABS( IINFO ) 664 IF( IINFO.LT.0 ) THEN 665 RETURN 666 ELSE 667 RESULT( 1 ) = ULPINV 668 GO TO 120 669 END IF 670 END IF 671* 672* Test 1: Check the decomposition A := Q * B * P' 673* 2: Check the orthogonality of Q 674* 3: Check the orthogonality of P 675* 4: Check the computation of Q' * C 676* 677 CALL ZBDT01( M, N, -1, A, LDA, Q, LDQ, BD, BE, P, LDP, 678 $ WORK, RWORK, RESULT( 1 ) ) 679 CALL ZUNT01( 'Columns', M, M, Q, LDQ, WORK, LWORK, RWORK, 680 $ RESULT( 2 ) ) 681 CALL ZUNT01( 'Rows', N, N, P, LDP, WORK, LWORK, RWORK, 682 $ RESULT( 3 ) ) 683 CALL ZBDT02( M, NRHS, C, LDC, CC, LDC, Q, LDQ, WORK, 684 $ RWORK, RESULT( 4 ) ) 685* 686* End of Loop -- Check for RESULT(j) > THRESH 687* 688 NTEST = 4 689 120 CONTINUE 690 NTESTT = NTESTT + NTEST 691* 692* Print out tests which fail. 693* 694 DO 130 JR = 1, NTEST 695 IF( RESULT( JR ).GE.THRESH ) THEN 696 IF( NERRS.EQ.0 ) 697 $ CALL DLAHD2( NOUNIT, 'ZBB' ) 698 NERRS = NERRS + 1 699 WRITE( NOUNIT, FMT = 9998 )M, N, K, IOLDSD, JTYPE, 700 $ JR, RESULT( JR ) 701 END IF 702 130 CONTINUE 703* 704 140 CONTINUE 705 150 CONTINUE 706 160 CONTINUE 707* 708* Summary 709* 710 CALL DLASUM( 'ZBB', NOUNIT, NERRS, NTESTT ) 711 RETURN 712* 713 9999 FORMAT( ' ZCHKBB: ', A, ' returned INFO=', I5, '.', / 9X, 'M=', 714 $ I5, ' N=', I5, ' K=', I5, ', JTYPE=', I5, ', ISEED=(', 715 $ 3( I5, ',' ), I5, ')' ) 716 9998 FORMAT( ' M =', I4, ' N=', I4, ', K=', I3, ', seed=', 717 $ 4( I4, ',' ), ' type ', I2, ', test(', I2, ')=', G10.3 ) 718* 719* End of ZCHKBB 720* 721 END 722