Annotation of rpl/lapack/lapack/dgeevx.f, revision 1.1
1.1 ! bertrand 1: SUBROUTINE DGEEVX( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, WR, WI,
! 2: $ VL, LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM,
! 3: $ RCONDE, RCONDV, WORK, LWORK, IWORK, INFO )
! 4: *
! 5: * -- LAPACK driver routine (version 3.2) --
! 6: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 7: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 8: * November 2006
! 9: *
! 10: * .. Scalar Arguments ..
! 11: CHARACTER BALANC, JOBVL, JOBVR, SENSE
! 12: INTEGER IHI, ILO, INFO, LDA, LDVL, LDVR, LWORK, N
! 13: DOUBLE PRECISION ABNRM
! 14: * ..
! 15: * .. Array Arguments ..
! 16: INTEGER IWORK( * )
! 17: DOUBLE PRECISION A( LDA, * ), RCONDE( * ), RCONDV( * ),
! 18: $ SCALE( * ), VL( LDVL, * ), VR( LDVR, * ),
! 19: $ WI( * ), WORK( * ), WR( * )
! 20: * ..
! 21: *
! 22: * Purpose
! 23: * =======
! 24: *
! 25: * DGEEVX computes for an N-by-N real nonsymmetric matrix A, the
! 26: * eigenvalues and, optionally, the left and/or right eigenvectors.
! 27: *
! 28: * Optionally also, it computes a balancing transformation to improve
! 29: * the conditioning of the eigenvalues and eigenvectors (ILO, IHI,
! 30: * SCALE, and ABNRM), reciprocal condition numbers for the eigenvalues
! 31: * (RCONDE), and reciprocal condition numbers for the right
! 32: * eigenvectors (RCONDV).
! 33: *
! 34: * The right eigenvector v(j) of A satisfies
! 35: * A * v(j) = lambda(j) * v(j)
! 36: * where lambda(j) is its eigenvalue.
! 37: * The left eigenvector u(j) of A satisfies
! 38: * u(j)**H * A = lambda(j) * u(j)**H
! 39: * where u(j)**H denotes the conjugate transpose of u(j).
! 40: *
! 41: * The computed eigenvectors are normalized to have Euclidean norm
! 42: * equal to 1 and largest component real.
! 43: *
! 44: * Balancing a matrix means permuting the rows and columns to make it
! 45: * more nearly upper triangular, and applying a diagonal similarity
! 46: * transformation D * A * D**(-1), where D is a diagonal matrix, to
! 47: * make its rows and columns closer in norm and the condition numbers
! 48: * of its eigenvalues and eigenvectors smaller. The computed
! 49: * reciprocal condition numbers correspond to the balanced matrix.
! 50: * Permuting rows and columns will not change the condition numbers
! 51: * (in exact arithmetic) but diagonal scaling will. For further
! 52: * explanation of balancing, see section 4.10.2 of the LAPACK
! 53: * Users' Guide.
! 54: *
! 55: * Arguments
! 56: * =========
! 57: *
! 58: * BALANC (input) CHARACTER*1
! 59: * Indicates how the input matrix should be diagonally scaled
! 60: * and/or permuted to improve the conditioning of its
! 61: * eigenvalues.
! 62: * = 'N': Do not diagonally scale or permute;
! 63: * = 'P': Perform permutations to make the matrix more nearly
! 64: * upper triangular. Do not diagonally scale;
! 65: * = 'S': Diagonally scale the matrix, i.e. replace A by
! 66: * D*A*D**(-1), where D is a diagonal matrix chosen
! 67: * to make the rows and columns of A more equal in
! 68: * norm. Do not permute;
! 69: * = 'B': Both diagonally scale and permute A.
! 70: *
! 71: * Computed reciprocal condition numbers will be for the matrix
! 72: * after balancing and/or permuting. Permuting does not change
! 73: * condition numbers (in exact arithmetic), but balancing does.
! 74: *
! 75: * JOBVL (input) CHARACTER*1
! 76: * = 'N': left eigenvectors of A are not computed;
! 77: * = 'V': left eigenvectors of A are computed.
! 78: * If SENSE = 'E' or 'B', JOBVL must = 'V'.
! 79: *
! 80: * JOBVR (input) CHARACTER*1
! 81: * = 'N': right eigenvectors of A are not computed;
! 82: * = 'V': right eigenvectors of A are computed.
! 83: * If SENSE = 'E' or 'B', JOBVR must = 'V'.
! 84: *
! 85: * SENSE (input) CHARACTER*1
! 86: * Determines which reciprocal condition numbers are computed.
! 87: * = 'N': None are computed;
! 88: * = 'E': Computed for eigenvalues only;
! 89: * = 'V': Computed for right eigenvectors only;
! 90: * = 'B': Computed for eigenvalues and right eigenvectors.
! 91: *
! 92: * If SENSE = 'E' or 'B', both left and right eigenvectors
! 93: * must also be computed (JOBVL = 'V' and JOBVR = 'V').
! 94: *
! 95: * N (input) INTEGER
! 96: * The order of the matrix A. N >= 0.
! 97: *
! 98: * A (input/output) DOUBLE PRECISION array, dimension (LDA,N)
! 99: * On entry, the N-by-N matrix A.
! 100: * On exit, A has been overwritten. If JOBVL = 'V' or
! 101: * JOBVR = 'V', A contains the real Schur form of the balanced
! 102: * version of the input matrix A.
! 103: *
! 104: * LDA (input) INTEGER
! 105: * The leading dimension of the array A. LDA >= max(1,N).
! 106: *
! 107: * WR (output) DOUBLE PRECISION array, dimension (N)
! 108: * WI (output) DOUBLE PRECISION array, dimension (N)
! 109: * WR and WI contain the real and imaginary parts,
! 110: * respectively, of the computed eigenvalues. Complex
! 111: * conjugate pairs of eigenvalues will appear consecutively
! 112: * with the eigenvalue having the positive imaginary part
! 113: * first.
! 114: *
! 115: * VL (output) DOUBLE PRECISION array, dimension (LDVL,N)
! 116: * If JOBVL = 'V', the left eigenvectors u(j) are stored one
! 117: * after another in the columns of VL, in the same order
! 118: * as their eigenvalues.
! 119: * If JOBVL = 'N', VL is not referenced.
! 120: * If the j-th eigenvalue is real, then u(j) = VL(:,j),
! 121: * the j-th column of VL.
! 122: * If the j-th and (j+1)-st eigenvalues form a complex
! 123: * conjugate pair, then u(j) = VL(:,j) + i*VL(:,j+1) and
! 124: * u(j+1) = VL(:,j) - i*VL(:,j+1).
! 125: *
! 126: * LDVL (input) INTEGER
! 127: * The leading dimension of the array VL. LDVL >= 1; if
! 128: * JOBVL = 'V', LDVL >= N.
! 129: *
! 130: * VR (output) DOUBLE PRECISION array, dimension (LDVR,N)
! 131: * If JOBVR = 'V', the right eigenvectors v(j) are stored one
! 132: * after another in the columns of VR, in the same order
! 133: * as their eigenvalues.
! 134: * If JOBVR = 'N', VR is not referenced.
! 135: * If the j-th eigenvalue is real, then v(j) = VR(:,j),
! 136: * the j-th column of VR.
! 137: * If the j-th and (j+1)-st eigenvalues form a complex
! 138: * conjugate pair, then v(j) = VR(:,j) + i*VR(:,j+1) and
! 139: * v(j+1) = VR(:,j) - i*VR(:,j+1).
! 140: *
! 141: * LDVR (input) INTEGER
! 142: * The leading dimension of the array VR. LDVR >= 1, and if
! 143: * JOBVR = 'V', LDVR >= N.
! 144: *
! 145: * ILO (output) INTEGER
! 146: * IHI (output) INTEGER
! 147: * ILO and IHI are integer values determined when A was
! 148: * balanced. The balanced A(i,j) = 0 if I > J and
! 149: * J = 1,...,ILO-1 or I = IHI+1,...,N.
! 150: *
! 151: * SCALE (output) DOUBLE PRECISION array, dimension (N)
! 152: * Details of the permutations and scaling factors applied
! 153: * when balancing A. If P(j) is the index of the row and column
! 154: * interchanged with row and column j, and D(j) is the scaling
! 155: * factor applied to row and column j, then
! 156: * SCALE(J) = P(J), for J = 1,...,ILO-1
! 157: * = D(J), for J = ILO,...,IHI
! 158: * = P(J) for J = IHI+1,...,N.
! 159: * The order in which the interchanges are made is N to IHI+1,
! 160: * then 1 to ILO-1.
! 161: *
! 162: * ABNRM (output) DOUBLE PRECISION
! 163: * The one-norm of the balanced matrix (the maximum
! 164: * of the sum of absolute values of elements of any column).
! 165: *
! 166: * RCONDE (output) DOUBLE PRECISION array, dimension (N)
! 167: * RCONDE(j) is the reciprocal condition number of the j-th
! 168: * eigenvalue.
! 169: *
! 170: * RCONDV (output) DOUBLE PRECISION array, dimension (N)
! 171: * RCONDV(j) is the reciprocal condition number of the j-th
! 172: * right eigenvector.
! 173: *
! 174: * WORK (workspace/output) DOUBLE PRECISION array, dimension (MAX(1,LWORK))
! 175: * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
! 176: *
! 177: * LWORK (input) INTEGER
! 178: * The dimension of the array WORK. If SENSE = 'N' or 'E',
! 179: * LWORK >= max(1,2*N), and if JOBVL = 'V' or JOBVR = 'V',
! 180: * LWORK >= 3*N. If SENSE = 'V' or 'B', LWORK >= N*(N+6).
! 181: * For good performance, LWORK must generally be larger.
! 182: *
! 183: * If LWORK = -1, then a workspace query is assumed; the routine
! 184: * only calculates the optimal size of the WORK array, returns
! 185: * this value as the first entry of the WORK array, and no error
! 186: * message related to LWORK is issued by XERBLA.
! 187: *
! 188: * IWORK (workspace) INTEGER array, dimension (2*N-2)
! 189: * If SENSE = 'N' or 'E', not referenced.
! 190: *
! 191: * INFO (output) INTEGER
! 192: * = 0: successful exit
! 193: * < 0: if INFO = -i, the i-th argument had an illegal value.
! 194: * > 0: if INFO = i, the QR algorithm failed to compute all the
! 195: * eigenvalues, and no eigenvectors or condition numbers
! 196: * have been computed; elements 1:ILO-1 and i+1:N of WR
! 197: * and WI contain eigenvalues which have converged.
! 198: *
! 199: * =====================================================================
! 200: *
! 201: * .. Parameters ..
! 202: DOUBLE PRECISION ZERO, ONE
! 203: PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
! 204: * ..
! 205: * .. Local Scalars ..
! 206: LOGICAL LQUERY, SCALEA, WANTVL, WANTVR, WNTSNB, WNTSNE,
! 207: $ WNTSNN, WNTSNV
! 208: CHARACTER JOB, SIDE
! 209: INTEGER HSWORK, I, ICOND, IERR, ITAU, IWRK, K, MAXWRK,
! 210: $ MINWRK, NOUT
! 211: DOUBLE PRECISION ANRM, BIGNUM, CS, CSCALE, EPS, R, SCL, SMLNUM,
! 212: $ SN
! 213: * ..
! 214: * .. Local Arrays ..
! 215: LOGICAL SELECT( 1 )
! 216: DOUBLE PRECISION DUM( 1 )
! 217: * ..
! 218: * .. External Subroutines ..
! 219: EXTERNAL DGEBAK, DGEBAL, DGEHRD, DHSEQR, DLABAD, DLACPY,
! 220: $ DLARTG, DLASCL, DORGHR, DROT, DSCAL, DTREVC,
! 221: $ DTRSNA, XERBLA
! 222: * ..
! 223: * .. External Functions ..
! 224: LOGICAL LSAME
! 225: INTEGER IDAMAX, ILAENV
! 226: DOUBLE PRECISION DLAMCH, DLANGE, DLAPY2, DNRM2
! 227: EXTERNAL LSAME, IDAMAX, ILAENV, DLAMCH, DLANGE, DLAPY2,
! 228: $ DNRM2
! 229: * ..
! 230: * .. Intrinsic Functions ..
! 231: INTRINSIC MAX, SQRT
! 232: * ..
! 233: * .. Executable Statements ..
! 234: *
! 235: * Test the input arguments
! 236: *
! 237: INFO = 0
! 238: LQUERY = ( LWORK.EQ.-1 )
! 239: WANTVL = LSAME( JOBVL, 'V' )
! 240: WANTVR = LSAME( JOBVR, 'V' )
! 241: WNTSNN = LSAME( SENSE, 'N' )
! 242: WNTSNE = LSAME( SENSE, 'E' )
! 243: WNTSNV = LSAME( SENSE, 'V' )
! 244: WNTSNB = LSAME( SENSE, 'B' )
! 245: IF( .NOT.( LSAME( BALANC, 'N' ) .OR. LSAME( BALANC,
! 246: $ 'S' ) .OR. LSAME( BALANC, 'P' ) .OR. LSAME( BALANC, 'B' ) ) )
! 247: $ THEN
! 248: INFO = -1
! 249: ELSE IF( ( .NOT.WANTVL ) .AND. ( .NOT.LSAME( JOBVL, 'N' ) ) ) THEN
! 250: INFO = -2
! 251: ELSE IF( ( .NOT.WANTVR ) .AND. ( .NOT.LSAME( JOBVR, 'N' ) ) ) THEN
! 252: INFO = -3
! 253: ELSE IF( .NOT.( WNTSNN .OR. WNTSNE .OR. WNTSNB .OR. WNTSNV ) .OR.
! 254: $ ( ( WNTSNE .OR. WNTSNB ) .AND. .NOT.( WANTVL .AND.
! 255: $ WANTVR ) ) ) THEN
! 256: INFO = -4
! 257: ELSE IF( N.LT.0 ) THEN
! 258: INFO = -5
! 259: ELSE IF( LDA.LT.MAX( 1, N ) ) THEN
! 260: INFO = -7
! 261: ELSE IF( LDVL.LT.1 .OR. ( WANTVL .AND. LDVL.LT.N ) ) THEN
! 262: INFO = -11
! 263: ELSE IF( LDVR.LT.1 .OR. ( WANTVR .AND. LDVR.LT.N ) ) THEN
! 264: INFO = -13
! 265: END IF
! 266: *
! 267: * Compute workspace
! 268: * (Note: Comments in the code beginning "Workspace:" describe the
! 269: * minimal amount of workspace needed at that point in the code,
! 270: * as well as the preferred amount for good performance.
! 271: * NB refers to the optimal block size for the immediately
! 272: * following subroutine, as returned by ILAENV.
! 273: * HSWORK refers to the workspace preferred by DHSEQR, as
! 274: * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
! 275: * the worst case.)
! 276: *
! 277: IF( INFO.EQ.0 ) THEN
! 278: IF( N.EQ.0 ) THEN
! 279: MINWRK = 1
! 280: MAXWRK = 1
! 281: ELSE
! 282: MAXWRK = N + N*ILAENV( 1, 'DGEHRD', ' ', N, 1, N, 0 )
! 283: *
! 284: IF( WANTVL ) THEN
! 285: CALL DHSEQR( 'S', 'V', N, 1, N, A, LDA, WR, WI, VL, LDVL,
! 286: $ WORK, -1, INFO )
! 287: ELSE IF( WANTVR ) THEN
! 288: CALL DHSEQR( 'S', 'V', N, 1, N, A, LDA, WR, WI, VR, LDVR,
! 289: $ WORK, -1, INFO )
! 290: ELSE
! 291: IF( WNTSNN ) THEN
! 292: CALL DHSEQR( 'E', 'N', N, 1, N, A, LDA, WR, WI, VR,
! 293: $ LDVR, WORK, -1, INFO )
! 294: ELSE
! 295: CALL DHSEQR( 'S', 'N', N, 1, N, A, LDA, WR, WI, VR,
! 296: $ LDVR, WORK, -1, INFO )
! 297: END IF
! 298: END IF
! 299: HSWORK = WORK( 1 )
! 300: *
! 301: IF( ( .NOT.WANTVL ) .AND. ( .NOT.WANTVR ) ) THEN
! 302: MINWRK = 2*N
! 303: IF( .NOT.WNTSNN )
! 304: $ MINWRK = MAX( MINWRK, N*N+6*N )
! 305: MAXWRK = MAX( MAXWRK, HSWORK )
! 306: IF( .NOT.WNTSNN )
! 307: $ MAXWRK = MAX( MAXWRK, N*N + 6*N )
! 308: ELSE
! 309: MINWRK = 3*N
! 310: IF( ( .NOT.WNTSNN ) .AND. ( .NOT.WNTSNE ) )
! 311: $ MINWRK = MAX( MINWRK, N*N + 6*N )
! 312: MAXWRK = MAX( MAXWRK, HSWORK )
! 313: MAXWRK = MAX( MAXWRK, N + ( N - 1 )*ILAENV( 1, 'DORGHR',
! 314: $ ' ', N, 1, N, -1 ) )
! 315: IF( ( .NOT.WNTSNN ) .AND. ( .NOT.WNTSNE ) )
! 316: $ MAXWRK = MAX( MAXWRK, N*N + 6*N )
! 317: MAXWRK = MAX( MAXWRK, 3*N )
! 318: END IF
! 319: MAXWRK = MAX( MAXWRK, MINWRK )
! 320: END IF
! 321: WORK( 1 ) = MAXWRK
! 322: *
! 323: IF( LWORK.LT.MINWRK .AND. .NOT.LQUERY ) THEN
! 324: INFO = -21
! 325: END IF
! 326: END IF
! 327: *
! 328: IF( INFO.NE.0 ) THEN
! 329: CALL XERBLA( 'DGEEVX', -INFO )
! 330: RETURN
! 331: ELSE IF( LQUERY ) THEN
! 332: RETURN
! 333: END IF
! 334: *
! 335: * Quick return if possible
! 336: *
! 337: IF( N.EQ.0 )
! 338: $ RETURN
! 339: *
! 340: * Get machine constants
! 341: *
! 342: EPS = DLAMCH( 'P' )
! 343: SMLNUM = DLAMCH( 'S' )
! 344: BIGNUM = ONE / SMLNUM
! 345: CALL DLABAD( SMLNUM, BIGNUM )
! 346: SMLNUM = SQRT( SMLNUM ) / EPS
! 347: BIGNUM = ONE / SMLNUM
! 348: *
! 349: * Scale A if max element outside range [SMLNUM,BIGNUM]
! 350: *
! 351: ICOND = 0
! 352: ANRM = DLANGE( 'M', N, N, A, LDA, DUM )
! 353: SCALEA = .FALSE.
! 354: IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN
! 355: SCALEA = .TRUE.
! 356: CSCALE = SMLNUM
! 357: ELSE IF( ANRM.GT.BIGNUM ) THEN
! 358: SCALEA = .TRUE.
! 359: CSCALE = BIGNUM
! 360: END IF
! 361: IF( SCALEA )
! 362: $ CALL DLASCL( 'G', 0, 0, ANRM, CSCALE, N, N, A, LDA, IERR )
! 363: *
! 364: * Balance the matrix and compute ABNRM
! 365: *
! 366: CALL DGEBAL( BALANC, N, A, LDA, ILO, IHI, SCALE, IERR )
! 367: ABNRM = DLANGE( '1', N, N, A, LDA, DUM )
! 368: IF( SCALEA ) THEN
! 369: DUM( 1 ) = ABNRM
! 370: CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, 1, 1, DUM, 1, IERR )
! 371: ABNRM = DUM( 1 )
! 372: END IF
! 373: *
! 374: * Reduce to upper Hessenberg form
! 375: * (Workspace: need 2*N, prefer N+N*NB)
! 376: *
! 377: ITAU = 1
! 378: IWRK = ITAU + N
! 379: CALL DGEHRD( N, ILO, IHI, A, LDA, WORK( ITAU ), WORK( IWRK ),
! 380: $ LWORK-IWRK+1, IERR )
! 381: *
! 382: IF( WANTVL ) THEN
! 383: *
! 384: * Want left eigenvectors
! 385: * Copy Householder vectors to VL
! 386: *
! 387: SIDE = 'L'
! 388: CALL DLACPY( 'L', N, N, A, LDA, VL, LDVL )
! 389: *
! 390: * Generate orthogonal matrix in VL
! 391: * (Workspace: need 2*N-1, prefer N+(N-1)*NB)
! 392: *
! 393: CALL DORGHR( N, ILO, IHI, VL, LDVL, WORK( ITAU ), WORK( IWRK ),
! 394: $ LWORK-IWRK+1, IERR )
! 395: *
! 396: * Perform QR iteration, accumulating Schur vectors in VL
! 397: * (Workspace: need 1, prefer HSWORK (see comments) )
! 398: *
! 399: IWRK = ITAU
! 400: CALL DHSEQR( 'S', 'V', N, ILO, IHI, A, LDA, WR, WI, VL, LDVL,
! 401: $ WORK( IWRK ), LWORK-IWRK+1, INFO )
! 402: *
! 403: IF( WANTVR ) THEN
! 404: *
! 405: * Want left and right eigenvectors
! 406: * Copy Schur vectors to VR
! 407: *
! 408: SIDE = 'B'
! 409: CALL DLACPY( 'F', N, N, VL, LDVL, VR, LDVR )
! 410: END IF
! 411: *
! 412: ELSE IF( WANTVR ) THEN
! 413: *
! 414: * Want right eigenvectors
! 415: * Copy Householder vectors to VR
! 416: *
! 417: SIDE = 'R'
! 418: CALL DLACPY( 'L', N, N, A, LDA, VR, LDVR )
! 419: *
! 420: * Generate orthogonal matrix in VR
! 421: * (Workspace: need 2*N-1, prefer N+(N-1)*NB)
! 422: *
! 423: CALL DORGHR( N, ILO, IHI, VR, LDVR, WORK( ITAU ), WORK( IWRK ),
! 424: $ LWORK-IWRK+1, IERR )
! 425: *
! 426: * Perform QR iteration, accumulating Schur vectors in VR
! 427: * (Workspace: need 1, prefer HSWORK (see comments) )
! 428: *
! 429: IWRK = ITAU
! 430: CALL DHSEQR( 'S', 'V', N, ILO, IHI, A, LDA, WR, WI, VR, LDVR,
! 431: $ WORK( IWRK ), LWORK-IWRK+1, INFO )
! 432: *
! 433: ELSE
! 434: *
! 435: * Compute eigenvalues only
! 436: * If condition numbers desired, compute Schur form
! 437: *
! 438: IF( WNTSNN ) THEN
! 439: JOB = 'E'
! 440: ELSE
! 441: JOB = 'S'
! 442: END IF
! 443: *
! 444: * (Workspace: need 1, prefer HSWORK (see comments) )
! 445: *
! 446: IWRK = ITAU
! 447: CALL DHSEQR( JOB, 'N', N, ILO, IHI, A, LDA, WR, WI, VR, LDVR,
! 448: $ WORK( IWRK ), LWORK-IWRK+1, INFO )
! 449: END IF
! 450: *
! 451: * If INFO > 0 from DHSEQR, then quit
! 452: *
! 453: IF( INFO.GT.0 )
! 454: $ GO TO 50
! 455: *
! 456: IF( WANTVL .OR. WANTVR ) THEN
! 457: *
! 458: * Compute left and/or right eigenvectors
! 459: * (Workspace: need 3*N)
! 460: *
! 461: CALL DTREVC( SIDE, 'B', SELECT, N, A, LDA, VL, LDVL, VR, LDVR,
! 462: $ N, NOUT, WORK( IWRK ), IERR )
! 463: END IF
! 464: *
! 465: * Compute condition numbers if desired
! 466: * (Workspace: need N*N+6*N unless SENSE = 'E')
! 467: *
! 468: IF( .NOT.WNTSNN ) THEN
! 469: CALL DTRSNA( SENSE, 'A', SELECT, N, A, LDA, VL, LDVL, VR, LDVR,
! 470: $ RCONDE, RCONDV, N, NOUT, WORK( IWRK ), N, IWORK,
! 471: $ ICOND )
! 472: END IF
! 473: *
! 474: IF( WANTVL ) THEN
! 475: *
! 476: * Undo balancing of left eigenvectors
! 477: *
! 478: CALL DGEBAK( BALANC, 'L', N, ILO, IHI, SCALE, N, VL, LDVL,
! 479: $ IERR )
! 480: *
! 481: * Normalize left eigenvectors and make largest component real
! 482: *
! 483: DO 20 I = 1, N
! 484: IF( WI( I ).EQ.ZERO ) THEN
! 485: SCL = ONE / DNRM2( N, VL( 1, I ), 1 )
! 486: CALL DSCAL( N, SCL, VL( 1, I ), 1 )
! 487: ELSE IF( WI( I ).GT.ZERO ) THEN
! 488: SCL = ONE / DLAPY2( DNRM2( N, VL( 1, I ), 1 ),
! 489: $ DNRM2( N, VL( 1, I+1 ), 1 ) )
! 490: CALL DSCAL( N, SCL, VL( 1, I ), 1 )
! 491: CALL DSCAL( N, SCL, VL( 1, I+1 ), 1 )
! 492: DO 10 K = 1, N
! 493: WORK( K ) = VL( K, I )**2 + VL( K, I+1 )**2
! 494: 10 CONTINUE
! 495: K = IDAMAX( N, WORK, 1 )
! 496: CALL DLARTG( VL( K, I ), VL( K, I+1 ), CS, SN, R )
! 497: CALL DROT( N, VL( 1, I ), 1, VL( 1, I+1 ), 1, CS, SN )
! 498: VL( K, I+1 ) = ZERO
! 499: END IF
! 500: 20 CONTINUE
! 501: END IF
! 502: *
! 503: IF( WANTVR ) THEN
! 504: *
! 505: * Undo balancing of right eigenvectors
! 506: *
! 507: CALL DGEBAK( BALANC, 'R', N, ILO, IHI, SCALE, N, VR, LDVR,
! 508: $ IERR )
! 509: *
! 510: * Normalize right eigenvectors and make largest component real
! 511: *
! 512: DO 40 I = 1, N
! 513: IF( WI( I ).EQ.ZERO ) THEN
! 514: SCL = ONE / DNRM2( N, VR( 1, I ), 1 )
! 515: CALL DSCAL( N, SCL, VR( 1, I ), 1 )
! 516: ELSE IF( WI( I ).GT.ZERO ) THEN
! 517: SCL = ONE / DLAPY2( DNRM2( N, VR( 1, I ), 1 ),
! 518: $ DNRM2( N, VR( 1, I+1 ), 1 ) )
! 519: CALL DSCAL( N, SCL, VR( 1, I ), 1 )
! 520: CALL DSCAL( N, SCL, VR( 1, I+1 ), 1 )
! 521: DO 30 K = 1, N
! 522: WORK( K ) = VR( K, I )**2 + VR( K, I+1 )**2
! 523: 30 CONTINUE
! 524: K = IDAMAX( N, WORK, 1 )
! 525: CALL DLARTG( VR( K, I ), VR( K, I+1 ), CS, SN, R )
! 526: CALL DROT( N, VR( 1, I ), 1, VR( 1, I+1 ), 1, CS, SN )
! 527: VR( K, I+1 ) = ZERO
! 528: END IF
! 529: 40 CONTINUE
! 530: END IF
! 531: *
! 532: * Undo scaling if necessary
! 533: *
! 534: 50 CONTINUE
! 535: IF( SCALEA ) THEN
! 536: CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, N-INFO, 1, WR( INFO+1 ),
! 537: $ MAX( N-INFO, 1 ), IERR )
! 538: CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, N-INFO, 1, WI( INFO+1 ),
! 539: $ MAX( N-INFO, 1 ), IERR )
! 540: IF( INFO.EQ.0 ) THEN
! 541: IF( ( WNTSNV .OR. WNTSNB ) .AND. ICOND.EQ.0 )
! 542: $ CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, N, 1, RCONDV, N,
! 543: $ IERR )
! 544: ELSE
! 545: CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, ILO-1, 1, WR, N,
! 546: $ IERR )
! 547: CALL DLASCL( 'G', 0, 0, CSCALE, ANRM, ILO-1, 1, WI, N,
! 548: $ IERR )
! 549: END IF
! 550: END IF
! 551: *
! 552: WORK( 1 ) = MAXWRK
! 553: RETURN
! 554: *
! 555: * End of DGEEVX
! 556: *
! 557: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dgeevx.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack