LAPACK  3.10.0 LAPACK: Linear Algebra PACKage

◆ 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.

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
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:151
Here is the call graph for this function:
Here is the caller graph for this function: