![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to ztgevc.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack|
Return to ztgevc.f CVS log | Up to [local] / rpl / lapack / lapack |
1.1 bertrand 1: SUBROUTINE ZTGEVC( SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL,
2: $ LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO )
3: *
4: * -- LAPACK routine (version 3.2) --
5: * -- LAPACK is a software package provided by Univ. of Tennessee, --
6: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7: * November 2006
8: *
9: * .. Scalar Arguments ..
10: CHARACTER HOWMNY, SIDE
11: INTEGER INFO, LDP, LDS, LDVL, LDVR, M, MM, N
12: * ..
13: * .. Array Arguments ..
14: LOGICAL SELECT( * )
15: DOUBLE PRECISION RWORK( * )
16: COMPLEX*16 P( LDP, * ), S( LDS, * ), VL( LDVL, * ),
17: $ VR( LDVR, * ), WORK( * )
18: * ..
19: *
20: *
21: * Purpose
22: * =======
23: *
24: * ZTGEVC computes some or all of the right and/or left eigenvectors of
25: * a pair of complex matrices (S,P), where S and P are upper triangular.
26: * Matrix pairs of this type are produced by the generalized Schur
27: * factorization of a complex matrix pair (A,B):
28: *
29: * A = Q*S*Z**H, B = Q*P*Z**H
30: *
31: * as computed by ZGGHRD + ZHGEQZ.
32: *
33: * The right eigenvector x and the left eigenvector y of (S,P)
34: * corresponding to an eigenvalue w are defined by:
35: *
36: * S*x = w*P*x, (y**H)*S = w*(y**H)*P,
37: *
38: * where y**H denotes the conjugate tranpose of y.
39: * The eigenvalues are not input to this routine, but are computed
40: * directly from the diagonal elements of S and P.
41: *
42: * This routine returns the matrices X and/or Y of right and left
43: * eigenvectors of (S,P), or the products Z*X and/or Q*Y,
44: * where Z and Q are input matrices.
45: * If Q and Z are the unitary factors from the generalized Schur
46: * factorization of a matrix pair (A,B), then Z*X and Q*Y
47: * are the matrices of right and left eigenvectors of (A,B).
48: *
49: * Arguments
50: * =========
51: *
52: * SIDE (input) CHARACTER*1
53: * = 'R': compute right eigenvectors only;
54: * = 'L': compute left eigenvectors only;
55: * = 'B': compute both right and left eigenvectors.
56: *
57: * HOWMNY (input) CHARACTER*1
58: * = 'A': compute all right and/or left eigenvectors;
59: * = 'B': compute all right and/or left eigenvectors,
60: * backtransformed by the matrices in VR and/or VL;
61: * = 'S': compute selected right and/or left eigenvectors,
62: * specified by the logical array SELECT.
63: *
64: * SELECT (input) LOGICAL array, dimension (N)
65: * If HOWMNY='S', SELECT specifies the eigenvectors to be
66: * computed. The eigenvector corresponding to the j-th
67: * eigenvalue is computed if SELECT(j) = .TRUE..
68: * Not referenced if HOWMNY = 'A' or 'B'.
69: *
70: * N (input) INTEGER
71: * The order of the matrices S and P. N >= 0.
72: *
73: * S (input) COMPLEX*16 array, dimension (LDS,N)
74: * The upper triangular matrix S from a generalized Schur
75: * factorization, as computed by ZHGEQZ.
76: *
77: * LDS (input) INTEGER
78: * The leading dimension of array S. LDS >= max(1,N).
79: *
80: * P (input) COMPLEX*16 array, dimension (LDP,N)
81: * The upper triangular matrix P from a generalized Schur
82: * factorization, as computed by ZHGEQZ. P must have real
83: * diagonal elements.
84: *
85: * LDP (input) INTEGER
86: * The leading dimension of array P. LDP >= max(1,N).
87: *
88: * VL (input/output) COMPLEX*16 array, dimension (LDVL,MM)
89: * On entry, if SIDE = 'L' or 'B' and HOWMNY = 'B', VL must
90: * contain an N-by-N matrix Q (usually the unitary matrix Q
91: * of left Schur vectors returned by ZHGEQZ).
92: * On exit, if SIDE = 'L' or 'B', VL contains:
93: * if HOWMNY = 'A', the matrix Y of left eigenvectors of (S,P);
94: * if HOWMNY = 'B', the matrix Q*Y;
95: * if HOWMNY = 'S', the left eigenvectors of (S,P) specified by
96: * SELECT, stored consecutively in the columns of
97: * VL, in the same order as their eigenvalues.
98: * Not referenced if SIDE = 'R'.
99: *
100: * LDVL (input) INTEGER
101: * The leading dimension of array VL. LDVL >= 1, and if
102: * SIDE = 'L' or 'l' or 'B' or 'b', LDVL >= N.
103: *
104: * VR (input/output) COMPLEX*16 array, dimension (LDVR,MM)
105: * On entry, if SIDE = 'R' or 'B' and HOWMNY = 'B', VR must
106: * contain an N-by-N matrix Q (usually the unitary matrix Z
107: * of right Schur vectors returned by ZHGEQZ).
108: * On exit, if SIDE = 'R' or 'B', VR contains:
109: * if HOWMNY = 'A', the matrix X of right eigenvectors of (S,P);
110: * if HOWMNY = 'B', the matrix Z*X;
111: * if HOWMNY = 'S', the right eigenvectors of (S,P) specified by
112: * SELECT, stored consecutively in the columns of
113: * VR, in the same order as their eigenvalues.
114: * Not referenced if SIDE = 'L'.
115: *
116: * LDVR (input) INTEGER
117: * The leading dimension of the array VR. LDVR >= 1, and if
118: * SIDE = 'R' or 'B', LDVR >= N.
119: *
120: * MM (input) INTEGER
121: * The number of columns in the arrays VL and/or VR. MM >= M.
122: *
123: * M (output) INTEGER
124: * The number of columns in the arrays VL and/or VR actually
125: * used to store the eigenvectors. If HOWMNY = 'A' or 'B', M
126: * is set to N. Each selected eigenvector occupies one column.
127: *
128: * WORK (workspace) COMPLEX*16 array, dimension (2*N)
129: *
130: * RWORK (workspace) DOUBLE PRECISION array, dimension (2*N)
131: *
132: * INFO (output) INTEGER
133: * = 0: successful exit.
134: * < 0: if INFO = -i, the i-th argument had an illegal value.
135: *
136: * =====================================================================
137: *
138: * .. Parameters ..
139: DOUBLE PRECISION ZERO, ONE
140: PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0 )
141: COMPLEX*16 CZERO, CONE
142: PARAMETER ( CZERO = ( 0.0D+0, 0.0D+0 ),
143: $ CONE = ( 1.0D+0, 0.0D+0 ) )
144: * ..
145: * .. Local Scalars ..
146: LOGICAL COMPL, COMPR, ILALL, ILBACK, ILBBAD, ILCOMP,
147: $ LSA, LSB
148: INTEGER I, IBEG, IEIG, IEND, IHWMNY, IM, ISIDE, ISRC,
149: $ J, JE, JR
150: DOUBLE PRECISION ACOEFA, ACOEFF, ANORM, ASCALE, BCOEFA, BIG,
151: $ BIGNUM, BNORM, BSCALE, DMIN, SAFMIN, SBETA,
152: $ SCALE, SMALL, TEMP, ULP, XMAX
153: COMPLEX*16 BCOEFF, CA, CB, D, SALPHA, SUM, SUMA, SUMB, X
154: * ..
155: * .. External Functions ..
156: LOGICAL LSAME
157: DOUBLE PRECISION DLAMCH
158: COMPLEX*16 ZLADIV
159: EXTERNAL LSAME, DLAMCH, ZLADIV
160: * ..
161: * .. External Subroutines ..
162: EXTERNAL DLABAD, XERBLA, ZGEMV
163: * ..
164: * .. Intrinsic Functions ..
165: INTRINSIC ABS, DBLE, DCMPLX, DCONJG, DIMAG, MAX, MIN
166: * ..
167: * .. Statement Functions ..
168: DOUBLE PRECISION ABS1
169: * ..
170: * .. Statement Function definitions ..
171: ABS1( X ) = ABS( DBLE( X ) ) + ABS( DIMAG( X ) )
172: * ..
173: * .. Executable Statements ..
174: *
175: * Decode and Test the input parameters
176: *
177: IF( LSAME( HOWMNY, 'A' ) ) THEN
178: IHWMNY = 1
179: ILALL = .TRUE.
180: ILBACK = .FALSE.
181: ELSE IF( LSAME( HOWMNY, 'S' ) ) THEN
182: IHWMNY = 2
183: ILALL = .FALSE.
184: ILBACK = .FALSE.
185: ELSE IF( LSAME( HOWMNY, 'B' ) ) THEN
186: IHWMNY = 3
187: ILALL = .TRUE.
188: ILBACK = .TRUE.
189: ELSE
190: IHWMNY = -1
191: END IF
192: *
193: IF( LSAME( SIDE, 'R' ) ) THEN
194: ISIDE = 1
195: COMPL = .FALSE.
196: COMPR = .TRUE.
197: ELSE IF( LSAME( SIDE, 'L' ) ) THEN
198: ISIDE = 2
199: COMPL = .TRUE.
200: COMPR = .FALSE.
201: ELSE IF( LSAME( SIDE, 'B' ) ) THEN
202: ISIDE = 3
203: COMPL = .TRUE.
204: COMPR = .TRUE.
205: ELSE
206: ISIDE = -1
207: END IF
208: *
209: INFO = 0
210: IF( ISIDE.LT.0 ) THEN
211: INFO = -1
212: ELSE IF( IHWMNY.LT.0 ) THEN
213: INFO = -2
214: ELSE IF( N.LT.0 ) THEN
215: INFO = -4
216: ELSE IF( LDS.LT.MAX( 1, N ) ) THEN
217: INFO = -6
218: ELSE IF( LDP.LT.MAX( 1, N ) ) THEN
219: INFO = -8
220: END IF
221: IF( INFO.NE.0 ) THEN
222: CALL XERBLA( 'ZTGEVC', -INFO )
223: RETURN
224: END IF
225: *
226: * Count the number of eigenvectors
227: *
228: IF( .NOT.ILALL ) THEN
229: IM = 0
230: DO 10 J = 1, N
231: IF( SELECT( J ) )
232: $ IM = IM + 1
233: 10 CONTINUE
234: ELSE
235: IM = N
236: END IF
237: *
238: * Check diagonal of B
239: *
240: ILBBAD = .FALSE.
241: DO 20 J = 1, N
242: IF( DIMAG( P( J, J ) ).NE.ZERO )
243: $ ILBBAD = .TRUE.
244: 20 CONTINUE
245: *
246: IF( ILBBAD ) THEN
247: INFO = -7
248: ELSE IF( COMPL .AND. LDVL.LT.N .OR. LDVL.LT.1 ) THEN
249: INFO = -10
250: ELSE IF( COMPR .AND. LDVR.LT.N .OR. LDVR.LT.1 ) THEN
251: INFO = -12
252: ELSE IF( MM.LT.IM ) THEN
253: INFO = -13
254: END IF
255: IF( INFO.NE.0 ) THEN
256: CALL XERBLA( 'ZTGEVC', -INFO )
257: RETURN
258: END IF
259: *
260: * Quick return if possible
261: *
262: M = IM
263: IF( N.EQ.0 )
264: $ RETURN
265: *
266: * Machine Constants
267: *
268: SAFMIN = DLAMCH( 'Safe minimum' )
269: BIG = ONE / SAFMIN
270: CALL DLABAD( SAFMIN, BIG )
271: ULP = DLAMCH( 'Epsilon' )*DLAMCH( 'Base' )
272: SMALL = SAFMIN*N / ULP
273: BIG = ONE / SMALL
274: BIGNUM = ONE / ( SAFMIN*N )
275: *
276: * Compute the 1-norm of each column of the strictly upper triangular
277: * part of A and B to check for possible overflow in the triangular
278: * solver.
279: *
280: ANORM = ABS1( S( 1, 1 ) )
281: BNORM = ABS1( P( 1, 1 ) )
282: RWORK( 1 ) = ZERO
283: RWORK( N+1 ) = ZERO
284: DO 40 J = 2, N
285: RWORK( J ) = ZERO
286: RWORK( N+J ) = ZERO
287: DO 30 I = 1, J - 1
288: RWORK( J ) = RWORK( J ) + ABS1( S( I, J ) )
289: RWORK( N+J ) = RWORK( N+J ) + ABS1( P( I, J ) )
290: 30 CONTINUE
291: ANORM = MAX( ANORM, RWORK( J )+ABS1( S( J, J ) ) )
292: BNORM = MAX( BNORM, RWORK( N+J )+ABS1( P( J, J ) ) )
293: 40 CONTINUE
294: *
295: ASCALE = ONE / MAX( ANORM, SAFMIN )
296: BSCALE = ONE / MAX( BNORM, SAFMIN )
297: *
298: * Left eigenvectors
299: *
300: IF( COMPL ) THEN
301: IEIG = 0
302: *
303: * Main loop over eigenvalues
304: *
305: DO 140 JE = 1, N
306: IF( ILALL ) THEN
307: ILCOMP = .TRUE.
308: ELSE
309: ILCOMP = SELECT( JE )
310: END IF
311: IF( ILCOMP ) THEN
312: IEIG = IEIG + 1
313: *
314: IF( ABS1( S( JE, JE ) ).LE.SAFMIN .AND.
315: $ ABS( DBLE( P( JE, JE ) ) ).LE.SAFMIN ) THEN
316: *
317: * Singular matrix pencil -- return unit eigenvector
318: *
319: DO 50 JR = 1, N
320: VL( JR, IEIG ) = CZERO
321: 50 CONTINUE
322: VL( IEIG, IEIG ) = CONE
323: GO TO 140
324: END IF
325: *
326: * Non-singular eigenvalue:
327: * Compute coefficients a and b in
328: * H
329: * y ( a A - b B ) = 0
330: *
331: TEMP = ONE / MAX( ABS1( S( JE, JE ) )*ASCALE,
332: $ ABS( DBLE( P( JE, JE ) ) )*BSCALE, SAFMIN )
333: SALPHA = ( TEMP*S( JE, JE ) )*ASCALE
334: SBETA = ( TEMP*DBLE( P( JE, JE ) ) )*BSCALE
335: ACOEFF = SBETA*ASCALE
336: BCOEFF = SALPHA*BSCALE
337: *
338: * Scale to avoid underflow
339: *
340: LSA = ABS( SBETA ).GE.SAFMIN .AND. ABS( ACOEFF ).LT.SMALL
341: LSB = ABS1( SALPHA ).GE.SAFMIN .AND. ABS1( BCOEFF ).LT.
342: $ SMALL
343: *
344: SCALE = ONE
345: IF( LSA )
346: $ SCALE = ( SMALL / ABS( SBETA ) )*MIN( ANORM, BIG )
347: IF( LSB )
348: $ SCALE = MAX( SCALE, ( SMALL / ABS1( SALPHA ) )*
349: $ MIN( BNORM, BIG ) )
350: IF( LSA .OR. LSB ) THEN
351: SCALE = MIN( SCALE, ONE /
352: $ ( SAFMIN*MAX( ONE, ABS( ACOEFF ),
353: $ ABS1( BCOEFF ) ) ) )
354: IF( LSA ) THEN
355: ACOEFF = ASCALE*( SCALE*SBETA )
356: ELSE
357: ACOEFF = SCALE*ACOEFF
358: END IF
359: IF( LSB ) THEN
360: BCOEFF = BSCALE*( SCALE*SALPHA )
361: ELSE
362: BCOEFF = SCALE*BCOEFF
363: END IF
364: END IF
365: *
366: ACOEFA = ABS( ACOEFF )
367: BCOEFA = ABS1( BCOEFF )
368: XMAX = ONE
369: DO 60 JR = 1, N
370: WORK( JR ) = CZERO
371: 60 CONTINUE
372: WORK( JE ) = CONE
373: DMIN = MAX( ULP*ACOEFA*ANORM, ULP*BCOEFA*BNORM, SAFMIN )
374: *
375: * H
376: * Triangular solve of (a A - b B) y = 0
377: *
378: * H
379: * (rowwise in (a A - b B) , or columnwise in a A - b B)
380: *
381: DO 100 J = JE + 1, N
382: *
383: * Compute
384: * j-1
385: * SUM = sum conjg( a*S(k,j) - b*P(k,j) )*x(k)
386: * k=je
387: * (Scale if necessary)
388: *
389: TEMP = ONE / XMAX
390: IF( ACOEFA*RWORK( J )+BCOEFA*RWORK( N+J ).GT.BIGNUM*
391: $ TEMP ) THEN
392: DO 70 JR = JE, J - 1
393: WORK( JR ) = TEMP*WORK( JR )
394: 70 CONTINUE
395: XMAX = ONE
396: END IF
397: SUMA = CZERO
398: SUMB = CZERO
399: *
400: DO 80 JR = JE, J - 1
401: SUMA = SUMA + DCONJG( S( JR, J ) )*WORK( JR )
402: SUMB = SUMB + DCONJG( P( JR, J ) )*WORK( JR )
403: 80 CONTINUE
404: SUM = ACOEFF*SUMA - DCONJG( BCOEFF )*SUMB
405: *
406: * Form x(j) = - SUM / conjg( a*S(j,j) - b*P(j,j) )
407: *
408: * with scaling and perturbation of the denominator
409: *
410: D = DCONJG( ACOEFF*S( J, J )-BCOEFF*P( J, J ) )
411: IF( ABS1( D ).LE.DMIN )
412: $ D = DCMPLX( DMIN )
413: *
414: IF( ABS1( D ).LT.ONE ) THEN
415: IF( ABS1( SUM ).GE.BIGNUM*ABS1( D ) ) THEN
416: TEMP = ONE / ABS1( SUM )
417: DO 90 JR = JE, J - 1
418: WORK( JR ) = TEMP*WORK( JR )
419: 90 CONTINUE
420: XMAX = TEMP*XMAX
421: SUM = TEMP*SUM
422: END IF
423: END IF
424: WORK( J ) = ZLADIV( -SUM, D )
425: XMAX = MAX( XMAX, ABS1( WORK( J ) ) )
426: 100 CONTINUE
427: *
428: * Back transform eigenvector if HOWMNY='B'.
429: *
430: IF( ILBACK ) THEN
431: CALL ZGEMV( 'N', N, N+1-JE, CONE, VL( 1, JE ), LDVL,
432: $ WORK( JE ), 1, CZERO, WORK( N+1 ), 1 )
433: ISRC = 2
434: IBEG = 1
435: ELSE
436: ISRC = 1
437: IBEG = JE
438: END IF
439: *
440: * Copy and scale eigenvector into column of VL
441: *
442: XMAX = ZERO
443: DO 110 JR = IBEG, N
444: XMAX = MAX( XMAX, ABS1( WORK( ( ISRC-1 )*N+JR ) ) )
445: 110 CONTINUE
446: *
447: IF( XMAX.GT.SAFMIN ) THEN
448: TEMP = ONE / XMAX
449: DO 120 JR = IBEG, N
450: VL( JR, IEIG ) = TEMP*WORK( ( ISRC-1 )*N+JR )
451: 120 CONTINUE
452: ELSE
453: IBEG = N + 1
454: END IF
455: *
456: DO 130 JR = 1, IBEG - 1
457: VL( JR, IEIG ) = CZERO
458: 130 CONTINUE
459: *
460: END IF
461: 140 CONTINUE
462: END IF
463: *
464: * Right eigenvectors
465: *
466: IF( COMPR ) THEN
467: IEIG = IM + 1
468: *
469: * Main loop over eigenvalues
470: *
471: DO 250 JE = N, 1, -1
472: IF( ILALL ) THEN
473: ILCOMP = .TRUE.
474: ELSE
475: ILCOMP = SELECT( JE )
476: END IF
477: IF( ILCOMP ) THEN
478: IEIG = IEIG - 1
479: *
480: IF( ABS1( S( JE, JE ) ).LE.SAFMIN .AND.
481: $ ABS( DBLE( P( JE, JE ) ) ).LE.SAFMIN ) THEN
482: *
483: * Singular matrix pencil -- return unit eigenvector
484: *
485: DO 150 JR = 1, N
486: VR( JR, IEIG ) = CZERO
487: 150 CONTINUE
488: VR( IEIG, IEIG ) = CONE
489: GO TO 250
490: END IF
491: *
492: * Non-singular eigenvalue:
493: * Compute coefficients a and b in
494: *
495: * ( a A - b B ) x = 0
496: *
497: TEMP = ONE / MAX( ABS1( S( JE, JE ) )*ASCALE,
498: $ ABS( DBLE( P( JE, JE ) ) )*BSCALE, SAFMIN )
499: SALPHA = ( TEMP*S( JE, JE ) )*ASCALE
500: SBETA = ( TEMP*DBLE( P( JE, JE ) ) )*BSCALE
501: ACOEFF = SBETA*ASCALE
502: BCOEFF = SALPHA*BSCALE
503: *
504: * Scale to avoid underflow
505: *
506: LSA = ABS( SBETA ).GE.SAFMIN .AND. ABS( ACOEFF ).LT.SMALL
507: LSB = ABS1( SALPHA ).GE.SAFMIN .AND. ABS1( BCOEFF ).LT.
508: $ SMALL
509: *
510: SCALE = ONE
511: IF( LSA )
512: $ SCALE = ( SMALL / ABS( SBETA ) )*MIN( ANORM, BIG )
513: IF( LSB )
514: $ SCALE = MAX( SCALE, ( SMALL / ABS1( SALPHA ) )*
515: $ MIN( BNORM, BIG ) )
516: IF( LSA .OR. LSB ) THEN
517: SCALE = MIN( SCALE, ONE /
518: $ ( SAFMIN*MAX( ONE, ABS( ACOEFF ),
519: $ ABS1( BCOEFF ) ) ) )
520: IF( LSA ) THEN
521: ACOEFF = ASCALE*( SCALE*SBETA )
522: ELSE
523: ACOEFF = SCALE*ACOEFF
524: END IF
525: IF( LSB ) THEN
526: BCOEFF = BSCALE*( SCALE*SALPHA )
527: ELSE
528: BCOEFF = SCALE*BCOEFF
529: END IF
530: END IF
531: *
532: ACOEFA = ABS( ACOEFF )
533: BCOEFA = ABS1( BCOEFF )
534: XMAX = ONE
535: DO 160 JR = 1, N
536: WORK( JR ) = CZERO
537: 160 CONTINUE
538: WORK( JE ) = CONE
539: DMIN = MAX( ULP*ACOEFA*ANORM, ULP*BCOEFA*BNORM, SAFMIN )
540: *
541: * Triangular solve of (a A - b B) x = 0 (columnwise)
542: *
543: * WORK(1:j-1) contains sums w,
544: * WORK(j+1:JE) contains x
545: *
546: DO 170 JR = 1, JE - 1
547: WORK( JR ) = ACOEFF*S( JR, JE ) - BCOEFF*P( JR, JE )
548: 170 CONTINUE
549: WORK( JE ) = CONE
550: *
551: DO 210 J = JE - 1, 1, -1
552: *
553: * Form x(j) := - w(j) / d
554: * with scaling and perturbation of the denominator
555: *
556: D = ACOEFF*S( J, J ) - BCOEFF*P( J, J )
557: IF( ABS1( D ).LE.DMIN )
558: $ D = DCMPLX( DMIN )
559: *
560: IF( ABS1( D ).LT.ONE ) THEN
561: IF( ABS1( WORK( J ) ).GE.BIGNUM*ABS1( D ) ) THEN
562: TEMP = ONE / ABS1( WORK( J ) )
563: DO 180 JR = 1, JE
564: WORK( JR ) = TEMP*WORK( JR )
565: 180 CONTINUE
566: END IF
567: END IF
568: *
569: WORK( J ) = ZLADIV( -WORK( J ), D )
570: *
571: IF( J.GT.1 ) THEN
572: *
573: * w = w + x(j)*(a S(*,j) - b P(*,j) ) with scaling
574: *
575: IF( ABS1( WORK( J ) ).GT.ONE ) THEN
576: TEMP = ONE / ABS1( WORK( J ) )
577: IF( ACOEFA*RWORK( J )+BCOEFA*RWORK( N+J ).GE.
578: $ BIGNUM*TEMP ) THEN
579: DO 190 JR = 1, JE
580: WORK( JR ) = TEMP*WORK( JR )
581: 190 CONTINUE
582: END IF
583: END IF
584: *
585: CA = ACOEFF*WORK( J )
586: CB = BCOEFF*WORK( J )
587: DO 200 JR = 1, J - 1
588: WORK( JR ) = WORK( JR ) + CA*S( JR, J ) -
589: $ CB*P( JR, J )
590: 200 CONTINUE
591: END IF
592: 210 CONTINUE
593: *
594: * Back transform eigenvector if HOWMNY='B'.
595: *
596: IF( ILBACK ) THEN
597: CALL ZGEMV( 'N', N, JE, CONE, VR, LDVR, WORK, 1,
598: $ CZERO, WORK( N+1 ), 1 )
599: ISRC = 2
600: IEND = N
601: ELSE
602: ISRC = 1
603: IEND = JE
604: END IF
605: *
606: * Copy and scale eigenvector into column of VR
607: *
608: XMAX = ZERO
609: DO 220 JR = 1, IEND
610: XMAX = MAX( XMAX, ABS1( WORK( ( ISRC-1 )*N+JR ) ) )
611: 220 CONTINUE
612: *
613: IF( XMAX.GT.SAFMIN ) THEN
614: TEMP = ONE / XMAX
615: DO 230 JR = 1, IEND
616: VR( JR, IEIG ) = TEMP*WORK( ( ISRC-1 )*N+JR )
617: 230 CONTINUE
618: ELSE
619: IEND = 0
620: END IF
621: *
622: DO 240 JR = IEND + 1, N
623: VR( JR, IEIG ) = CZERO
624: 240 CONTINUE
625: *
626: END IF
627: 250 CONTINUE
628: END IF
629: *
630: RETURN
631: *
632: * End of ZTGEVC
633: *
634: END