LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine dgegv ( character  JOBVL,
character  JOBVR,
integer  N,
double precision, dimension( lda, * )  A,
integer  LDA,
double precision, dimension( ldb, * )  B,
integer  LDB,
double precision, dimension( * )  ALPHAR,
double precision, dimension( * )  ALPHAI,
double precision, dimension( * )  BETA,
double precision, dimension( ldvl, * )  VL,
integer  LDVL,
double precision, dimension( ldvr, * )  VR,
integer  LDVR,
double precision, dimension( * )  WORK,
integer  LWORK,
integer  INFO 
)

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

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

Purpose:
 This routine is deprecated and has been replaced by routine DGGEV.

 DGEGV computes the eigenvalues and, optionally, the left and/or right
 eigenvectors of a real matrix pair (A,B).
 Given two square matrices A and B,
 the generalized nonsymmetric eigenvalue problem (GNEP) is to find the
 eigenvalues lambda and corresponding (non-zero) eigenvectors x such
 that

    A*x = lambda*B*x.

 An alternate form is to find the eigenvalues mu and corresponding
 eigenvectors y such that

    mu*A*y = B*y.

 These two forms are equivalent with mu = 1/lambda and x = y if
 neither lambda nor mu is zero.  In order to deal with the case that
 lambda or mu is zero or small, two values alpha and beta are returned
 for each eigenvalue, such that lambda = alpha/beta and
 mu = beta/alpha.

 The vectors x and y in the above equations are right eigenvectors of
 the matrix pair (A,B).  Vectors u and v satisfying

    u**H*A = lambda*u**H*B  or  mu*v**H*A = v**H*B

 are left eigenvectors of (A,B).

 Note: this routine performs "full balancing" on A and B
Parameters
[in]JOBVL
          JOBVL is CHARACTER*1
          = 'N':  do not compute the left generalized eigenvectors;
          = 'V':  compute the left generalized eigenvectors (returned
                  in VL).
[in]JOBVR
          JOBVR is CHARACTER*1
          = 'N':  do not compute the right generalized eigenvectors;
          = 'V':  compute the right generalized eigenvectors (returned
                  in VR).
[in]N
          N is INTEGER
          The order of the matrices A, B, VL, and VR.  N >= 0.
[in,out]A
          A is DOUBLE PRECISION array, dimension (LDA, N)
          On entry, the matrix A.
          If JOBVL = 'V' or JOBVR = 'V', then on exit A
          contains the real Schur form of A from the generalized Schur
          factorization of the pair (A,B) after balancing.
          If no eigenvectors were computed, then only the diagonal
          blocks from the Schur form will be correct.  See DGGHRD and
          DHGEQZ for details.
[in]LDA
          LDA is INTEGER
          The leading dimension of A.  LDA >= max(1,N).
[in,out]B
          B is DOUBLE PRECISION array, dimension (LDB, N)
          On entry, the matrix B.
          If JOBVL = 'V' or JOBVR = 'V', then on exit B contains the
          upper triangular matrix obtained from B in the generalized
          Schur factorization of the pair (A,B) after balancing.
          If no eigenvectors were computed, then only those elements of
          B corresponding to the diagonal blocks from the Schur form of
          A will be correct.  See DGGHRD and DHGEQZ for details.
[in]LDB
          LDB is INTEGER
          The leading dimension of B.  LDB >= max(1,N).
[out]ALPHAR
          ALPHAR is DOUBLE PRECISION array, dimension (N)
          The real parts of each scalar alpha defining an eigenvalue of
          GNEP.
[out]ALPHAI
          ALPHAI is DOUBLE PRECISION array, dimension (N)
          The imaginary parts of each scalar alpha defining an
          eigenvalue of GNEP.  If ALPHAI(j) is zero, then the j-th
          eigenvalue is real; if positive, then the j-th and
          (j+1)-st eigenvalues are a complex conjugate pair, with
          ALPHAI(j+1) = -ALPHAI(j).
[out]BETA
          BETA is DOUBLE PRECISION array, dimension (N)
          The scalars beta that define the eigenvalues of GNEP.
          
          Together, the quantities alpha = (ALPHAR(j),ALPHAI(j)) and
          beta = BETA(j) represent the j-th eigenvalue of the matrix
          pair (A,B), in one of the forms lambda = alpha/beta or
          mu = beta/alpha.  Since either lambda or mu may overflow,
          they should not, in general, be computed.
[out]VL
          VL is DOUBLE PRECISION array, dimension (LDVL,N)
          If JOBVL = 'V', the left eigenvectors u(j) are stored
          in the columns of VL, in the same order as their eigenvalues.
          If the j-th eigenvalue is real, then u(j) = VL(:,j).
          If the j-th and (j+1)-st eigenvalues form a complex conjugate
          pair, then
             u(j) = VL(:,j) + i*VL(:,j+1)
          and
            u(j+1) = VL(:,j) - i*VL(:,j+1).

          Each eigenvector is scaled so that its largest component has
          abs(real part) + abs(imag. part) = 1, except for eigenvectors
          corresponding to an eigenvalue with alpha = beta = 0, which
          are set to zero.
          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 DOUBLE PRECISION array, dimension (LDVR,N)
          If JOBVR = 'V', the right eigenvectors x(j) are stored
          in the columns of VR, in the same order as their eigenvalues.
          If the j-th eigenvalue is real, then x(j) = VR(:,j).
          If the j-th and (j+1)-st eigenvalues form a complex conjugate
          pair, then
            x(j) = VR(:,j) + i*VR(:,j+1)
          and
            x(j+1) = VR(:,j) - i*VR(:,j+1).

          Each eigenvector is scaled so that its largest component has
          abs(real part) + abs(imag. part) = 1, except for eigenvalues
          corresponding to an eigenvalue with alpha = beta = 0, which
          are set to zero.
          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 DOUBLE PRECISION 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,8*N).
          For good performance, LWORK must generally be larger.
          To compute the optimal value of LWORK, call ILAENV to get
          blocksizes (for DGEQRF, DORMQR, and DORGQR.)  Then compute:
          NB  -- MAX of the blocksizes for DGEQRF, DORMQR, and DORGQR;
          The optimal LWORK is:
              2*N + MAX( 6*N, N*(NB+1) ).

          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]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 ALPHAR(j), ALPHAI(j), and BETA(j)
                should be correct for j=INFO+1,...,N.
          > N:  errors that usually indicate LAPACK problems:
                =N+1: error return from DGGBAL
                =N+2: error return from DGEQRF
                =N+3: error return from DORMQR
                =N+4: error return from DORGQR
                =N+5: error return from DGGHRD
                =N+6: error return from DHGEQZ (other than failed
                                                iteration)
                =N+7: error return from DTGEVC
                =N+8: error return from DGGBAK (computing VL)
                =N+9: error return from DGGBAK (computing VR)
                =N+10: error return from DLASCL (various calls)
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Date
November 2011
Further Details:
  Balancing
  ---------

  This driver calls DGGBAL to both permute and scale rows and columns
  of A and B.  The permutations PL and PR are chosen so that PL*A*PR
  and PL*B*R will be upper triangular except for the diagonal blocks
  A(i:j,i:j) and B(i:j,i:j), with i and j as close together as
  possible.  The diagonal scaling matrices DL and DR are chosen so
  that the pair  DL*PL*A*PR*DR, DL*PL*B*PR*DR have elements close to
  one (except for the elements that start out zero.)

  After the eigenvalues and eigenvectors of the balanced matrices
  have been computed, DGGBAK transforms the eigenvectors back to what
  they would have been (in perfect arithmetic) if they had not been
  balanced.

  Contents of A and B on Exit
  -------- -- - --- - -- ----

  If any eigenvectors are computed (either JOBVL='V' or JOBVR='V' or
  both), then on exit the arrays A and B will contain the real Schur
  form[*] of the "balanced" versions of A and B.  If no eigenvectors
  are computed, then only the diagonal blocks will be correct.

  [*] See DHGEQZ, DGEGS, or read the book "Matrix Computations",
      by Golub & van Loan, pub. by Johns Hopkins U. Press.

Definition at line 308 of file dgegv.f.

308 *
309 * -- LAPACK driver routine (version 3.4.0) --
310 * -- LAPACK is a software package provided by Univ. of Tennessee, --
311 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
312 * November 2011
313 *
314 * .. Scalar Arguments ..
315  CHARACTER jobvl, jobvr
316  INTEGER info, lda, ldb, ldvl, ldvr, lwork, n
317 * ..
318 * .. Array Arguments ..
319  DOUBLE PRECISION a( lda, * ), alphai( * ), alphar( * ),
320  $ b( ldb, * ), beta( * ), vl( ldvl, * ),
321  $ vr( ldvr, * ), work( * )
322 * ..
323 *
324 * =====================================================================
325 *
326 * .. Parameters ..
327  DOUBLE PRECISION zero, one
328  parameter ( zero = 0.0d0, one = 1.0d0 )
329 * ..
330 * .. Local Scalars ..
331  LOGICAL ilimit, ilv, ilvl, ilvr, lquery
332  CHARACTER chtemp
333  INTEGER icols, ihi, iinfo, ijobvl, ijobvr, ileft, ilo,
334  $ in, iright, irows, itau, iwork, jc, jr, lopt,
335  $ lwkmin, lwkopt, nb, nb1, nb2, nb3
336  DOUBLE PRECISION absai, absar, absb, anrm, anrm1, anrm2, bnrm,
337  $ bnrm1, bnrm2, eps, onepls, safmax, safmin,
338  $ salfai, salfar, sbeta, scale, temp
339 * ..
340 * .. Local Arrays ..
341  LOGICAL ldumma( 1 )
342 * ..
343 * .. External Subroutines ..
344  EXTERNAL dgeqrf, dggbak, dggbal, dgghrd, dhgeqz, dlacpy,
346 * ..
347 * .. External Functions ..
348  LOGICAL lsame
349  INTEGER ilaenv
350  DOUBLE PRECISION dlamch, dlange
351  EXTERNAL lsame, ilaenv, dlamch, dlange
352 * ..
353 * .. Intrinsic Functions ..
354  INTRINSIC abs, int, max
355 * ..
356 * .. Executable Statements ..
357 *
358 * Decode the input arguments
359 *
360  IF( lsame( jobvl, 'N' ) ) THEN
361  ijobvl = 1
362  ilvl = .false.
363  ELSE IF( lsame( jobvl, 'V' ) ) THEN
364  ijobvl = 2
365  ilvl = .true.
366  ELSE
367  ijobvl = -1
368  ilvl = .false.
369  END IF
370 *
371  IF( lsame( jobvr, 'N' ) ) THEN
372  ijobvr = 1
373  ilvr = .false.
374  ELSE IF( lsame( jobvr, 'V' ) ) THEN
375  ijobvr = 2
376  ilvr = .true.
377  ELSE
378  ijobvr = -1
379  ilvr = .false.
380  END IF
381  ilv = ilvl .OR. ilvr
382 *
383 * Test the input arguments
384 *
385  lwkmin = max( 8*n, 1 )
386  lwkopt = lwkmin
387  work( 1 ) = lwkopt
388  lquery = ( lwork.EQ.-1 )
389  info = 0
390  IF( ijobvl.LE.0 ) THEN
391  info = -1
392  ELSE IF( ijobvr.LE.0 ) THEN
393  info = -2
394  ELSE IF( n.LT.0 ) THEN
395  info = -3
396  ELSE IF( lda.LT.max( 1, n ) ) THEN
397  info = -5
398  ELSE IF( ldb.LT.max( 1, n ) ) THEN
399  info = -7
400  ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
401  info = -12
402  ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
403  info = -14
404  ELSE IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
405  info = -16
406  END IF
407 *
408  IF( info.EQ.0 ) THEN
409  nb1 = ilaenv( 1, 'DGEQRF', ' ', n, n, -1, -1 )
410  nb2 = ilaenv( 1, 'DORMQR', ' ', n, n, n, -1 )
411  nb3 = ilaenv( 1, 'DORGQR', ' ', n, n, n, -1 )
412  nb = max( nb1, nb2, nb3 )
413  lopt = 2*n + max( 6*n, n*( nb+1 ) )
414  work( 1 ) = lopt
415  END IF
416 *
417  IF( info.NE.0 ) THEN
418  CALL xerbla( 'DGEGV ', -info )
419  RETURN
420  ELSE IF( lquery ) THEN
421  RETURN
422  END IF
423 *
424 * Quick return if possible
425 *
426  IF( n.EQ.0 )
427  $ RETURN
428 *
429 * Get machine constants
430 *
431  eps = dlamch( 'E' )*dlamch( 'B' )
432  safmin = dlamch( 'S' )
433  safmin = safmin + safmin
434  safmax = one / safmin
435  onepls = one + ( 4*eps )
436 *
437 * Scale A
438 *
439  anrm = dlange( 'M', n, n, a, lda, work )
440  anrm1 = anrm
441  anrm2 = one
442  IF( anrm.LT.one ) THEN
443  IF( safmax*anrm.LT.one ) THEN
444  anrm1 = safmin
445  anrm2 = safmax*anrm
446  END IF
447  END IF
448 *
449  IF( anrm.GT.zero ) THEN
450  CALL dlascl( 'G', -1, -1, anrm, one, n, n, a, lda, iinfo )
451  IF( iinfo.NE.0 ) THEN
452  info = n + 10
453  RETURN
454  END IF
455  END IF
456 *
457 * Scale B
458 *
459  bnrm = dlange( 'M', n, n, b, ldb, work )
460  bnrm1 = bnrm
461  bnrm2 = one
462  IF( bnrm.LT.one ) THEN
463  IF( safmax*bnrm.LT.one ) THEN
464  bnrm1 = safmin
465  bnrm2 = safmax*bnrm
466  END IF
467  END IF
468 *
469  IF( bnrm.GT.zero ) THEN
470  CALL dlascl( 'G', -1, -1, bnrm, one, n, n, b, ldb, iinfo )
471  IF( iinfo.NE.0 ) THEN
472  info = n + 10
473  RETURN
474  END IF
475  END IF
476 *
477 * Permute the matrix to make it more nearly triangular
478 * Workspace layout: (8*N words -- "work" requires 6*N words)
479 * left_permutation, right_permutation, work...
480 *
481  ileft = 1
482  iright = n + 1
483  iwork = iright + n
484  CALL dggbal( 'P', n, a, lda, b, ldb, ilo, ihi, work( ileft ),
485  $ work( iright ), work( iwork ), iinfo )
486  IF( iinfo.NE.0 ) THEN
487  info = n + 1
488  GO TO 120
489  END IF
490 *
491 * Reduce B to triangular form, and initialize VL and/or VR
492 * Workspace layout: ("work..." must have at least N words)
493 * left_permutation, right_permutation, tau, work...
494 *
495  irows = ihi + 1 - ilo
496  IF( ilv ) THEN
497  icols = n + 1 - ilo
498  ELSE
499  icols = irows
500  END IF
501  itau = iwork
502  iwork = itau + irows
503  CALL dgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
504  $ work( iwork ), lwork+1-iwork, iinfo )
505  IF( iinfo.GE.0 )
506  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
507  IF( iinfo.NE.0 ) THEN
508  info = n + 2
509  GO TO 120
510  END IF
511 *
512  CALL dormqr( 'L', 'T', irows, icols, irows, b( ilo, ilo ), ldb,
513  $ work( itau ), a( ilo, ilo ), lda, work( iwork ),
514  $ lwork+1-iwork, iinfo )
515  IF( iinfo.GE.0 )
516  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
517  IF( iinfo.NE.0 ) THEN
518  info = n + 3
519  GO TO 120
520  END IF
521 *
522  IF( ilvl ) THEN
523  CALL dlaset( 'Full', n, n, zero, one, vl, ldvl )
524  CALL dlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
525  $ vl( ilo+1, ilo ), ldvl )
526  CALL dorgqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
527  $ work( itau ), work( iwork ), lwork+1-iwork,
528  $ iinfo )
529  IF( iinfo.GE.0 )
530  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
531  IF( iinfo.NE.0 ) THEN
532  info = n + 4
533  GO TO 120
534  END IF
535  END IF
536 *
537  IF( ilvr )
538  $ CALL dlaset( 'Full', n, n, zero, one, vr, ldvr )
539 *
540 * Reduce to generalized Hessenberg form
541 *
542  IF( ilv ) THEN
543 *
544 * Eigenvectors requested -- work on whole matrix.
545 *
546  CALL dgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
547  $ ldvl, vr, ldvr, iinfo )
548  ELSE
549  CALL dgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
550  $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, iinfo )
551  END IF
552  IF( iinfo.NE.0 ) THEN
553  info = n + 5
554  GO TO 120
555  END IF
556 *
557 * Perform QZ algorithm
558 * Workspace layout: ("work..." must have at least 1 word)
559 * left_permutation, right_permutation, work...
560 *
561  iwork = itau
562  IF( ilv ) THEN
563  chtemp = 'S'
564  ELSE
565  chtemp = 'E'
566  END IF
567  CALL dhgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
568  $ alphar, alphai, beta, vl, ldvl, vr, ldvr,
569  $ work( iwork ), lwork+1-iwork, iinfo )
570  IF( iinfo.GE.0 )
571  $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
572  IF( iinfo.NE.0 ) THEN
573  IF( iinfo.GT.0 .AND. iinfo.LE.n ) THEN
574  info = iinfo
575  ELSE IF( iinfo.GT.n .AND. iinfo.LE.2*n ) THEN
576  info = iinfo - n
577  ELSE
578  info = n + 6
579  END IF
580  GO TO 120
581  END IF
582 *
583  IF( ilv ) THEN
584 *
585 * Compute Eigenvectors (DTGEVC requires 6*N words of workspace)
586 *
587  IF( ilvl ) THEN
588  IF( ilvr ) THEN
589  chtemp = 'B'
590  ELSE
591  chtemp = 'L'
592  END IF
593  ELSE
594  chtemp = 'R'
595  END IF
596 *
597  CALL dtgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
598  $ vr, ldvr, n, in, work( iwork ), iinfo )
599  IF( iinfo.NE.0 ) THEN
600  info = n + 7
601  GO TO 120
602  END IF
603 *
604 * Undo balancing on VL and VR, rescale
605 *
606  IF( ilvl ) THEN
607  CALL dggbak( 'P', 'L', n, ilo, ihi, work( ileft ),
608  $ work( iright ), n, vl, ldvl, iinfo )
609  IF( iinfo.NE.0 ) THEN
610  info = n + 8
611  GO TO 120
612  END IF
613  DO 50 jc = 1, n
614  IF( alphai( jc ).LT.zero )
615  $ GO TO 50
616  temp = zero
617  IF( alphai( jc ).EQ.zero ) THEN
618  DO 10 jr = 1, n
619  temp = max( temp, abs( vl( jr, jc ) ) )
620  10 CONTINUE
621  ELSE
622  DO 20 jr = 1, n
623  temp = max( temp, abs( vl( jr, jc ) )+
624  $ abs( vl( jr, jc+1 ) ) )
625  20 CONTINUE
626  END IF
627  IF( temp.LT.safmin )
628  $ GO TO 50
629  temp = one / temp
630  IF( alphai( jc ).EQ.zero ) THEN
631  DO 30 jr = 1, n
632  vl( jr, jc ) = vl( jr, jc )*temp
633  30 CONTINUE
634  ELSE
635  DO 40 jr = 1, n
636  vl( jr, jc ) = vl( jr, jc )*temp
637  vl( jr, jc+1 ) = vl( jr, jc+1 )*temp
638  40 CONTINUE
639  END IF
640  50 CONTINUE
641  END IF
642  IF( ilvr ) THEN
643  CALL dggbak( 'P', 'R', n, ilo, ihi, work( ileft ),
644  $ work( iright ), n, vr, ldvr, iinfo )
645  IF( iinfo.NE.0 ) THEN
646  info = n + 9
647  GO TO 120
648  END IF
649  DO 100 jc = 1, n
650  IF( alphai( jc ).LT.zero )
651  $ GO TO 100
652  temp = zero
653  IF( alphai( jc ).EQ.zero ) THEN
654  DO 60 jr = 1, n
655  temp = max( temp, abs( vr( jr, jc ) ) )
656  60 CONTINUE
657  ELSE
658  DO 70 jr = 1, n
659  temp = max( temp, abs( vr( jr, jc ) )+
660  $ abs( vr( jr, jc+1 ) ) )
661  70 CONTINUE
662  END IF
663  IF( temp.LT.safmin )
664  $ GO TO 100
665  temp = one / temp
666  IF( alphai( jc ).EQ.zero ) THEN
667  DO 80 jr = 1, n
668  vr( jr, jc ) = vr( jr, jc )*temp
669  80 CONTINUE
670  ELSE
671  DO 90 jr = 1, n
672  vr( jr, jc ) = vr( jr, jc )*temp
673  vr( jr, jc+1 ) = vr( jr, jc+1 )*temp
674  90 CONTINUE
675  END IF
676  100 CONTINUE
677  END IF
678 *
679 * End of eigenvector calculation
680 *
681  END IF
682 *
683 * Undo scaling in alpha, beta
684 *
685 * Note: this does not give the alpha and beta for the unscaled
686 * problem.
687 *
688 * Un-scaling is limited to avoid underflow in alpha and beta
689 * if they are significant.
690 *
691  DO 110 jc = 1, n
692  absar = abs( alphar( jc ) )
693  absai = abs( alphai( jc ) )
694  absb = abs( beta( jc ) )
695  salfar = anrm*alphar( jc )
696  salfai = anrm*alphai( jc )
697  sbeta = bnrm*beta( jc )
698  ilimit = .false.
699  scale = one
700 *
701 * Check for significant underflow in ALPHAI
702 *
703  IF( abs( salfai ).LT.safmin .AND. absai.GE.
704  $ max( safmin, eps*absar, eps*absb ) ) THEN
705  ilimit = .true.
706  scale = ( onepls*safmin / anrm1 ) /
707  $ max( onepls*safmin, anrm2*absai )
708 *
709  ELSE IF( salfai.EQ.zero ) THEN
710 *
711 * If insignificant underflow in ALPHAI, then make the
712 * conjugate eigenvalue real.
713 *
714  IF( alphai( jc ).LT.zero .AND. jc.GT.1 ) THEN
715  alphai( jc-1 ) = zero
716  ELSE IF( alphai( jc ).GT.zero .AND. jc.LT.n ) THEN
717  alphai( jc+1 ) = zero
718  END IF
719  END IF
720 *
721 * Check for significant underflow in ALPHAR
722 *
723  IF( abs( salfar ).LT.safmin .AND. absar.GE.
724  $ max( safmin, eps*absai, eps*absb ) ) THEN
725  ilimit = .true.
726  scale = max( scale, ( onepls*safmin / anrm1 ) /
727  $ max( onepls*safmin, anrm2*absar ) )
728  END IF
729 *
730 * Check for significant underflow in BETA
731 *
732  IF( abs( sbeta ).LT.safmin .AND. absb.GE.
733  $ max( safmin, eps*absar, eps*absai ) ) THEN
734  ilimit = .true.
735  scale = max( scale, ( onepls*safmin / bnrm1 ) /
736  $ max( onepls*safmin, bnrm2*absb ) )
737  END IF
738 *
739 * Check for possible overflow when limiting scaling
740 *
741  IF( ilimit ) THEN
742  temp = ( scale*safmin )*max( abs( salfar ), abs( salfai ),
743  $ abs( sbeta ) )
744  IF( temp.GT.one )
745  $ scale = scale / temp
746  IF( scale.LT.one )
747  $ ilimit = .false.
748  END IF
749 *
750 * Recompute un-scaled ALPHAR, ALPHAI, BETA if necessary.
751 *
752  IF( ilimit ) THEN
753  salfar = ( scale*alphar( jc ) )*anrm
754  salfai = ( scale*alphai( jc ) )*anrm
755  sbeta = ( scale*beta( jc ) )*bnrm
756  END IF
757  alphar( jc ) = salfar
758  alphai( jc ) = salfai
759  beta( jc ) = sbeta
760  110 CONTINUE
761 *
762  120 CONTINUE
763  work( 1 ) = lwkopt
764 *
765  RETURN
766 *
767 * End of DGEGV
768 *
subroutine dgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
DGGHRD
Definition: dgghrd.f:209
subroutine dlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
DLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: dlaset.f:112
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65
subroutine dhgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
DHGEQZ
Definition: dhgeqz.f:306
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:105
subroutine dlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
DLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: dlascl.f:145
subroutine dormqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
DORMQR
Definition: dormqr.f:169
subroutine dtgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, INFO)
DTGEVC
Definition: dtgevc.f:297
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
double precision function dlange(NORM, M, N, A, LDA, WORK)
DLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: dlange.f:116
subroutine dggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
DGGBAL
Definition: dggbal.f:179
subroutine dgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
DGEQRF
Definition: dgeqrf.f:138
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine dggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
DGGBAK
Definition: dggbak.f:149
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine dorgqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
DORGQR
Definition: dorgqr.f:130

Here is the call graph for this function: