LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
cgegv.f
Go to the documentation of this file.
1 *> \brief <b> CGEEVX 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 CGEGV + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cgegv.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cgegv.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cgegv.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE CGEGV( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA,
22 * VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO )
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER JOBVL, JOBVR
26 * INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
27 * ..
28 * .. Array Arguments ..
29 * REAL RWORK( * )
30 * COMPLEX A( LDA, * ), ALPHA( * ), B( LDB, * ),
31 * $ BETA( * ), VL( LDVL, * ), VR( LDVR, * ),
32 * $ WORK( * )
33 * ..
34 *
35 *
36 *> \par Purpose:
37 * =============
38 *>
39 *> \verbatim
40 *>
41 *> This routine is deprecated and has been replaced by routine CGGEV.
42 *>
43 *> CGEGV computes the eigenvalues and, optionally, the left and/or right
44 *> eigenvectors of a complex matrix pair (A,B).
45 *> Given two square matrices A and B,
46 *> the generalized nonsymmetric eigenvalue problem (GNEP) is to find the
47 *> eigenvalues lambda and corresponding (non-zero) eigenvectors x such
48 *> that
49 *> A*x = lambda*B*x.
50 *>
51 *> An alternate form is to find the eigenvalues mu and corresponding
52 *> eigenvectors y such that
53 *> mu*A*y = B*y.
54 *>
55 *> These two forms are equivalent with mu = 1/lambda and x = y if
56 *> neither lambda nor mu is zero. In order to deal with the case that
57 *> lambda or mu is zero or small, two values alpha and beta are returned
58 *> for each eigenvalue, such that lambda = alpha/beta and
59 *> mu = beta/alpha.
60 *>
61 *> The vectors x and y in the above equations are right eigenvectors of
62 *> the matrix pair (A,B). Vectors u and v satisfying
63 *> u**H*A = lambda*u**H*B or mu*v**H*A = v**H*B
64 *> are left eigenvectors of (A,B).
65 *>
66 *> Note: this routine performs "full balancing" on A and B
67 *> \endverbatim
68 *
69 * Arguments:
70 * ==========
71 *
72 *> \param[in] JOBVL
73 *> \verbatim
74 *> JOBVL is CHARACTER*1
75 *> = 'N': do not compute the left generalized eigenvectors;
76 *> = 'V': compute the left generalized eigenvectors (returned
77 *> in VL).
78 *> \endverbatim
79 *>
80 *> \param[in] JOBVR
81 *> \verbatim
82 *> JOBVR is CHARACTER*1
83 *> = 'N': do not compute the right generalized eigenvectors;
84 *> = 'V': compute the right generalized eigenvectors (returned
85 *> in VR).
86 *> \endverbatim
87 *>
88 *> \param[in] N
89 *> \verbatim
90 *> N is INTEGER
91 *> The order of the matrices A, B, VL, and VR. N >= 0.
92 *> \endverbatim
93 *>
94 *> \param[in,out] A
95 *> \verbatim
96 *> A is COMPLEX array, dimension (LDA, N)
97 *> On entry, the matrix A.
98 *> If JOBVL = 'V' or JOBVR = 'V', then on exit A
99 *> contains the Schur form of A from the generalized Schur
100 *> factorization of the pair (A,B) after balancing. If no
101 *> eigenvectors were computed, then only the diagonal elements
102 *> of the Schur form will be correct. See CGGHRD and CHGEQZ
103 *> for details.
104 *> \endverbatim
105 *>
106 *> \param[in] LDA
107 *> \verbatim
108 *> LDA is INTEGER
109 *> The leading dimension of A. LDA >= max(1,N).
110 *> \endverbatim
111 *>
112 *> \param[in,out] B
113 *> \verbatim
114 *> B is COMPLEX array, dimension (LDB, N)
115 *> On entry, the matrix B.
116 *> If JOBVL = 'V' or JOBVR = 'V', then on exit B contains the
117 *> upper triangular matrix obtained from B in the generalized
118 *> Schur factorization of the pair (A,B) after balancing.
119 *> If no eigenvectors were computed, then only the diagonal
120 *> elements of B will be correct. See CGGHRD and CHGEQZ for
121 *> details.
122 *> \endverbatim
123 *>
124 *> \param[in] LDB
125 *> \verbatim
126 *> LDB is INTEGER
127 *> The leading dimension of B. LDB >= max(1,N).
128 *> \endverbatim
129 *>
130 *> \param[out] ALPHA
131 *> \verbatim
132 *> ALPHA is COMPLEX array, dimension (N)
133 *> The complex scalars alpha that define the eigenvalues of
134 *> GNEP.
135 *> \endverbatim
136 *>
137 *> \param[out] BETA
138 *> \verbatim
139 *> BETA is COMPLEX array, dimension (N)
140 *> The complex scalars beta that define the eigenvalues of GNEP.
141 *>
142 *> Together, the quantities alpha = ALPHA(j) and beta = BETA(j)
143 *> represent the j-th eigenvalue of the matrix pair (A,B), in
144 *> one of the forms lambda = alpha/beta or mu = beta/alpha.
145 *> Since either lambda or mu may overflow, they should not,
146 *> in general, be computed.
147 *> \endverbatim
148 *>
149 *> \param[out] VL
150 *> \verbatim
151 *> VL is COMPLEX array, dimension (LDVL,N)
152 *> If JOBVL = 'V', the left eigenvectors u(j) are stored
153 *> in the columns of VL, in the same order as their eigenvalues.
154 *> Each eigenvector is scaled so that its largest component has
155 *> abs(real part) + abs(imag. part) = 1, except for eigenvectors
156 *> corresponding to an eigenvalue with alpha = beta = 0, which
157 *> are set to zero.
158 *> Not referenced if JOBVL = 'N'.
159 *> \endverbatim
160 *>
161 *> \param[in] LDVL
162 *> \verbatim
163 *> LDVL is INTEGER
164 *> The leading dimension of the matrix VL. LDVL >= 1, and
165 *> if JOBVL = 'V', LDVL >= N.
166 *> \endverbatim
167 *>
168 *> \param[out] VR
169 *> \verbatim
170 *> VR is COMPLEX array, dimension (LDVR,N)
171 *> If JOBVR = 'V', the right eigenvectors x(j) are stored
172 *> in the columns of VR, in the same order as their eigenvalues.
173 *> Each eigenvector is scaled so that its largest component has
174 *> abs(real part) + abs(imag. part) = 1, except for eigenvectors
175 *> corresponding to an eigenvalue with alpha = beta = 0, which
176 *> are set to zero.
177 *> Not referenced if JOBVR = 'N'.
178 *> \endverbatim
179 *>
180 *> \param[in] LDVR
181 *> \verbatim
182 *> LDVR is INTEGER
183 *> The leading dimension of the matrix VR. LDVR >= 1, and
184 *> if JOBVR = 'V', LDVR >= N.
185 *> \endverbatim
186 *>
187 *> \param[out] WORK
188 *> \verbatim
189 *> WORK is COMPLEX array, dimension (MAX(1,LWORK))
190 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
191 *> \endverbatim
192 *>
193 *> \param[in] LWORK
194 *> \verbatim
195 *> LWORK is INTEGER
196 *> The dimension of the array WORK. LWORK >= max(1,2*N).
197 *> For good performance, LWORK must generally be larger.
198 *> To compute the optimal value of LWORK, call ILAENV to get
199 *> blocksizes (for CGEQRF, CUNMQR, and CUNGQR.) Then compute:
200 *> NB -- MAX of the blocksizes for CGEQRF, CUNMQR, and CUNGQR;
201 *> The optimal LWORK is MAX( 2*N, N*(NB+1) ).
202 *>
203 *> If LWORK = -1, then a workspace query is assumed; the routine
204 *> only calculates the optimal size of the WORK array, returns
205 *> this value as the first entry of the WORK array, and no error
206 *> message related to LWORK is issued by XERBLA.
207 *> \endverbatim
208 *>
209 *> \param[out] RWORK
210 *> \verbatim
211 *> RWORK is REAL array, dimension (8*N)
212 *> \endverbatim
213 *>
214 *> \param[out] INFO
215 *> \verbatim
216 *> INFO is INTEGER
217 *> = 0: successful exit
218 *> < 0: if INFO = -i, the i-th argument had an illegal value.
219 *> =1,...,N:
220 *> The QZ iteration failed. No eigenvectors have been
221 *> calculated, but ALPHA(j) and BETA(j) should be
222 *> correct for j=INFO+1,...,N.
223 *> > N: errors that usually indicate LAPACK problems:
224 *> =N+1: error return from CGGBAL
225 *> =N+2: error return from CGEQRF
226 *> =N+3: error return from CUNMQR
227 *> =N+4: error return from CUNGQR
228 *> =N+5: error return from CGGHRD
229 *> =N+6: error return from CHGEQZ (other than failed
230 *> iteration)
231 *> =N+7: error return from CTGEVC
232 *> =N+8: error return from CGGBAK (computing VL)
233 *> =N+9: error return from CGGBAK (computing VR)
234 *> =N+10: error return from CLASCL (various calls)
235 *> \endverbatim
236 *
237 * Authors:
238 * ========
239 *
240 *> \author Univ. of Tennessee
241 *> \author Univ. of California Berkeley
242 *> \author Univ. of Colorado Denver
243 *> \author NAG Ltd.
244 *
245 *> \ingroup complexGEeigen
246 *
247 *> \par Further Details:
248 * =====================
249 *>
250 *> \verbatim
251 *>
252 *> Balancing
253 *> ---------
254 *>
255 *> This driver calls CGGBAL to both permute and scale rows and columns
256 *> of A and B. The permutations PL and PR are chosen so that PL*A*PR
257 *> and PL*B*R will be upper triangular except for the diagonal blocks
258 *> A(i:j,i:j) and B(i:j,i:j), with i and j as close together as
259 *> possible. The diagonal scaling matrices DL and DR are chosen so
260 *> that the pair DL*PL*A*PR*DR, DL*PL*B*PR*DR have elements close to
261 *> one (except for the elements that start out zero.)
262 *>
263 *> After the eigenvalues and eigenvectors of the balanced matrices
264 *> have been computed, CGGBAK transforms the eigenvectors back to what
265 *> they would have been (in perfect arithmetic) if they had not been
266 *> balanced.
267 *>
268 *> Contents of A and B on Exit
269 *> -------- -- - --- - -- ----
270 *>
271 *> If any eigenvectors are computed (either JOBVL='V' or JOBVR='V' or
272 *> both), then on exit the arrays A and B will contain the complex Schur
273 *> form[*] of the "balanced" versions of A and B. If no eigenvectors
274 *> are computed, then only the diagonal blocks will be correct.
275 *>
276 *> [*] In other words, upper triangular form.
277 *> \endverbatim
278 *>
279 * =====================================================================
280  SUBROUTINE cgegv( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA,
281  $ VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO )
282 *
283 * -- LAPACK driver routine --
284 * -- LAPACK is a software package provided by Univ. of Tennessee, --
285 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
286 *
287 * .. Scalar Arguments ..
288  CHARACTER JOBVL, JOBVR
289  INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
290 * ..
291 * .. Array Arguments ..
292  REAL RWORK( * )
293  COMPLEX A( LDA, * ), ALPHA( * ), B( LDB, * ),
294  $ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
295  $ work( * )
296 * ..
297 *
298 * =====================================================================
299 *
300 * .. Parameters ..
301  REAL ZERO, ONE
302  parameter( zero = 0.0e0, one = 1.0e0 )
303  COMPLEX CZERO, CONE
304  parameter( czero = ( 0.0e0, 0.0e0 ),
305  $ cone = ( 1.0e0, 0.0e0 ) )
306 * ..
307 * .. Local Scalars ..
308  LOGICAL ILIMIT, ILV, ILVL, ILVR, LQUERY
309  CHARACTER CHTEMP
310  INTEGER ICOLS, IHI, IINFO, IJOBVL, IJOBVR, ILEFT, ILO,
311  $ in, iright, irows, irwork, itau, iwork, jc, jr,
312  $ lopt, lwkmin, lwkopt, nb, nb1, nb2, nb3
313  REAL ABSAI, ABSAR, ABSB, ANRM, ANRM1, ANRM2, BNRM,
314  $ bnrm1, bnrm2, eps, safmax, safmin, salfai,
315  $ salfar, sbeta, scale, temp
316  COMPLEX X
317 * ..
318 * .. Local Arrays ..
319  LOGICAL LDUMMA( 1 )
320 * ..
321 * .. External Subroutines ..
322  EXTERNAL cgeqrf, cggbak, cggbal, cgghrd, chgeqz, clacpy,
324 * ..
325 * .. External Functions ..
326  LOGICAL LSAME
327  INTEGER ILAENV
328  REAL CLANGE, SLAMCH
329  EXTERNAL ilaenv, lsame, clange, slamch
330 * ..
331 * .. Intrinsic Functions ..
332  INTRINSIC abs, aimag, cmplx, int, max, real
333 * ..
334 * .. Statement Functions ..
335  REAL ABS1
336 * ..
337 * .. Statement Function definitions ..
338  abs1( x ) = abs( real( x ) ) + abs( aimag( x ) )
339 * ..
340 * .. Executable Statements ..
341 *
342 * Decode the input arguments
343 *
344  IF( lsame( jobvl, 'N' ) ) THEN
345  ijobvl = 1
346  ilvl = .false.
347  ELSE IF( lsame( jobvl, 'V' ) ) THEN
348  ijobvl = 2
349  ilvl = .true.
350  ELSE
351  ijobvl = -1
352  ilvl = .false.
353  END IF
354 *
355  IF( lsame( jobvr, 'N' ) ) THEN
356  ijobvr = 1
357  ilvr = .false.
358  ELSE IF( lsame( jobvr, 'V' ) ) THEN
359  ijobvr = 2
360  ilvr = .true.
361  ELSE
362  ijobvr = -1
363  ilvr = .false.
364  END IF
365  ilv = ilvl .OR. ilvr
366 *
367 * Test the input arguments
368 *
369  lwkmin = max( 2*n, 1 )
370  lwkopt = lwkmin
371  work( 1 ) = lwkopt
372  lquery = ( lwork.EQ.-1 )
373  info = 0
374  IF( ijobvl.LE.0 ) THEN
375  info = -1
376  ELSE IF( ijobvr.LE.0 ) THEN
377  info = -2
378  ELSE IF( n.LT.0 ) THEN
379  info = -3
380  ELSE IF( lda.LT.max( 1, n ) ) THEN
381  info = -5
382  ELSE IF( ldb.LT.max( 1, n ) ) THEN
383  info = -7
384  ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
385  info = -11
386  ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
387  info = -13
388  ELSE IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
389  info = -15
390  END IF
391 *
392  IF( info.EQ.0 ) THEN
393  nb1 = ilaenv( 1, 'CGEQRF', ' ', n, n, -1, -1 )
394  nb2 = ilaenv( 1, 'CUNMQR', ' ', n, n, n, -1 )
395  nb3 = ilaenv( 1, 'CUNGQR', ' ', n, n, n, -1 )
396  nb = max( nb1, nb2, nb3 )
397  lopt = max( 2*n, n*(nb+1) )
398  work( 1 ) = lopt
399  END IF
400 *
401  IF( info.NE.0 ) THEN
402  CALL xerbla( 'CGEGV ', -info )
403  RETURN
404  ELSE IF( lquery ) THEN
405  RETURN
406  END IF
407 *
408 * Quick return if possible
409 *
410  IF( n.EQ.0 )
411  $ RETURN
412 *
413 * Get machine constants
414 *
415  eps = slamch( 'E' )*slamch( 'B' )
416  safmin = slamch( 'S' )
417  safmin = safmin + safmin
418  safmax = one / safmin
419 *
420 * Scale A
421 *
422  anrm = clange( 'M', n, n, a, lda, rwork )
423  anrm1 = anrm
424  anrm2 = one
425  IF( anrm.LT.one ) THEN
426  IF( safmax*anrm.LT.one ) THEN
427  anrm1 = safmin
428  anrm2 = safmax*anrm
429  END IF
430  END IF
431 *
432  IF( anrm.GT.zero ) THEN
433  CALL clascl( 'G', -1, -1, anrm, one, n, n, a, lda, iinfo )
434  IF( iinfo.NE.0 ) THEN
435  info = n + 10
436  RETURN
437  END IF
438  END IF
439 *
440 * Scale B
441 *
442  bnrm = clange( 'M', n, n, b, ldb, rwork )
443  bnrm1 = bnrm
444  bnrm2 = one
445  IF( bnrm.LT.one ) THEN
446  IF( safmax*bnrm.LT.one ) THEN
447  bnrm1 = safmin
448  bnrm2 = safmax*bnrm
449  END IF
450  END IF
451 *
452  IF( bnrm.GT.zero ) THEN
453  CALL clascl( 'G', -1, -1, bnrm, one, n, n, b, ldb, iinfo )
454  IF( iinfo.NE.0 ) THEN
455  info = n + 10
456  RETURN
457  END IF
458  END IF
459 *
460 * Permute the matrix to make it more nearly triangular
461 * Also "balance" the matrix.
462 *
463  ileft = 1
464  iright = n + 1
465  irwork = iright + n
466  CALL cggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
467  $ rwork( iright ), rwork( irwork ), iinfo )
468  IF( iinfo.NE.0 ) THEN
469  info = n + 1
470  GO TO 80
471  END IF
472 *
473 * Reduce B to triangular form, and initialize VL and/or VR
474 *
475  irows = ihi + 1 - ilo
476  IF( ilv ) THEN
477  icols = n + 1 - ilo
478  ELSE
479  icols = irows
480  END IF
481  itau = 1
482  iwork = itau + irows
483  CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
484  $ work( iwork ), lwork+1-iwork, iinfo )
485  IF( iinfo.GE.0 )
486  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
487  IF( iinfo.NE.0 ) THEN
488  info = n + 2
489  GO TO 80
490  END IF
491 *
492  CALL cunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
493  $ work( itau ), a( ilo, ilo ), lda, work( iwork ),
494  $ lwork+1-iwork, iinfo )
495  IF( iinfo.GE.0 )
496  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
497  IF( iinfo.NE.0 ) THEN
498  info = n + 3
499  GO TO 80
500  END IF
501 *
502  IF( ilvl ) THEN
503  CALL claset( 'Full', n, n, czero, cone, vl, ldvl )
504  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
505  $ vl( ilo+1, ilo ), ldvl )
506  CALL cungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
507  $ work( itau ), work( iwork ), lwork+1-iwork,
508  $ iinfo )
509  IF( iinfo.GE.0 )
510  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
511  IF( iinfo.NE.0 ) THEN
512  info = n + 4
513  GO TO 80
514  END IF
515  END IF
516 *
517  IF( ilvr )
518  $ CALL claset( 'Full', n, n, czero, cone, vr, ldvr )
519 *
520 * Reduce to generalized Hessenberg form
521 *
522  IF( ilv ) THEN
523 *
524 * Eigenvectors requested -- work on whole matrix.
525 *
526  CALL cgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
527  $ ldvl, vr, ldvr, iinfo )
528  ELSE
529  CALL cgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
530  $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, iinfo )
531  END IF
532  IF( iinfo.NE.0 ) THEN
533  info = n + 5
534  GO TO 80
535  END IF
536 *
537 * Perform QZ algorithm
538 *
539  iwork = itau
540  IF( ilv ) THEN
541  chtemp = 'S'
542  ELSE
543  chtemp = 'E'
544  END IF
545  CALL chgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
546  $ alpha, beta, vl, ldvl, vr, ldvr, work( iwork ),
547  $ lwork+1-iwork, rwork( irwork ), iinfo )
548  IF( iinfo.GE.0 )
549  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
550  IF( iinfo.NE.0 ) THEN
551  IF( iinfo.GT.0 .AND. iinfo.LE.n ) THEN
552  info = iinfo
553  ELSE IF( iinfo.GT.n .AND. iinfo.LE.2*n ) THEN
554  info = iinfo - n
555  ELSE
556  info = n + 6
557  END IF
558  GO TO 80
559  END IF
560 *
561  IF( ilv ) THEN
562 *
563 * Compute Eigenvectors
564 *
565  IF( ilvl ) THEN
566  IF( ilvr ) THEN
567  chtemp = 'B'
568  ELSE
569  chtemp = 'L'
570  END IF
571  ELSE
572  chtemp = 'R'
573  END IF
574 *
575  CALL ctgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
576  $ vr, ldvr, n, in, work( iwork ), rwork( irwork ),
577  $ iinfo )
578  IF( iinfo.NE.0 ) THEN
579  info = n + 7
580  GO TO 80
581  END IF
582 *
583 * Undo balancing on VL and VR, rescale
584 *
585  IF( ilvl ) THEN
586  CALL cggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
587  $ rwork( iright ), n, vl, ldvl, iinfo )
588  IF( iinfo.NE.0 ) THEN
589  info = n + 8
590  GO TO 80
591  END IF
592  DO 30 jc = 1, n
593  temp = zero
594  DO 10 jr = 1, n
595  temp = max( temp, abs1( vl( jr, jc ) ) )
596  10 CONTINUE
597  IF( temp.LT.safmin )
598  $ GO TO 30
599  temp = one / temp
600  DO 20 jr = 1, n
601  vl( jr, jc ) = vl( jr, jc )*temp
602  20 CONTINUE
603  30 CONTINUE
604  END IF
605  IF( ilvr ) THEN
606  CALL cggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
607  $ rwork( iright ), n, vr, ldvr, iinfo )
608  IF( iinfo.NE.0 ) THEN
609  info = n + 9
610  GO TO 80
611  END IF
612  DO 60 jc = 1, n
613  temp = zero
614  DO 40 jr = 1, n
615  temp = max( temp, abs1( vr( jr, jc ) ) )
616  40 CONTINUE
617  IF( temp.LT.safmin )
618  $ GO TO 60
619  temp = one / temp
620  DO 50 jr = 1, n
621  vr( jr, jc ) = vr( jr, jc )*temp
622  50 CONTINUE
623  60 CONTINUE
624  END IF
625 *
626 * End of eigenvector calculation
627 *
628  END IF
629 *
630 * Undo scaling in alpha, beta
631 *
632 * Note: this does not give the alpha and beta for the unscaled
633 * problem.
634 *
635 * Un-scaling is limited to avoid underflow in alpha and beta
636 * if they are significant.
637 *
638  DO 70 jc = 1, n
639  absar = abs( real( alpha( jc ) ) )
640  absai = abs( aimag( alpha( jc ) ) )
641  absb = abs( real( beta( jc ) ) )
642  salfar = anrm*real( alpha( jc ) )
643  salfai = anrm*aimag( alpha( jc ) )
644  sbeta = bnrm*real( beta( jc ) )
645  ilimit = .false.
646  scale = one
647 *
648 * Check for significant underflow in imaginary part of ALPHA
649 *
650  IF( abs( salfai ).LT.safmin .AND. absai.GE.
651  $ max( safmin, eps*absar, eps*absb ) ) THEN
652  ilimit = .true.
653  scale = ( safmin / anrm1 ) / max( safmin, anrm2*absai )
654  END IF
655 *
656 * Check for significant underflow in real part of ALPHA
657 *
658  IF( abs( salfar ).LT.safmin .AND. absar.GE.
659  $ max( safmin, eps*absai, eps*absb ) ) THEN
660  ilimit = .true.
661  scale = max( scale, ( safmin / anrm1 ) /
662  $ max( safmin, anrm2*absar ) )
663  END IF
664 *
665 * Check for significant underflow in BETA
666 *
667  IF( abs( sbeta ).LT.safmin .AND. absb.GE.
668  $ max( safmin, eps*absar, eps*absai ) ) THEN
669  ilimit = .true.
670  scale = max( scale, ( safmin / bnrm1 ) /
671  $ max( safmin, bnrm2*absb ) )
672  END IF
673 *
674 * Check for possible overflow when limiting scaling
675 *
676  IF( ilimit ) THEN
677  temp = ( scale*safmin )*max( abs( salfar ), abs( salfai ),
678  $ abs( sbeta ) )
679  IF( temp.GT.one )
680  $ scale = scale / temp
681  IF( scale.LT.one )
682  $ ilimit = .false.
683  END IF
684 *
685 * Recompute un-scaled ALPHA, BETA if necessary.
686 *
687  IF( ilimit ) THEN
688  salfar = ( scale*real( alpha( jc ) ) )*anrm
689  salfai = ( scale*aimag( alpha( jc ) ) )*anrm
690  sbeta = ( scale*beta( jc ) )*bnrm
691  END IF
692  alpha( jc ) = cmplx( salfar, salfai )
693  beta( jc ) = sbeta
694  70 CONTINUE
695 *
696  80 CONTINUE
697  work( 1 ) = lwkopt
698 *
699  RETURN
700 *
701 * End of CGEGV
702 *
703  END
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:177
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:148
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:284
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:145
subroutine ctgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTGEVC
Definition: ctgevc.f:219
subroutine cgegv(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO)
CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
Definition: cgegv.f:282
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: claset.f:106
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:143
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:103
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:168
subroutine cgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
CGGHRD
Definition: cgghrd.f:204
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:128