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>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dggev3.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack