303 SUBROUTINE dgeevx( BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, WR, WI,
304 $ VL, LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM,
305 $ RCONDE, RCONDV, WORK, LWORK, IWORK, INFO )
313 CHARACTER BALANC, JOBVL, JOBVR, SENSE
314 INTEGER IHI, ILO, INFO, LDA, LDVL, LDVR, LWORK, N
315 DOUBLE PRECISION ABNRM
319 DOUBLE PRECISION A( LDA, * ), RCONDE( * ), RCONDV( * ),
320 $ scale( * ), vl( ldvl, * ), vr( ldvr, * ),
321 $ wi( * ), work( * ), wr( * )
327 DOUBLE PRECISION ZERO, ONE
328 PARAMETER ( ZERO = 0.0d0, one = 1.0d0 )
331 LOGICAL LQUERY, SCALEA, WANTVL, WANTVR, WNTSNB, WNTSNE,
334 INTEGER HSWORK, I, ICOND, IERR, ITAU, IWRK, K,
335 $ lwork_trevc, maxwrk, minwrk, nout
336 DOUBLE PRECISION ANRM, BIGNUM, CS, CSCALE, EPS, R, SCL, SMLNUM,
341 DOUBLE PRECISION DUM( 1 )
350 INTEGER IDAMAX, ILAENV
351 DOUBLE PRECISION DLAMCH, DLANGE, DLAPY2, DNRM2
352 EXTERNAL lsame, idamax, ilaenv, dlamch, dlange, dlapy2,
363 lquery = ( lwork.EQ.-1 )
364 wantvl = lsame( jobvl,
'V' )
365 wantvr = lsame( jobvr,
'V' )
366 wntsnn = lsame( sense,
'N' )
367 wntsne = lsame( sense,
'E' )
368 wntsnv = lsame( sense,
'V' )
369 wntsnb = lsame( sense,
'B' )
370 IF( .NOT.( lsame( balanc,
'N' ) .OR. lsame( balanc,
'S' )
371 $ .OR. lsame( balanc,
'P' ) .OR. lsame( balanc,
'B' ) ) )
374 ELSE IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl,
'N' ) ) )
THEN
376 ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr,
'N' ) ) )
THEN
378 ELSE IF( .NOT.( wntsnn .OR. wntsne .OR. wntsnb .OR. wntsnv ) .OR.
379 $ ( ( wntsne .OR. wntsnb ) .AND. .NOT.( wantvl .AND.
382 ELSE IF( n.LT.0 )
THEN
384 ELSE IF( lda.LT.max( 1, n ) )
THEN
386 ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) )
THEN
388 ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) )
THEN
407 maxwrk = n + n*ilaenv( 1,
'DGEHRD',
' ', n, 1, n, 0 )
410 CALL dtrevc3(
'L',
'B',
SELECT, n, a, lda,
411 $ vl, ldvl, vr, ldvr,
412 $ n, nout, work, -1, ierr )
413 lwork_trevc = int( work(1) )
414 maxwrk = max( maxwrk, n + lwork_trevc )
415 CALL dhseqr(
'S',
'V', n, 1, n, a, lda, wr, wi, vl, ldvl,
417 ELSE IF( wantvr )
THEN
418 CALL dtrevc3(
'R',
'B',
SELECT, n, a, lda,
419 $ vl, ldvl, vr, ldvr,
420 $ n, nout, work, -1, ierr )
421 lwork_trevc = int( work(1) )
422 maxwrk = max( maxwrk, n + lwork_trevc )
423 CALL dhseqr(
'S',
'V', n, 1, n, a, lda, wr, wi, vr, ldvr,
427 CALL dhseqr(
'E',
'N', n, 1, n, a, lda, wr, wi, vr,
428 $ ldvr, work, -1, info )
430 CALL dhseqr(
'S',
'N', n, 1, n, a, lda, wr, wi, vr,
431 $ ldvr, work, -1, info )
434 hswork = int( work(1) )
436 IF( ( .NOT.wantvl ) .AND. ( .NOT.wantvr ) )
THEN
439 $ minwrk = max( minwrk, n*n+6*n )
440 maxwrk = max( maxwrk, hswork )
442 $ maxwrk = max( maxwrk, n*n + 6*n )
445 IF( ( .NOT.wntsnn ) .AND. ( .NOT.wntsne ) )
446 $ minwrk = max( minwrk, n*n + 6*n )
447 maxwrk = max( maxwrk, hswork )
448 maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1,
'DORGHR',
449 $
' ', n, 1, n, -1 ) )
450 IF( ( .NOT.wntsnn ) .AND. ( .NOT.wntsne ) )
451 $ maxwrk = max( maxwrk, n*n + 6*n )
452 maxwrk = max( maxwrk, 3*n )
454 maxwrk = max( maxwrk, minwrk )
458 IF( lwork.LT.minwrk .AND. .NOT.lquery )
THEN
464 CALL xerbla(
'DGEEVX', -info )
466 ELSE IF( lquery )
THEN
478 smlnum = dlamch(
'S' )
479 bignum = one / smlnum
480 smlnum = sqrt( smlnum ) / eps
481 bignum = one / smlnum
486 anrm = dlange(
'M', n, n, a, lda, dum )
488 IF( anrm.GT.zero .AND. anrm.LT.smlnum )
THEN
491 ELSE IF( anrm.GT.bignum )
THEN
496 $
CALL dlascl(
'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
500 CALL dgebal( balanc, n, a, lda, ilo, ihi, scale, ierr )
501 abnrm = dlange(
'1', n, n, a, lda, dum )
504 CALL dlascl(
'G', 0, 0, cscale, anrm, 1, 1, dum, 1, ierr )
513 CALL dgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
514 $ lwork-iwrk+1, ierr )
522 CALL dlacpy(
'L', n, n, a, lda, vl, ldvl )
527 CALL dorghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
528 $ lwork-iwrk+1, ierr )
534 CALL dhseqr(
'S',
'V', n, ilo, ihi, a, lda, wr, wi, vl, ldvl,
535 $ work( iwrk ), lwork-iwrk+1, info )
543 CALL dlacpy(
'F', n, n, vl, ldvl, vr, ldvr )
546 ELSE IF( wantvr )
THEN
552 CALL dlacpy(
'L', n, n, a, lda, vr, ldvr )
557 CALL dorghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
558 $ lwork-iwrk+1, ierr )
564 CALL dhseqr(
'S',
'V', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,
565 $ work( iwrk ), lwork-iwrk+1, info )
581 CALL dhseqr( job,
'N', n, ilo, ihi, a, lda, wr, wi, vr, ldvr,
582 $ work( iwrk ), lwork-iwrk+1, info )
590 IF( wantvl .OR. wantvr )
THEN
595 CALL dtrevc3( side,
'B',
SELECT, n, a, lda, vl, ldvl, vr, ldvr,
596 $ n, nout, work( iwrk ), lwork-iwrk+1, ierr )
602 IF( .NOT.wntsnn )
THEN
603 CALL dtrsna( sense,
'A',
SELECT, n, a, lda, vl, ldvl, vr, ldvr,
604 $ rconde, rcondv, n, nout, work( iwrk ), n, iwork,
612 CALL dgebak( balanc,
'L', n, ilo, ihi, scale, n, vl, ldvl,
618 IF( wi( i ).EQ.zero )
THEN
619 scl = one / dnrm2( n, vl( 1, i ), 1 )
620 CALL dscal( n, scl, vl( 1, i ), 1 )
621 ELSE IF( wi( i ).GT.zero )
THEN
622 scl = one / dlapy2( dnrm2( n, vl( 1, i ), 1 ),
623 $ dnrm2( n, vl( 1, i+1 ), 1 ) )
624 CALL dscal( n, scl, vl( 1, i ), 1 )
625 CALL dscal( n, scl, vl( 1, i+1 ), 1 )
627 work( k ) = vl( k, i )**2 + vl( k, i+1 )**2
629 k = idamax( n, work, 1 )
630 CALL dlartg( vl( k, i ), vl( k, i+1 ), cs, sn, r )
631 CALL drot( n, vl( 1, i ), 1, vl( 1, i+1 ), 1, cs, sn )
641 CALL dgebak( balanc,
'R', n, ilo, ihi, scale, n, vr, ldvr,
647 IF( wi( i ).EQ.zero )
THEN
648 scl = one / dnrm2( n, vr( 1, i ), 1 )
649 CALL dscal( n, scl, vr( 1, i ), 1 )
650 ELSE IF( wi( i ).GT.zero )
THEN
651 scl = one / dlapy2( dnrm2( n, vr( 1, i ), 1 ),
652 $ dnrm2( n, vr( 1, i+1 ), 1 ) )
653 CALL dscal( n, scl, vr( 1, i ), 1 )
654 CALL dscal( n, scl, vr( 1, i+1 ), 1 )
656 work( k ) = vr( k, i )**2 + vr( k, i+1 )**2
658 k = idamax( n, work, 1 )
659 CALL dlartg( vr( k, i ), vr( k, i+1 ), cs, sn, r )
660 CALL drot( n, vr( 1, i ), 1, vr( 1, i+1 ), 1, cs, sn )
670 CALL dlascl(
'G', 0, 0, cscale, anrm, n-info, 1, wr( info+1 ),
671 $ max( n-info, 1 ), ierr )
672 CALL dlascl(
'G', 0, 0, cscale, anrm, n-info, 1, wi( info+1 ),
673 $ max( n-info, 1 ), ierr )
675 IF( ( wntsnv .OR. wntsnb ) .AND. icond.EQ.0 )
676 $
CALL dlascl(
'G', 0, 0, cscale, anrm, n, 1, rcondv, n,
679 CALL dlascl(
'G', 0, 0, cscale, anrm, ilo-1, 1, wr, n,
681 CALL dlascl(
'G', 0, 0, cscale, anrm, ilo-1, 1, wi, n,
subroutine xerbla(srname, info)
subroutine dgebak(job, side, n, ilo, ihi, scale, m, v, ldv, info)
DGEBAK
subroutine dgebal(job, n, a, lda, ilo, ihi, scale, info)
DGEBAL
subroutine dgeevx(balanc, jobvl, jobvr, sense, n, a, lda, wr, wi, vl, ldvl, vr, ldvr, ilo, ihi, scale, abnrm, rconde, rcondv, work, lwork, iwork, info)
DGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
subroutine dgehrd(n, ilo, ihi, a, lda, tau, work, lwork, info)
DGEHRD
subroutine dhseqr(job, compz, n, ilo, ihi, h, ldh, wr, wi, z, ldz, work, lwork, info)
DHSEQR
subroutine dlacpy(uplo, m, n, a, lda, b, ldb)
DLACPY copies all or part of one two-dimensional array to another.
subroutine dlartg(f, g, c, s, r)
DLARTG generates a plane rotation with real cosine and real sine.
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 drot(n, dx, incx, dy, incy, c, s)
DROT
subroutine dscal(n, da, dx, incx)
DSCAL
subroutine dtrevc3(side, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, mm, m, work, lwork, info)
DTREVC3
subroutine dtrsna(job, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, s, sep, mm, m, work, ldwork, iwork, info)
DTRSNA
subroutine dorghr(n, ilo, ihi, a, lda, tau, work, lwork, info)
DORGHR