1*> \brief \b CGET23 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 CGET23( COMP, ISRT, BALANC, JTYPE, THRESH, ISEED, 12* NOUNIT, N, A, LDA, H, W, W1, VL, LDVL, VR, 13* LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN, 14* RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT, 15* WORK, LWORK, RWORK, INFO ) 16* 17* .. Scalar Arguments .. 18* LOGICAL COMP 19* CHARACTER BALANC 20* INTEGER INFO, ISRT, JTYPE, LDA, LDLRE, LDVL, LDVR, 21* $ LWORK, N, NOUNIT 22* REAL THRESH 23* .. 24* .. Array Arguments .. 25* INTEGER ISEED( 4 ) 26* REAL RCDEIN( * ), RCDVIN( * ), RCNDE1( * ), 27* $ RCNDV1( * ), RCONDE( * ), RCONDV( * ), 28* $ RESULT( 11 ), RWORK( * ), SCALE( * ), 29* $ SCALE1( * ) 30* COMPLEX A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ), 31* $ VL( LDVL, * ), VR( LDVR, * ), W( * ), W1( * ), 32* $ WORK( * ) 33* .. 34* 35* 36*> \par Purpose: 37* ============= 38*> 39*> \verbatim 40*> 41*> CGET23 checks the nonsymmetric eigenvalue problem driver CGEEVX. 42*> If COMP = .FALSE., the first 8 of the following tests will be 43*> performed on the input matrix A, and also test 9 if LWORK is 44*> sufficiently large. 45*> if COMP is .TRUE. all 11 tests will be performed. 46*> 47*> (1) | A * VR - VR * W | / ( n |A| ulp ) 48*> 49*> Here VR is the matrix of unit right eigenvectors. 50*> W is a diagonal matrix with diagonal entries W(j). 51*> 52*> (2) | A**H * VL - VL * W**H | / ( n |A| ulp ) 53*> 54*> Here VL is the matrix of unit left eigenvectors, A**H is the 55*> conjugate transpose of A, and W is as above. 56*> 57*> (3) | |VR(i)| - 1 | / ulp and largest component real 58*> 59*> VR(i) denotes the i-th column of VR. 60*> 61*> (4) | |VL(i)| - 1 | / ulp and largest component real 62*> 63*> VL(i) denotes the i-th column of VL. 64*> 65*> (5) 0 if W(full) = W(partial), 1/ulp otherwise 66*> 67*> W(full) denotes the eigenvalues computed when VR, VL, RCONDV 68*> and RCONDE are also computed, and W(partial) denotes the 69*> eigenvalues computed when only some of VR, VL, RCONDV, and 70*> RCONDE are computed. 71*> 72*> (6) 0 if VR(full) = VR(partial), 1/ulp otherwise 73*> 74*> VR(full) denotes the right eigenvectors computed when VL, RCONDV 75*> and RCONDE are computed, and VR(partial) denotes the result 76*> when only some of VL and RCONDV are computed. 77*> 78*> (7) 0 if VL(full) = VL(partial), 1/ulp otherwise 79*> 80*> VL(full) denotes the left eigenvectors computed when VR, RCONDV 81*> and RCONDE are computed, and VL(partial) denotes the result 82*> when only some of VR and RCONDV are computed. 83*> 84*> (8) 0 if SCALE, ILO, IHI, ABNRM (full) = 85*> SCALE, ILO, IHI, ABNRM (partial) 86*> 1/ulp otherwise 87*> 88*> SCALE, ILO, IHI and ABNRM describe how the matrix is balanced. 89*> (full) is when VR, VL, RCONDE and RCONDV are also computed, and 90*> (partial) is when some are not computed. 91*> 92*> (9) 0 if RCONDV(full) = RCONDV(partial), 1/ulp otherwise 93*> 94*> RCONDV(full) denotes the reciprocal condition numbers of the 95*> right eigenvectors computed when VR, VL and RCONDE are also 96*> computed. RCONDV(partial) denotes the reciprocal condition 97*> numbers when only some of VR, VL and RCONDE are computed. 98*> 99*> (10) |RCONDV - RCDVIN| / cond(RCONDV) 100*> 101*> RCONDV is the reciprocal right eigenvector condition number 102*> computed by CGEEVX and RCDVIN (the precomputed true value) 103*> is supplied as input. cond(RCONDV) is the condition number of 104*> RCONDV, and takes errors in computing RCONDV into account, so 105*> that the resulting quantity should be O(ULP). cond(RCONDV) is 106*> essentially given by norm(A)/RCONDE. 107*> 108*> (11) |RCONDE - RCDEIN| / cond(RCONDE) 109*> 110*> RCONDE is the reciprocal eigenvalue condition number 111*> computed by CGEEVX and RCDEIN (the precomputed true value) 112*> is supplied as input. cond(RCONDE) is the condition number 113*> of RCONDE, and takes errors in computing RCONDE into account, 114*> so that the resulting quantity should be O(ULP). cond(RCONDE) 115*> is essentially given by norm(A)/RCONDV. 116*> \endverbatim 117* 118* Arguments: 119* ========== 120* 121*> \param[in] COMP 122*> \verbatim 123*> COMP is LOGICAL 124*> COMP describes which input tests to perform: 125*> = .FALSE. if the computed condition numbers are not to 126*> be tested against RCDVIN and RCDEIN 127*> = .TRUE. if they are to be compared 128*> \endverbatim 129*> 130*> \param[in] ISRT 131*> \verbatim 132*> ISRT is INTEGER 133*> If COMP = .TRUE., ISRT indicates in how the eigenvalues 134*> corresponding to values in RCDVIN and RCDEIN are ordered: 135*> = 0 means the eigenvalues are sorted by 136*> increasing real part 137*> = 1 means the eigenvalues are sorted by 138*> increasing imaginary part 139*> If COMP = .FALSE., ISRT is not referenced. 140*> \endverbatim 141*> 142*> \param[in] BALANC 143*> \verbatim 144*> BALANC is CHARACTER 145*> Describes the balancing option to be tested. 146*> = 'N' for no permuting or diagonal scaling 147*> = 'P' for permuting but no diagonal scaling 148*> = 'S' for no permuting but diagonal scaling 149*> = 'B' for permuting and diagonal scaling 150*> \endverbatim 151*> 152*> \param[in] JTYPE 153*> \verbatim 154*> JTYPE is INTEGER 155*> Type of input matrix. Used to label output if error occurs. 156*> \endverbatim 157*> 158*> \param[in] THRESH 159*> \verbatim 160*> THRESH is REAL 161*> A test will count as "failed" if the "error", computed as 162*> described above, exceeds THRESH. Note that the error 163*> is scaled to be O(1), so THRESH should be a reasonably 164*> small multiple of 1, e.g., 10 or 100. In particular, 165*> it should not depend on the precision (single vs. double) 166*> or the size of the matrix. It must be at least zero. 167*> \endverbatim 168*> 169*> \param[in] ISEED 170*> \verbatim 171*> ISEED is INTEGER array, dimension (4) 172*> If COMP = .FALSE., the random number generator seed 173*> used to produce matrix. 174*> If COMP = .TRUE., ISEED(1) = the number of the example. 175*> Used to label output if error occurs. 176*> \endverbatim 177*> 178*> \param[in] NOUNIT 179*> \verbatim 180*> NOUNIT is INTEGER 181*> The FORTRAN unit number for printing out error messages 182*> (e.g., if a routine returns INFO not equal to 0.) 183*> \endverbatim 184*> 185*> \param[in] N 186*> \verbatim 187*> N is INTEGER 188*> The dimension of A. N must be at least 0. 189*> \endverbatim 190*> 191*> \param[in,out] A 192*> \verbatim 193*> A is COMPLEX array, dimension (LDA,N) 194*> Used to hold the matrix whose eigenvalues are to be 195*> computed. 196*> \endverbatim 197*> 198*> \param[in] LDA 199*> \verbatim 200*> LDA is INTEGER 201*> The leading dimension of A, and H. LDA must be at 202*> least 1 and at least N. 203*> \endverbatim 204*> 205*> \param[out] H 206*> \verbatim 207*> H is COMPLEX array, dimension (LDA,N) 208*> Another copy of the test matrix A, modified by CGEEVX. 209*> \endverbatim 210*> 211*> \param[out] W 212*> \verbatim 213*> W is COMPLEX array, dimension (N) 214*> Contains the eigenvalues of A. 215*> \endverbatim 216*> 217*> \param[out] W1 218*> \verbatim 219*> W1 is COMPLEX array, dimension (N) 220*> Like W, this array contains the eigenvalues of A, 221*> but those computed when CGEEVX only computes a partial 222*> eigendecomposition, i.e. not the eigenvalues and left 223*> and right eigenvectors. 224*> \endverbatim 225*> 226*> \param[out] VL 227*> \verbatim 228*> VL is COMPLEX array, dimension (LDVL,N) 229*> VL holds the computed left eigenvectors. 230*> \endverbatim 231*> 232*> \param[in] LDVL 233*> \verbatim 234*> LDVL is INTEGER 235*> Leading dimension of VL. Must be at least max(1,N). 236*> \endverbatim 237*> 238*> \param[out] VR 239*> \verbatim 240*> VR is COMPLEX array, dimension (LDVR,N) 241*> VR holds the computed right eigenvectors. 242*> \endverbatim 243*> 244*> \param[in] LDVR 245*> \verbatim 246*> LDVR is INTEGER 247*> Leading dimension of VR. Must be at least max(1,N). 248*> \endverbatim 249*> 250*> \param[out] LRE 251*> \verbatim 252*> LRE is COMPLEX array, dimension (LDLRE,N) 253*> LRE holds the computed right or left eigenvectors. 254*> \endverbatim 255*> 256*> \param[in] LDLRE 257*> \verbatim 258*> LDLRE is INTEGER 259*> Leading dimension of LRE. Must be at least max(1,N). 260*> \endverbatim 261*> 262*> \param[out] RCONDV 263*> \verbatim 264*> RCONDV is REAL array, dimension (N) 265*> RCONDV holds the computed reciprocal condition numbers 266*> for eigenvectors. 267*> \endverbatim 268*> 269*> \param[out] RCNDV1 270*> \verbatim 271*> RCNDV1 is REAL array, dimension (N) 272*> RCNDV1 holds more computed reciprocal condition numbers 273*> for eigenvectors. 274*> \endverbatim 275*> 276*> \param[in] RCDVIN 277*> \verbatim 278*> RCDVIN is REAL array, dimension (N) 279*> When COMP = .TRUE. RCDVIN holds the precomputed reciprocal 280*> condition numbers for eigenvectors to be compared with 281*> RCONDV. 282*> \endverbatim 283*> 284*> \param[out] RCONDE 285*> \verbatim 286*> RCONDE is REAL array, dimension (N) 287*> RCONDE holds the computed reciprocal condition numbers 288*> for eigenvalues. 289*> \endverbatim 290*> 291*> \param[out] RCNDE1 292*> \verbatim 293*> RCNDE1 is REAL array, dimension (N) 294*> RCNDE1 holds more computed reciprocal condition numbers 295*> for eigenvalues. 296*> \endverbatim 297*> 298*> \param[in] RCDEIN 299*> \verbatim 300*> RCDEIN is REAL array, dimension (N) 301*> When COMP = .TRUE. RCDEIN holds the precomputed reciprocal 302*> condition numbers for eigenvalues to be compared with 303*> RCONDE. 304*> \endverbatim 305*> 306*> \param[out] SCALE 307*> \verbatim 308*> SCALE is REAL array, dimension (N) 309*> Holds information describing balancing of matrix. 310*> \endverbatim 311*> 312*> \param[out] SCALE1 313*> \verbatim 314*> SCALE1 is REAL array, dimension (N) 315*> Holds information describing balancing of matrix. 316*> \endverbatim 317*> 318*> \param[out] RESULT 319*> \verbatim 320*> RESULT is REAL array, dimension (11) 321*> The values computed by the 11 tests described above. 322*> The values are currently limited to 1/ulp, to avoid 323*> overflow. 324*> \endverbatim 325*> 326*> \param[out] WORK 327*> \verbatim 328*> WORK is COMPLEX array, dimension (LWORK) 329*> \endverbatim 330*> 331*> \param[in] LWORK 332*> \verbatim 333*> LWORK is INTEGER 334*> The number of entries in WORK. This must be at least 335*> 2*N, and 2*N+N**2 if tests 9, 10 or 11 are to be performed. 336*> \endverbatim 337*> 338*> \param[out] RWORK 339*> \verbatim 340*> RWORK is REAL array, dimension (2*N) 341*> \endverbatim 342*> 343*> \param[out] INFO 344*> \verbatim 345*> INFO is INTEGER 346*> If 0, successful exit. 347*> If <0, input parameter -INFO had an incorrect value. 348*> If >0, CGEEVX returned an error code, the absolute 349*> value of which is returned. 350*> \endverbatim 351* 352* Authors: 353* ======== 354* 355*> \author Univ. of Tennessee 356*> \author Univ. of California Berkeley 357*> \author Univ. of Colorado Denver 358*> \author NAG Ltd. 359* 360*> \ingroup complex_eig 361* 362* ===================================================================== 363 SUBROUTINE CGET23( COMP, ISRT, BALANC, JTYPE, THRESH, ISEED, 364 $ NOUNIT, N, A, LDA, H, W, W1, VL, LDVL, VR, 365 $ LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN, 366 $ RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT, 367 $ WORK, LWORK, RWORK, INFO ) 368* 369* -- LAPACK test routine -- 370* -- LAPACK is a software package provided by Univ. of Tennessee, -- 371* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 372* 373* .. Scalar Arguments .. 374 LOGICAL COMP 375 CHARACTER BALANC 376 INTEGER INFO, ISRT, JTYPE, LDA, LDLRE, LDVL, LDVR, 377 $ LWORK, N, NOUNIT 378 REAL THRESH 379* .. 380* .. Array Arguments .. 381 INTEGER ISEED( 4 ) 382 REAL RCDEIN( * ), RCDVIN( * ), RCNDE1( * ), 383 $ RCNDV1( * ), RCONDE( * ), RCONDV( * ), 384 $ RESULT( 11 ), RWORK( * ), SCALE( * ), 385 $ SCALE1( * ) 386 COMPLEX A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ), 387 $ VL( LDVL, * ), VR( LDVR, * ), W( * ), W1( * ), 388 $ WORK( * ) 389* .. 390* 391* ===================================================================== 392* 393* .. Parameters .. 394 REAL ZERO, ONE, TWO 395 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0 ) 396 REAL EPSIN 397 PARAMETER ( EPSIN = 5.9605E-8 ) 398* .. 399* .. Local Scalars .. 400 LOGICAL BALOK, NOBAL 401 CHARACTER SENSE 402 INTEGER I, IHI, IHI1, IINFO, ILO, ILO1, ISENS, ISENSM, 403 $ J, JJ, KMIN 404 REAL ABNRM, ABNRM1, EPS, SMLNUM, TNRM, TOL, TOLIN, 405 $ ULP, ULPINV, V, VMAX, VMX, VRICMP, VRIMIN, 406 $ VRMX, VTST 407 COMPLEX CTMP 408* .. 409* .. Local Arrays .. 410 CHARACTER SENS( 2 ) 411 REAL RES( 2 ) 412 COMPLEX CDUM( 1 ) 413* .. 414* .. External Functions .. 415 LOGICAL LSAME 416 REAL SCNRM2, SLAMCH 417 EXTERNAL LSAME, SCNRM2, SLAMCH 418* .. 419* .. External Subroutines .. 420 EXTERNAL CGEEVX, CGET22, CLACPY, XERBLA 421* .. 422* .. Intrinsic Functions .. 423 INTRINSIC ABS, AIMAG, MAX, MIN, REAL 424* .. 425* .. Data statements .. 426 DATA SENS / 'N', 'V' / 427* .. 428* .. Executable Statements .. 429* 430* Check for errors 431* 432 NOBAL = LSAME( BALANC, 'N' ) 433 BALOK = NOBAL .OR. LSAME( BALANC, 'P' ) .OR. 434 $ LSAME( BALANC, 'S' ) .OR. LSAME( BALANC, 'B' ) 435 INFO = 0 436 IF( ISRT.NE.0 .AND. ISRT.NE.1 ) THEN 437 INFO = -2 438 ELSE IF( .NOT.BALOK ) THEN 439 INFO = -3 440 ELSE IF( THRESH.LT.ZERO ) THEN 441 INFO = -5 442 ELSE IF( NOUNIT.LE.0 ) THEN 443 INFO = -7 444 ELSE IF( N.LT.0 ) THEN 445 INFO = -8 446 ELSE IF( LDA.LT.1 .OR. LDA.LT.N ) THEN 447 INFO = -10 448 ELSE IF( LDVL.LT.1 .OR. LDVL.LT.N ) THEN 449 INFO = -15 450 ELSE IF( LDVR.LT.1 .OR. LDVR.LT.N ) THEN 451 INFO = -17 452 ELSE IF( LDLRE.LT.1 .OR. LDLRE.LT.N ) THEN 453 INFO = -19 454 ELSE IF( LWORK.LT.2*N .OR. ( COMP .AND. LWORK.LT.2*N+N*N ) ) THEN 455 INFO = -30 456 END IF 457* 458 IF( INFO.NE.0 ) THEN 459 CALL XERBLA( 'CGET23', -INFO ) 460 RETURN 461 END IF 462* 463* Quick return if nothing to do 464* 465 DO 10 I = 1, 11 466 RESULT( I ) = -ONE 467 10 CONTINUE 468* 469 IF( N.EQ.0 ) 470 $ RETURN 471* 472* More Important constants 473* 474 ULP = SLAMCH( 'Precision' ) 475 SMLNUM = SLAMCH( 'S' ) 476 ULPINV = ONE / ULP 477* 478* Compute eigenvalues and eigenvectors, and test them 479* 480 IF( LWORK.GE.2*N+N*N ) THEN 481 SENSE = 'B' 482 ISENSM = 2 483 ELSE 484 SENSE = 'E' 485 ISENSM = 1 486 END IF 487 CALL CLACPY( 'F', N, N, A, LDA, H, LDA ) 488 CALL CGEEVX( BALANC, 'V', 'V', SENSE, N, H, LDA, W, VL, LDVL, VR, 489 $ LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV, WORK, 490 $ LWORK, RWORK, IINFO ) 491 IF( IINFO.NE.0 ) THEN 492 RESULT( 1 ) = ULPINV 493 IF( JTYPE.NE.22 ) THEN 494 WRITE( NOUNIT, FMT = 9998 )'CGEEVX1', IINFO, N, JTYPE, 495 $ BALANC, ISEED 496 ELSE 497 WRITE( NOUNIT, FMT = 9999 )'CGEEVX1', IINFO, N, ISEED( 1 ) 498 END IF 499 INFO = ABS( IINFO ) 500 RETURN 501 END IF 502* 503* Do Test (1) 504* 505 CALL CGET22( 'N', 'N', 'N', N, A, LDA, VR, LDVR, W, WORK, RWORK, 506 $ RES ) 507 RESULT( 1 ) = RES( 1 ) 508* 509* Do Test (2) 510* 511 CALL CGET22( 'C', 'N', 'C', N, A, LDA, VL, LDVL, W, WORK, RWORK, 512 $ RES ) 513 RESULT( 2 ) = RES( 1 ) 514* 515* Do Test (3) 516* 517 DO 30 J = 1, N 518 TNRM = SCNRM2( N, VR( 1, J ), 1 ) 519 RESULT( 3 ) = MAX( RESULT( 3 ), 520 $ MIN( ULPINV, ABS( TNRM-ONE ) / ULP ) ) 521 VMX = ZERO 522 VRMX = ZERO 523 DO 20 JJ = 1, N 524 VTST = ABS( VR( JJ, J ) ) 525 IF( VTST.GT.VMX ) 526 $ VMX = VTST 527 IF( AIMAG( VR( JJ, J ) ).EQ.ZERO .AND. 528 $ ABS( REAL( VR( JJ, J ) ) ).GT.VRMX ) 529 $ VRMX = ABS( REAL( VR( JJ, J ) ) ) 530 20 CONTINUE 531 IF( VRMX / VMX.LT.ONE-TWO*ULP ) 532 $ RESULT( 3 ) = ULPINV 533 30 CONTINUE 534* 535* Do Test (4) 536* 537 DO 50 J = 1, N 538 TNRM = SCNRM2( N, VL( 1, J ), 1 ) 539 RESULT( 4 ) = MAX( RESULT( 4 ), 540 $ MIN( ULPINV, ABS( TNRM-ONE ) / ULP ) ) 541 VMX = ZERO 542 VRMX = ZERO 543 DO 40 JJ = 1, N 544 VTST = ABS( VL( JJ, J ) ) 545 IF( VTST.GT.VMX ) 546 $ VMX = VTST 547 IF( AIMAG( VL( JJ, J ) ).EQ.ZERO .AND. 548 $ ABS( REAL( VL( JJ, J ) ) ).GT.VRMX ) 549 $ VRMX = ABS( REAL( VL( JJ, J ) ) ) 550 40 CONTINUE 551 IF( VRMX / VMX.LT.ONE-TWO*ULP ) 552 $ RESULT( 4 ) = ULPINV 553 50 CONTINUE 554* 555* Test for all options of computing condition numbers 556* 557 DO 200 ISENS = 1, ISENSM 558* 559 SENSE = SENS( ISENS ) 560* 561* Compute eigenvalues only, and test them 562* 563 CALL CLACPY( 'F', N, N, A, LDA, H, LDA ) 564 CALL CGEEVX( BALANC, 'N', 'N', SENSE, N, H, LDA, W1, CDUM, 1, 565 $ CDUM, 1, ILO1, IHI1, SCALE1, ABNRM1, RCNDE1, 566 $ RCNDV1, WORK, LWORK, RWORK, IINFO ) 567 IF( IINFO.NE.0 ) THEN 568 RESULT( 1 ) = ULPINV 569 IF( JTYPE.NE.22 ) THEN 570 WRITE( NOUNIT, FMT = 9998 )'CGEEVX2', IINFO, N, JTYPE, 571 $ BALANC, ISEED 572 ELSE 573 WRITE( NOUNIT, FMT = 9999 )'CGEEVX2', IINFO, N, 574 $ ISEED( 1 ) 575 END IF 576 INFO = ABS( IINFO ) 577 GO TO 190 578 END IF 579* 580* Do Test (5) 581* 582 DO 60 J = 1, N 583 IF( W( J ).NE.W1( J ) ) 584 $ RESULT( 5 ) = ULPINV 585 60 CONTINUE 586* 587* Do Test (8) 588* 589 IF( .NOT.NOBAL ) THEN 590 DO 70 J = 1, N 591 IF( SCALE( J ).NE.SCALE1( J ) ) 592 $ RESULT( 8 ) = ULPINV 593 70 CONTINUE 594 IF( ILO.NE.ILO1 ) 595 $ RESULT( 8 ) = ULPINV 596 IF( IHI.NE.IHI1 ) 597 $ RESULT( 8 ) = ULPINV 598 IF( ABNRM.NE.ABNRM1 ) 599 $ RESULT( 8 ) = ULPINV 600 END IF 601* 602* Do Test (9) 603* 604 IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN 605 DO 80 J = 1, N 606 IF( RCONDV( J ).NE.RCNDV1( J ) ) 607 $ RESULT( 9 ) = ULPINV 608 80 CONTINUE 609 END IF 610* 611* Compute eigenvalues and right eigenvectors, and test them 612* 613 CALL CLACPY( 'F', N, N, A, LDA, H, LDA ) 614 CALL CGEEVX( BALANC, 'N', 'V', SENSE, N, H, LDA, W1, CDUM, 1, 615 $ LRE, LDLRE, ILO1, IHI1, SCALE1, ABNRM1, RCNDE1, 616 $ RCNDV1, WORK, LWORK, RWORK, IINFO ) 617 IF( IINFO.NE.0 ) THEN 618 RESULT( 1 ) = ULPINV 619 IF( JTYPE.NE.22 ) THEN 620 WRITE( NOUNIT, FMT = 9998 )'CGEEVX3', IINFO, N, JTYPE, 621 $ BALANC, ISEED 622 ELSE 623 WRITE( NOUNIT, FMT = 9999 )'CGEEVX3', IINFO, N, 624 $ ISEED( 1 ) 625 END IF 626 INFO = ABS( IINFO ) 627 GO TO 190 628 END IF 629* 630* Do Test (5) again 631* 632 DO 90 J = 1, N 633 IF( W( J ).NE.W1( J ) ) 634 $ RESULT( 5 ) = ULPINV 635 90 CONTINUE 636* 637* Do Test (6) 638* 639 DO 110 J = 1, N 640 DO 100 JJ = 1, N 641 IF( VR( J, JJ ).NE.LRE( J, JJ ) ) 642 $ RESULT( 6 ) = ULPINV 643 100 CONTINUE 644 110 CONTINUE 645* 646* Do Test (8) again 647* 648 IF( .NOT.NOBAL ) THEN 649 DO 120 J = 1, N 650 IF( SCALE( J ).NE.SCALE1( J ) ) 651 $ RESULT( 8 ) = ULPINV 652 120 CONTINUE 653 IF( ILO.NE.ILO1 ) 654 $ RESULT( 8 ) = ULPINV 655 IF( IHI.NE.IHI1 ) 656 $ RESULT( 8 ) = ULPINV 657 IF( ABNRM.NE.ABNRM1 ) 658 $ RESULT( 8 ) = ULPINV 659 END IF 660* 661* Do Test (9) again 662* 663 IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN 664 DO 130 J = 1, N 665 IF( RCONDV( J ).NE.RCNDV1( J ) ) 666 $ RESULT( 9 ) = ULPINV 667 130 CONTINUE 668 END IF 669* 670* Compute eigenvalues and left eigenvectors, and test them 671* 672 CALL CLACPY( 'F', N, N, A, LDA, H, LDA ) 673 CALL CGEEVX( BALANC, 'V', 'N', SENSE, N, H, LDA, W1, LRE, 674 $ LDLRE, CDUM, 1, ILO1, IHI1, SCALE1, ABNRM1, 675 $ RCNDE1, RCNDV1, WORK, LWORK, RWORK, IINFO ) 676 IF( IINFO.NE.0 ) THEN 677 RESULT( 1 ) = ULPINV 678 IF( JTYPE.NE.22 ) THEN 679 WRITE( NOUNIT, FMT = 9998 )'CGEEVX4', IINFO, N, JTYPE, 680 $ BALANC, ISEED 681 ELSE 682 WRITE( NOUNIT, FMT = 9999 )'CGEEVX4', IINFO, N, 683 $ ISEED( 1 ) 684 END IF 685 INFO = ABS( IINFO ) 686 GO TO 190 687 END IF 688* 689* Do Test (5) again 690* 691 DO 140 J = 1, N 692 IF( W( J ).NE.W1( J ) ) 693 $ RESULT( 5 ) = ULPINV 694 140 CONTINUE 695* 696* Do Test (7) 697* 698 DO 160 J = 1, N 699 DO 150 JJ = 1, N 700 IF( VL( J, JJ ).NE.LRE( J, JJ ) ) 701 $ RESULT( 7 ) = ULPINV 702 150 CONTINUE 703 160 CONTINUE 704* 705* Do Test (8) again 706* 707 IF( .NOT.NOBAL ) THEN 708 DO 170 J = 1, N 709 IF( SCALE( J ).NE.SCALE1( J ) ) 710 $ RESULT( 8 ) = ULPINV 711 170 CONTINUE 712 IF( ILO.NE.ILO1 ) 713 $ RESULT( 8 ) = ULPINV 714 IF( IHI.NE.IHI1 ) 715 $ RESULT( 8 ) = ULPINV 716 IF( ABNRM.NE.ABNRM1 ) 717 $ RESULT( 8 ) = ULPINV 718 END IF 719* 720* Do Test (9) again 721* 722 IF( ISENS.EQ.2 .AND. N.GT.1 ) THEN 723 DO 180 J = 1, N 724 IF( RCONDV( J ).NE.RCNDV1( J ) ) 725 $ RESULT( 9 ) = ULPINV 726 180 CONTINUE 727 END IF 728* 729 190 CONTINUE 730* 731 200 CONTINUE 732* 733* If COMP, compare condition numbers to precomputed ones 734* 735 IF( COMP ) THEN 736 CALL CLACPY( 'F', N, N, A, LDA, H, LDA ) 737 CALL CGEEVX( 'N', 'V', 'V', 'B', N, H, LDA, W, VL, LDVL, VR, 738 $ LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV, 739 $ WORK, LWORK, RWORK, IINFO ) 740 IF( IINFO.NE.0 ) THEN 741 RESULT( 1 ) = ULPINV 742 WRITE( NOUNIT, FMT = 9999 )'CGEEVX5', IINFO, N, ISEED( 1 ) 743 INFO = ABS( IINFO ) 744 GO TO 250 745 END IF 746* 747* Sort eigenvalues and condition numbers lexicographically 748* to compare with inputs 749* 750 DO 220 I = 1, N - 1 751 KMIN = I 752 IF( ISRT.EQ.0 ) THEN 753 VRIMIN = REAL( W( I ) ) 754 ELSE 755 VRIMIN = AIMAG( W( I ) ) 756 END IF 757 DO 210 J = I + 1, N 758 IF( ISRT.EQ.0 ) THEN 759 VRICMP = REAL( W( J ) ) 760 ELSE 761 VRICMP = AIMAG( W( J ) ) 762 END IF 763 IF( VRICMP.LT.VRIMIN ) THEN 764 KMIN = J 765 VRIMIN = VRICMP 766 END IF 767 210 CONTINUE 768 CTMP = W( KMIN ) 769 W( KMIN ) = W( I ) 770 W( I ) = CTMP 771 VRIMIN = RCONDE( KMIN ) 772 RCONDE( KMIN ) = RCONDE( I ) 773 RCONDE( I ) = VRIMIN 774 VRIMIN = RCONDV( KMIN ) 775 RCONDV( KMIN ) = RCONDV( I ) 776 RCONDV( I ) = VRIMIN 777 220 CONTINUE 778* 779* Compare condition numbers for eigenvectors 780* taking their condition numbers into account 781* 782 RESULT( 10 ) = ZERO 783 EPS = MAX( EPSIN, ULP ) 784 V = MAX( REAL( N )*EPS*ABNRM, SMLNUM ) 785 IF( ABNRM.EQ.ZERO ) 786 $ V = ONE 787 DO 230 I = 1, N 788 IF( V.GT.RCONDV( I )*RCONDE( I ) ) THEN 789 TOL = RCONDV( I ) 790 ELSE 791 TOL = V / RCONDE( I ) 792 END IF 793 IF( V.GT.RCDVIN( I )*RCDEIN( I ) ) THEN 794 TOLIN = RCDVIN( I ) 795 ELSE 796 TOLIN = V / RCDEIN( I ) 797 END IF 798 TOL = MAX( TOL, SMLNUM / EPS ) 799 TOLIN = MAX( TOLIN, SMLNUM / EPS ) 800 IF( EPS*( RCDVIN( I )-TOLIN ).GT.RCONDV( I )+TOL ) THEN 801 VMAX = ONE / EPS 802 ELSE IF( RCDVIN( I )-TOLIN.GT.RCONDV( I )+TOL ) THEN 803 VMAX = ( RCDVIN( I )-TOLIN ) / ( RCONDV( I )+TOL ) 804 ELSE IF( RCDVIN( I )+TOLIN.LT.EPS*( RCONDV( I )-TOL ) ) THEN 805 VMAX = ONE / EPS 806 ELSE IF( RCDVIN( I )+TOLIN.LT.RCONDV( I )-TOL ) THEN 807 VMAX = ( RCONDV( I )-TOL ) / ( RCDVIN( I )+TOLIN ) 808 ELSE 809 VMAX = ONE 810 END IF 811 RESULT( 10 ) = MAX( RESULT( 10 ), VMAX ) 812 230 CONTINUE 813* 814* Compare condition numbers for eigenvalues 815* taking their condition numbers into account 816* 817 RESULT( 11 ) = ZERO 818 DO 240 I = 1, N 819 IF( V.GT.RCONDV( I ) ) THEN 820 TOL = ONE 821 ELSE 822 TOL = V / RCONDV( I ) 823 END IF 824 IF( V.GT.RCDVIN( I ) ) THEN 825 TOLIN = ONE 826 ELSE 827 TOLIN = V / RCDVIN( I ) 828 END IF 829 TOL = MAX( TOL, SMLNUM / EPS ) 830 TOLIN = MAX( TOLIN, SMLNUM / EPS ) 831 IF( EPS*( RCDEIN( I )-TOLIN ).GT.RCONDE( I )+TOL ) THEN 832 VMAX = ONE / EPS 833 ELSE IF( RCDEIN( I )-TOLIN.GT.RCONDE( I )+TOL ) THEN 834 VMAX = ( RCDEIN( I )-TOLIN ) / ( RCONDE( I )+TOL ) 835 ELSE IF( RCDEIN( I )+TOLIN.LT.EPS*( RCONDE( I )-TOL ) ) THEN 836 VMAX = ONE / EPS 837 ELSE IF( RCDEIN( I )+TOLIN.LT.RCONDE( I )-TOL ) THEN 838 VMAX = ( RCONDE( I )-TOL ) / ( RCDEIN( I )+TOLIN ) 839 ELSE 840 VMAX = ONE 841 END IF 842 RESULT( 11 ) = MAX( RESULT( 11 ), VMAX ) 843 240 CONTINUE 844 250 CONTINUE 845* 846 END IF 847* 848 9999 FORMAT( ' CGET23: ', A, ' returned INFO=', I6, '.', / 9X, 'N=', 849 $ I6, ', INPUT EXAMPLE NUMBER = ', I4 ) 850 9998 FORMAT( ' CGET23: ', A, ' returned INFO=', I6, '.', / 9X, 'N=', 851 $ I6, ', JTYPE=', I6, ', BALANC = ', A, ', ISEED=(', 852 $ 3( I5, ',' ), I5, ')' ) 853* 854 RETURN 855* 856* End of CGET23 857* 858 END 859