LAPACK 3.11.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 dlabad, xerbla, zgeqrf, zggbak, zggbal, zgghrd,
258 $ zunmqr
259* ..
260* .. External Functions ..
261 LOGICAL LSAME
262 INTEGER ILAENV
263 DOUBLE PRECISION DLAMCH, ZLANGE
264 EXTERNAL lsame, ilaenv, dlamch, zlange
265* ..
266* .. Intrinsic Functions ..
267 INTRINSIC abs, dble, dimag, max, sqrt
268* ..
269* .. Statement Functions ..
270 DOUBLE PRECISION ABS1
271* ..
272* .. Statement Function definitions ..
273 abs1( x ) = abs( dble( x ) ) + abs( dimag( x ) )
274* ..
275* .. Executable Statements ..
276*
277* Decode the input arguments
278*
279 IF( lsame( jobvl, 'N' ) ) THEN
280 ijobvl = 1
281 ilvl = .false.
282 ELSE IF( lsame( jobvl, 'V' ) ) THEN
283 ijobvl = 2
284 ilvl = .true.
285 ELSE
286 ijobvl = -1
287 ilvl = .false.
288 END IF
289*
290 IF( lsame( jobvr, 'N' ) ) THEN
291 ijobvr = 1
292 ilvr = .false.
293 ELSE IF( lsame( jobvr, 'V' ) ) THEN
294 ijobvr = 2
295 ilvr = .true.
296 ELSE
297 ijobvr = -1
298 ilvr = .false.
299 END IF
300 ilv = ilvl .OR. ilvr
301*
302* Test the input arguments
303*
304 info = 0
305 lquery = ( lwork.EQ.-1 )
306 IF( ijobvl.LE.0 ) THEN
307 info = -1
308 ELSE IF( ijobvr.LE.0 ) THEN
309 info = -2
310 ELSE IF( n.LT.0 ) THEN
311 info = -3
312 ELSE IF( lda.LT.max( 1, n ) ) THEN
313 info = -5
314 ELSE IF( ldb.LT.max( 1, n ) ) THEN
315 info = -7
316 ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
317 info = -11
318 ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
319 info = -13
320 END IF
321*
322* Compute workspace
323* (Note: Comments in the code beginning "Workspace:" describe the
324* minimal amount of workspace needed at that point in the code,
325* as well as the preferred amount for good performance.
326* NB refers to the optimal block size for the immediately
327* following subroutine, as returned by ILAENV. The workspace is
328* computed assuming ILO = 1 and IHI = N, the worst case.)
329*
330 IF( info.EQ.0 ) THEN
331 lwkmin = max( 1, 2*n )
332 lwkopt = max( 1, n + n*ilaenv( 1, 'ZGEQRF', ' ', n, 1, n, 0 ) )
333 lwkopt = max( lwkopt, n +
334 $ n*ilaenv( 1, 'ZUNMQR', ' ', n, 1, n, 0 ) )
335 IF( ilvl ) THEN
336 lwkopt = max( lwkopt, n +
337 $ n*ilaenv( 1, 'ZUNGQR', ' ', n, 1, n, -1 ) )
338 END IF
339 work( 1 ) = lwkopt
340*
341 IF( lwork.LT.lwkmin .AND. .NOT.lquery )
342 $ info = -15
343 END IF
344*
345 IF( info.NE.0 ) THEN
346 CALL xerbla( 'ZGGEV ', -info )
347 RETURN
348 ELSE IF( lquery ) THEN
349 RETURN
350 END IF
351*
352* Quick return if possible
353*
354 IF( n.EQ.0 )
355 $ RETURN
356*
357* Get machine constants
358*
359 eps = dlamch( 'E' )*dlamch( 'B' )
360 smlnum = dlamch( 'S' )
361 bignum = one / smlnum
362 CALL dlabad( smlnum, bignum )
363 smlnum = sqrt( smlnum ) / eps
364 bignum = one / smlnum
365*
366* Scale A if max element outside range [SMLNUM,BIGNUM]
367*
368 anrm = zlange( 'M', n, n, a, lda, rwork )
369 ilascl = .false.
370 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
371 anrmto = smlnum
372 ilascl = .true.
373 ELSE IF( anrm.GT.bignum ) THEN
374 anrmto = bignum
375 ilascl = .true.
376 END IF
377 IF( ilascl )
378 $ CALL zlascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
379*
380* Scale B if max element outside range [SMLNUM,BIGNUM]
381*
382 bnrm = zlange( 'M', n, n, b, ldb, rwork )
383 ilbscl = .false.
384 IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
385 bnrmto = smlnum
386 ilbscl = .true.
387 ELSE IF( bnrm.GT.bignum ) THEN
388 bnrmto = bignum
389 ilbscl = .true.
390 END IF
391 IF( ilbscl )
392 $ CALL zlascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
393*
394* Permute the matrices A, B to isolate eigenvalues if possible
395* (Real Workspace: need 6*N)
396*
397 ileft = 1
398 iright = n + 1
399 irwrk = iright + n
400 CALL zggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
401 $ rwork( iright ), rwork( irwrk ), ierr )
402*
403* Reduce B to triangular form (QR decomposition of B)
404* (Complex Workspace: need N, prefer N*NB)
405*
406 irows = ihi + 1 - ilo
407 IF( ilv ) THEN
408 icols = n + 1 - ilo
409 ELSE
410 icols = irows
411 END IF
412 itau = 1
413 iwrk = itau + irows
414 CALL zgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
415 $ work( iwrk ), lwork+1-iwrk, ierr )
416*
417* Apply the orthogonal transformation to matrix A
418* (Complex Workspace: need N, prefer N*NB)
419*
420 CALL zunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
421 $ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
422 $ lwork+1-iwrk, ierr )
423*
424* Initialize VL
425* (Complex Workspace: need N, prefer N*NB)
426*
427 IF( ilvl ) THEN
428 CALL zlaset( 'Full', n, n, czero, cone, vl, ldvl )
429 IF( irows.GT.1 ) THEN
430 CALL zlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
431 $ vl( ilo+1, ilo ), ldvl )
432 END IF
433 CALL zungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
434 $ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
435 END IF
436*
437* Initialize VR
438*
439 IF( ilvr )
440 $ CALL zlaset( 'Full', n, n, czero, cone, vr, ldvr )
441*
442* Reduce to generalized Hessenberg form
443*
444 IF( ilv ) THEN
445*
446* Eigenvectors requested -- work on whole matrix.
447*
448 CALL zgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
449 $ ldvl, vr, ldvr, ierr )
450 ELSE
451 CALL zgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
452 $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, ierr )
453 END IF
454*
455* Perform QZ algorithm (Compute eigenvalues, and optionally, the
456* Schur form and Schur vectors)
457* (Complex Workspace: need N)
458* (Real Workspace: need N)
459*
460 iwrk = itau
461 IF( ilv ) THEN
462 chtemp = 'S'
463 ELSE
464 chtemp = 'E'
465 END IF
466 CALL zhgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
467 $ alpha, beta, vl, ldvl, vr, ldvr, work( iwrk ),
468 $ lwork+1-iwrk, rwork( irwrk ), ierr )
469 IF( ierr.NE.0 ) THEN
470 IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
471 info = ierr
472 ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
473 info = ierr - n
474 ELSE
475 info = n + 1
476 END IF
477 GO TO 70
478 END IF
479*
480* Compute Eigenvectors
481* (Real Workspace: need 2*N)
482* (Complex Workspace: need 2*N)
483*
484 IF( ilv ) THEN
485 IF( ilvl ) THEN
486 IF( ilvr ) THEN
487 chtemp = 'B'
488 ELSE
489 chtemp = 'L'
490 END IF
491 ELSE
492 chtemp = 'R'
493 END IF
494*
495 CALL ztgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
496 $ vr, ldvr, n, in, work( iwrk ), rwork( irwrk ),
497 $ ierr )
498 IF( ierr.NE.0 ) THEN
499 info = n + 2
500 GO TO 70
501 END IF
502*
503* Undo balancing on VL and VR and normalization
504* (Workspace: none needed)
505*
506 IF( ilvl ) THEN
507 CALL zggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
508 $ rwork( iright ), n, vl, ldvl, ierr )
509 DO 30 jc = 1, n
510 temp = zero
511 DO 10 jr = 1, n
512 temp = max( temp, abs1( vl( jr, jc ) ) )
513 10 CONTINUE
514 IF( temp.LT.smlnum )
515 $ GO TO 30
516 temp = one / temp
517 DO 20 jr = 1, n
518 vl( jr, jc ) = vl( jr, jc )*temp
519 20 CONTINUE
520 30 CONTINUE
521 END IF
522 IF( ilvr ) THEN
523 CALL zggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
524 $ rwork( iright ), n, vr, ldvr, ierr )
525 DO 60 jc = 1, n
526 temp = zero
527 DO 40 jr = 1, n
528 temp = max( temp, abs1( vr( jr, jc ) ) )
529 40 CONTINUE
530 IF( temp.LT.smlnum )
531 $ GO TO 60
532 temp = one / temp
533 DO 50 jr = 1, n
534 vr( jr, jc ) = vr( jr, jc )*temp
535 50 CONTINUE
536 60 CONTINUE
537 END IF
538 END IF
539*
540* Undo scaling if necessary
541*
542 70 CONTINUE
543*
544 IF( ilascl )
545 $ CALL zlascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
546*
547 IF( ilbscl )
548 $ CALL zlascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
549*
550 work( 1 ) = lwkopt
551 RETURN
552*
553* End of ZGGEV
554*
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:53
subroutine zggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
ZGGBAL
Definition: zggbal.f:177
subroutine zggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
ZGGBAK
Definition: zggbak.f:148
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 ztgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
ZTGEVC
Definition: ztgevc.f:219
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 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 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 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 zungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
ZUNGQR
Definition: zungqr.f:128
subroutine zgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
ZGGHRD
Definition: zgghrd.f:204
subroutine zunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
ZUNMQR
Definition: zunmqr.f:167
subroutine zgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
ZGEQRF VARIANT: left-looking Level 3 BLAS of the algorithm.
Definition: zgeqrf.f:152
Here is the call graph for this function:
Here is the caller graph for this function: