LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ zggev()

subroutine zggev ( character  jobvl,
character  jobvr,
integer  n,
complex*16, dimension( lda, * )  a,
integer  lda,
complex*16, dimension( ldb, * )  b,
integer  ldb,
complex*16, dimension( * )  alpha,
complex*16, dimension( * )  beta,
complex*16, dimension( ldvl, * )  vl,
integer  ldvl,
complex*16, dimension( ldvr, * )  vr,
integer  ldvr,
complex*16, dimension( * )  work,
integer  lwork,
double precision, dimension( * )  rwork,
integer  info 
)

ZGGEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

Download ZGGEV + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 ZGGEV computes for a pair of N-by-N complex nonsymmetric matrices
 (A,B), the generalized eigenvalues, and optionally, the left and/or
 right generalized eigenvectors.

 A generalized eigenvalue for a pair of matrices (A,B) is a scalar
 lambda or a ratio alpha/beta = lambda, such that A - lambda*B is
 singular. It is usually represented as the pair (alpha,beta), as
 there is a reasonable interpretation for beta=0, and even for both
 being zero.

 The right generalized eigenvector v(j) corresponding to the
 generalized eigenvalue lambda(j) of (A,B) satisfies

              A * v(j) = lambda(j) * B * v(j).

 The left generalized eigenvector u(j) corresponding to the
 generalized eigenvalues lambda(j) of (A,B) satisfies

              u(j)**H * A = lambda(j) * u(j)**H * B

 where u(j)**H is the conjugate-transpose of u(j).
Parameters
[in]JOBVL
          JOBVL is CHARACTER*1
          = 'N':  do not compute the left generalized eigenvectors;
          = 'V':  compute the left generalized eigenvectors.
[in]JOBVR
          JOBVR is CHARACTER*1
          = 'N':  do not compute the right generalized eigenvectors;
          = 'V':  compute the right generalized eigenvectors.
[in]N
          N is INTEGER
          The order of the matrices A, B, VL, and VR.  N >= 0.
[in,out]A
          A is COMPLEX*16 array, dimension (LDA, N)
          On entry, the matrix A in the pair (A,B).
          On exit, A has been overwritten.
[in]LDA
          LDA is INTEGER
          The leading dimension of A.  LDA >= max(1,N).
[in,out]B
          B is COMPLEX*16 array, dimension (LDB, N)
          On entry, the matrix B in the pair (A,B).
          On exit, B has been overwritten.
[in]LDB
          LDB is INTEGER
          The leading dimension of B.  LDB >= max(1,N).
[out]ALPHA
          ALPHA is COMPLEX*16 array, dimension (N)
[out]BETA
          BETA is COMPLEX*16 array, dimension (N)
          On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the
          generalized eigenvalues.

          Note: the quotients ALPHA(j)/BETA(j) may easily over- or
          underflow, and BETA(j) may even be zero.  Thus, the user
          should avoid naively computing the ratio alpha/beta.
          However, ALPHA will be always less than and usually
          comparable with norm(A) in magnitude, and BETA always less
          than and usually comparable with norm(B).
[out]VL
          VL is COMPLEX*16 array, dimension (LDVL,N)
          If JOBVL = 'V', the left generalized eigenvectors u(j) are
          stored one after another in the columns of VL, in the same
          order as their eigenvalues.
          Each eigenvector is scaled so the largest component has
          abs(real part) + abs(imag. part) = 1.
          Not referenced if JOBVL = 'N'.
[in]LDVL
          LDVL is INTEGER
          The leading dimension of the matrix VL. LDVL >= 1, and
          if JOBVL = 'V', LDVL >= N.
[out]VR
          VR is COMPLEX*16 array, dimension (LDVR,N)
          If JOBVR = 'V', the right generalized eigenvectors v(j) are
          stored one after another in the columns of VR, in the same
          order as their eigenvalues.
          Each eigenvector is scaled so the largest component has
          abs(real part) + abs(imag. part) = 1.
          Not referenced if JOBVR = 'N'.
[in]LDVR
          LDVR is INTEGER
          The leading dimension of the matrix VR. LDVR >= 1, and
          if JOBVR = 'V', LDVR >= N.
[out]WORK
          WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
[in]LWORK
          LWORK is INTEGER
          The dimension of the array WORK.  LWORK >= max(1,2*N).
          For good performance, LWORK must generally be larger.

          If LWORK = -1, then a workspace query is assumed; the routine
          only calculates the optimal size of the WORK array, returns
          this value as the first entry of the WORK array, and no error
          message related to LWORK is issued by XERBLA.
[out]RWORK
          RWORK is DOUBLE PRECISION array, dimension (8*N)
[out]INFO
          INFO is INTEGER
          = 0:  successful exit
          < 0:  if INFO = -i, the i-th argument had an illegal value.
          =1,...,N:
                The QZ iteration failed.  No eigenvectors have been
                calculated, but ALPHA(j) and BETA(j) should be
                correct for j=INFO+1,...,N.
          > N:  =N+1: other then QZ iteration failed in ZHGEQZ,
                =N+2: error return from ZTGEVC.
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 215 of file zggev.f.

217*
218* -- LAPACK driver routine --
219* -- LAPACK is a software package provided by Univ. of Tennessee, --
220* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
221*
222* .. Scalar Arguments ..
223 CHARACTER JOBVL, JOBVR
224 INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
225* ..
226* .. Array Arguments ..
227 DOUBLE PRECISION RWORK( * )
228 COMPLEX*16 A( LDA, * ), ALPHA( * ), B( LDB, * ),
229 $ BETA( * ), VL( LDVL, * ), VR( LDVR, * ),
230 $ WORK( * )
231* ..
232*
233* =====================================================================
234*
235* .. Parameters ..
236 DOUBLE PRECISION ZERO, ONE
237 parameter( zero = 0.0d0, one = 1.0d0 )
238 COMPLEX*16 CZERO, CONE
239 parameter( czero = ( 0.0d0, 0.0d0 ),
240 $ cone = ( 1.0d0, 0.0d0 ) )
241* ..
242* .. Local Scalars ..
243 LOGICAL ILASCL, ILBSCL, ILV, ILVL, ILVR, LQUERY
244 CHARACTER CHTEMP
245 INTEGER ICOLS, IERR, IHI, IJOBVL, IJOBVR, ILEFT, ILO,
246 $ IN, IRIGHT, IROWS, IRWRK, ITAU, IWRK, JC, JR,
247 $ LWKMIN, LWKOPT
248 DOUBLE PRECISION ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS,
249 $ SMLNUM, TEMP
250 COMPLEX*16 X
251* ..
252* .. Local Arrays ..
253 LOGICAL LDUMMA( 1 )
254* ..
255* .. External Subroutines ..
256 EXTERNAL xerbla, zgeqrf, zggbak, zggbal, zgghrd, zhgeqz,
258* ..
259* .. External Functions ..
260 LOGICAL LSAME
261 INTEGER ILAENV
262 DOUBLE PRECISION DLAMCH, ZLANGE
263 EXTERNAL lsame, ilaenv, dlamch, zlange
264* ..
265* .. Intrinsic Functions ..
266 INTRINSIC abs, dble, dimag, max, sqrt
267* ..
268* .. Statement Functions ..
269 DOUBLE PRECISION ABS1
270* ..
271* .. Statement Function definitions ..
272 abs1( x ) = abs( dble( x ) ) + abs( dimag( x ) )
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 = -11
317 ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
318 info = -13
319 END IF
320*
321* Compute workspace
322* (Note: Comments in the code beginning "Workspace:" describe the
323* minimal amount of workspace needed at that point in the code,
324* as well as the preferred amount for good performance.
325* NB refers to the optimal block size for the immediately
326* following subroutine, as returned by ILAENV. The workspace is
327* computed assuming ILO = 1 and IHI = N, the worst case.)
328*
329 IF( info.EQ.0 ) THEN
330 lwkmin = max( 1, 2*n )
331 lwkopt = max( 1, n + n*ilaenv( 1, 'ZGEQRF', ' ', n, 1, n, 0 ) )
332 lwkopt = max( lwkopt, n +
333 $ n*ilaenv( 1, 'ZUNMQR', ' ', n, 1, n, 0 ) )
334 IF( ilvl ) THEN
335 lwkopt = max( lwkopt, n +
336 $ n*ilaenv( 1, 'ZUNGQR', ' ', n, 1, n, -1 ) )
337 END IF
338 work( 1 ) = lwkopt
339*
340 IF( lwork.LT.lwkmin .AND. .NOT.lquery )
341 $ info = -15
342 END IF
343*
344 IF( info.NE.0 ) THEN
345 CALL xerbla( 'ZGGEV ', -info )
346 RETURN
347 ELSE IF( lquery ) THEN
348 RETURN
349 END IF
350*
351* Quick return if possible
352*
353 IF( n.EQ.0 )
354 $ RETURN
355*
356* Get machine constants
357*
358 eps = dlamch( 'E' )*dlamch( 'B' )
359 smlnum = dlamch( 'S' )
360 bignum = one / smlnum
361 smlnum = sqrt( smlnum ) / eps
362 bignum = one / smlnum
363*
364* Scale A if max element outside range [SMLNUM,BIGNUM]
365*
366 anrm = zlange( 'M', n, n, a, lda, rwork )
367 ilascl = .false.
368 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
369 anrmto = smlnum
370 ilascl = .true.
371 ELSE IF( anrm.GT.bignum ) THEN
372 anrmto = bignum
373 ilascl = .true.
374 END IF
375 IF( ilascl )
376 $ CALL zlascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
377*
378* Scale B if max element outside range [SMLNUM,BIGNUM]
379*
380 bnrm = zlange( 'M', n, n, b, ldb, rwork )
381 ilbscl = .false.
382 IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
383 bnrmto = smlnum
384 ilbscl = .true.
385 ELSE IF( bnrm.GT.bignum ) THEN
386 bnrmto = bignum
387 ilbscl = .true.
388 END IF
389 IF( ilbscl )
390 $ CALL zlascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
391*
392* Permute the matrices A, B to isolate eigenvalues if possible
393* (Real Workspace: need 6*N)
394*
395 ileft = 1
396 iright = n + 1
397 irwrk = iright + n
398 CALL zggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
399 $ rwork( iright ), rwork( irwrk ), ierr )
400*
401* Reduce B to triangular form (QR decomposition of B)
402* (Complex Workspace: need N, prefer N*NB)
403*
404 irows = ihi + 1 - ilo
405 IF( ilv ) THEN
406 icols = n + 1 - ilo
407 ELSE
408 icols = irows
409 END IF
410 itau = 1
411 iwrk = itau + irows
412 CALL zgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
413 $ work( iwrk ), lwork+1-iwrk, ierr )
414*
415* Apply the orthogonal transformation to matrix A
416* (Complex Workspace: need N, prefer N*NB)
417*
418 CALL zunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
419 $ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
420 $ lwork+1-iwrk, ierr )
421*
422* Initialize VL
423* (Complex Workspace: need N, prefer N*NB)
424*
425 IF( ilvl ) THEN
426 CALL zlaset( 'Full', n, n, czero, cone, vl, ldvl )
427 IF( irows.GT.1 ) THEN
428 CALL zlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
429 $ vl( ilo+1, ilo ), ldvl )
430 END IF
431 CALL zungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
432 $ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
433 END IF
434*
435* Initialize VR
436*
437 IF( ilvr )
438 $ CALL zlaset( 'Full', n, n, czero, cone, vr, ldvr )
439*
440* Reduce to generalized Hessenberg form
441*
442 IF( ilv ) THEN
443*
444* Eigenvectors requested -- work on whole matrix.
445*
446 CALL zgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
447 $ ldvl, vr, ldvr, ierr )
448 ELSE
449 CALL zgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
450 $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, ierr )
451 END IF
452*
453* Perform QZ algorithm (Compute eigenvalues, and optionally, the
454* Schur form and Schur vectors)
455* (Complex Workspace: need N)
456* (Real Workspace: need N)
457*
458 iwrk = itau
459 IF( ilv ) THEN
460 chtemp = 'S'
461 ELSE
462 chtemp = 'E'
463 END IF
464 CALL zhgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
465 $ alpha, beta, vl, ldvl, vr, ldvr, work( iwrk ),
466 $ lwork+1-iwrk, rwork( irwrk ), ierr )
467 IF( ierr.NE.0 ) THEN
468 IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
469 info = ierr
470 ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
471 info = ierr - n
472 ELSE
473 info = n + 1
474 END IF
475 GO TO 70
476 END IF
477*
478* Compute Eigenvectors
479* (Real Workspace: need 2*N)
480* (Complex Workspace: need 2*N)
481*
482 IF( ilv ) THEN
483 IF( ilvl ) THEN
484 IF( ilvr ) THEN
485 chtemp = 'B'
486 ELSE
487 chtemp = 'L'
488 END IF
489 ELSE
490 chtemp = 'R'
491 END IF
492*
493 CALL ztgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
494 $ vr, ldvr, n, in, work( iwrk ), rwork( irwrk ),
495 $ ierr )
496 IF( ierr.NE.0 ) THEN
497 info = n + 2
498 GO TO 70
499 END IF
500*
501* Undo balancing on VL and VR and normalization
502* (Workspace: none needed)
503*
504 IF( ilvl ) THEN
505 CALL zggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
506 $ rwork( iright ), n, vl, ldvl, ierr )
507 DO 30 jc = 1, n
508 temp = zero
509 DO 10 jr = 1, n
510 temp = max( temp, abs1( vl( jr, jc ) ) )
511 10 CONTINUE
512 IF( temp.LT.smlnum )
513 $ GO TO 30
514 temp = one / temp
515 DO 20 jr = 1, n
516 vl( jr, jc ) = vl( jr, jc )*temp
517 20 CONTINUE
518 30 CONTINUE
519 END IF
520 IF( ilvr ) THEN
521 CALL zggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
522 $ rwork( iright ), n, vr, ldvr, ierr )
523 DO 60 jc = 1, n
524 temp = zero
525 DO 40 jr = 1, n
526 temp = max( temp, abs1( vr( jr, jc ) ) )
527 40 CONTINUE
528 IF( temp.LT.smlnum )
529 $ GO TO 60
530 temp = one / temp
531 DO 50 jr = 1, n
532 vr( jr, jc ) = vr( jr, jc )*temp
533 50 CONTINUE
534 60 CONTINUE
535 END IF
536 END IF
537*
538* Undo scaling if necessary
539*
540 70 CONTINUE
541*
542 IF( ilascl )
543 $ CALL zlascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
544*
545 IF( ilbscl )
546 $ CALL zlascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
547*
548 work( 1 ) = lwkopt
549 RETURN
550*
551* End of ZGGEV
552*
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
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
Definition ilaenv.f:162
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
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
double precision function zlange(norm, m, n, a, lda, work)
ZLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition zlange.f:115
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
logical function lsame(ca, cb)
LSAME
Definition lsame.f:48
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
Here is the call graph for this function:
Here is the caller graph for this function: