LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
zgeevx.f
Go to the documentation of this file.
1 *> \brief <b> ZGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors 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 ZGEEVX + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgeevx.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgeevx.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgeevx.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE ZGEEVX( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, W, VL,
22 * LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE,
23 * RCONDV, WORK, LWORK, RWORK, INFO )
24 *
25 * .. Scalar Arguments ..
26 * CHARACTER BALANC, JOBVL, JOBVR, SENSE
27 * INTEGER IHI, ILO, INFO, LDA, LDVL, LDVR, LWORK, N
28 * DOUBLE PRECISION ABNRM
29 * ..
30 * .. Array Arguments ..
31 * DOUBLE PRECISION RCONDE( * ), RCONDV( * ), RWORK( * ),
32 * $ SCALE( * )
33 * COMPLEX*16 A( LDA, * ), VL( LDVL, * ), VR( LDVR, * ),
34 * $ W( * ), WORK( * )
35 * ..
36 *
37 *
38 *> \par Purpose:
39 * =============
40 *>
41 *> \verbatim
42 *>
43 *> ZGEEVX computes for an N-by-N complex nonsymmetric matrix A, the
44 *> eigenvalues and, optionally, the left and/or right eigenvectors.
45 *>
46 *> Optionally also, it computes a balancing transformation to improve
47 *> the conditioning of the eigenvalues and eigenvectors (ILO, IHI,
48 *> SCALE, and ABNRM), reciprocal condition numbers for the eigenvalues
49 *> (RCONDE), and reciprocal condition numbers for the right
50 *> eigenvectors (RCONDV).
51 *>
52 *> The right eigenvector v(j) of A satisfies
53 *> A * v(j) = lambda(j) * v(j)
54 *> where lambda(j) is its eigenvalue.
55 *> The left eigenvector u(j) of A satisfies
56 *> u(j)**H * A = lambda(j) * u(j)**H
57 *> where u(j)**H denotes the conjugate transpose of u(j).
58 *>
59 *> The computed eigenvectors are normalized to have Euclidean norm
60 *> equal to 1 and largest component real.
61 *>
62 *> Balancing a matrix means permuting the rows and columns to make it
63 *> more nearly upper triangular, and applying a diagonal similarity
64 *> transformation D * A * D**(-1), where D is a diagonal matrix, to
65 *> make its rows and columns closer in norm and the condition numbers
66 *> of its eigenvalues and eigenvectors smaller. The computed
67 *> reciprocal condition numbers correspond to the balanced matrix.
68 *> Permuting rows and columns will not change the condition numbers
69 *> (in exact arithmetic) but diagonal scaling will. For further
70 *> explanation of balancing, see section 4.10.2 of the LAPACK
71 *> Users' Guide.
72 *> \endverbatim
73 *
74 * Arguments:
75 * ==========
76 *
77 *> \param[in] BALANC
78 *> \verbatim
79 *> BALANC is CHARACTER*1
80 *> Indicates how the input matrix should be diagonally scaled
81 *> and/or permuted to improve the conditioning of its
82 *> eigenvalues.
83 *> = 'N': Do not diagonally scale or permute;
84 *> = 'P': Perform permutations to make the matrix more nearly
85 *> upper triangular. Do not diagonally scale;
86 *> = 'S': Diagonally scale the matrix, ie. replace A by
87 *> D*A*D**(-1), where D is a diagonal matrix chosen
88 *> to make the rows and columns of A more equal in
89 *> norm. Do not permute;
90 *> = 'B': Both diagonally scale and permute A.
91 *>
92 *> Computed reciprocal condition numbers will be for the matrix
93 *> after balancing and/or permuting. Permuting does not change
94 *> condition numbers (in exact arithmetic), but balancing does.
95 *> \endverbatim
96 *>
97 *> \param[in] JOBVL
98 *> \verbatim
99 *> JOBVL is CHARACTER*1
100 *> = 'N': left eigenvectors of A are not computed;
101 *> = 'V': left eigenvectors of A are computed.
102 *> If SENSE = 'E' or 'B', JOBVL must = 'V'.
103 *> \endverbatim
104 *>
105 *> \param[in] JOBVR
106 *> \verbatim
107 *> JOBVR is CHARACTER*1
108 *> = 'N': right eigenvectors of A are not computed;
109 *> = 'V': right eigenvectors of A are computed.
110 *> If SENSE = 'E' or 'B', JOBVR must = 'V'.
111 *> \endverbatim
112 *>
113 *> \param[in] SENSE
114 *> \verbatim
115 *> SENSE is CHARACTER*1
116 *> Determines which reciprocal condition numbers are computed.
117 *> = 'N': None are computed;
118 *> = 'E': Computed for eigenvalues only;
119 *> = 'V': Computed for right eigenvectors only;
120 *> = 'B': Computed for eigenvalues and right eigenvectors.
121 *>
122 *> If SENSE = 'E' or 'B', both left and right eigenvectors
123 *> must also be computed (JOBVL = 'V' and JOBVR = 'V').
124 *> \endverbatim
125 *>
126 *> \param[in] N
127 *> \verbatim
128 *> N is INTEGER
129 *> The order of the matrix A. N >= 0.
130 *> \endverbatim
131 *>
132 *> \param[in,out] A
133 *> \verbatim
134 *> A is COMPLEX*16 array, dimension (LDA,N)
135 *> On entry, the N-by-N matrix A.
136 *> On exit, A has been overwritten. If JOBVL = 'V' or
137 *> JOBVR = 'V', A contains the Schur form of the balanced
138 *> version of the matrix A.
139 *> \endverbatim
140 *>
141 *> \param[in] LDA
142 *> \verbatim
143 *> LDA is INTEGER
144 *> The leading dimension of the array A. LDA >= max(1,N).
145 *> \endverbatim
146 *>
147 *> \param[out] W
148 *> \verbatim
149 *> W is COMPLEX*16 array, dimension (N)
150 *> W contains the computed eigenvalues.
151 *> \endverbatim
152 *>
153 *> \param[out] VL
154 *> \verbatim
155 *> VL is COMPLEX*16 array, dimension (LDVL,N)
156 *> If JOBVL = 'V', the left eigenvectors u(j) are stored one
157 *> after another in the columns of VL, in the same order
158 *> as their eigenvalues.
159 *> If JOBVL = 'N', VL is not referenced.
160 *> u(j) = VL(:,j), the j-th column of VL.
161 *> \endverbatim
162 *>
163 *> \param[in] LDVL
164 *> \verbatim
165 *> LDVL is INTEGER
166 *> The leading dimension of the array VL. LDVL >= 1; if
167 *> JOBVL = 'V', LDVL >= N.
168 *> \endverbatim
169 *>
170 *> \param[out] VR
171 *> \verbatim
172 *> VR is COMPLEX*16 array, dimension (LDVR,N)
173 *> If JOBVR = 'V', the right eigenvectors v(j) are stored one
174 *> after another in the columns of VR, in the same order
175 *> as their eigenvalues.
176 *> If JOBVR = 'N', VR is not referenced.
177 *> v(j) = VR(:,j), the j-th column of VR.
178 *> \endverbatim
179 *>
180 *> \param[in] LDVR
181 *> \verbatim
182 *> LDVR is INTEGER
183 *> The leading dimension of the array VR. LDVR >= 1; if
184 *> JOBVR = 'V', LDVR >= N.
185 *> \endverbatim
186 *>
187 *> \param[out] ILO
188 *> \verbatim
189 *> ILO is INTEGER
190 *> \endverbatim
191 *>
192 *> \param[out] IHI
193 *> \verbatim
194 *> IHI is INTEGER
195 *> ILO and IHI are integer values determined when A was
196 *> balanced. The balanced A(i,j) = 0 if I > J and
197 *> J = 1,...,ILO-1 or I = IHI+1,...,N.
198 *> \endverbatim
199 *>
200 *> \param[out] SCALE
201 *> \verbatim
202 *> SCALE is DOUBLE PRECISION array, dimension (N)
203 *> Details of the permutations and scaling factors applied
204 *> when balancing A. If P(j) is the index of the row and column
205 *> interchanged with row and column j, and D(j) is the scaling
206 *> factor applied to row and column j, then
207 *> SCALE(J) = P(J), for J = 1,...,ILO-1
208 *> = D(J), for J = ILO,...,IHI
209 *> = P(J) for J = IHI+1,...,N.
210 *> The order in which the interchanges are made is N to IHI+1,
211 *> then 1 to ILO-1.
212 *> \endverbatim
213 *>
214 *> \param[out] ABNRM
215 *> \verbatim
216 *> ABNRM is DOUBLE PRECISION
217 *> The one-norm of the balanced matrix (the maximum
218 *> of the sum of absolute values of elements of any column).
219 *> \endverbatim
220 *>
221 *> \param[out] RCONDE
222 *> \verbatim
223 *> RCONDE is DOUBLE PRECISION array, dimension (N)
224 *> RCONDE(j) is the reciprocal condition number of the j-th
225 *> eigenvalue.
226 *> \endverbatim
227 *>
228 *> \param[out] RCONDV
229 *> \verbatim
230 *> RCONDV is DOUBLE PRECISION array, dimension (N)
231 *> RCONDV(j) is the reciprocal condition number of the j-th
232 *> right eigenvector.
233 *> \endverbatim
234 *>
235 *> \param[out] WORK
236 *> \verbatim
237 *> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
238 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
239 *> \endverbatim
240 *>
241 *> \param[in] LWORK
242 *> \verbatim
243 *> LWORK is INTEGER
244 *> The dimension of the array WORK. If SENSE = 'N' or 'E',
245 *> LWORK >= max(1,2*N), and if SENSE = 'V' or 'B',
246 *> LWORK >= N*N+2*N.
247 *> For good performance, LWORK must generally be larger.
248 *>
249 *> If LWORK = -1, then a workspace query is assumed; the routine
250 *> only calculates the optimal size of the WORK array, returns
251 *> this value as the first entry of the WORK array, and no error
252 *> message related to LWORK is issued by XERBLA.
253 *> \endverbatim
254 *>
255 *> \param[out] RWORK
256 *> \verbatim
257 *> RWORK is DOUBLE PRECISION array, dimension (2*N)
258 *> \endverbatim
259 *>
260 *> \param[out] INFO
261 *> \verbatim
262 *> INFO is INTEGER
263 *> = 0: successful exit
264 *> < 0: if INFO = -i, the i-th argument had an illegal value.
265 *> > 0: if INFO = i, the QR algorithm failed to compute all the
266 *> eigenvalues, and no eigenvectors or condition numbers
267 *> have been computed; elements 1:ILO-1 and i+1:N of W
268 *> contain eigenvalues which have converged.
269 *> \endverbatim
270 *
271 * Authors:
272 * ========
273 *
274 *> \author Univ. of Tennessee
275 *> \author Univ. of California Berkeley
276 *> \author Univ. of Colorado Denver
277 *> \author NAG Ltd.
278 *
279 *> \date November 2011
280 *
281 *> \ingroup complex16GEeigen
282 *
283 * =====================================================================
284  SUBROUTINE zgeevx( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, W, VL,
285  $ ldvl, vr, ldvr, ilo, ihi, scale, abnrm, rconde,
286  $ rcondv, work, lwork, rwork, info )
287 *
288 * -- LAPACK driver routine (version 3.4.0) --
289 * -- LAPACK is a software package provided by Univ. of Tennessee, --
290 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
291 * November 2011
292 *
293 * .. Scalar Arguments ..
294  CHARACTER balanc, jobvl, jobvr, sense
295  INTEGER ihi, ilo, info, lda, ldvl, ldvr, lwork, n
296  DOUBLE PRECISION abnrm
297 * ..
298 * .. Array Arguments ..
299  DOUBLE PRECISION rconde( * ), rcondv( * ), rwork( * ),
300  $ scale( * )
301  COMPLEX*16 a( lda, * ), vl( ldvl, * ), vr( ldvr, * ),
302  $ w( * ), work( * )
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 lquery, scalea, wantvl, wantvr, wntsnb, wntsne,
313  $ wntsnn, wntsnv
314  CHARACTER job, side
315  INTEGER hswork, i, icond, ierr, itau, iwrk, k, maxwrk,
316  $ minwrk, nout
317  DOUBLE PRECISION anrm, bignum, cscale, eps, scl, smlnum
318  COMPLEX*16 tmp
319 * ..
320 * .. Local Arrays ..
321  LOGICAL select( 1 )
322  DOUBLE PRECISION dum( 1 )
323 * ..
324 * .. External Subroutines ..
325  EXTERNAL dlabad, dlascl, xerbla, zdscal, zgebak, zgebal,
327  $ ztrsna, zunghr
328 * ..
329 * .. External Functions ..
330  LOGICAL lsame
331  INTEGER idamax, ilaenv
332  DOUBLE PRECISION dlamch, dznrm2, zlange
333  EXTERNAL lsame, idamax, ilaenv, dlamch, dznrm2, zlange
334 * ..
335 * .. Intrinsic Functions ..
336  INTRINSIC dble, dcmplx, dconjg, dimag, max, sqrt
337 * ..
338 * .. Executable Statements ..
339 *
340 * Test the input arguments
341 *
342  info = 0
343  lquery = ( lwork.EQ.-1 )
344  wantvl = lsame( jobvl, 'V' )
345  wantvr = lsame( jobvr, 'V' )
346  wntsnn = lsame( sense, 'N' )
347  wntsne = lsame( sense, 'E' )
348  wntsnv = lsame( sense, 'V' )
349  wntsnb = lsame( sense, 'B' )
350  IF( .NOT.( lsame( balanc, 'N' ) .OR. lsame( balanc, 'S' ) .OR.
351  $ lsame( balanc, 'P' ) .OR. lsame( balanc, 'B' ) ) ) THEN
352  info = -1
353  ELSE IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl, 'N' ) ) ) THEN
354  info = -2
355  ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr, 'N' ) ) ) THEN
356  info = -3
357  ELSE IF( .NOT.( wntsnn .OR. wntsne .OR. wntsnb .OR. wntsnv ) .OR.
358  $ ( ( wntsne .OR. wntsnb ) .AND. .NOT.( wantvl .AND.
359  $ wantvr ) ) ) THEN
360  info = -4
361  ELSE IF( n.LT.0 ) THEN
362  info = -5
363  ELSE IF( lda.LT.max( 1, n ) ) THEN
364  info = -7
365  ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) ) THEN
366  info = -10
367  ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) ) THEN
368  info = -12
369  END IF
370 *
371 * Compute workspace
372 * (Note: Comments in the code beginning "Workspace:" describe the
373 * minimal amount of workspace needed at that point in the code,
374 * as well as the preferred amount for good performance.
375 * CWorkspace refers to complex workspace, and RWorkspace to real
376 * workspace. NB refers to the optimal block size for the
377 * immediately following subroutine, as returned by ILAENV.
378 * HSWORK refers to the workspace preferred by ZHSEQR, as
379 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
380 * the worst case.)
381 *
382  IF( info.EQ.0 ) THEN
383  IF( n.EQ.0 ) THEN
384  minwrk = 1
385  maxwrk = 1
386  ELSE
387  maxwrk = n + n*ilaenv( 1, 'ZGEHRD', ' ', n, 1, n, 0 )
388 *
389  IF( wantvl ) THEN
390  CALL zhseqr( 'S', 'V', n, 1, n, a, lda, w, vl, ldvl,
391  $ work, -1, info )
392  ELSE IF( wantvr ) THEN
393  CALL zhseqr( 'S', 'V', n, 1, n, a, lda, w, vr, ldvr,
394  $ work, -1, info )
395  ELSE
396  IF( wntsnn ) THEN
397  CALL zhseqr( 'E', 'N', n, 1, n, a, lda, w, vr, ldvr,
398  $ work, -1, info )
399  ELSE
400  CALL zhseqr( 'S', 'N', n, 1, n, a, lda, w, vr, ldvr,
401  $ work, -1, info )
402  END IF
403  END IF
404  hswork = work( 1 )
405 *
406  IF( ( .NOT.wantvl ) .AND. ( .NOT.wantvr ) ) THEN
407  minwrk = 2*n
408  IF( .NOT.( wntsnn .OR. wntsne ) )
409  $ minwrk = max( minwrk, n*n + 2*n )
410  maxwrk = max( maxwrk, hswork )
411  IF( .NOT.( wntsnn .OR. wntsne ) )
412  $ maxwrk = max( maxwrk, n*n + 2*n )
413  ELSE
414  minwrk = 2*n
415  IF( .NOT.( wntsnn .OR. wntsne ) )
416  $ minwrk = max( minwrk, n*n + 2*n )
417  maxwrk = max( maxwrk, hswork )
418  maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'ZUNGHR',
419  $ ' ', n, 1, n, -1 ) )
420  IF( .NOT.( wntsnn .OR. wntsne ) )
421  $ maxwrk = max( maxwrk, n*n + 2*n )
422  maxwrk = max( maxwrk, 2*n )
423  END IF
424  maxwrk = max( maxwrk, minwrk )
425  END IF
426  work( 1 ) = maxwrk
427 *
428  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
429  info = -20
430  END IF
431  END IF
432 *
433  IF( info.NE.0 ) THEN
434  CALL xerbla( 'ZGEEVX', -info )
435  return
436  ELSE IF( lquery ) THEN
437  return
438  END IF
439 *
440 * Quick return if possible
441 *
442  IF( n.EQ.0 )
443  $ return
444 *
445 * Get machine constants
446 *
447  eps = dlamch( 'P' )
448  smlnum = dlamch( 'S' )
449  bignum = one / smlnum
450  CALL dlabad( smlnum, bignum )
451  smlnum = sqrt( smlnum ) / eps
452  bignum = one / smlnum
453 *
454 * Scale A if max element outside range [SMLNUM,BIGNUM]
455 *
456  icond = 0
457  anrm = zlange( 'M', n, n, a, lda, dum )
458  scalea = .false.
459  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
460  scalea = .true.
461  cscale = smlnum
462  ELSE IF( anrm.GT.bignum ) THEN
463  scalea = .true.
464  cscale = bignum
465  END IF
466  IF( scalea )
467  $ CALL zlascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
468 *
469 * Balance the matrix and compute ABNRM
470 *
471  CALL zgebal( balanc, n, a, lda, ilo, ihi, scale, ierr )
472  abnrm = zlange( '1', n, n, a, lda, dum )
473  IF( scalea ) THEN
474  dum( 1 ) = abnrm
475  CALL dlascl( 'G', 0, 0, cscale, anrm, 1, 1, dum, 1, ierr )
476  abnrm = dum( 1 )
477  END IF
478 *
479 * Reduce to upper Hessenberg form
480 * (CWorkspace: need 2*N, prefer N+N*NB)
481 * (RWorkspace: none)
482 *
483  itau = 1
484  iwrk = itau + n
485  CALL zgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
486  $ lwork-iwrk+1, ierr )
487 *
488  IF( wantvl ) THEN
489 *
490 * Want left eigenvectors
491 * Copy Householder vectors to VL
492 *
493  side = 'L'
494  CALL zlacpy( 'L', n, n, a, lda, vl, ldvl )
495 *
496 * Generate unitary matrix in VL
497 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
498 * (RWorkspace: none)
499 *
500  CALL zunghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
501  $ lwork-iwrk+1, ierr )
502 *
503 * Perform QR iteration, accumulating Schur vectors in VL
504 * (CWorkspace: need 1, prefer HSWORK (see comments) )
505 * (RWorkspace: none)
506 *
507  iwrk = itau
508  CALL zhseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,
509  $ work( iwrk ), lwork-iwrk+1, info )
510 *
511  IF( wantvr ) THEN
512 *
513 * Want left and right eigenvectors
514 * Copy Schur vectors to VR
515 *
516  side = 'B'
517  CALL zlacpy( 'F', n, n, vl, ldvl, vr, ldvr )
518  END IF
519 *
520  ELSE IF( wantvr ) THEN
521 *
522 * Want right eigenvectors
523 * Copy Householder vectors to VR
524 *
525  side = 'R'
526  CALL zlacpy( 'L', n, n, a, lda, vr, ldvr )
527 *
528 * Generate unitary matrix in VR
529 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
530 * (RWorkspace: none)
531 *
532  CALL zunghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
533  $ lwork-iwrk+1, ierr )
534 *
535 * Perform QR iteration, accumulating Schur vectors in VR
536 * (CWorkspace: need 1, prefer HSWORK (see comments) )
537 * (RWorkspace: none)
538 *
539  iwrk = itau
540  CALL zhseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,
541  $ work( iwrk ), lwork-iwrk+1, info )
542 *
543  ELSE
544 *
545 * Compute eigenvalues only
546 * If condition numbers desired, compute Schur form
547 *
548  IF( wntsnn ) THEN
549  job = 'E'
550  ELSE
551  job = 'S'
552  END IF
553 *
554 * (CWorkspace: need 1, prefer HSWORK (see comments) )
555 * (RWorkspace: none)
556 *
557  iwrk = itau
558  CALL zhseqr( job, 'N', n, ilo, ihi, a, lda, w, vr, ldvr,
559  $ work( iwrk ), lwork-iwrk+1, info )
560  END IF
561 *
562 * If INFO > 0 from ZHSEQR, then quit
563 *
564  IF( info.GT.0 )
565  $ go to 50
566 *
567  IF( wantvl .OR. wantvr ) THEN
568 *
569 * Compute left and/or right eigenvectors
570 * (CWorkspace: need 2*N)
571 * (RWorkspace: need N)
572 *
573  CALL ztrevc( side, 'B', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
574  $ n, nout, work( iwrk ), rwork, ierr )
575  END IF
576 *
577 * Compute condition numbers if desired
578 * (CWorkspace: need N*N+2*N unless SENSE = 'E')
579 * (RWorkspace: need 2*N unless SENSE = 'E')
580 *
581  IF( .NOT.wntsnn ) THEN
582  CALL ztrsna( sense, 'A', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
583  $ rconde, rcondv, n, nout, work( iwrk ), n, rwork,
584  $ icond )
585  END IF
586 *
587  IF( wantvl ) THEN
588 *
589 * Undo balancing of left eigenvectors
590 *
591  CALL zgebak( balanc, 'L', n, ilo, ihi, scale, n, vl, ldvl,
592  $ ierr )
593 *
594 * Normalize left eigenvectors and make largest component real
595 *
596  DO 20 i = 1, n
597  scl = one / dznrm2( n, vl( 1, i ), 1 )
598  CALL zdscal( n, scl, vl( 1, i ), 1 )
599  DO 10 k = 1, n
600  rwork( k ) = dble( vl( k, i ) )**2 +
601  $ dimag( vl( k, i ) )**2
602  10 continue
603  k = idamax( n, rwork, 1 )
604  tmp = dconjg( vl( k, i ) ) / sqrt( rwork( k ) )
605  CALL zscal( n, tmp, vl( 1, i ), 1 )
606  vl( k, i ) = dcmplx( dble( vl( k, i ) ), zero )
607  20 continue
608  END IF
609 *
610  IF( wantvr ) THEN
611 *
612 * Undo balancing of right eigenvectors
613 *
614  CALL zgebak( balanc, 'R', n, ilo, ihi, scale, n, vr, ldvr,
615  $ ierr )
616 *
617 * Normalize right eigenvectors and make largest component real
618 *
619  DO 40 i = 1, n
620  scl = one / dznrm2( n, vr( 1, i ), 1 )
621  CALL zdscal( n, scl, vr( 1, i ), 1 )
622  DO 30 k = 1, n
623  rwork( k ) = dble( vr( k, i ) )**2 +
624  $ dimag( vr( k, i ) )**2
625  30 continue
626  k = idamax( n, rwork, 1 )
627  tmp = dconjg( vr( k, i ) ) / sqrt( rwork( k ) )
628  CALL zscal( n, tmp, vr( 1, i ), 1 )
629  vr( k, i ) = dcmplx( dble( vr( k, i ) ), zero )
630  40 continue
631  END IF
632 *
633 * Undo scaling if necessary
634 *
635  50 continue
636  IF( scalea ) THEN
637  CALL zlascl( 'G', 0, 0, cscale, anrm, n-info, 1, w( info+1 ),
638  $ max( n-info, 1 ), ierr )
639  IF( info.EQ.0 ) THEN
640  IF( ( wntsnv .OR. wntsnb ) .AND. icond.EQ.0 )
641  $ CALL dlascl( 'G', 0, 0, cscale, anrm, n, 1, rcondv, n,
642  $ ierr )
643  ELSE
644  CALL zlascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, w, n, ierr )
645  END IF
646  END IF
647 *
648  work( 1 ) = maxwrk
649  return
650 *
651 * End of ZGEEVX
652 *
653  END