LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zgegv.f
Go to the documentation of this file.
1*> \brief <b> ZGEGV computes the eigenvalues and, optionally, the left and/or right eigenvectors of a complex matrix pair (A,B).</b>
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download ZGEGV + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgegv.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgegv.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgegv.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE ZGEGV( 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* DOUBLE PRECISION RWORK( * )
30* COMPLEX*16 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 ZGGEV.
42*>
43*> ZGEGV 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*16 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 ZGGHRD and ZHGEQZ
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*16 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 ZGGHRD and ZHGEQZ 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*16 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*16 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*16 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*16 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*16 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 ZGEQRF, ZUNMQR, and ZUNGQR.) Then compute:
200*> NB -- MAX of the blocksizes for ZGEQRF, ZUNMQR, and ZUNGQR;
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 DOUBLE PRECISION 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 ZGGBAL
225*> =N+2: error return from ZGEQRF
226*> =N+3: error return from ZUNMQR
227*> =N+4: error return from ZUNGQR
228*> =N+5: error return from ZGGHRD
229*> =N+6: error return from ZHGEQZ (other than failed
230*> iteration)
231*> =N+7: error return from ZTGEVC
232*> =N+8: error return from ZGGBAK (computing VL)
233*> =N+9: error return from ZGGBAK (computing VR)
234*> =N+10: error return from ZLASCL (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 complex16GEeigen
246*
247*> \par Further Details:
248* =====================
249*>
250*> \verbatim
251*>
252*> Balancing
253*> ---------
254*>
255*> This driver calls ZGGBAL 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, ZGGBAK 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 zgegv( 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 DOUBLE PRECISION RWORK( * )
293 COMPLEX*16 A( LDA, * ), ALPHA( * ), B( LDB, * ),
294 $ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
295 $ work( * )
296* ..
297*
298* =====================================================================
299*
300* .. Parameters ..
301 DOUBLE PRECISION ZERO, ONE
302 parameter( zero = 0.0d0, one = 1.0d0 )
303 COMPLEX*16 CZERO, CONE
304 parameter( czero = ( 0.0d0, 0.0d0 ),
305 $ cone = ( 1.0d0, 0.0d0 ) )
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 DOUBLE PRECISION ABSAI, ABSAR, ABSB, ANRM, ANRM1, ANRM2, BNRM,
314 $ bnrm1, bnrm2, eps, safmax, safmin, salfai,
315 $ salfar, sbeta, scale, temp
316 COMPLEX*16 X
317* ..
318* .. Local Arrays ..
319 LOGICAL LDUMMA( 1 )
320* ..
321* .. External Subroutines ..
322 EXTERNAL xerbla, zgeqrf, zggbak, zggbal, zgghrd, zhgeqz,
324* ..
325* .. External Functions ..
326 LOGICAL LSAME
327 INTEGER ILAENV
328 DOUBLE PRECISION DLAMCH, ZLANGE
329 EXTERNAL lsame, ilaenv, dlamch, zlange
330* ..
331* .. Intrinsic Functions ..
332 INTRINSIC abs, dble, dcmplx, dimag, int, max
333* ..
334* .. Statement Functions ..
335 DOUBLE PRECISION ABS1
336* ..
337* .. Statement Function definitions ..
338 abs1( x ) = abs( dble( x ) ) + abs( dimag( 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, 'ZGEQRF', ' ', n, n, -1, -1 )
394 nb2 = ilaenv( 1, 'ZUNMQR', ' ', n, n, n, -1 )
395 nb3 = ilaenv( 1, 'ZUNGQR', ' ', 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( 'ZGEGV ', -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 = dlamch( 'E' )*dlamch( 'B' )
416 safmin = dlamch( 'S' )
417 safmin = safmin + safmin
418 safmax = one / safmin
419*
420* Scale A
421*
422 anrm = zlange( '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 zlascl( '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 = zlange( '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 zlascl( '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 zggbal( '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 zgeqrf( 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 zunmqr( '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 zlaset( 'Full', n, n, czero, cone, vl, ldvl )
504 CALL zlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
505 $ vl( ilo+1, ilo ), ldvl )
506 CALL zungqr( 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 zlaset( '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 zgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
527 $ ldvl, vr, ldvr, iinfo )
528 ELSE
529 CALL zgghrd( '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 zhgeqz( 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 ztgevc( 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 zggbak( '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 zggbak( '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( dble( alpha( jc ) ) )
640 absai = abs( dimag( alpha( jc ) ) )
641 absb = abs( dble( beta( jc ) ) )
642 salfar = anrm*dble( alpha( jc ) )
643 salfai = anrm*dimag( alpha( jc ) )
644 sbeta = bnrm*dble( 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*dble( alpha( jc ) ) )*anrm
689 salfai = ( scale*dimag( alpha( jc ) ) )*anrm
690 sbeta = ( scale*beta( jc ) )*bnrm
691 END IF
692 alpha( jc ) = dcmplx( salfar, salfai )
693 beta( jc ) = sbeta
694 70 CONTINUE
695*
696 80 CONTINUE
697 work( 1 ) = lwkopt
698*
699 RETURN
700*
701* End of ZGEGV
702*
703 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine zgeqrf(m, n, a, lda, tau, work, lwork, info)
ZGEQRF
Definition zgeqrf.f:146
subroutine zggbak(job, side, n, ilo, ihi, lscale, rscale, m, v, ldv, info)
ZGGBAK
Definition zggbak.f:148
subroutine zggbal(job, n, a, lda, b, ldb, ilo, ihi, lscale, rscale, work, info)
ZGGBAL
Definition zggbal.f:177
subroutine zgghrd(compq, compz, n, ilo, ihi, a, lda, b, ldb, q, ldq, z, ldz, info)
ZGGHRD
Definition zgghrd.f:204
subroutine zhgeqz(job, compq, compz, n, ilo, ihi, h, ldh, t, ldt, alpha, beta, q, ldq, z, ldz, work, lwork, rwork, info)
ZHGEQZ
Definition zhgeqz.f:284
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:103
subroutine zlascl(type, kl, ku, cfrom, cto, m, n, a, lda, info)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition zlascl.f:143
subroutine zlaset(uplo, m, n, alpha, beta, a, lda)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition zlaset.f:106
subroutine ztgevc(side, howmny, select, n, s, lds, p, ldp, vl, ldvl, vr, ldvr, mm, m, work, rwork, info)
ZTGEVC
Definition ztgevc.f:219
subroutine zungqr(m, n, k, a, lda, tau, work, lwork, info)
ZUNGQR
Definition zungqr.f:128
subroutine zunmqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
ZUNMQR
Definition zunmqr.f:167
subroutine zgegv(jobvl, jobvr, n, a, lda, b, ldb, alpha, beta, vl, ldvl, vr, ldvr, work, lwork, rwork, info)
ZGEGV computes the eigenvalues and, optionally, the left and/or right eigenvectors of a complex matri...
Definition zgegv.f:282