387 SUBROUTINE sggevx( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, B, LDB,
388 $ ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR, ILO,
389 $ IHI, LSCALE, RSCALE, ABNRM, BBNRM, RCONDE,
390 $ RCONDV, WORK, LWORK, IWORK, BWORK, INFO )
397 CHARACTER BALANC, JOBVL, JOBVR, SENSE
398 INTEGER IHI, ILO, INFO, LDA, LDB, LDVL, LDVR, LWORK, N
404 REAL A( LDA, * ), ALPHAI( * ), ALPHAR( * ),
405 $ b( ldb, * ), beta( * ), lscale( * ),
406 $ rconde( * ), rcondv( * ), rscale( * ),
407 $ vl( ldvl, * ), vr( ldvr, * ), work( * )
414 PARAMETER ( ZERO = 0.0e+0, one = 1.0e+0 )
417 LOGICAL ILASCL, ILBSCL, ILV, ILVL, ILVR, LQUERY, NOSCL,
418 $ PAIR, WANTSB, WANTSE, WANTSN, WANTSV
420 INTEGER I, ICOLS, IERR, IJOBVL, IJOBVR, IN, IROWS,
421 $ ITAU, IWRK, IWRK1, J, JC, JR, M, MAXWRK,
423 REAL ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS,
438 EXTERNAL lsame, ilaenv, slamch, slange
441 INTRINSIC abs, max, sqrt
447 IF( lsame( jobvl,
'N' ) )
THEN
450 ELSE IF( lsame( jobvl,
'V' ) )
THEN
458 IF( lsame( jobvr,
'N' ) )
THEN
461 ELSE IF( lsame( jobvr,
'V' ) )
THEN
470 noscl = lsame( balanc,
'N' ) .OR. lsame( balanc,
'P' )
471 wantsn = lsame( sense,
'N' )
472 wantse = lsame( sense,
'E' )
473 wantsv = lsame( sense,
'V' )
474 wantsb = lsame( sense,
'B' )
479 lquery = ( lwork.EQ.-1 )
480 IF( .NOT.( noscl .OR. lsame( balanc,
'S' ) .OR.
481 $ lsame( balanc,
'B' ) ) )
THEN
483 ELSE IF( ijobvl.LE.0 )
THEN
485 ELSE IF( ijobvr.LE.0 )
THEN
487 ELSE IF( .NOT.( wantsn .OR. wantse .OR. wantsb .OR. wantsv ) )
490 ELSE IF( n.LT.0 )
THEN
492 ELSE IF( lda.LT.max( 1, n ) )
THEN
494 ELSE IF( ldb.LT.max( 1, n ) )
THEN
496 ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) )
THEN
498 ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) )
THEN
515 IF( noscl .AND. .NOT.ilv )
THEN
522 ELSE IF( wantsv .OR. wantsb )
THEN
523 minwrk = 2*n*( n + 4 ) + 16
526 maxwrk = max( maxwrk,
527 $ n + n*ilaenv( 1,
'SGEQRF',
' ', n, 1, n, 0 ) )
528 maxwrk = max( maxwrk,
529 $ n + n*ilaenv( 1,
'SORMQR',
' ', n, 1, n, 0 ) )
531 maxwrk = max( maxwrk, n +
532 $ n*ilaenv( 1,
'SORGQR',
' ', n, 1, n, 0 ) )
537 IF( lwork.LT.minwrk .AND. .NOT.lquery )
THEN
543 CALL xerbla(
'SGGEVX', -info )
545 ELSE IF( lquery )
THEN
558 smlnum = slamch(
'S' )
559 bignum = one / smlnum
560 CALL slabad( smlnum, bignum )
561 smlnum = sqrt( smlnum ) / eps
562 bignum = one / smlnum
566 anrm = slange(
'M', n, n, a, lda, work )
568 IF( anrm.GT.zero .AND. anrm.LT.smlnum )
THEN
571 ELSE IF( anrm.GT.bignum )
THEN
576 $
CALL slascl(
'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
580 bnrm = slange(
'M', n, n, b, ldb, work )
582 IF( bnrm.GT.zero .AND. bnrm.LT.smlnum )
THEN
585 ELSE IF( bnrm.GT.bignum )
THEN
590 $
CALL slascl(
'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
595 CALL sggbal( balanc, n, a, lda, b, ldb, ilo, ihi, lscale, rscale,
600 abnrm = slange(
'1', n, n, a, lda, work( 1 ) )
603 CALL slascl(
'G', 0, 0, anrmto, anrm, 1, 1, work( 1 ), 1,
608 bbnrm = slange(
'1', n, n, b, ldb, work( 1 ) )
611 CALL slascl(
'G', 0, 0, bnrmto, bnrm, 1, 1, work( 1 ), 1,
619 irows = ihi + 1 - ilo
620 IF( ilv .OR. .NOT.wantsn )
THEN
627 CALL sgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
628 $ work( iwrk ), lwork+1-iwrk, ierr )
633 CALL sormqr(
'L',
'T', irows, icols, irows, b( ilo, ilo ), ldb,
634 $ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
635 $ lwork+1-iwrk, ierr )
641 CALL slaset(
'Full', n, n, zero, one, vl, ldvl )
642 IF( irows.GT.1 )
THEN
643 CALL slacpy(
'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
644 $ vl( ilo+1, ilo ), ldvl )
646 CALL sorgqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
647 $ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
651 $
CALL slaset(
'Full', n, n, zero, one, vr, ldvr )
656 IF( ilv .OR. .NOT.wantsn )
THEN
660 CALL sgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
661 $ ldvl, vr, ldvr, ierr )
663 CALL sgghrd(
'N',
'N', irows, 1, irows, a( ilo, ilo ), lda,
664 $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, ierr )
671 IF( ilv .OR. .NOT.wantsn )
THEN
677 CALL shgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
678 $ alphar, alphai, beta, vl, ldvl, vr, ldvr, work,
681 IF( ierr.GT.0 .AND. ierr.LE.n )
THEN
683 ELSE IF( ierr.GT.n .AND. ierr.LE.2*n )
THEN
696 IF( ilv .OR. .NOT.wantsn )
THEN
708 CALL stgevc( chtemp,
'B', ldumma, n, a, lda, b, ldb, vl,
709 $ ldvl, vr, ldvr, n, in, work, ierr )
716 IF( .NOT.wantsn )
THEN
735 IF( a( i+1, i ).NE.zero )
THEN
746 ELSE IF( mm.EQ.2 )
THEN
748 bwork( i+1 ) = .true.
757 IF( wantse .OR. wantsb )
THEN
758 CALL stgevc(
'B',
'S', bwork, n, a, lda, b, ldb,
759 $ work( 1 ), n, work( iwrk ), n, mm, m,
760 $ work( iwrk1 ), ierr )
767 CALL stgsna( sense,
'S', bwork, n, a, lda, b, ldb,
768 $ work( 1 ), n, work( iwrk ), n, rconde( i ),
769 $ rcondv( i ), mm, m, work( iwrk1 ),
770 $ lwork-iwrk1+1, iwork, ierr )
780 CALL sggbak( balanc,
'L', n, ilo, ihi, lscale, rscale, n, vl,
784 IF( alphai( jc ).LT.zero )
787 IF( alphai( jc ).EQ.zero )
THEN
789 temp = max( temp, abs( vl( jr, jc ) ) )
793 temp = max( temp, abs( vl( jr, jc ) )+
794 $ abs( vl( jr, jc+1 ) ) )
800 IF( alphai( jc ).EQ.zero )
THEN
802 vl( jr, jc ) = vl( jr, jc )*temp
806 vl( jr, jc ) = vl( jr, jc )*temp
807 vl( jr, jc+1 ) = vl( jr, jc+1 )*temp
813 CALL sggbak( balanc,
'R', n, ilo, ihi, lscale, rscale, n, vr,
816 IF( alphai( jc ).LT.zero )
819 IF( alphai( jc ).EQ.zero )
THEN
821 temp = max( temp, abs( vr( jr, jc ) ) )
825 temp = max( temp, abs( vr( jr, jc ) )+
826 $ abs( vr( jr, jc+1 ) ) )
832 IF( alphai( jc ).EQ.zero )
THEN
834 vr( jr, jc ) = vr( jr, jc )*temp
838 vr( jr, jc ) = vr( jr, jc )*temp
839 vr( jr, jc+1 ) = vr( jr, jc+1 )*temp
850 CALL slascl(
'G', 0, 0, anrmto, anrm, n, 1, alphar, n, ierr )
851 CALL slascl(
'G', 0, 0, anrmto, anrm, n, 1, alphai, n, ierr )
855 CALL slascl(
'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
subroutine slabad(SMALL, LARGE)
SLABAD
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
subroutine slaset(UPLO, M, N, ALPHA, BETA, A, LDA)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
subroutine xerbla(SRNAME, INFO)
XERBLA
subroutine sggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
SGGBAK
subroutine sggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
SGGBAL
subroutine stgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, INFO)
STGEVC
subroutine sgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
SGEQRF
subroutine shgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHAR, ALPHAI, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
SHGEQZ
subroutine sggevx(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)
SGGEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
subroutine stgsna(JOB, HOWMNY, SELECT, N, A, LDA, B, LDB, VL, LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK, INFO)
STGSNA
subroutine sorgqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
SORGQR
subroutine sormqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
SORMQR
subroutine sgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
SGGHRD