404 SUBROUTINE ddrvev( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
405 $ nounit, a, lda, h, wr, wi, wr1, wi1, vl, ldvl,
406 $ vr, ldvr, lre, ldlre, result, work, nwork,
415 INTEGER INFO, LDA, LDLRE, LDVL, LDVR, NOUNIT, NSIZES,
417 DOUBLE PRECISION THRESH
421 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
422 DOUBLE PRECISION A( lda, * ), H( lda, * ), LRE( ldlre, * ),
423 $ result( 7 ), vl( ldvl, * ), vr( ldvr, * ),
424 $ wi( * ), wi1( * ), work( * ), wr( * ), wr1( * )
430 DOUBLE PRECISION ZERO, ONE
431 parameter ( zero = 0.0d0, one = 1.0d0 )
433 parameter ( two = 2.0d0 )
435 parameter ( maxtyp = 21 )
440 INTEGER IINFO, IMODE, ITYPE, IWK, J, JCOL, JJ, JSIZE,
441 $ jtype, mtypes, n, nerrs, nfail, nmax, nnwork,
442 $ ntest, ntestf, ntestt
443 DOUBLE PRECISION ANORM, COND, CONDS, OVFL, RTULP, RTULPI, TNRM,
444 $ ulp, ulpinv, unfl, vmx, vrmx, vtst
447 CHARACTER ADUMMA( 1 )
448 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KCONDS( maxtyp ),
449 $ kmagn( maxtyp ), kmode( maxtyp ),
451 DOUBLE PRECISION DUM( 1 ), RES( 2 )
454 DOUBLE PRECISION DLAMCH, DLAPY2, DNRM2
455 EXTERNAL dlamch, dlapy2, dnrm2
462 INTRINSIC abs, max, min, sqrt
465 DATA ktype / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
466 DATA kmagn / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
468 DATA kmode / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
469 $ 1, 5, 5, 5, 4, 3, 1 /
470 DATA kconds / 3*0, 5*0, 4*1, 6*2, 3*0 /
474 path( 1: 1 ) =
'Double precision'
488 nmax = max( nmax, nn( j ) )
495 IF( nsizes.LT.0 )
THEN
497 ELSE IF( badnn )
THEN
499 ELSE IF( ntypes.LT.0 )
THEN
501 ELSE IF( thresh.LT.zero )
THEN
503 ELSE IF( nounit.LE.0 )
THEN
505 ELSE IF( lda.LT.1 .OR. lda.LT.nmax )
THEN
507 ELSE IF( ldvl.LT.1 .OR. ldvl.LT.nmax )
THEN
509 ELSE IF( ldvr.LT.1 .OR. ldvr.LT.nmax )
THEN
511 ELSE IF( ldlre.LT.1 .OR. ldlre.LT.nmax )
THEN
513 ELSE IF( 5*nmax+2*nmax**2.GT.nwork )
THEN
518 CALL xerbla(
'DDRVEV', -info )
524 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
529 unfl = dlamch(
'Safe minimum' )
532 ulp = dlamch(
'Precision' )
541 DO 270 jsize = 1, nsizes
543 IF( nsizes.NE.1 )
THEN
544 mtypes = min( maxtyp, ntypes )
546 mtypes = min( maxtyp+1, ntypes )
549 DO 260 jtype = 1, mtypes
550 IF( .NOT.dotype( jtype ) )
556 ioldsd( j ) = iseed( j )
575 IF( mtypes.GT.maxtyp )
578 itype = ktype( jtype )
579 imode = kmode( jtype )
583 GO TO ( 30, 40, 50 )kmagn( jtype )
599 CALL dlaset(
'Full', lda, n, zero, zero, a, lda )
607 IF( itype.EQ.1 )
THEN
610 ELSE IF( itype.EQ.2 )
THEN
615 a( jcol, jcol ) = anorm
618 ELSE IF( itype.EQ.3 )
THEN
623 a( jcol, jcol ) = anorm
625 $ a( jcol, jcol-1 ) = one
628 ELSE IF( itype.EQ.4 )
THEN
632 CALL dlatms( n, n,
'S', iseed,
'S', work, imode, cond,
633 $ anorm, 0, 0,
'N', a, lda, work( n+1 ),
636 ELSE IF( itype.EQ.5 )
THEN
640 CALL dlatms( n, n,
'S', iseed,
'S', work, imode, cond,
641 $ anorm, n, n,
'N', a, lda, work( n+1 ),
644 ELSE IF( itype.EQ.6 )
THEN
648 IF( kconds( jtype ).EQ.1 )
THEN
650 ELSE IF( kconds( jtype ).EQ.2 )
THEN
657 CALL dlatme( n,
'S', iseed, work, imode, cond, one,
658 $ adumma,
'T',
'T',
'T', work( n+1 ), 4,
659 $ conds, n, n, anorm, a, lda, work( 2*n+1 ),
662 ELSE IF( itype.EQ.7 )
THEN
666 CALL dlatmr( n, n,
'S', iseed,
'S', work, 6, one, one,
667 $
'T',
'N', work( n+1 ), 1, one,
668 $ work( 2*n+1 ), 1, one,
'N', idumma, 0, 0,
669 $ zero, anorm,
'NO', a, lda, iwork, iinfo )
671 ELSE IF( itype.EQ.8 )
THEN
675 CALL dlatmr( n, n,
'S', iseed,
'S', work, 6, one, one,
676 $
'T',
'N', work( n+1 ), 1, one,
677 $ work( 2*n+1 ), 1, one,
'N', idumma, n, n,
678 $ zero, anorm,
'NO', a, lda, iwork, iinfo )
680 ELSE IF( itype.EQ.9 )
THEN
684 CALL dlatmr( n, n,
'S', iseed,
'N', work, 6, one, one,
685 $
'T',
'N', work( n+1 ), 1, one,
686 $ work( 2*n+1 ), 1, one,
'N', idumma, n, n,
687 $ zero, anorm,
'NO', a, lda, iwork, iinfo )
689 CALL dlaset(
'Full', 2, n, zero, zero, a, lda )
690 CALL dlaset(
'Full', n-3, 1, zero, zero, a( 3, 1 ),
692 CALL dlaset(
'Full', n-3, 2, zero, zero, a( 3, n-1 ),
694 CALL dlaset(
'Full', 1, n, zero, zero, a( n, 1 ),
698 ELSE IF( itype.EQ.10 )
THEN
702 CALL dlatmr( n, n,
'S', iseed,
'N', work, 6, one, one,
703 $
'T',
'N', work( n+1 ), 1, one,
704 $ work( 2*n+1 ), 1, one,
'N', idumma, n, 0,
705 $ zero, anorm,
'NO', a, lda, iwork, iinfo )
712 IF( iinfo.NE.0 )
THEN
713 WRITE( nounit, fmt = 9993 )
'Generator', iinfo, n, jtype,
727 nnwork = 5*n + 2*n**2
729 nnwork = max( nnwork, 1 )
739 CALL dlacpy(
'F', n, n, a, lda, h, lda )
740 CALL dgeev(
'V',
'V', n, h, lda, wr, wi, vl, ldvl, vr,
741 $ ldvr, work, nnwork, iinfo )
742 IF( iinfo.NE.0 )
THEN
744 WRITE( nounit, fmt = 9993 )
'DGEEV1', iinfo, n, jtype,
752 CALL dget22(
'N',
'N',
'N', n, a, lda, vr, ldvr, wr, wi,
754 result( 1 ) = res( 1 )
758 CALL dget22(
'T',
'N',
'T', n, a, lda, vl, ldvl, wr, wi,
760 result( 2 ) = res( 1 )
766 IF( wi( j ).EQ.zero )
THEN
767 tnrm = dnrm2( n, vr( 1, j ), 1 )
768 ELSE IF( wi( j ).GT.zero )
THEN
769 tnrm = dlapy2( dnrm2( n, vr( 1, j ), 1 ),
770 $ dnrm2( n, vr( 1, j+1 ), 1 ) )
772 result( 3 ) = max( result( 3 ),
773 $ min( ulpinv, abs( tnrm-one ) / ulp ) )
774 IF( wi( j ).GT.zero )
THEN
778 vtst = dlapy2( vr( jj, j ), vr( jj, j+1 ) )
781 IF( vr( jj, j+1 ).EQ.zero .AND.
782 $ abs( vr( jj, j ) ).GT.vrmx )
783 $ vrmx = abs( vr( jj, j ) )
785 IF( vrmx / vmx.LT.one-two*ulp )
786 $ result( 3 ) = ulpinv
794 IF( wi( j ).EQ.zero )
THEN
795 tnrm = dnrm2( n, vl( 1, j ), 1 )
796 ELSE IF( wi( j ).GT.zero )
THEN
797 tnrm = dlapy2( dnrm2( n, vl( 1, j ), 1 ),
798 $ dnrm2( n, vl( 1, j+1 ), 1 ) )
800 result( 4 ) = max( result( 4 ),
801 $ min( ulpinv, abs( tnrm-one ) / ulp ) )
802 IF( wi( j ).GT.zero )
THEN
806 vtst = dlapy2( vl( jj, j ), vl( jj, j+1 ) )
809 IF( vl( jj, j+1 ).EQ.zero .AND.
810 $ abs( vl( jj, j ) ).GT.vrmx )
811 $ vrmx = abs( vl( jj, j ) )
813 IF( vrmx / vmx.LT.one-two*ulp )
814 $ result( 4 ) = ulpinv
820 CALL dlacpy(
'F', n, n, a, lda, h, lda )
821 CALL dgeev(
'N',
'N', n, h, lda, wr1, wi1, dum, 1, dum,
822 $ 1, work, nnwork, iinfo )
823 IF( iinfo.NE.0 )
THEN
825 WRITE( nounit, fmt = 9993 )
'DGEEV2', iinfo, n, jtype,
834 IF( wr( j ).NE.wr1( j ) .OR. wi( j ).NE.wi1( j ) )
835 $ result( 5 ) = ulpinv
840 CALL dlacpy(
'F', n, n, a, lda, h, lda )
841 CALL dgeev(
'N',
'V', n, h, lda, wr1, wi1, dum, 1, lre,
842 $ ldlre, work, nnwork, iinfo )
843 IF( iinfo.NE.0 )
THEN
845 WRITE( nounit, fmt = 9993 )
'DGEEV3', iinfo, n, jtype,
854 IF( wr( j ).NE.wr1( j ) .OR. wi( j ).NE.wi1( j ) )
855 $ result( 5 ) = ulpinv
862 IF( vr( j, jj ).NE.lre( j, jj ) )
863 $ result( 6 ) = ulpinv
869 CALL dlacpy(
'F', n, n, a, lda, h, lda )
870 CALL dgeev(
'V',
'N', n, h, lda, wr1, wi1, lre, ldlre,
871 $ dum, 1, work, nnwork, iinfo )
872 IF( iinfo.NE.0 )
THEN
874 WRITE( nounit, fmt = 9993 )
'DGEEV4', iinfo, n, jtype,
883 IF( wr( j ).NE.wr1( j ) .OR. wi( j ).NE.wi1( j ) )
884 $ result( 5 ) = ulpinv
891 IF( vl( j, jj ).NE.lre( j, jj ) )
892 $ result( 7 ) = ulpinv
903 IF( result( j ).GE.zero )
905 IF( result( j ).GE.thresh )
910 $ ntestf = ntestf + 1
911 IF( ntestf.EQ.1 )
THEN
912 WRITE( nounit, fmt = 9999 )path
913 WRITE( nounit, fmt = 9998 )
914 WRITE( nounit, fmt = 9997 )
915 WRITE( nounit, fmt = 9996 )
916 WRITE( nounit, fmt = 9995 )thresh
921 IF( result( j ).GE.thresh )
THEN
922 WRITE( nounit, fmt = 9994 )n, iwk, ioldsd, jtype,
927 nerrs = nerrs + nfail
928 ntestt = ntestt + ntest
936 CALL dlasum( path, nounit, nerrs, ntestt )
938 9999
FORMAT( / 1x, a3,
' -- Real Eigenvalue-Eigenvector Decomposition',
939 $
' Driver', /
' Matrix types (see DDRVEV for details): ' )
941 9998
FORMAT( /
' Special Matrices:', /
' 1=Zero matrix. ',
942 $
' ',
' 5=Diagonal: geometr. spaced entries.',
943 $ /
' 2=Identity matrix. ',
' 6=Diagona',
944 $
'l: clustered entries.', /
' 3=Transposed Jordan block. ',
945 $
' ',
' 7=Diagonal: large, evenly spaced.', /
' ',
946 $
'4=Diagonal: evenly spaced entries. ',
' 8=Diagonal: s',
947 $
'mall, evenly spaced.' )
948 9997
FORMAT(
' Dense, Non-Symmetric Matrices:', /
' 9=Well-cond., ev',
949 $
'enly spaced eigenvals.',
' 14=Ill-cond., geomet. spaced e',
950 $
'igenals.', /
' 10=Well-cond., geom. spaced eigenvals. ',
951 $
' 15=Ill-conditioned, clustered e.vals.', /
' 11=Well-cond',
952 $
'itioned, clustered e.vals. ',
' 16=Ill-cond., random comp',
953 $
'lex ', /
' 12=Well-cond., random complex ', 6x,
' ',
954 $
' 17=Ill-cond., large rand. complx ', /
' 13=Ill-condi',
955 $
'tioned, evenly spaced. ',
' 18=Ill-cond., small rand.',
957 9996
FORMAT(
' 19=Matrix with random O(1) entries. ',
' 21=Matrix ',
958 $
'with small random entries.', /
' 20=Matrix with large ran',
959 $
'dom entries. ', / )
960 9995
FORMAT(
' Tests performed with test threshold =', f8.2,
961 $ / /
' 1 = | A VR - VR W | / ( n |A| ulp ) ',
962 $ /
' 2 = | transpose(A) VL - VL W | / ( n |A| ulp ) ',
963 $ /
' 3 = | |VR(i)| - 1 | / ulp ',
964 $ /
' 4 = | |VL(i)| - 1 | / ulp ',
965 $ /
' 5 = 0 if W same no matter if VR or VL computed,',
966 $
' 1/ulp otherwise', /
967 $
' 6 = 0 if VR same no matter if VL computed,',
968 $
' 1/ulp otherwise', /
969 $
' 7 = 0 if VL same no matter if VR computed,',
970 $
' 1/ulp otherwise', / )
971 9994
FORMAT(
' N=', i5,
', IWK=', i2,
', seed=', 4( i4,
',' ),
972 $
' type ', i2,
', test(', i2,
')=', g10.3 )
973 9993
FORMAT(
' DDRVEV: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
974 $ i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
subroutine dget22(TRANSA, TRANSE, TRANSW, N, A, LDA, E, LDE, WR, WI, WORK, RESULT)
DGET22
subroutine dlatmr(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, RSIGN, GRADE, DL, MODEL, CONDL, DR, MODER, CONDR, PIVTNG, IPIVOT, KL, KU, SPARSE, ANORM, PACK, A, LDA, IWORK, INFO)
DLATMR
subroutine dlatme(N, DIST, ISEED, D, MODE, COND, DMAX, EI, RSIGN, UPPER, SIM, DS, MODES, CONDS, KL, KU, ANORM, A, LDA, WORK, INFO)
DLATME
subroutine ddrvev(NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, NOUNIT, A, LDA, H, WR, WI, WR1, WI1, VL, LDVL, VR, LDVR, LRE, LDLRE, RESULT, WORK, NWORK, IWORK, INFO)
DDRVEV
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 xerbla(SRNAME, INFO)
XERBLA
subroutine dlabad(SMALL, LARGE)
DLABAD
subroutine dlasum(TYPE, IOUNIT, IE, NRUN)
DLASUM
subroutine dgeev(JOBVL, JOBVR, N, A, LDA, WR, WI, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)
DGEEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices ...
subroutine dlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
DLATMS