389 SUBROUTINE dggevx( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, B, LDB,
390 $ alphar, alphai, beta, vl, ldvl, vr, ldvr, ilo,
391 $ ihi, lscale, rscale, abnrm, bbnrm, rconde,
392 $ rcondv, work, lwork, iwork, bwork, info )
400 CHARACTER BALANC, JOBVL, JOBVR, SENSE
401 INTEGER IHI, ILO, INFO, LDA, LDB, LDVL, LDVR, LWORK, N
402 DOUBLE PRECISION ABNRM, BBNRM
407 DOUBLE PRECISION A( lda, * ), ALPHAI( * ), ALPHAR( * ),
408 $ b( ldb, * ), beta( * ), lscale( * ),
409 $ rconde( * ), rcondv( * ), rscale( * ),
410 $ vl( ldvl, * ), vr( ldvr, * ), work( * )
416 DOUBLE PRECISION ZERO, ONE
417 parameter( zero = 0.0d+0, one = 1.0d+0 )
420 LOGICAL ILASCL, ILBSCL, ILV, ILVL, ILVR, LQUERY, NOSCL,
421 $ pair, wantsb, wantse, wantsn, wantsv
423 INTEGER I, ICOLS, IERR, IJOBVL, IJOBVR, IN, IROWS,
424 $ itau, iwrk, iwrk1, j, jc, jr, m, maxwrk,
426 DOUBLE PRECISION ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS,
440 DOUBLE PRECISION DLAMCH, DLANGE
441 EXTERNAL lsame, ilaenv, dlamch, dlange
444 INTRINSIC abs, max, sqrt
450 IF( lsame( jobvl,
'N' ) )
THEN
453 ELSE IF( lsame( jobvl,
'V' ) )
THEN
461 IF( lsame( jobvr,
'N' ) )
THEN
464 ELSE IF( lsame( jobvr,
'V' ) )
THEN
473 noscl = lsame( balanc,
'N' ) .OR. lsame( balanc,
'P' )
474 wantsn = lsame( sense,
'N' )
475 wantse = lsame( sense,
'E' )
476 wantsv = lsame( sense,
'V' )
477 wantsb = lsame( sense,
'B' )
482 lquery = ( lwork.EQ.-1 )
483 IF( .NOT.( lsame( balanc,
'N' ) .OR. lsame( balanc,
484 $
'S' ) .OR. lsame( balanc,
'P' ) .OR. lsame( balanc,
'B' ) ) )
487 ELSE IF( ijobvl.LE.0 )
THEN
489 ELSE IF( ijobvr.LE.0 )
THEN
491 ELSE IF( .NOT.( wantsn .OR. wantse .OR. wantsb .OR. wantsv ) )
494 ELSE IF( n.LT.0 )
THEN
496 ELSE IF( lda.LT.max( 1, n ) )
THEN
498 ELSE IF( ldb.LT.max( 1, n ) )
THEN
500 ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) )
THEN
502 ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) )
THEN
519 IF( noscl .AND. .NOT.ilv )
THEN
524 IF( wantse .OR. wantsb )
THEN
527 IF( wantsv .OR. wantsb )
THEN
528 minwrk = max( minwrk, 2*n*( n + 4 ) + 16 )
531 maxwrk = max( maxwrk,
532 $ n + n*ilaenv( 1,
'DGEQRF',
' ', n, 1, n, 0 ) )
533 maxwrk = max( maxwrk,
534 $ n + n*ilaenv( 1,
'DORMQR',
' ', n, 1, n, 0 ) )
536 maxwrk = max( maxwrk, n +
537 $ n*ilaenv( 1,
'DORGQR',
' ', n, 1, n, 0 ) )
542 IF( lwork.LT.minwrk .AND. .NOT.lquery )
THEN
548 CALL
xerbla(
'DGGEVX', -info )
550 ELSE IF( lquery )
THEN
563 smlnum = dlamch(
'S' )
564 bignum = one / smlnum
565 CALL
dlabad( smlnum, bignum )
566 smlnum = sqrt( smlnum ) / eps
567 bignum = one / smlnum
571 anrm = dlange(
'M', n, n, a, lda, work )
573 IF( anrm.GT.zero .AND. anrm.LT.smlnum )
THEN
576 ELSE IF( anrm.GT.bignum )
THEN
581 $ CALL
dlascl(
'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
585 bnrm = dlange(
'M', n, n, b, ldb, work )
587 IF( bnrm.GT.zero .AND. bnrm.LT.smlnum )
THEN
590 ELSE IF( bnrm.GT.bignum )
THEN
595 $ CALL
dlascl(
'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
600 CALL
dggbal( balanc, n, a, lda, b, ldb, ilo, ihi, lscale, rscale,
605 abnrm = dlange(
'1', n, n, a, lda, work( 1 ) )
608 CALL
dlascl(
'G', 0, 0, anrmto, anrm, 1, 1, work( 1 ), 1,
613 bbnrm = dlange(
'1', n, n, b, ldb, work( 1 ) )
616 CALL
dlascl(
'G', 0, 0, bnrmto, bnrm, 1, 1, work( 1 ), 1,
624 irows = ihi + 1 - ilo
625 IF( ilv .OR. .NOT.wantsn )
THEN
632 CALL
dgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
633 $ work( iwrk ), lwork+1-iwrk, ierr )
638 CALL
dormqr(
'L',
'T', irows, icols, irows, b( ilo, ilo ), ldb,
639 $ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
640 $ lwork+1-iwrk, ierr )
646 CALL
dlaset(
'Full', n, n, zero, one, vl, ldvl )
647 IF( irows.GT.1 )
THEN
648 CALL
dlacpy(
'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
649 $ vl( ilo+1, ilo ), ldvl )
651 CALL
dorgqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
652 $ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
656 $ CALL
dlaset(
'Full', n, n, zero, one, vr, ldvr )
661 IF( ilv .OR. .NOT.wantsn )
THEN
665 CALL
dgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
666 $ ldvl, vr, ldvr, ierr )
668 CALL
dgghrd(
'N',
'N', irows, 1, irows, a( ilo, ilo ), lda,
669 $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, ierr )
676 IF( ilv .OR. .NOT.wantsn )
THEN
682 CALL
dhgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
683 $ alphar, alphai, beta, vl, ldvl, vr, ldvr, work,
686 IF( ierr.GT.0 .AND. ierr.LE.n )
THEN
688 ELSE IF( ierr.GT.n .AND. ierr.LE.2*n )
THEN
701 IF( ilv .OR. .NOT.wantsn )
THEN
713 CALL
dtgevc( chtemp,
'B', ldumma, n, a, lda, b, ldb, vl,
714 $ ldvl, vr, ldvr, n, in, work, ierr )
721 IF( .NOT.wantsn )
THEN
740 IF( a( i+1, i ).NE.zero )
THEN
751 ELSE IF( mm.EQ.2 )
THEN
753 bwork( i+1 ) = .true.
762 IF( wantse .OR. wantsb )
THEN
763 CALL
dtgevc(
'B',
'S', bwork, n, a, lda, b, ldb,
764 $ work( 1 ), n, work( iwrk ), n, mm, m,
765 $ work( iwrk1 ), ierr )
772 CALL
dtgsna( sense,
'S', bwork, n, a, lda, b, ldb,
773 $ work( 1 ), n, work( iwrk ), n, rconde( i ),
774 $ rcondv( i ), mm, m, work( iwrk1 ),
775 $ lwork-iwrk1+1, iwork, ierr )
785 CALL
dggbak( balanc,
'L', n, ilo, ihi, lscale, rscale, n, vl,
789 IF( alphai( jc ).LT.zero )
792 IF( alphai( jc ).EQ.zero )
THEN
794 temp = max( temp, abs( vl( jr, jc ) ) )
798 temp = max( temp, abs( vl( jr, jc ) )+
799 $ abs( vl( jr, jc+1 ) ) )
805 IF( alphai( jc ).EQ.zero )
THEN
807 vl( jr, jc ) = vl( jr, jc )*temp
811 vl( jr, jc ) = vl( jr, jc )*temp
812 vl( jr, jc+1 ) = vl( jr, jc+1 )*temp
818 CALL
dggbak( balanc,
'R', n, ilo, ihi, lscale, rscale, n, vr,
821 IF( alphai( jc ).LT.zero )
824 IF( alphai( jc ).EQ.zero )
THEN
826 temp = max( temp, abs( vr( jr, jc ) ) )
830 temp = max( temp, abs( vr( jr, jc ) )+
831 $ abs( vr( jr, jc+1 ) ) )
837 IF( alphai( jc ).EQ.zero )
THEN
839 vr( jr, jc ) = vr( jr, jc )*temp
843 vr( jr, jc ) = vr( jr, jc )*temp
844 vr( jr, jc+1 ) = vr( jr, jc+1 )*temp
855 CALL
dlascl(
'G', 0, 0, anrmto, anrm, n, 1, alphar, n, ierr )
856 CALL
dlascl(
'G', 0, 0, anrmto, anrm, n, 1, alphai, n, ierr )
860 CALL
dlascl(
'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
subroutine dtgsna(JOB, HOWMNY, SELECT, N, A, LDA, B, LDB, VL, LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK, INFO)
DTGSNA
subroutine dgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
DGGHRD
subroutine dgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
DGEQRF
subroutine dormqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
DORMQR
subroutine xerbla(SRNAME, INFO)
XERBLA
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.
subroutine dlabad(SMALL, LARGE)
DLABAD
subroutine dorgqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
DORGQR
subroutine dhgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
DHGEQZ
subroutine dggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
DGGBAK
subroutine dtgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, INFO)
DTGEVC
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
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...
subroutine dggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
DGGBAL
subroutine dggevx(BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, B, LDB, ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR, ILO, IHI, LSCALE, RSCALE, ABNRM, BBNRM, RCONDE, RCONDV, WORK, LWORK, IWORK, BWORK, INFO)
DGGEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices ...