297 SUBROUTINE sdrgvx( NSIZE, THRESH, NIN, NOUT, A, LDA, B, AI, BI,
298 $ ALPHAR, ALPHAI, BETA, VL, VR, ILO, IHI, LSCALE,
299 $ RSCALE, S, STRU, DIF, DIFTRU, WORK, LWORK,
300 $ IWORK, LIWORK, RESULT, BWORK, INFO )
307 INTEGER IHI, ILO, INFO, LDA, LIWORK, LWORK, NIN, NOUT,
314 REAL A( LDA, * ), AI( LDA, * ), ALPHAI( * ),
315 $ alphar( * ), b( lda, * ), beta( * ),
316 $ bi( lda, * ), dif( * ), diftru( * ),
317 $ lscale( * ), result( 4 ), rscale( * ), s( * ),
318 $ stru( * ), vl( lda, * ), vr( lda, * ),
325 REAL , ONE, TEN, TNTH, HALF
326 PARAMETER ( ZERO = 0.0e+0, one = 1.0e+0, ten = 1.0e+1,
327 $ tnth = 1.0e-1, half = 0.5d+0 )
330 INTEGER I, IPTYPE, IWA, IWB
332REAL ABNORM, ANORM, BNORM, RATIO1, RATIO2, THRSH2,
341 EXTERNAL ILAENV, SLAMCH, SLANGE
347 INTRINSIC abs,
max, sqrt
357 IF( nsize.LT.0 )
THEN
359 ELSE IF( thresh.LT.zero )
THEN
361 ELSE IF( nin.LE.0 )
THEN
363 ELSE IF( nout.LE.0 )
THEN
365 ELSE IF( lda.LT.1 .OR. lda.LT.nmax )
THEN
367 ELSE IF( liwork.LT.nmax+6 )
THEN
379 IF( info.EQ.0 .AND. lwork.GE.1 )
THEN
380 minwrk = 2*nmax*nmax + 12*nmax + 16
381 maxwrk = 6*nmax + nmax*ilaenv( 1,
'SGEQRF',
' ', nmax, 1, nmax,
383 maxwrk =
max( maxwrk, 2*nmax*nmax+12*nmax+16 )
387 IF( lwork.LT.minwrk )
411 WEIGHT( 4 ) = ONE / WEIGHT( 2 )
412 WEIGHT( 5 ) = ONE / WEIGHT( 1 )
422 CALL SLATM6( IPTYPE, 5, A, LDA, B, VR, LDA, VL,
423 $ LDA, WEIGHT( IWA ), WEIGHT( IWB ),
424 $ WEIGHT( IWX ), WEIGHT( IWY ), STRU,
431 CALL SLACPY( 'f
', N, N, A, LDA, AI, LDA )
432 CALL SLACPY( 'f
', N, N, B, LDA, BI, LDA )
434 CALL SGGEVX( 'n
', 'v
', 'v
', 'b
', N, AI, LDA, BI,
435 $ LDA, ALPHAR, ALPHAI, BETA, VL, LDA,
436 $ VR, LDA, ILO, IHI, LSCALE, RSCALE,
437 $ ANORM, BNORM, S, DIF, WORK, LWORK,
438 $ IWORK, BWORK, LINFO )
439.NE.
IF( LINFO0 ) THEN
441 WRITE( NOUT, FMT = 9999 )'sggevx', LINFO, N,
448 CALL SLACPY( 'full
', N, N, AI, LDA, WORK, N )
449 CALL SLACPY( 'full
', N, N, BI, LDA, WORK( N*N+1 ),
451 ABNORM = SLANGE( 'fro
', N, 2*N, WORK, N, WORK )
456 CALL SGET52( .TRUE., N, A, LDA, B, LDA, VL, LDA,
457 $ ALPHAR, ALPHAI, BETA, WORK,
459.GT.
IF( RESULT( 2 )THRESH ) THEN
460 WRITE( NOUT, FMT = 9998 )'left
', 'sggevx',
461 $ RESULT( 2 ), N, IPTYPE, IWA, IWB, IWX, IWY
465 CALL SGET52( .FALSE., N, A, LDA, B, LDA, VR, LDA,
466 $ ALPHAR, ALPHAI, BETA, WORK,
468.GT.
IF( RESULT( 3 )THRESH ) THEN
469 WRITE( NOUT, FMT = 9998 )'right
', 'sggevx',
470 $ RESULT( 3 ), N, IPTYPE, IWA, IWB, IWX, IWY
477.EQ.
IF( S( I )ZERO ) THEN
478.GT.
IF( STRU( I )ABNORM*ULP )
479 $ RESULT( 3 ) = ULPINV
480.EQ.
ELSE IF( STRU( I )ZERO ) THEN
481.GT.
IF( S( I )ABNORM*ULP )
482 $ RESULT( 3 ) = ULPINV
484 WORK( I ) = MAX( ABS( STRU( I ) / S( I ) ),
485 $ ABS( S( I ) / STRU( I ) ) )
486 RESULT( 3 ) = MAX( RESULT( 3 ), WORK( I ) )
493.EQ.
IF( DIF( 1 )ZERO ) THEN
494.GT.
IF( DIFTRU( 1 )ABNORM*ULP )
495 $ RESULT( 4 ) = ULPINV
496.EQ.
ELSE IF( DIFTRU( 1 )ZERO ) THEN
497.GT.
IF( DIF( 1 )ABNORM*ULP )
498 $ RESULT( 4 ) = ULPINV
499.EQ.
ELSE IF( DIF( 5 )ZERO ) THEN
500.GT.
IF( DIFTRU( 5 )ABNORM*ULP )
501 $ RESULT( 4 ) = ULPINV
502.EQ.
ELSE IF( DIFTRU( 5 )ZERO ) THEN
503.GT.
IF( DIF( 5 )ABNORM*ULP )
504 $ RESULT( 4 ) = ULPINV
506 RATIO1 = MAX( ABS( DIFTRU( 1 ) / DIF( 1 ) ),
507 $ ABS( DIF( 1 ) / DIFTRU( 1 ) ) )
508 RATIO2 = MAX( ABS( DIFTRU( 5 ) / DIF( 5 ) ),
509 $ ABS( DIF( 5 ) / DIFTRU( 5 ) ) )
510 RESULT( 4 ) = MAX( RATIO1, RATIO2 )
518.GE..AND..GE..OR.
IF( ( RESULT( J )THRSH2 J4 )
519.GE..AND..LE.
$ ( RESULT( J )THRESH J3 ) )
525.EQ.
IF( NERRS0 ) THEN
526 WRITE( NOUT, FMT = 9997 )'sxv
'
532 WRITE( NOUT, FMT = 9995 )
533 WRITE( NOUT, FMT = 9994 )
534 WRITE( NOUT, FMT = 9993 )
538 WRITE( NOUT, FMT = 9992 )'''',
543.LT.
IF( RESULT( J )10000.0 ) THEN
544 WRITE( NOUT, FMT = 9991 )IPTYPE, IWA,
545 $ IWB, IWX, IWY, J, RESULT( J )
547 WRITE( NOUT, FMT = 9990 )IPTYPE, IWA,
548 $ IWB, IWX, IWY, J, RESULT( J )
568 READ( NIN, FMT = *, END = 150 )N
572 READ( NIN, FMT = * )( A( I, J ), J = 1, N )
575 READ( NIN, FMT = * )( B( I, J ), J = 1, N )
577 READ( NIN, FMT = * )( STRU( I ), I = 1, N )
578 READ( NIN, FMT = * )( DIFTRU( I ), I = 1, N )
586 CALL SLACPY( 'f
', N, N, A, LDA, AI, LDA )
587 CALL SLACPY( 'f
', N, N, B, LDA, BI, LDA )
589 CALL SGGEVX( 'n
', 'v
', 'v
', 'b
', N, AI, LDA, BI, LDA, ALPHAR,
590 $ ALPHAI, BETA, VL, LDA, VR, LDA, ILO, IHI, LSCALE,
591 $ RSCALE, ANORM, BNORM, S, DIF, WORK, LWORK, IWORK,
594.NE.
IF( LINFO0 ) THEN
596 WRITE( NOUT, FMT = 9987 )'sggevx', LINFO, N, NPTKNT
602 CALL SLACPY( 'full
', N, N, AI, LDA, WORK, N )
603 CALL SLACPY( 'full
', N, N, BI, LDA, WORK( N*N+1 ), N )
604 ABNORM = SLANGE( 'fro
', N, 2*N, WORK, N, WORK )
609 CALL SGET52( .TRUE., N, A, LDA, B, LDA, VL, LDA, ALPHAR, ALPHAI,
610 $ BETA, WORK, RESULT( 1 ) )
611.GT.
IF( RESULT( 2 )THRESH ) THEN
612 WRITE( NOUT, FMT = 9986 )'left
', 'sggevx', RESULT( 2 ), N,
617 CALL SGET52( .FALSE., N, A, LDA, B, LDA, VR, LDA, ALPHAR, ALPHAI,
618 $ BETA, WORK, RESULT( 2 ) )
619.GT.
IF( RESULT( 3 )THRESH ) THEN
620 WRITE( NOUT, FMT = 9986 )'right
', 'sggevx', RESULT( 3 ), N,
628.EQ.
IF( S( I )ZERO ) THEN
629.GT.
IF( STRU( I )ABNORM*ULP )
630 $ RESULT( 3 ) = ULPINV
631.EQ.
ELSE IF( STRU( I )ZERO ) THEN
632.GT.
IF( S( I )ABNORM*ULP )
633 $ RESULT( 3 ) = ULPINV
635 WORK( I ) = MAX( ABS( STRU( I ) / S( I ) ),
636 $ ABS( S( I ) / STRU( I ) ) )
637 RESULT( 3 ) = MAX( RESULT( 3 ), WORK( I ) )
644.EQ.
IF( DIF( 1 )ZERO ) THEN
645.GT.
IF( DIFTRU( 1 )ABNORM*ULP )
646 $ RESULT( 4 ) = ULPINV
647.EQ.
ELSE IF( DIFTRU( 1 )ZERO ) THEN
648.GT.
IF( DIF( 1 )ABNORM*ULP )
649 $ RESULT( 4 ) = ULPINV
650.EQ.
ELSE IF( DIF( 5 )ZERO ) THEN
651.GT.
IF( DIFTRU( 5 )ABNORM*ULP )
652 $ RESULT( 4 ) = ULPINV
653.EQ.
ELSE IF( DIFTRU( 5 )ZERO ) THEN
654.GT.
IF( DIF( 5 )ABNORM*ULP )
655 $ RESULT( 4 ) = ULPINV
657 RATIO1 = MAX( ABS( DIFTRU( 1 ) / DIF( 1 ) ),
658 $ ABS( DIF( 1 ) / DIFTRU( 1 ) ) )
659 RATIO2 = MAX( ABS( DIFTRU( 5 ) / DIF( 5 ) ),
660 $ ABS( DIF( 5 ) / DIFTRU( 5 ) ) )
661 RESULT( 4 ) = MAX( RATIO1, RATIO2 )
669.GE.
IF( RESULT( J )THRSH2 ) THEN
674.EQ.
IF( NERRS0 ) THEN
675 WRITE( NOUT, FMT = 9997 )'sxv
'
681 WRITE( NOUT, FMT = 9996 )
685 WRITE( NOUT, FMT = 9992 )'''', 'transpose
', ''''
689.LT.
IF( RESULT( J )10000.0 ) THEN
690 WRITE( NOUT, FMT = 9989 )NPTKNT, N, J, RESULT( J )
692 WRITE( NOUT, FMT = 9988 )NPTKNT, N, J, RESULT( J )
704 CALL ALASVM( 'sxv
', NOUT, NERRS, NTESTT, 0 )
710 9999 FORMAT( ' sdrgvx:
', A, ' returned info=
', I6, '.
', / 9X, 'n=
',
711 $ I6, ', jtype=
', I6, ')
' )
713 9998 FORMAT( ' sdrgvx:
', A, ' eigenvectors from
', A, ' incorrectly
',
714 $ 'normalized.
', / ' bits of error=
', 0P, G10.3, ',
', 9X,
715 $ 'n=
', I6, ', jtype=
', I6, ', iwa=
', I5, ', iwb=
', I5,
716 $ ', iwx=
', I5, ', iwy=
', I5 )
718 9997 FORMAT( / 1X, A3, ' -- real expert eigenvalue/vector
',
719 $ ' problem driver
' )
721 9996 FORMAT( ' input
' )
723 9995 FORMAT( ' matrix types:
', / )
725 9994 FORMAT( ' TYPE 1: da is diagonal, db is identity,
',
726 $ / ' a = y^(-h) da x^(-1), b = y^(-h) db x^(-1)
',
727 $ / ' yh and x are left and right eigenvectors.
', / )
729 9993 FORMAT( ' TYPE 2: da is quasi-diagonal, db is identity,
',
730 $ / ' a = y^(-h) da x^(-1), b = y^(-h) db x^(-1)
',
731 $ / ' yh and x are left and right eigenvectors.
', / )
733 9992 FORMAT( / ' tests performed:
', / 4X,
734 $ ' a is
alpha, b is beta, l is a left eigenvector, ', / 4x,
735 $
' r is a right eigenvector and ', a,
' means ', a,
'.',
736 $ /
' 1 = max | ( b A - a B )', a,
' l | / const.',
737 $ /
' 2 = max | ( b A - a B ) r | / const.',
738 $ /
' 3 = max ( Sest/Stru, Stru/Sest ) ',
739 $
' over all eigenvalues', /
740 $
' 4 = max( DIFest/DIFtru, DIFtru/DIFest ) ',
741 $
' over the 1st and 5th eigenvectors'
743 9991
FORMAT(
' Type='','' IWA=', i2,
', IWB=', i2,
', IWX=',
744 $ i2,
', IWY=', i2,
', result ', i2,
' is', 0p, f8.2 )
745 9990
FORMAT(
' Type=', i2,
',',
' IWA=', i2,
', IWB='', IWX=',
746 $ i2,
', IWY=', i2,
', result ', i2,
' is', 1p, e10.3 )
747 9989
FORMAT(
' Input example #', i2,
', matrix order=', i4,
',',
748 $
' result ', i2,
' is', 0p, f8.2 )
749 9988
FORMAT(
' Input example #', i2,
', matrix order=', i4,
',',
750 $
' result ', i2,
' is', 1p, e10.3 )
751 9
FORMAT' SDRGVX: '' returned INFO=''.', / 9x,
'N=',
752 $ i6,
', Input example #', i2,
')' )
754 9986
FORMAT(
' SDRGVX: ', a,
' Eigenvectors from ', a,
' incorrectly ',
755 $
'normalized.', /
' Bits of error=', 0p, g10.3,
',', 9x,
756 $
'N=', i6,
', Input Example #', i2,
')' )
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 slatm6(type, n, a, lda, b, x, ldx, y, ldy, alpha, beta, wx, wy, s, dif)
SLATM6
subroutine sdrgvx(nsize, thresh, nin, nout, a, lda, b, ai, bi, alphar, alphai, beta, vl, vr, ilo, ihi, lscale, rscale, s, stru, dif, diftru, work, lwork, iwork, liwork, result, bwork, info)
SDRGVX
subroutine sget52(left, n, a, lda, b, ldb, e, lde, alphar, alphai, beta, work, result)
SGET52