1*> \brief <b> DGEESX computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices</b> 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8*> \htmlonly 9*> Download DGEESX + dependencies 10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgeesx.f"> 11*> [TGZ]</a> 12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgeesx.f"> 13*> [ZIP]</a> 14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgeesx.f"> 15*> [TXT]</a> 16*> \endhtmlonly 17* 18* Definition: 19* =========== 20* 21* SUBROUTINE DGEESX( JOBVS, SORT, SELECT, SENSE, N, A, LDA, SDIM, 22* WR, WI, VS, LDVS, RCONDE, RCONDV, WORK, LWORK, 23* IWORK, LIWORK, BWORK, INFO ) 24* 25* .. Scalar Arguments .. 26* CHARACTER JOBVS, SENSE, SORT 27* INTEGER INFO, LDA, LDVS, LIWORK, LWORK, N, SDIM 28* DOUBLE PRECISION RCONDE, RCONDV 29* .. 30* .. Array Arguments .. 31* LOGICAL BWORK( * ) 32* INTEGER IWORK( * ) 33* DOUBLE PRECISION A( LDA, * ), VS( LDVS, * ), WI( * ), WORK( * ), 34* $ WR( * ) 35* .. 36* .. Function Arguments .. 37* LOGICAL SELECT 38* EXTERNAL SELECT 39* .. 40* 41* 42*> \par Purpose: 43* ============= 44*> 45*> \verbatim 46*> 47*> DGEESX computes for an N-by-N real nonsymmetric matrix A, the 48*> eigenvalues, the real Schur form T, and, optionally, the matrix of 49*> Schur vectors Z. This gives the Schur factorization A = Z*T*(Z**T). 50*> 51*> Optionally, it also orders the eigenvalues on the diagonal of the 52*> real Schur form so that selected eigenvalues are at the top left; 53*> computes a reciprocal condition number for the average of the 54*> selected eigenvalues (RCONDE); and computes a reciprocal condition 55*> number for the right invariant subspace corresponding to the 56*> selected eigenvalues (RCONDV). The leading columns of Z form an 57*> orthonormal basis for this invariant subspace. 58*> 59*> For further explanation of the reciprocal condition numbers RCONDE 60*> and RCONDV, see Section 4.10 of the LAPACK Users' Guide (where 61*> these quantities are called s and sep respectively). 62*> 63*> A real matrix is in real Schur form if it is upper quasi-triangular 64*> with 1-by-1 and 2-by-2 blocks. 2-by-2 blocks will be standardized in 65*> the form 66*> [ a b ] 67*> [ c a ] 68*> 69*> where b*c < 0. The eigenvalues of such a block are a +- sqrt(bc). 70*> \endverbatim 71* 72* Arguments: 73* ========== 74* 75*> \param[in] JOBVS 76*> \verbatim 77*> JOBVS is CHARACTER*1 78*> = 'N': Schur vectors are not computed; 79*> = 'V': Schur vectors are computed. 80*> \endverbatim 81*> 82*> \param[in] SORT 83*> \verbatim 84*> SORT is CHARACTER*1 85*> Specifies whether or not to order the eigenvalues on the 86*> diagonal of the Schur form. 87*> = 'N': Eigenvalues are not ordered; 88*> = 'S': Eigenvalues are ordered (see SELECT). 89*> \endverbatim 90*> 91*> \param[in] SELECT 92*> \verbatim 93*> SELECT is procedure) LOGICAL FUNCTION of two DOUBLE PRECISION arguments 94*> SELECT must be declared EXTERNAL in the calling subroutine. 95*> If SORT = 'S', SELECT is used to select eigenvalues to sort 96*> to the top left of the Schur form. 97*> If SORT = 'N', SELECT is not referenced. 98*> An eigenvalue WR(j)+sqrt(-1)*WI(j) is selected if 99*> SELECT(WR(j),WI(j)) is true; i.e., if either one of a 100*> complex conjugate pair of eigenvalues is selected, then both 101*> are. Note that a selected complex eigenvalue may no longer 102*> satisfy SELECT(WR(j),WI(j)) = .TRUE. after ordering, since 103*> ordering may change the value of complex eigenvalues 104*> (especially if the eigenvalue is ill-conditioned); in this 105*> case INFO may be set to N+3 (see INFO below). 106*> \endverbatim 107*> 108*> \param[in] SENSE 109*> \verbatim 110*> SENSE is CHARACTER*1 111*> Determines which reciprocal condition numbers are computed. 112*> = 'N': None are computed; 113*> = 'E': Computed for average of selected eigenvalues only; 114*> = 'V': Computed for selected right invariant subspace only; 115*> = 'B': Computed for both. 116*> If SENSE = 'E', 'V' or 'B', SORT must equal 'S'. 117*> \endverbatim 118*> 119*> \param[in] N 120*> \verbatim 121*> N is INTEGER 122*> The order of the matrix A. N >= 0. 123*> \endverbatim 124*> 125*> \param[in,out] A 126*> \verbatim 127*> A is DOUBLE PRECISION array, dimension (LDA, N) 128*> On entry, the N-by-N matrix A. 129*> On exit, A is overwritten by its real Schur form T. 130*> \endverbatim 131*> 132*> \param[in] LDA 133*> \verbatim 134*> LDA is INTEGER 135*> The leading dimension of the array A. LDA >= max(1,N). 136*> \endverbatim 137*> 138*> \param[out] SDIM 139*> \verbatim 140*> SDIM is INTEGER 141*> If SORT = 'N', SDIM = 0. 142*> If SORT = 'S', SDIM = number of eigenvalues (after sorting) 143*> for which SELECT is true. (Complex conjugate 144*> pairs for which SELECT is true for either 145*> eigenvalue count as 2.) 146*> \endverbatim 147*> 148*> \param[out] WR 149*> \verbatim 150*> WR is DOUBLE PRECISION array, dimension (N) 151*> \endverbatim 152*> 153*> \param[out] WI 154*> \verbatim 155*> WI is DOUBLE PRECISION array, dimension (N) 156*> WR and WI contain the real and imaginary parts, respectively, 157*> of the computed eigenvalues, in the same order that they 158*> appear on the diagonal of the output Schur form T. Complex 159*> conjugate pairs of eigenvalues appear consecutively with the 160*> eigenvalue having the positive imaginary part first. 161*> \endverbatim 162*> 163*> \param[out] VS 164*> \verbatim 165*> VS is DOUBLE PRECISION array, dimension (LDVS,N) 166*> If JOBVS = 'V', VS contains the orthogonal matrix Z of Schur 167*> vectors. 168*> If JOBVS = 'N', VS is not referenced. 169*> \endverbatim 170*> 171*> \param[in] LDVS 172*> \verbatim 173*> LDVS is INTEGER 174*> The leading dimension of the array VS. LDVS >= 1, and if 175*> JOBVS = 'V', LDVS >= N. 176*> \endverbatim 177*> 178*> \param[out] RCONDE 179*> \verbatim 180*> RCONDE is DOUBLE PRECISION 181*> If SENSE = 'E' or 'B', RCONDE contains the reciprocal 182*> condition number for the average of the selected eigenvalues. 183*> Not referenced if SENSE = 'N' or 'V'. 184*> \endverbatim 185*> 186*> \param[out] RCONDV 187*> \verbatim 188*> RCONDV is DOUBLE PRECISION 189*> If SENSE = 'V' or 'B', RCONDV contains the reciprocal 190*> condition number for the selected right invariant subspace. 191*> Not referenced if SENSE = 'N' or 'E'. 192*> \endverbatim 193*> 194*> \param[out] WORK 195*> \verbatim 196*> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK)) 197*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK. 198*> \endverbatim 199*> 200*> \param[in] LWORK 201*> \verbatim 202*> LWORK is INTEGER 203*> The dimension of the array WORK. LWORK >= max(1,3*N). 204*> Also, if SENSE = 'E' or 'V' or 'B', 205*> LWORK >= N+2*SDIM*(N-SDIM), where SDIM is the number of 206*> selected eigenvalues computed by this routine. Note that 207*> N+2*SDIM*(N-SDIM) <= N+N*N/2. Note also that an error is only 208*> returned if LWORK < max(1,3*N), but if SENSE = 'E' or 'V' or 209*> 'B' this may not be large enough. 210*> For good performance, LWORK must generally be larger. 211*> 212*> If LWORK = -1, then a workspace query is assumed; the routine 213*> only calculates upper bounds on the optimal sizes of the 214*> arrays WORK and IWORK, returns these values as the first 215*> entries of the WORK and IWORK arrays, and no error messages 216*> related to LWORK or LIWORK are issued by XERBLA. 217*> \endverbatim 218*> 219*> \param[out] IWORK 220*> \verbatim 221*> IWORK is INTEGER array, dimension (MAX(1,LIWORK)) 222*> On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK. 223*> \endverbatim 224*> 225*> \param[in] LIWORK 226*> \verbatim 227*> LIWORK is INTEGER 228*> The dimension of the array IWORK. 229*> LIWORK >= 1; if SENSE = 'V' or 'B', LIWORK >= SDIM*(N-SDIM). 230*> Note that SDIM*(N-SDIM) <= N*N/4. Note also that an error is 231*> only returned if LIWORK < 1, but if SENSE = 'V' or 'B' this 232*> may not be large enough. 233*> 234*> If LIWORK = -1, then a workspace query is assumed; the 235*> routine only calculates upper bounds on the optimal sizes of 236*> the arrays WORK and IWORK, returns these values as the first 237*> entries of the WORK and IWORK arrays, and no error messages 238*> related to LWORK or LIWORK are issued by XERBLA. 239*> \endverbatim 240*> 241*> \param[out] BWORK 242*> \verbatim 243*> BWORK is LOGICAL array, dimension (N) 244*> Not referenced if SORT = 'N'. 245*> \endverbatim 246*> 247*> \param[out] INFO 248*> \verbatim 249*> INFO is INTEGER 250*> = 0: successful exit 251*> < 0: if INFO = -i, the i-th argument had an illegal value. 252*> > 0: if INFO = i, and i is 253*> <= N: the QR algorithm failed to compute all the 254*> eigenvalues; elements 1:ILO-1 and i+1:N of WR and WI 255*> contain those eigenvalues which have converged; if 256*> JOBVS = 'V', VS contains the transformation which 257*> reduces A to its partially converged Schur form. 258*> = N+1: the eigenvalues could not be reordered because some 259*> eigenvalues were too close to separate (the problem 260*> is very ill-conditioned); 261*> = N+2: after reordering, roundoff changed values of some 262*> complex eigenvalues so that leading eigenvalues in 263*> the Schur form no longer satisfy SELECT=.TRUE. This 264*> could also be caused by underflow due to scaling. 265*> \endverbatim 266* 267* Authors: 268* ======== 269* 270*> \author Univ. of Tennessee 271*> \author Univ. of California Berkeley 272*> \author Univ. of Colorado Denver 273*> \author NAG Ltd. 274* 275*> \date November 2011 276* 277*> \ingroup doubleGEeigen 278* 279* ===================================================================== 280 SUBROUTINE DGEESX( JOBVS, SORT, SELECT, SENSE, N, A, LDA, SDIM, 281 $ WR, WI, VS, LDVS, RCONDE, RCONDV, WORK, LWORK, 282 $ IWORK, LIWORK, BWORK, INFO ) 283* 284* -- LAPACK driver routine (version 3.4.0) -- 285* -- LAPACK is a software package provided by Univ. of Tennessee, -- 286* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 287* November 2011 288* 289* .. Scalar Arguments .. 290 CHARACTER JOBVS, SENSE, SORT 291 INTEGER INFO, LDA, LDVS, LIWORK, LWORK, N, SDIM 292 DOUBLE PRECISION RCONDE, RCONDV 293* .. 294* .. Array Arguments .. 295 LOGICAL BWORK( * ) 296 INTEGER IWORK( * ) 297 DOUBLE PRECISION A( LDA, * ), VS( LDVS, * ), WI( * ), WORK( * ), 298 $ WR( * ) 299* .. 300* .. Function Arguments .. 301 LOGICAL SELECT 302 EXTERNAL SELECT 303* .. 304* 305* ===================================================================== 306* 307* .. Parameters .. 308 DOUBLE PRECISION ZERO, ONE 309 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 ) 310* .. 311* .. Local Scalars .. 312 LOGICAL CURSL, LASTSL, LQUERY, LST2SL, SCALEA, WANTSB, 313 $ WANTSE, WANTSN, WANTST, WANTSV, WANTVS 314 INTEGER HSWORK, I, I1, I2, IBAL, ICOND, IERR, IEVAL, 315 $ IHI, ILO, INXT, IP, ITAU, IWRK, LIWRK, LWRK, 316 $ MAXWRK, MINWRK 317 DOUBLE PRECISION ANRM, BIGNUM, CSCALE, EPS, SMLNUM 318* .. 319* .. Local Arrays .. 320 DOUBLE PRECISION DUM( 1 ) 321* .. 322* .. External Subroutines .. 323 EXTERNAL DCOPY, DGEBAK, DGEBAL, DGEHRD, DHSEQR, DLACPY, 324 $ DLASCL, DORGHR, DSWAP, DTRSEN, XERBLA 325* .. 326* .. External Functions .. 327 LOGICAL LSAME 328 INTEGER ILAENV 329 DOUBLE PRECISION DLAMCH, DLANGE 330 EXTERNAL LSAME, ILAENV, DLABAD, DLAMCH, DLANGE 331* .. 332* .. Intrinsic Functions .. 333 INTRINSIC MAX, SQRT 334* .. 335* .. Executable Statements .. 336* 337* Test the input arguments 338* 339 INFO = 0 340 WANTVS = LSAME( JOBVS, 'V' ) 341 WANTST = LSAME( SORT, 'S' ) 342 WANTSN = LSAME( SENSE, 'N' ) 343 WANTSE = LSAME( SENSE, 'E' ) 344 WANTSV = LSAME( SENSE, 'V' ) 345 WANTSB = LSAME( SENSE, 'B' ) 346 LQUERY = ( LWORK.EQ.-1 .OR. LIWORK.EQ.-1 ) 347* 348 IF( ( .NOT.WANTVS ) .AND. ( .NOT.LSAME( JOBVS, 'N' ) ) ) THEN 349 INFO = -1 350 ELSE IF( ( .NOT.WANTST ) .AND. ( .NOT.LSAME( SORT, 'N' ) ) ) THEN 351 INFO = -2 352 ELSE IF( .NOT.( WANTSN .OR. WANTSE .OR. WANTSV .OR. WANTSB ) .OR. 353 $ ( .NOT.WANTST .AND. .NOT.WANTSN ) ) THEN 354 INFO = -4 355 ELSE IF( N.LT.0 ) THEN 356 INFO = -5 357 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 358 INFO = -7 359 ELSE IF( LDVS.LT.1 .OR. ( WANTVS .AND. LDVS.LT.N ) ) THEN 360 INFO = -12 361 END IF 362* 363* Compute workspace 364* (Note: Comments in the code beginning "RWorkspace:" describe the 365* minimal amount of real workspace needed at that point in the 366* code, as well as the preferred amount for good performance. 367* IWorkspace refers to integer workspace. 368* NB refers to the optimal block size for the immediately 369* following subroutine, as returned by ILAENV. 370* HSWORK refers to the workspace preferred by DHSEQR, as 371* calculated below. HSWORK is computed assuming ILO=1 and IHI=N, 372* the worst case. 373* If SENSE = 'E', 'V' or 'B', then the amount of workspace needed 374* depends on SDIM, which is computed by the routine DTRSEN later 375* in the code.) 376* 377 IF( INFO.EQ.0 ) THEN 378 LIWRK = 1 379 IF( N.EQ.0 ) THEN 380 MINWRK = 1 381 LWRK = 1 382 ELSE 383 MAXWRK = 2*N + N*ILAENV( 1, 'DGEHRD', ' ', N, 1, N, 0 ) 384 MINWRK = 3*N 385* 386 CALL DHSEQR( 'S', JOBVS, N, 1, N, A, LDA, WR, WI, VS, LDVS, 387 $ WORK, -1, IEVAL ) 388 HSWORK = WORK( 1 ) 389* 390 IF( .NOT.WANTVS ) THEN 391 MAXWRK = MAX( MAXWRK, N + HSWORK ) 392 ELSE 393 MAXWRK = MAX( MAXWRK, 2*N + ( N - 1 )*ILAENV( 1, 394 $ 'DORGHR', ' ', N, 1, N, -1 ) ) 395 MAXWRK = MAX( MAXWRK, N + HSWORK ) 396 END IF 397 LWRK = MAXWRK 398 IF( .NOT.WANTSN ) 399 $ LWRK = MAX( LWRK, N + ( N*N )/2 ) 400 IF( WANTSV .OR. WANTSB ) 401 $ LIWRK = ( N*N )/4 402 END IF 403 IWORK( 1 ) = LIWRK 404 WORK( 1 ) = LWRK 405* 406 IF( LWORK.LT.MINWRK .AND. .NOT.LQUERY ) THEN 407 INFO = -16 408 ELSE IF( LIWORK.LT.1 .AND. .NOT.LQUERY ) THEN 409 INFO = -18 410 END IF 411 END IF 412* 413 IF( INFO.NE.0 ) THEN 414 CALL XERBLA( 'DGEESX', -INFO ) 415 RETURN 416 ELSE IF( LQUERY ) THEN 417 RETURN 418 END IF 419* 420* Quick return if possible 421* 422 IF( N.EQ.0 ) THEN 423 SDIM = 0 424 RETURN 425 END IF 426* 427* Get machine constants 428* 429 EPS = DLAMCH( 'P' ) 430 SMLNUM = DLAMCH( 'S' ) 431 BIGNUM = ONE / SMLNUM 432 CALL DLABAD( SMLNUM, BIGNUM ) 433 SMLNUM = SQRT( SMLNUM ) / EPS 434 BIGNUM = ONE / SMLNUM 435* 436* Scale A if max element outside range [SMLNUM,BIGNUM] 437* 438 ANRM = DLANGE( 'M', N, N, A, LDA, DUM ) 439 SCALEA = .FALSE. 440 IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN 441 SCALEA = .TRUE. 442 CSCALE = SMLNUM 443 ELSE IF( ANRM.GT.BIGNUM ) THEN 444 SCALEA = .TRUE. 445 CSCALE = BIGNUM 446 END IF 447 IF( SCALEA ) 448 $ CALL DLASCL( 'G', 0, 0, ANRM, CSCALE, N, N, A, LDA, IERR ) 449* 450* Permute the matrix to make it more nearly triangular 451* (RWorkspace: need N) 452* 453 IBAL = 1 454 CALL DGEBAL( 'P', N, A, LDA, ILO, IHI, WORK( IBAL ), IERR ) 455* 456* Reduce to upper Hessenberg form 457* (RWorkspace: need 3*N, prefer 2*N+N*NB) 458* 459 ITAU = N + IBAL 460 IWRK = N + ITAU 461 CALL DGEHRD( N, ILO, IHI, A, LDA, WORK( ITAU ), WORK( IWRK ), 462 $ LWORK-IWRK+1, IERR ) 463* 464 IF( WANTVS ) THEN 465* 466* Copy Householder vectors to VS 467* 468 CALL DLACPY( 'L', N, N, A, LDA, VS, LDVS ) 469* 470* Generate orthogonal matrix in VS 471* (RWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) 472* 473 CALL DORGHR( N, ILO, IHI, VS, LDVS, WORK( ITAU ), WORK( IWRK ), 474 $ LWORK-IWRK+1, IERR ) 475 END IF 476* 477 SDIM = 0 478* 479* Perform QR iteration, accumulating Schur vectors in VS if desired 480* (RWorkspace: need N+1, prefer N+HSWORK (see comments) ) 481* 482 IWRK = ITAU 483 CALL DHSEQR( 'S', JOBVS, N, ILO, IHI, A, LDA, WR, WI, VS, LDVS, 484 $ WORK( IWRK ), LWORK-IWRK+1, IEVAL ) 485 IF( IEVAL.GT.0 ) 486 $ INFO = IEVAL 487* 488* Sort eigenvalues if desired 489* 490 IF( WANTST .AND. INFO.EQ.0 ) THEN 491 IF( SCALEA ) THEN 492 CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, N, 1, WR, N, IERR ) 493 CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, N, 1, WI, N, IERR ) 494 END IF 495 DO 10 I = 1, N 496 BWORK( I ) = SELECT( WR( I ), WI( I ) ) 497 10 CONTINUE 498* 499* Reorder eigenvalues, transform Schur vectors, and compute 500* reciprocal condition numbers 501* (RWorkspace: if SENSE is not 'N', need N+2*SDIM*(N-SDIM) 502* otherwise, need N ) 503* (IWorkspace: if SENSE is 'V' or 'B', need SDIM*(N-SDIM) 504* otherwise, need 0 ) 505* 506 CALL DTRSEN( SENSE, JOBVS, BWORK, N, A, LDA, VS, LDVS, WR, WI, 507 $ SDIM, RCONDE, RCONDV, WORK( IWRK ), LWORK-IWRK+1, 508 $ IWORK, LIWORK, ICOND ) 509 IF( .NOT.WANTSN ) 510 $ MAXWRK = MAX( MAXWRK, N+2*SDIM*( N-SDIM ) ) 511 IF( ICOND.EQ.-15 ) THEN 512* 513* Not enough real workspace 514* 515 INFO = -16 516 ELSE IF( ICOND.EQ.-17 ) THEN 517* 518* Not enough integer workspace 519* 520 INFO = -18 521 ELSE IF( ICOND.GT.0 ) THEN 522* 523* DTRSEN failed to reorder or to restore standard Schur form 524* 525 INFO = ICOND + N 526 END IF 527 END IF 528* 529 IF( WANTVS ) THEN 530* 531* Undo balancing 532* (RWorkspace: need N) 533* 534 CALL DGEBAK( 'P', 'R', N, ILO, IHI, WORK( IBAL ), N, VS, LDVS, 535 $ IERR ) 536 END IF 537* 538 IF( SCALEA ) THEN 539* 540* Undo scaling for the Schur form of A 541* 542 CALL DLASCL( 'H', 0, 0, CSCALE, ANRM, N, N, A, LDA, IERR ) 543 CALL DCOPY( N, A, LDA+1, WR, 1 ) 544 IF( ( WANTSV .OR. WANTSB ) .AND. INFO.EQ.0 ) THEN 545 DUM( 1 ) = RCONDV 546 CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, 1, 1, DUM, 1, IERR ) 547 RCONDV = DUM( 1 ) 548 END IF 549 IF( CSCALE.EQ.SMLNUM ) THEN 550* 551* If scaling back towards underflow, adjust WI if an 552* offdiagonal element of a 2-by-2 block in the Schur form 553* underflows. 554* 555 IF( IEVAL.GT.0 ) THEN 556 I1 = IEVAL + 1 557 I2 = IHI - 1 558 CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, ILO-1, 1, WI, N, 559 $ IERR ) 560 ELSE IF( WANTST ) THEN 561 I1 = 1 562 I2 = N - 1 563 ELSE 564 I1 = ILO 565 I2 = IHI - 1 566 END IF 567 INXT = I1 - 1 568 DO 20 I = I1, I2 569 IF( I.LT.INXT ) 570 $ GO TO 20 571 IF( WI( I ).EQ.ZERO ) THEN 572 INXT = I + 1 573 ELSE 574 IF( A( I+1, I ).EQ.ZERO ) THEN 575 WI( I ) = ZERO 576 WI( I+1 ) = ZERO 577 ELSE IF( A( I+1, I ).NE.ZERO .AND. A( I, I+1 ).EQ. 578 $ ZERO ) THEN 579 WI( I ) = ZERO 580 WI( I+1 ) = ZERO 581 IF( I.GT.1 ) 582 $ CALL DSWAP( I-1, A( 1, I ), 1, A( 1, I+1 ), 1 ) 583 IF( N.GT.I+1 ) 584 $ CALL DSWAP( N-I-1, A( I, I+2 ), LDA, 585 $ A( I+1, I+2 ), LDA ) 586 CALL DSWAP( N, VS( 1, I ), 1, VS( 1, I+1 ), 1 ) 587 A( I, I+1 ) = A( I+1, I ) 588 A( I+1, I ) = ZERO 589 END IF 590 INXT = I + 2 591 END IF 592 20 CONTINUE 593 END IF 594 CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, N-IEVAL, 1, 595 $ WI( IEVAL+1 ), MAX( N-IEVAL, 1 ), IERR ) 596 END IF 597* 598 IF( WANTST .AND. INFO.EQ.0 ) THEN 599* 600* Check if reordering successful 601* 602 LASTSL = .TRUE. 603 LST2SL = .TRUE. 604 SDIM = 0 605 IP = 0 606 DO 30 I = 1, N 607 CURSL = SELECT( WR( I ), WI( I ) ) 608 IF( WI( I ).EQ.ZERO ) THEN 609 IF( CURSL ) 610 $ SDIM = SDIM + 1 611 IP = 0 612 IF( CURSL .AND. .NOT.LASTSL ) 613 $ INFO = N + 2 614 ELSE 615 IF( IP.EQ.1 ) THEN 616* 617* Last eigenvalue of conjugate pair 618* 619 CURSL = CURSL .OR. LASTSL 620 LASTSL = CURSL 621 IF( CURSL ) 622 $ SDIM = SDIM + 2 623 IP = -1 624 IF( CURSL .AND. .NOT.LST2SL ) 625 $ INFO = N + 2 626 ELSE 627* 628* First eigenvalue of conjugate pair 629* 630 IP = 1 631 END IF 632 END IF 633 LST2SL = LASTSL 634 LASTSL = CURSL 635 30 CONTINUE 636 END IF 637* 638 WORK( 1 ) = MAXWRK 639 IF( WANTSV .OR. WANTSB ) THEN 640 IWORK( 1 ) = MAX( 1, SDIM*( N-SDIM ) ) 641 ELSE 642 IWORK( 1 ) = 1 643 END IF 644* 645 RETURN 646* 647* End of DGEESX 648* 649 END 650