Annotation of rpl/lapack/lapack/dggev3.f, revision 1.1
1.1 ! bertrand 1: *> \brief <b> DGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices (blocked algorithm)</b>
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: *> \htmlonly
! 9: *> Download DGGEV3 + dependencies
! 10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dggev3.f">
! 11: *> [TGZ]</a>
! 12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dggev3.f">
! 13: *> [ZIP]</a>
! 14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dggev3.f">
! 15: *> [TXT]</a>
! 16: *> \endhtmlonly
! 17: *
! 18: * Definition:
! 19: * ===========
! 20: *
! 21: * SUBROUTINE DGGEV3( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR,
! 22: * $ ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK,
! 23: * $ INFO )
! 24: *
! 25: * .. Scalar Arguments ..
! 26: * CHARACTER JOBVL, JOBVR
! 27: * INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
! 28: * ..
! 29: * .. Array Arguments ..
! 30: * DOUBLE PRECISION A( LDA, * ), ALPHAI( * ), ALPHAR( * ),
! 31: * $ B( LDB, * ), BETA( * ), VL( LDVL, * ),
! 32: * $ VR( LDVR, * ), WORK( * )
! 33: * ..
! 34: *
! 35: *
! 36: *> \par Purpose:
! 37: * =============
! 38: *>
! 39: *> \verbatim
! 40: *>
! 41: *> DGGEV3 computes for a pair of N-by-N real nonsymmetric matrices (A,B)
! 42: *> the generalized eigenvalues, and optionally, the left and/or right
! 43: *> generalized eigenvectors.
! 44: *>
! 45: *> A generalized eigenvalue for a pair of matrices (A,B) is a scalar
! 46: *> lambda or a ratio alpha/beta = lambda, such that A - lambda*B is
! 47: *> singular. It is usually represented as the pair (alpha,beta), as
! 48: *> there is a reasonable interpretation for beta=0, and even for both
! 49: *> being zero.
! 50: *>
! 51: *> The right eigenvector v(j) corresponding to the eigenvalue lambda(j)
! 52: *> of (A,B) satisfies
! 53: *>
! 54: *> A * v(j) = lambda(j) * B * v(j).
! 55: *>
! 56: *> The left eigenvector u(j) corresponding to the eigenvalue lambda(j)
! 57: *> of (A,B) satisfies
! 58: *>
! 59: *> u(j)**H * A = lambda(j) * u(j)**H * B .
! 60: *>
! 61: *> where u(j)**H is the conjugate-transpose of u(j).
! 62: *>
! 63: *> \endverbatim
! 64: *
! 65: * Arguments:
! 66: * ==========
! 67: *
! 68: *> \param[in] JOBVL
! 69: *> \verbatim
! 70: *> JOBVL is CHARACTER*1
! 71: *> = 'N': do not compute the left generalized eigenvectors;
! 72: *> = 'V': compute the left generalized eigenvectors.
! 73: *> \endverbatim
! 74: *>
! 75: *> \param[in] JOBVR
! 76: *> \verbatim
! 77: *> JOBVR is CHARACTER*1
! 78: *> = 'N': do not compute the right generalized eigenvectors;
! 79: *> = 'V': compute the right generalized eigenvectors.
! 80: *> \endverbatim
! 81: *>
! 82: *> \param[in] N
! 83: *> \verbatim
! 84: *> N is INTEGER
! 85: *> The order of the matrices A, B, VL, and VR. N >= 0.
! 86: *> \endverbatim
! 87: *>
! 88: *> \param[in,out] A
! 89: *> \verbatim
! 90: *> A is DOUBLE PRECISION array, dimension (LDA, N)
! 91: *> On entry, the matrix A in the pair (A,B).
! 92: *> On exit, A has been overwritten.
! 93: *> \endverbatim
! 94: *>
! 95: *> \param[in] LDA
! 96: *> \verbatim
! 97: *> LDA is INTEGER
! 98: *> The leading dimension of A. LDA >= max(1,N).
! 99: *> \endverbatim
! 100: *>
! 101: *> \param[in,out] B
! 102: *> \verbatim
! 103: *> B is DOUBLE PRECISION array, dimension (LDB, N)
! 104: *> On entry, the matrix B in the pair (A,B).
! 105: *> On exit, B has been overwritten.
! 106: *> \endverbatim
! 107: *>
! 108: *> \param[in] LDB
! 109: *> \verbatim
! 110: *> LDB is INTEGER
! 111: *> The leading dimension of B. LDB >= max(1,N).
! 112: *> \endverbatim
! 113: *>
! 114: *> \param[out] ALPHAR
! 115: *> \verbatim
! 116: *> ALPHAR is DOUBLE PRECISION array, dimension (N)
! 117: *> \endverbatim
! 118: *>
! 119: *> \param[out] ALPHAI
! 120: *> \verbatim
! 121: *> ALPHAI is DOUBLE PRECISION array, dimension (N)
! 122: *> \endverbatim
! 123: *>
! 124: *> \param[out] BETA
! 125: *> \verbatim
! 126: *> BETA is DOUBLE PRECISION array, dimension (N)
! 127: *> On exit, (ALPHAR(j) + ALPHAI(j)*i)/BETA(j), j=1,...,N, will
! 128: *> be the generalized eigenvalues. If ALPHAI(j) is zero, then
! 129: *> the j-th eigenvalue is real; if positive, then the j-th and
! 130: *> (j+1)-st eigenvalues are a complex conjugate pair, with
! 131: *> ALPHAI(j+1) negative.
! 132: *>
! 133: *> Note: the quotients ALPHAR(j)/BETA(j) and ALPHAI(j)/BETA(j)
! 134: *> may easily over- or underflow, and BETA(j) may even be zero.
! 135: *> Thus, the user should avoid naively computing the ratio
! 136: *> alpha/beta. However, ALPHAR and ALPHAI will be always less
! 137: *> than and usually comparable with norm(A) in magnitude, and
! 138: *> BETA always less than and usually comparable with norm(B).
! 139: *> \endverbatim
! 140: *>
! 141: *> \param[out] VL
! 142: *> \verbatim
! 143: *> VL is DOUBLE PRECISION array, dimension (LDVL,N)
! 144: *> If JOBVL = 'V', the left eigenvectors u(j) are stored one
! 145: *> after another in the columns of VL, in the same order as
! 146: *> their eigenvalues. If the j-th eigenvalue is real, then
! 147: *> u(j) = VL(:,j), the j-th column of VL. If the j-th and
! 148: *> (j+1)-th eigenvalues form a complex conjugate pair, then
! 149: *> u(j) = VL(:,j)+i*VL(:,j+1) and u(j+1) = VL(:,j)-i*VL(:,j+1).
! 150: *> Each eigenvector is scaled so the largest component has
! 151: *> abs(real part)+abs(imag. part)=1.
! 152: *> Not referenced if JOBVL = 'N'.
! 153: *> \endverbatim
! 154: *>
! 155: *> \param[in] LDVL
! 156: *> \verbatim
! 157: *> LDVL is INTEGER
! 158: *> The leading dimension of the matrix VL. LDVL >= 1, and
! 159: *> if JOBVL = 'V', LDVL >= N.
! 160: *> \endverbatim
! 161: *>
! 162: *> \param[out] VR
! 163: *> \verbatim
! 164: *> VR is DOUBLE PRECISION array, dimension (LDVR,N)
! 165: *> If JOBVR = 'V', the right eigenvectors v(j) are stored one
! 166: *> after another in the columns of VR, in the same order as
! 167: *> their eigenvalues. If the j-th eigenvalue is real, then
! 168: *> v(j) = VR(:,j), the j-th column of VR. If the j-th and
! 169: *> (j+1)-th eigenvalues form a complex conjugate pair, then
! 170: *> v(j) = VR(:,j)+i*VR(:,j+1) and v(j+1) = VR(:,j)-i*VR(:,j+1).
! 171: *> Each eigenvector is scaled so the largest component has
! 172: *> abs(real part)+abs(imag. part)=1.
! 173: *> Not referenced if JOBVR = 'N'.
! 174: *> \endverbatim
! 175: *>
! 176: *> \param[in] LDVR
! 177: *> \verbatim
! 178: *> LDVR is INTEGER
! 179: *> The leading dimension of the matrix VR. LDVR >= 1, and
! 180: *> if JOBVR = 'V', LDVR >= N.
! 181: *> \endverbatim
! 182: *>
! 183: *> \param[out] WORK
! 184: *> \verbatim
! 185: *> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
! 186: *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
! 187: *> \endverbatim
! 188: *>
! 189: *> \param[in] LWORK
! 190: *> \verbatim
! 191: *> LWORK is INTEGER
! 192: *>
! 193: *> If LWORK = -1, then a workspace query is assumed; the routine
! 194: *> only calculates the optimal size of the WORK array, returns
! 195: *> this value as the first entry of the WORK array, and no error
! 196: *> message related to LWORK is issued by XERBLA.
! 197: *> \endverbatim
! 198: *>
! 199: *> \param[out] INFO
! 200: *> \verbatim
! 201: *> INFO is INTEGER
! 202: *> = 0: successful exit
! 203: *> < 0: if INFO = -i, the i-th argument had an illegal value.
! 204: *> = 1,...,N:
! 205: *> The QZ iteration failed. No eigenvectors have been
! 206: *> calculated, but ALPHAR(j), ALPHAI(j), and BETA(j)
! 207: *> should be correct for j=INFO+1,...,N.
! 208: *> > N: =N+1: other than QZ iteration failed in DHGEQZ.
! 209: *> =N+2: error return from DTGEVC.
! 210: *> \endverbatim
! 211: *
! 212: * Authors:
! 213: * ========
! 214: *
! 215: *> \author Univ. of Tennessee
! 216: *> \author Univ. of California Berkeley
! 217: *> \author Univ. of Colorado Denver
! 218: *> \author NAG Ltd.
! 219: *
! 220: *> \date January 2015
! 221: *
! 222: *> \ingroup doubleGEeigen
! 223: *
! 224: * =====================================================================
! 225: SUBROUTINE DGGEV3( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR,
! 226: $ ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK,
! 227: $ INFO )
! 228: *
! 229: * -- LAPACK driver routine (version 3.6.0) --
! 230: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 231: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 232: * January 2015
! 233: *
! 234: * .. Scalar Arguments ..
! 235: CHARACTER JOBVL, JOBVR
! 236: INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
! 237: * ..
! 238: * .. Array Arguments ..
! 239: DOUBLE PRECISION A( LDA, * ), ALPHAI( * ), ALPHAR( * ),
! 240: $ B( LDB, * ), BETA( * ), VL( LDVL, * ),
! 241: $ VR( LDVR, * ), WORK( * )
! 242: * ..
! 243: *
! 244: * =====================================================================
! 245: *
! 246: * .. Parameters ..
! 247: DOUBLE PRECISION ZERO, ONE
! 248: PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
! 249: * ..
! 250: * .. Local Scalars ..
! 251: LOGICAL ILASCL, ILBSCL, ILV, ILVL, ILVR, LQUERY
! 252: CHARACTER CHTEMP
! 253: INTEGER ICOLS, IERR, IHI, IJOBVL, IJOBVR, ILEFT, ILO,
! 254: $ IN, IRIGHT, IROWS, ITAU, IWRK, JC, JR, LWKOPT
! 255: DOUBLE PRECISION ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS,
! 256: $ SMLNUM, TEMP
! 257: * ..
! 258: * .. Local Arrays ..
! 259: LOGICAL LDUMMA( 1 )
! 260: * ..
! 261: * .. External Subroutines ..
! 262: EXTERNAL DGEQRF, DGGBAK, DGGBAL, DGGHD3, DHGEQZ, DLABAD,
! 263: $ DLACPY, DLASCL, DLASET, DORGQR, DORMQR, DTGEVC,
! 264: $ XERBLA
! 265: * ..
! 266: * .. External Functions ..
! 267: LOGICAL LSAME
! 268: DOUBLE PRECISION DLAMCH, DLANGE
! 269: EXTERNAL LSAME, DLAMCH, DLANGE
! 270: * ..
! 271: * .. Intrinsic Functions ..
! 272: INTRINSIC ABS, MAX, SQRT
! 273: * ..
! 274: * .. Executable Statements ..
! 275: *
! 276: * Decode the input arguments
! 277: *
! 278: IF( LSAME( JOBVL, 'N' ) ) THEN
! 279: IJOBVL = 1
! 280: ILVL = .FALSE.
! 281: ELSE IF( LSAME( JOBVL, 'V' ) ) THEN
! 282: IJOBVL = 2
! 283: ILVL = .TRUE.
! 284: ELSE
! 285: IJOBVL = -1
! 286: ILVL = .FALSE.
! 287: END IF
! 288: *
! 289: IF( LSAME( JOBVR, 'N' ) ) THEN
! 290: IJOBVR = 1
! 291: ILVR = .FALSE.
! 292: ELSE IF( LSAME( JOBVR, 'V' ) ) THEN
! 293: IJOBVR = 2
! 294: ILVR = .TRUE.
! 295: ELSE
! 296: IJOBVR = -1
! 297: ILVR = .FALSE.
! 298: END IF
! 299: ILV = ILVL .OR. ILVR
! 300: *
! 301: * Test the input arguments
! 302: *
! 303: INFO = 0
! 304: LQUERY = ( LWORK.EQ.-1 )
! 305: IF( IJOBVL.LE.0 ) THEN
! 306: INFO = -1
! 307: ELSE IF( IJOBVR.LE.0 ) THEN
! 308: INFO = -2
! 309: ELSE IF( N.LT.0 ) THEN
! 310: INFO = -3
! 311: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
! 312: INFO = -5
! 313: ELSE IF( LDB.LT.MAX( 1, N ) ) THEN
! 314: INFO = -7
! 315: ELSE IF( LDVL.LT.1 .OR. ( ILVL .AND. LDVL.LT.N ) ) THEN
! 316: INFO = -12
! 317: ELSE IF( LDVR.LT.1 .OR. ( ILVR .AND. LDVR.LT.N ) ) THEN
! 318: INFO = -14
! 319: ELSE IF( LWORK.LT.MAX( 1, 8*N ) .AND. .NOT.LQUERY ) THEN
! 320: INFO = -16
! 321: END IF
! 322: *
! 323: * Compute workspace
! 324: *
! 325: IF( INFO.EQ.0 ) THEN
! 326: CALL DGEQRF( N, N, B, LDB, WORK, WORK, -1, IERR )
! 327: LWKOPT = MAX(1, 8*N, 3*N+INT( WORK( 1 ) ) )
! 328: CALL DORMQR( 'L', 'T', N, N, N, B, LDB, WORK, A, LDA, WORK, -1,
! 329: $ IERR )
! 330: LWKOPT = MAX( LWKOPT, 3*N+INT( WORK ( 1 ) ) )
! 331: IF( ILVL ) THEN
! 332: CALL DORGQR( N, N, N, VL, LDVL, WORK, WORK, -1, IERR )
! 333: LWKOPT = MAX( LWKOPT, 3*N+INT( WORK ( 1 ) ) )
! 334: END IF
! 335: IF( ILV ) THEN
! 336: CALL DGGHD3( JOBVL, JOBVR, N, 1, N, A, LDA, B, LDB, VL,
! 337: $ LDVL, VR, LDVR, WORK, -1, IERR )
! 338: LWKOPT = MAX( LWKOPT, 3*N+INT( WORK ( 1 ) ) )
! 339: CALL DHGEQZ( 'S', JOBVL, JOBVR, N, 1, N, A, LDA, B, LDB,
! 340: $ ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR,
! 341: $ WORK, -1, IERR )
! 342: LWKOPT = MAX( LWKOPT, 2*N+INT( WORK ( 1 ) ) )
! 343: ELSE
! 344: CALL DGGHD3( 'N', 'N', N, 1, N, A, LDA, B, LDB, VL, LDVL,
! 345: $ VR, LDVR, WORK, -1, IERR )
! 346: LWKOPT = MAX( LWKOPT, 3*N+INT( WORK ( 1 ) ) )
! 347: CALL DHGEQZ( 'E', JOBVL, JOBVR, N, 1, N, A, LDA, B, LDB,
! 348: $ ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR,
! 349: $ WORK, -1, IERR )
! 350: LWKOPT = MAX( LWKOPT, 2*N+INT( WORK ( 1 ) ) )
! 351: END IF
! 352:
! 353: WORK( 1 ) = LWKOPT
! 354: END IF
! 355: *
! 356: IF( INFO.NE.0 ) THEN
! 357: CALL XERBLA( 'DGGEV3 ', -INFO )
! 358: RETURN
! 359: ELSE IF( LQUERY ) THEN
! 360: RETURN
! 361: END IF
! 362: *
! 363: * Quick return if possible
! 364: *
! 365: IF( N.EQ.0 )
! 366: $ RETURN
! 367: *
! 368: * Get machine constants
! 369: *
! 370: EPS = DLAMCH( 'P' )
! 371: SMLNUM = DLAMCH( 'S' )
! 372: BIGNUM = ONE / SMLNUM
! 373: CALL DLABAD( SMLNUM, BIGNUM )
! 374: SMLNUM = SQRT( SMLNUM ) / EPS
! 375: BIGNUM = ONE / SMLNUM
! 376: *
! 377: * Scale A if max element outside range [SMLNUM,BIGNUM]
! 378: *
! 379: ANRM = DLANGE( 'M', N, N, A, LDA, WORK )
! 380: ILASCL = .FALSE.
! 381: IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
! 382: ANRMTO = SMLNUM
! 383: ILASCL = .TRUE.
! 384: ELSE IF( ANRM.GT.BIGNUM ) THEN
! 385: ANRMTO = BIGNUM
! 386: ILASCL = .TRUE.
! 387: END IF
! 388: IF( ILASCL )
! 389: $ CALL DLASCL( 'G', 0, 0, ANRM, ANRMTO, N, N, A, LDA, IERR )
! 390: *
! 391: * Scale B if max element outside range [SMLNUM,BIGNUM]
! 392: *
! 393: BNRM = DLANGE( 'M', N, N, B, LDB, WORK )
! 394: ILBSCL = .FALSE.
! 395: IF( BNRM.GT.ZERO .AND. BNRM.LT.SMLNUM ) THEN
! 396: BNRMTO = SMLNUM
! 397: ILBSCL = .TRUE.
! 398: ELSE IF( BNRM.GT.BIGNUM ) THEN
! 399: BNRMTO = BIGNUM
! 400: ILBSCL = .TRUE.
! 401: END IF
! 402: IF( ILBSCL )
! 403: $ CALL DLASCL( 'G', 0, 0, BNRM, BNRMTO, N, N, B, LDB, IERR )
! 404: *
! 405: * Permute the matrices A, B to isolate eigenvalues if possible
! 406: *
! 407: ILEFT = 1
! 408: IRIGHT = N + 1
! 409: IWRK = IRIGHT + N
! 410: CALL DGGBAL( 'P', N, A, LDA, B, LDB, ILO, IHI, WORK( ILEFT ),
! 411: $ WORK( IRIGHT ), WORK( IWRK ), IERR )
! 412: *
! 413: * Reduce B to triangular form (QR decomposition of B)
! 414: *
! 415: IROWS = IHI + 1 - ILO
! 416: IF( ILV ) THEN
! 417: ICOLS = N + 1 - ILO
! 418: ELSE
! 419: ICOLS = IROWS
! 420: END IF
! 421: ITAU = IWRK
! 422: IWRK = ITAU + IROWS
! 423: CALL DGEQRF( IROWS, ICOLS, B( ILO, ILO ), LDB, WORK( ITAU ),
! 424: $ WORK( IWRK ), LWORK+1-IWRK, IERR )
! 425: *
! 426: * Apply the orthogonal transformation to matrix A
! 427: *
! 428: CALL DORMQR( 'L', 'T', IROWS, ICOLS, IROWS, B( ILO, ILO ), LDB,
! 429: $ WORK( ITAU ), A( ILO, ILO ), LDA, WORK( IWRK ),
! 430: $ LWORK+1-IWRK, IERR )
! 431: *
! 432: * Initialize VL
! 433: *
! 434: IF( ILVL ) THEN
! 435: CALL DLASET( 'Full', N, N, ZERO, ONE, VL, LDVL )
! 436: IF( IROWS.GT.1 ) THEN
! 437: CALL DLACPY( 'L', IROWS-1, IROWS-1, B( ILO+1, ILO ), LDB,
! 438: $ VL( ILO+1, ILO ), LDVL )
! 439: END IF
! 440: CALL DORGQR( IROWS, IROWS, IROWS, VL( ILO, ILO ), LDVL,
! 441: $ WORK( ITAU ), WORK( IWRK ), LWORK+1-IWRK, IERR )
! 442: END IF
! 443: *
! 444: * Initialize VR
! 445: *
! 446: IF( ILVR )
! 447: $ CALL DLASET( 'Full', N, N, ZERO, ONE, VR, LDVR )
! 448: *
! 449: * Reduce to generalized Hessenberg form
! 450: *
! 451: IF( ILV ) THEN
! 452: *
! 453: * Eigenvectors requested -- work on whole matrix.
! 454: *
! 455: CALL DGGHD3( JOBVL, JOBVR, N, ILO, IHI, A, LDA, B, LDB, VL,
! 456: $ LDVL, VR, LDVR, WORK( IWRK ), LWORK+1-IWRK, IERR )
! 457: ELSE
! 458: CALL DGGHD3( 'N', 'N', IROWS, 1, IROWS, A( ILO, ILO ), LDA,
! 459: $ B( ILO, ILO ), LDB, VL, LDVL, VR, LDVR,
! 460: $ WORK( IWRK ), LWORK+1-IWRK, IERR )
! 461: END IF
! 462: *
! 463: * Perform QZ algorithm (Compute eigenvalues, and optionally, the
! 464: * Schur forms and Schur vectors)
! 465: *
! 466: IWRK = ITAU
! 467: IF( ILV ) THEN
! 468: CHTEMP = 'S'
! 469: ELSE
! 470: CHTEMP = 'E'
! 471: END IF
! 472: CALL DHGEQZ( CHTEMP, JOBVL, JOBVR, N, ILO, IHI, A, LDA, B, LDB,
! 473: $ ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR,
! 474: $ WORK( IWRK ), LWORK+1-IWRK, IERR )
! 475: IF( IERR.NE.0 ) THEN
! 476: IF( IERR.GT.0 .AND. IERR.LE.N ) THEN
! 477: INFO = IERR
! 478: ELSE IF( IERR.GT.N .AND. IERR.LE.2*N ) THEN
! 479: INFO = IERR - N
! 480: ELSE
! 481: INFO = N + 1
! 482: END IF
! 483: GO TO 110
! 484: END IF
! 485: *
! 486: * Compute Eigenvectors
! 487: *
! 488: IF( ILV ) THEN
! 489: IF( ILVL ) THEN
! 490: IF( ILVR ) THEN
! 491: CHTEMP = 'B'
! 492: ELSE
! 493: CHTEMP = 'L'
! 494: END IF
! 495: ELSE
! 496: CHTEMP = 'R'
! 497: END IF
! 498: CALL DTGEVC( CHTEMP, 'B', LDUMMA, N, A, LDA, B, LDB, VL, LDVL,
! 499: $ VR, LDVR, N, IN, WORK( IWRK ), IERR )
! 500: IF( IERR.NE.0 ) THEN
! 501: INFO = N + 2
! 502: GO TO 110
! 503: END IF
! 504: *
! 505: * Undo balancing on VL and VR and normalization
! 506: *
! 507: IF( ILVL ) THEN
! 508: CALL DGGBAK( 'P', 'L', N, ILO, IHI, WORK( ILEFT ),
! 509: $ WORK( IRIGHT ), N, VL, LDVL, IERR )
! 510: DO 50 JC = 1, N
! 511: IF( ALPHAI( JC ).LT.ZERO )
! 512: $ GO TO 50
! 513: TEMP = ZERO
! 514: IF( ALPHAI( JC ).EQ.ZERO ) THEN
! 515: DO 10 JR = 1, N
! 516: TEMP = MAX( TEMP, ABS( VL( JR, JC ) ) )
! 517: 10 CONTINUE
! 518: ELSE
! 519: DO 20 JR = 1, N
! 520: TEMP = MAX( TEMP, ABS( VL( JR, JC ) )+
! 521: $ ABS( VL( JR, JC+1 ) ) )
! 522: 20 CONTINUE
! 523: END IF
! 524: IF( TEMP.LT.SMLNUM )
! 525: $ GO TO 50
! 526: TEMP = ONE / TEMP
! 527: IF( ALPHAI( JC ).EQ.ZERO ) THEN
! 528: DO 30 JR = 1, N
! 529: VL( JR, JC ) = VL( JR, JC )*TEMP
! 530: 30 CONTINUE
! 531: ELSE
! 532: DO 40 JR = 1, N
! 533: VL( JR, JC ) = VL( JR, JC )*TEMP
! 534: VL( JR, JC+1 ) = VL( JR, JC+1 )*TEMP
! 535: 40 CONTINUE
! 536: END IF
! 537: 50 CONTINUE
! 538: END IF
! 539: IF( ILVR ) THEN
! 540: CALL DGGBAK( 'P', 'R', N, ILO, IHI, WORK( ILEFT ),
! 541: $ WORK( IRIGHT ), N, VR, LDVR, IERR )
! 542: DO 100 JC = 1, N
! 543: IF( ALPHAI( JC ).LT.ZERO )
! 544: $ GO TO 100
! 545: TEMP = ZERO
! 546: IF( ALPHAI( JC ).EQ.ZERO ) THEN
! 547: DO 60 JR = 1, N
! 548: TEMP = MAX( TEMP, ABS( VR( JR, JC ) ) )
! 549: 60 CONTINUE
! 550: ELSE
! 551: DO 70 JR = 1, N
! 552: TEMP = MAX( TEMP, ABS( VR( JR, JC ) )+
! 553: $ ABS( VR( JR, JC+1 ) ) )
! 554: 70 CONTINUE
! 555: END IF
! 556: IF( TEMP.LT.SMLNUM )
! 557: $ GO TO 100
! 558: TEMP = ONE / TEMP
! 559: IF( ALPHAI( JC ).EQ.ZERO ) THEN
! 560: DO 80 JR = 1, N
! 561: VR( JR, JC ) = VR( JR, JC )*TEMP
! 562: 80 CONTINUE
! 563: ELSE
! 564: DO 90 JR = 1, N
! 565: VR( JR, JC ) = VR( JR, JC )*TEMP
! 566: VR( JR, JC+1 ) = VR( JR, JC+1 )*TEMP
! 567: 90 CONTINUE
! 568: END IF
! 569: 100 CONTINUE
! 570: END IF
! 571: *
! 572: * End of eigenvector calculation
! 573: *
! 574: END IF
! 575: *
! 576: * Undo scaling if necessary
! 577: *
! 578: 110 CONTINUE
! 579: *
! 580: IF( ILASCL ) THEN
! 581: CALL DLASCL( 'G', 0, 0, ANRMTO, ANRM, N, 1, ALPHAR, N, IERR )
! 582: CALL DLASCL( 'G', 0, 0, ANRMTO, ANRM, N, 1, ALPHAI, N, IERR )
! 583: END IF
! 584: *
! 585: IF( ILBSCL ) THEN
! 586: CALL DLASCL( 'G', 0, 0, BNRMTO, BNRM, N, 1, BETA, N, IERR )
! 587: END IF
! 588: *
! 589: WORK( 1 ) = LWKOPT
! 590: RETURN
! 591: *
! 592: * End of DGGEV3
! 593: *
! 594: END
CVSweb interface <joel.bertrand@systella.fr>