416 SUBROUTINE schkhs( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
417 $ NOUNIT, A, LDA, H, T1, T2, U, LDU, Z, UZ, WR1,
418 $ WI1, WR2, WI2, WR3, WI3, EVECTL, EVECTR,
419 $ EVECTY, EVECTX, UU, TAU, WORK, NWORK, IWORK,
420 $ SELECT, RESULT, INFO )
427 INTEGER INFO, LDA, LDU, NOUNIT, NSIZES, NTYPES, NWORK
431 LOGICAL DOTYPE( * ), SELECT( * )
432 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
433 REAL A( LDA, * ), EVECTL( LDU, * ),
434 $ EVECTR( LDU, * ), EVECTX( LDU, * ),
435 $ evecty( ldu, * ), h( lda, * ), result( 16 ),
436 $ t1( lda, * ), t2( lda, * ), tau( * ),
437 $ u( ldu, * ), uu( ldu, * ), uz( ldu, * ),
438 $ wi1( * ), wi2( * ), wi3( * ), work( * ),
439 $ wr1( * ), wr2( * ), wr3( * ), z( ldu, * )
446 PARAMETER ( ZERO = 0.0, one = 1.0 )
448 PARAMETER ( MAXTYP = 21 )
452 INTEGER I, IHI, IINFO, ILO, IMODE, IN, ITYPE, J, JCOL,
453 $ JJ, JSIZE, JTYPE, K, MTYPES, N, N1, NERRS,
454 $ NMATS, NMAX, NSELC, NSELR, NTEST, NTESTT
455 REAL ANINV, ANORM, COND, CONDS, OVFL, RTOVFL, RTULP,
456 $ rtulpi, rtunfl, temp1, temp2, ulp, ulpinv, unfl
459 CHARACTER ADUMMA( 1 )
460 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KCONDS( MAXTYP ),
461 $ KMAGN( MAXTYP ), KMODE( MAXTYP ),
476 INTRINSIC abs, max, min, real, sqrt
479 DATA ktype / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
480 DATA kmagn / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
482 DATA kmode / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
483 $ 1, 5, 5, 5, 4, 3, 1 /
484 DATA kconds / 3*0, 5*0, 4*1, 6*2, 3*0 /
496 nmax = max( nmax, nn( j ) )
503 IF( nsizes.LT.0 )
THEN
505 ELSE IF( badnn )
THEN
507 ELSE IF( ntypes.LT.0 )
THEN
509 ELSE IF( thresh.LT.zero )
THEN
511 ELSE IF( lda.LE.1 .OR. lda.LT.nmax )
THEN
513 ELSE IF( ldu.LE.1 .OR. ldu.LT.nmax )
THEN
515 ELSE IF( 4*nmax*nmax+2.GT.nwork )
THEN
520 CALL xerbla(
'SCHKHS', -info )
526 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
531 unfl = slamch(
'Safe minimum' )
532 ovfl = slamch(
'Overflow' )
533 ulp = slamch(
'Epsilon' )*slamch(
'Base' )
535 rtunfl = sqrt( unfl )
536 rtovfl = sqrt( ovfl )
545 DO 270 jsize = 1, nsizes
550 aninv = one / real( n1 )
552 IF( nsizes.NE.1 )
THEN
553 mtypes = min( maxtyp, ntypes )
555 mtypes = min( maxtyp+1, ntypes )
558 DO 260 jtype = 1, mtypes
559 IF( .NOT.dotype( jtype ) )
567 ioldsd( j ) = iseed( j )
592 IF( mtypes.GT.maxtyp )
595 itype = ktype( jtype )
596 imode = kmode( jtype )
600 GO TO ( 40, 50, 60 )kmagn( jtype )
607 anorm = ( rtovfl*ulp )*aninv
611 anorm = rtunfl*n*ulpinv
616 CALL slaset(
'Full', lda, n, zero, zero, a, lda )
622 IF( itype.EQ.1 )
THEN
628 ELSE IF( itype.EQ.2 )
THEN
633 a( jcol, jcol ) = anorm
636 ELSE IF( itype.EQ.3 )
THEN
641 a( jcol, jcol ) = anorm
643 $ a( jcol, jcol-1 ) = one
646 ELSE IF( itype.EQ.4 )
THEN
650 CALL slatms( n, n,
'S', iseed,
'S', work, imode, cond,
651 $ anorm, 0, 0,
'N', a, lda, work( n+1 ),
654 ELSE IF( itype.EQ.5 )
THEN
658 CALL slatms( n, n,
'S', iseed,
'S', work, imode, cond,
659 $ anorm, n, n,
'N', a, lda, work( n+1 ),
662 ELSE IF( itype.EQ.6 )
THEN
666 IF( kconds( jtype ).EQ.1 )
THEN
668 ELSE IF( kconds( jtype ).EQ.2 )
THEN
675 CALL slatme( n,
'S', iseed, work, imode, cond, one,
676 $ adumma,
'T',
'T',
'T', work( n+1 ), 4,
677 $ conds, n, n, anorm, a, lda, work( 2*n+1 ),
680 ELSE IF( itype.EQ.7 )
THEN
684 CALL slatmr( n, n,
'S', iseed,
'S', work, 6, one, one,
685 $
'T',
'N', work( n+1 ), 1, one,
686 $ work( 2*n+1 ), 1, one,
'N', idumma, 0, 0,
687 $ zero, anorm,
'NO', a, lda, iwork, iinfo )
689 ELSE IF( itype.EQ.8 )
THEN
693 CALL slatmr( n, n,
'S', iseed,
'S', work, 6, one, one,
694 $
'T',
'N', work( n+1 ), 1, one,
695 $ work( 2*n+1 ), 1, one,
'N', idumma, n, n,
696 $ zero, anorm,
'NO', a, lda, iwork, iinfo )
698 ELSE IF( itype.EQ.9 )
THEN
702 CALL slatmr( 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, n,
705 $ zero, anorm,
'NO', a, lda, iwork, iinfo )
707 ELSE IF( itype.EQ.10 )
THEN
711 CALL slatmr( n, n,
'S', iseed,
'N', work, 6, one, one,
712 $
'T',
'N', work( n+1 ), 1, one,
713 $ work( 2*n+1 ), 1, one,
'N', idumma, n, 0,
714 $ zero, anorm,
'NO', a, lda, iwork, iinfo )
721 IF( iinfo.NE.0 )
THEN
722 WRITE( nounit, fmt = 9999 )
'Generator', iinfo, n, jtype,
732 CALL slacpy(
' ', n, n, a, lda, h, lda )
739 CALL sgehrd( n, ilo, ihi, h, lda, work, work( n+1 ),
742 IF( iinfo.NE.0 )
THEN
744 WRITE( nounit, fmt = 9999 )
'SGEHRD', iinfo, n, jtype,
753 u( i, j ) = h( i, j )
754 uu( i, j ) = h( i, j )
758 CALL scopy( n-1, work, 1, tau, 1 )
759 CALL sorghr( n, ilo, ihi, u, ldu, work, work( n+1 ),
763 CALL shst01( n, ilo, ihi, a, lda, h, lda, u, ldu, work,
764 $ nwork, result( 1 ) )
770 CALL slacpy(
' ', n, n, h, lda, t2, lda )
774 CALL shseqr(
'E',
'N', n, ilo, ihi, t2, lda, wr3, wi3, uz,
775 $ ldu, work, nwork, iinfo )
776 IF( iinfo.NE.0 )
THEN
777 WRITE( nounit, fmt = 9999 )
'SHSEQR(E)', iinfo, n, jtype,
779 IF( iinfo.LE.n+2 )
THEN
787 CALL slacpy(
' ', n, n, h, lda, t2, lda )
789 CALL shseqr(
'S',
'N', n, ilo, ihi, t2, lda, wr2, wi2, uz,
790 $ ldu, work, nwork, iinfo )
791 IF( iinfo.NE.0 .AND. iinfo.LE.n+2 )
THEN
792 WRITE( nounit, fmt = 9999 )
'SHSEQR(S)', iinfo, n, jtype,
801 CALL slacpy(
' ', n, n, h, lda, t1, lda )
802 CALL slacpy(
' ', n, n, u, ldu, uz, ldu )
804 CALL shseqr(
'S',
'V', n, ilo, ihi, t1, lda, wr1, wi1, uz,
805 $ ldu, work, nwork, iinfo )
806 IF( iinfo.NE.0 .AND. iinfo.LE.n+2 )
THEN
807 WRITE( nounit, fmt = 9999 )
'SHSEQR(V)', iinfo, n, jtype,
815 CALL sgemm(
'T',
'N', n, n, n, one, u, ldu, uz, ldu, zero,
822 CALL shst01( n, ilo, ihi, h, lda, t1, lda, z, ldu, work,
823 $ nwork, result( 3 ) )
828 CALL shst01( n, ilo, ihi, a, lda, t1, lda, uz, ldu, work,
829 $ nwork, result( 5 ) )
833 CALL sget10( n, n, t2, lda, t1, lda, work, result( 7 ) )
840 temp1 = max( temp1, abs( wr1( j ) )+abs( wi1( j ) ),
841 $ abs( wr2( j ) )+abs( wi2( j ) ) )
842 temp2 = max( temp2, abs( wr1( j )-wr2( j ) )+
843 $ abs( wi1( j )-wi2( j ) ) )
846 result( 8 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
861 IF( wi1( j ).EQ.zero )
THEN
862 IF( nselr.LT.max( n / 4, 1 ) )
THEN
866 SELECT( j ) = .false.
870 IF( nselc.LT.max( n / 4, 1 ) )
THEN
873 SELECT( j-1 ) = .false.
875 SELECT( j ) = .false.
876 SELECT( j-1 ) = .false.
883 CALL strevc(
'Right',
'All',
SELECT, n, t1, lda, dumma, ldu,
884 $ evectr, ldu, n, in, work, iinfo )
885 IF( iinfo.NE.0 )
THEN
886 WRITE( nounit, fmt = 9999 )
'STREVC(R,A)', iinfo, n,
894 CALL sget22(
'N',
'N',
'N', n, t1, lda, evectr, ldu, wr1,
895 $ wi1, work, dumma( 1 ) )
896 result( 9 ) = dumma( 1 )
897 IF( dumma( 2 ).GT.thresh )
THEN
898 WRITE( nounit, fmt = 9998 )
'Right',
'STREVC',
899 $ dumma( 2 ), n, jtype, ioldsd
905 CALL strevc(
'Right',
'Some',
SELECT, n, t1, lda, dumma,
906 $ ldu, evectl, ldu, n, in, work, iinfo )
907 IF( iinfo.NE.0 )
THEN
908 WRITE( nounit, fmt = 9999 )
'STREVC(R,S)', iinfo, n,
917 IF(
SELECT( j ) .AND. wi1( j ).EQ.zero )
THEN
919 IF( evectr( jj, j ).NE.evectl( jj, k ) )
THEN
925 ELSE IF(
SELECT( j ) .AND. wi1( j ).NE.zero )
THEN
927 IF( evectr( jj, j ).NE.evectl( jj, k ) .OR.
928 $ evectr( jj, j+1 ).NE.evectl( jj, k+1 ) )
THEN
938 $
WRITE( nounit, fmt = 9997 )
'Right',
'STREVC', n, jtype,
944 result( 10 ) = ulpinv
945 CALL strevc(
'Left',
'All',
SELECT, n, t1, lda, evectl, ldu,
946 $ dumma, ldu, n, in, work, iinfo )
947 IF( iinfo.NE.0 )
THEN
948 WRITE( nounit, fmt = 9999 )
'STREVC(L,A)', iinfo, n,
956 CALL sget22(
'Trans',
'N',
'Conj', n, t1, lda, evectl, ldu,
957 $ wr1, wi1, work, dumma( 3 ) )
958 result( 10 ) = dumma( 3 )
959 IF( dumma( 4 ).GT.thresh )
THEN
960 WRITE( nounit, fmt = 9998 )
'Left',
'STREVC', dumma( 4 ),
967 CALL strevc(
'Left',
'Some',
SELECT, n, t1, lda, evectr,
968 $ ldu, dumma, ldu, n, in, work, iinfo )
969 IF( iinfo.NE.0 )
THEN
970 WRITE( nounit, fmt = 9999 )
'STREVC(L,S)', iinfo, n,
979 IF(
SELECT( j ) .AND. wi1( j ).EQ.zero )
THEN
981 IF( evectl( jj, j ).NE.evectr( jj, k ) )
THEN
987 ELSE IF(
SELECT( j ) .AND. wi1( j ).NE.zero )
THEN
989 IF( evectl( jj, j ).NE.evectr( jj, k ) .OR.
990 $ evectl( jj, j+1 ).NE.evectr( jj, k+1 ) )
THEN
1000 $
WRITE( nounit, fmt = 9997 )
'Left',
'STREVC', n, jtype,
1006 result( 11 ) = ulpinv
1008 SELECT( j ) = .true.
1011 CALL shsein(
'Right',
'Qr',
'Ninitv',
SELECT, n, h, lda,
1012 $ wr3, wi3, dumma, ldu, evectx, ldu, n1, in,
1013 $ work, iwork, iwork, iinfo )
1014 IF( iinfo.NE.0 )
THEN
1015 WRITE( nounit, fmt = 9999 )
'SHSEIN(R)', iinfo, n, jtype,
1026 CALL sget22(
'N',
'N',
'N', n, h, lda, evectx, ldu, wr3,
1027 $ wi3, work, dumma( 1 ) )
1028 IF( dumma( 1 ).LT.ulpinv )
1029 $ result( 11 ) = dumma( 1 )*aninv
1030 IF( dumma( 2 ).GT.thresh )
THEN
1031 WRITE( nounit, fmt = 9998 )
'Right',
'SHSEIN',
1032 $ dumma( 2 ), n, jtype, ioldsd
1039 result( 12 ) = ulpinv
1041 SELECT( j ) = .true.
1044 CALL shsein(
'Left',
'Qr',
'Ninitv',
SELECT, n, h, lda, wr3,
1045 $ wi3, evecty, ldu, dumma, ldu, n1, in, work,
1046 $ iwork, iwork, iinfo )
1047 IF( iinfo.NE.0 )
THEN
1048 WRITE( nounit, fmt = 9999 )
'SHSEIN(L)', iinfo, n, jtype,
1059 CALL sget22(
'C',
'N',
'C', n, h, lda, evecty, ldu, wr3,
1060 $ wi3, work, dumma( 3 ) )
1061 IF( dumma( 3 ).LT.ulpinv )
1062 $ result( 12 ) = dumma( 3 )*aninv
1063 IF( dumma( 4 ).GT.thresh )
THEN
1064 WRITE( nounit, fmt = 9998 )
'Left',
'SHSEIN',
1065 $ dumma( 4 ), n, jtype, ioldsd
1072 result( 13 ) = ulpinv
1074 CALL sormhr(
'Left',
'No transpose', n, n, ilo, ihi, uu,
1075 $ ldu, tau, evectx, ldu, work, nwork, iinfo )
1076 IF( iinfo.NE.0 )
THEN
1077 WRITE( nounit, fmt = 9999 )
'SORMHR(R)', iinfo, n, jtype,
1088 CALL sget22(
'N',
'N',
'N', n, a, lda, evectx, ldu, wr3,
1089 $ wi3, work, dumma( 1 ) )
1090 IF( dumma( 1 ).LT.ulpinv )
1091 $ result( 13 ) = dumma( 1 )*aninv
1097 result( 14 ) = ulpinv
1099 CALL sormhr(
'Left',
'No transpose', n, n, ilo, ihi, uu,
1100 $ ldu, tau, evecty, ldu, work, nwork, iinfo )
1101 IF( iinfo.NE.0 )
THEN
1102 WRITE( nounit, fmt = 9999 )
'SORMHR(L)', iinfo, n, jtype,
1113 CALL sget22(
'C',
'N',
'C', n, a, lda, evecty, ldu, wr3,
1114 $ wi3, work, dumma( 3 ) )
1115 IF( dumma( 3 ).LT.ulpinv )
1116 $ result( 14 ) = dumma( 3 )*aninv
1124 result( 15 ) = ulpinv
1126 CALL slacpy(
' ', n, n, uz, ldu, evectr, ldu )
1128 CALL strevc3(
'Right',
'Back',
SELECT, n, t1, lda, dumma,
1129 $ ldu, evectr, ldu, n, in, work, nwork, iinfo )
1130 IF( iinfo.NE.0 )
THEN
1131 WRITE( nounit, fmt = 9999 )
'STREVC3(R,B)', iinfo, n,
1141 CALL sget22(
'N',
'N',
'N', n, a, lda, evectr, ldu, wr1,
1142 $ wi1, work, dumma( 1 ) )
1143 result( 15 ) = dumma( 1 )
1144 IF( dumma( 2 ).GT.thresh )
THEN
1145 WRITE( nounit, fmt = 9998 )
'Right',
'STREVC3',
1146 $ dumma( 2 ), n, jtype, ioldsd
1152 result( 16 ) = ulpinv
1154 CALL slacpy(
' ', n, n, uz, ldu, evectl, ldu )
1156 CALL strevc3(
'Left',
'Back',
SELECT, n, t1, lda, evectl,
1157 $ ldu, dumma, ldu, n, in, work, nwork, iinfo )
1158 IF( iinfo.NE.0 )
THEN
1159 WRITE( nounit, fmt = 9999 )
'STREVC3(L,B)', iinfo, n,
1169 CALL sget22(
'Trans',
'N',
'Conj', n, a, lda, evectl, ldu,
1170 $ wr1, wi1, work, dumma( 3 ) )
1171 result( 16 ) = dumma( 3 )
1172 IF( dumma( 4 ).GT.thresh )
THEN
1173 WRITE( nounit, fmt = 9998 )
'Left',
'STREVC3', dumma( 4 ),
1181 ntestt = ntestt + ntest
1182 CALL slafts(
'SHS', n, n, jtype, ntest, result, ioldsd,
1183 $ thresh, nounit, nerrs )
1190 CALL slasum(
'SHS', nounit, nerrs, ntestt )
1194 9999
FORMAT(
' SCHKHS: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
1195 $ i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
1196 9998
FORMAT(
' SCHKHS: ', a,
' Eigenvectors from ', a,
' incorrectly ',
1197 $
'normalized.', /
' Bits of error=', 0p, g10.3,
',', 9x,
1198 $
'N=', i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
1200 9997
FORMAT(
' SCHKHS: Selected ', a,
' Eigenvectors from ', a,
1201 $
' do not match other eigenvectors ', 9x,
'N=', i6,
1202 $
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
subroutine xerbla(srname, info)
subroutine scopy(n, sx, incx, sy, incy)
SCOPY
subroutine sgehrd(n, ilo, ihi, a, lda, tau, work, lwork, info)
SGEHRD
subroutine sgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
SGEMM
subroutine shsein(side, eigsrc, initv, select, n, h, ldh, wr, wi, vl, ldvl, vr, ldvr, mm, m, work, ifaill, ifailr, info)
SHSEIN
subroutine shseqr(job, compz, n, ilo, ihi, h, ldh, wr, wi, z, ldz, work, lwork, info)
SHSEQR
subroutine slacpy(uplo, m, n, a, lda, b, ldb)
SLACPY copies all or part of one two-dimensional array to another.
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 strevc3(side, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, mm, m, work, lwork, info)
STREVC3
subroutine strevc(side, howmny, select, n, t, ldt, vl, ldvl, vr, ldvr, mm, m, work, info)
STREVC
subroutine sorghr(n, ilo, ihi, a, lda, tau, work, lwork, info)
SORGHR
subroutine sormhr(side, trans, m, n, ilo, ihi, a, lda, tau, c, ldc, work, lwork, info)
SORMHR
subroutine schkhs(nsizes, nn, ntypes, dotype, iseed, thresh, nounit, a, lda, h, t1, t2, u, ldu, z, uz, wr1, wi1, wr2, wi2, wr3, wi3, evectl, evectr, evecty, evectx, uu, tau, work, nwork, iwork, select, result, info)
SCHKHS
subroutine sget10(m, n, a, lda, b, ldb, work, result)
SGET10
subroutine sget22(transa, transe, transw, n, a, lda, e, lde, wr, wi, work, result)
SGET22
subroutine shst01(n, ilo, ihi, a, lda, h, ldh, q, ldq, work, lwork, result)
SHST01
subroutine slafts(type, m, n, imat, ntests, result, iseed, thresh, iounit, ie)
SLAFTS
subroutine slasum(type, iounit, ie, nrun)
SLASUM
subroutine slatme(n, dist, iseed, d, mode, cond, dmax, ei, rsign, upper, sim, ds, modes, conds, kl, ku, anorm, a, lda, work, info)
SLATME
subroutine slatmr(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)
SLATMR
subroutine slatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
SLATMS