404 SUBROUTINE sdrgev3( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
405 $ NOUNIT, A, LDA, B, S, T, Q, LDQ, Z, QE, LDQE,
406 $ ALPHAR, ALPHAI, BETA, ALPHR1, ALPHI1, BETA1,
407 $ WORK, LWORK, RESULT, INFO )
414 INTEGER INFO, LDA, LDQ, LDQE, LWORK, NOUNIT, NSIZES,
420 INTEGER ISEED( 4 ), NN( * )
421 REAL A( LDA, * ), ALPHAI( * ), ALPHI1( * ),
422 $ alphar( * ), alphr1( * ), b( lda, * ),
423 $ beta( * ), beta1( * ), q( ldq, * ),
424 $ qe( ldqe, * ), result( * ), s( lda, * ),
425 $ t( lda, * ), work( * ), z( ldq, * )
432 PARAMETER ( ZERO = 0.0e+0, one = 1.0e+0 )
434 parameter( maxtyp = 26 )
438 INTEGER I, IADD, IERR, IN, J, JC, JR, JSIZE, JTYPE,
439 $ MAXWRK, MINWRK, MTYPES, N, N1, NERRS, NMATS,
441 REAL SAFMAX, SAFMIN, ULP, ULPINV
444 INTEGER IASIGN( MAXTYP ), IBSIGN( MAXTYP ),
445 $ IOLDSD( 4 ), KADD( 6 ), KAMAGN( MAXTYP ),
446 $ KATYPE( MAXTYP ), KAZERO( MAXTYP ),
447 $ kbmagn( maxtyp ), kbtype( maxtyp ),
448 $ kbzero( maxtyp ), kclass( maxtyp ),
449 $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
455 EXTERNAL ILAENV, SLAMCH, SLARND
462 INTRINSIC abs, max, min, real, sign
465 DATA kclass / 15*1, 10*2, 1*3 /
466 DATA kz1 / 0, 1, 2, 1, 3, 3 /
467 DATA kz2 / 0, 0, 1, 2, 1, 1 /
468 DATA kadd / 0, 0, 0, 0, 3, 2 /
469 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
470 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
471 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
472 $ 1, 1, -4, 2, -4, 8*8, 0 /
473 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
475 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
477 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
479 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
481 DATA ktrian / 16*0, 10*1 /
482 DATA iasign / 6*0, 2, 0, 2*2, 2*0, 3*2, 0, 2, 3*0,
484 DATA ibsign / 7*0, 2, 2*0, 2*2, 2*0, 2, 0, 2, 9*0 /
495 nmax = max( nmax, nn( j ) )
500 IF( nsizes.LT.0 )
THEN
502 ELSE IF( badnn )
THEN
504 ELSE IF( ntypes.LT.0 )
THEN
506 ELSE IF( thresh.LT.zero )
THEN
508 ELSE IF( lda.LE.1 .OR. lda.LT.nmax )
THEN
510 ELSE IF( ldq.LE.1 .OR. ldq.LT.nmax )
THEN
512 ELSE IF( ldqe.LE.1 .OR. ldqe.LT.nmax )
THEN
524 IF( info.EQ.0 .AND. lwork.GE.1 )
THEN
525 minwrk = max( 1, 8*nmax, nmax*( nmax+1 ) )
526 maxwrk = 7*nmax + nmax*ilaenv( 1,
'SGEQRF',
' ', nmax, 1, nmax,
528 maxwrk = max( maxwrk, nmax*( nmax+1 ) )
532 IF( lwork.LT.minwrk )
536 CALL xerbla(
'SDRGEV3', -info )
542 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
545 safmin = slamch(
'Safe minimum' )
546 ulp = slamch(
'Epsilon' )*slamch(
'Base' )
547 safmin = safmin / ulp
548 safmax = one / safmin
562 DO 220 jsize = 1, nsizes
565 rmagn( 2 ) = safmax*ulp / real( n1 )
566 rmagn( 3 ) = safmin*ulpinv*n1
568 IF( nsizes.NE.1 )
THEN
569 mtypes = min( maxtyp, ntypes )
571 mtypes = min( maxtyp+1, ntypes )
574 DO 210 jtype = 1, mtypes
575 IF( .NOT.dotype( jtype ) )
582 ioldsd( j ) = iseed( j )
608 IF( mtypes.GT.maxtyp )
611 IF( kclass( jtype ).LT.3 )
THEN
615 IF( abs( katype( jtype ) ).EQ.3 )
THEN
616 in = 2*( ( n-1 ) / 2 ) + 1
618 $
CALL slaset(
'Full', n, n, zero, zero, a, lda )
622 CALL slatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
623 $ kz2( kazero( jtype ) ), iasign( jtype ),
624 $ rmagn( kamagn( jtype ) ), ulp,
625 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
627 iadd = kadd( kazero( jtype ) )
628 IF( iadd.GT.0 .AND. iadd.LE.n )
629 $ a( iadd, iadd ) = one
633 IF( abs( kbtype( jtype ) ).EQ.3 )
THEN
634 in = 2*( ( n-1 ) / 2 ) + 1
636 $
CALL slaset(
'Full', n, n, zero, zero, b, lda )
640 CALL slatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
641 $ kz2( kbzero( jtype ) ), ibsign( jtype ),
642 $ rmagn( kbmagn( jtype ) ), one,
643 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
645 iadd = kadd( kbzero( jtype ) )
646 IF( iadd.NE.0 .AND. iadd.LE.n )
647 $ b( iadd, iadd ) = one
649 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 )
THEN
658 q( jr, jc ) = slarnd( 3, iseed )
659 z( jr, jc ) = slarnd( 3, iseed )
661 CALL slarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
663 work( 2*n+jc ) = sign( one, q( jc, jc ) )
665 CALL slarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
667 work( 3*n+jc ) = sign( one, z( jc, jc ) )
672 work( 3*n ) = sign( one, slarnd( 2, iseed ) )
675 work( 4*n ) = sign( one, slarnd( 2, iseed ) )
681 a( jr, jc ) = work( 2*n+jr )*work( 3*n+jc )*
683 b( jr, jc ) = work( 2*n+jr )*work( 3*n+jc )*
687 CALL sorm2r(
'L',
'N', n, n, n-1, q, ldq, work, a,
688 $ lda, work( 2*n+1 ), ierr )
691 CALL sorm2r(
'R',
'T', n, n, n-1, z, ldq, work( n+1 ),
692 $ a, lda, work( 2*n+1 ), ierr )
695 CALL sorm2r(
'L',
'N', n, n, n-1, q, ldq, work, b,
696 $ lda, work( 2*n+1 ), ierr )
699 CALL sorm2r(
'R',
'T', n, n, n-1, z, ldq, work( n+1 ),
700 $ b, lda, work( 2*n+1 ), ierr )
710 a( jr, jc ) = rmagn( kamagn( jtype ) )*
712 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
721 WRITE( nounit, fmt = 9999 )
'Generator', ierr, n, jtype,
743 CALL slacpy(
' ', n, n, a, lda, s, lda )
744 CALL slacpy(
' ', n, n, b, lda, t, lda )
745 CALL sggev3(
'V',
'V', n, s, lda, t, lda, alphar, alphai,
746 $ beta, q, ldq, z, ldq, work, lwork, ierr )
747 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
749 WRITE( nounit, fmt = 9999 )
'SGGEV31', ierr, n, jtype,
757 CALL sget52( .true., n, a, lda, b, lda, q, ldq, alphar,
758 $ alphai, beta, work, result( 1 ) )
759 IF( result( 2 ).GT.thresh )
THEN
760 WRITE( nounit, fmt = 9998 )
'Left',
'SGGEV31',
761 $ result( 2 ), n, jtype, ioldsd
766 CALL sget52( .false., n, a, lda, b, lda, z, ldq, alphar,
767 $ alphai, beta, work, result( 3 ) )
768 IF( result( 4 ).GT.thresh )
THEN
769 WRITE( nounit, fmt = 9998 )
'Right',
'SGGEV31',
770 $ result( 4 ), n, jtype, ioldsd
775 CALL slacpy(
' ', n, n, a, lda, s, lda )
776 CALL slacpy(
' ', n, n, b, lda, t, lda )
777 CALL sggev3(
'N',
'N', n, s, lda, t, lda, alphr1, alphi1,
778 $ beta1, q, ldq, z, ldq, work, lwork, ierr )
779 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
781 WRITE( nounit, fmt = 9999 )
'SGGEV32', ierr, n, jtype,
788 IF( alphar( j ).NE.alphr1( j ) .OR.
789 $ beta( j ).NE. beta1( j ) )
THEN
797 CALL slacpy(
' ', n, n, a, lda, s, lda )
798 CALL slacpy(
' ', n, n, b, lda, t, lda )
799 CALL sggev3(
'V',
'N', n, s, lda, t, lda, alphr1, alphi1,
800 $ beta1, qe, ldqe, z, ldq, work, lwork, ierr )
801 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
803 WRITE( nounit, fmt = 9999 )
'SGGEV33', ierr, n, jtype,
810 IF( alphar( j ).NE.alphr1( j ) .OR. alphai( j ).NE.
811 $ alphi1( j ) .OR. beta( j ).NE.beta1( j ) )
812 $ result( 6 ) = ulpinv
817 IF( q( j, jc ).NE.qe( j, jc ) )
818 $ result( 6 ) = ulpinv
825 CALL slacpy(
' ', n, n, a, lda, s, lda )
826 CALL slacpy(
' ', n, n, b, lda, t, lda )
827 CALL sggev3(
'N',
'V', n, s, lda, t, lda, alphr1, alphi1,
828 $ beta1, q, ldq, qe, ldqe, work, lwork, ierr )
829 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
831 WRITE( nounit, fmt = 9999 )
'SGGEV34', ierr, n, jtype,
838 IF( alphar( j ).NE.alphr1( j ) .OR. alphai( j ).NE.
839 $ alphi1( j ) .OR. beta( j ).NE.beta1( j ) )
840 $ result( 7 ) = ulpinv
845 IF( z( j, jc ).NE.qe( j, jc ) )
846 $ result( 7 ) = ulpinv
859 IF( result( jr ).GE.thresh )
THEN
864 IF( nerrs.EQ.0 )
THEN
865 WRITE( nounit, fmt = 9997 )
'SGV'
869 WRITE( nounit, fmt = 9996 )
870 WRITE( nounit, fmt = 9995 )
871 WRITE( nounit, fmt = 9994 )
'Orthogonal'
875 WRITE( nounit, fmt = 9993 )
879 IF( result( jr ).LT.10000.0 )
THEN
880 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
883 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
894 CALL alasvm(
'SGV', nounit, nerrs, ntestt, 0 )
900 9999
FORMAT(
' SDRGEV3: ', a,
' returned INFO=', i6,
'.', / 3x,
'N=',
901 $ i6,
', JTYPE=', i6,
', ISEED=(', 4( i4,
',' ), i5,
')' )
903 9998
FORMAT(
' SDRGEV3: ', a,
' Eigenvectors from ', a,
904 $
' incorrectly normalized.', /
' Bits of error=', 0p, g10.3,
905 $
',', 3x,
'N=', i4,
', JTYPE=', i3,
', ISEED=(',
906 $ 4( i4,
',' ), i5,
')' )
908 9997
FORMAT( / 1x, a3,
' -- Real Generalized eigenvalue problem driver'
911 9996
FORMAT(
' Matrix types (see SDRGEV3 for details): ' )
913 9995
FORMAT(
' Special Matrices:', 23x,
914 $
'(J''=transposed Jordan block)',
915 $ /
' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
916 $
'6=(diag(J'',I), diag(I,J''))', /
' Diagonal Matrices: ( ',
917 $
'D=diag(0,1,2,...) )', /
' 7=(D,I) 9=(large*D, small*I',
918 $
') 11=(large*I, small*D) 13=(large*D, large*I)', /
919 $
' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
920 $
' 14=(small*D, small*I)', /
' 15=(D, reversed D)' )
921 9994
FORMAT(
' Matrices Rotated by Random ', a,
' Matrices U, V:',
922 $ /
' 16=Transposed Jordan Blocks 19=geometric ',
923 $
'alpha, beta=0,1', /
' 17=arithm. alpha&beta ',
924 $
' 20=arithmetic alpha, beta=0,1', /
' 18=clustered ',
925 $
'alpha, beta=0,1 21=random alpha, beta=0,1',
926 $ /
' Large & Small Matrices:', /
' 22=(large, small) ',
927 $
'23=(small,large) 24=(small,small) 25=(large,large)',
928 $ /
' 26=random O(1) matrices.' )
930 9993
FORMAT( /
' Tests performed: ',
931 $ /
' 1 = max | ( b A - a B )''*l | / const.,',
932 $ /
' 2 = | |VR(i)| - 1 | / ulp,',
933 $ /
' 3 = max | ( b A - a B )*r | / const.',
934 $ /
' 4 = | |VL(i)| - 1 | / ulp,',
935 $ /
' 5 = 0 if W same no matter if r or l computed,',
936 $ /
' 6 = 0 if l same no matter if l computed,',
937 $ /
' 7 = 0 if r same no matter if r computed,', / 1x )
938 9992
FORMAT(
' Matrix order=', i5,
', type=', i2,
', seed=',
939 $ 4( i4,
',' ),
' result ', i2,
' is', 0p, f8.2 )
940 9991
FORMAT(
' Matrix order=', i5,
', type=', i2,
', seed=',
941 $ 4( i4,
',' ),
' result ', i2,
' is', 1p, e10.3 )
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
subroutine xlaenv(ispec, nvalue)
XLAENV
subroutine xerbla(srname, info)
subroutine sggev3(jobvl, jobvr, n, a, lda, b, ldb, alphar, alphai, beta, vl, ldvl, vr, ldvr, work, lwork, info)
SGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices (...
subroutine slacpy(uplo, m, n, a, lda, b, ldb)
SLACPY copies all or part of one two-dimensional array to another.
subroutine slarfg(n, alpha, x, incx, tau)
SLARFG generates an elementary reflector (Householder matrix).
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 sorm2r(side, trans, m, n, k, a, lda, tau, c, ldc, work, info)
SORM2R multiplies a general matrix by the orthogonal matrix from a QR factorization determined by sge...
subroutine sdrgev3(nsizes, nn, ntypes, dotype, iseed, thresh, nounit, a, lda, b, s, t, q, ldq, z, qe, ldqe, alphar, alphai, beta, alphr1, alphi1, beta1, work, lwork, result, info)
SDRGEV3
subroutine sget52(left, n, a, lda, b, ldb, e, lde, alphar, alphai, beta, work, result)
SGET52
subroutine slatm4(itype, n, nz1, nz2, isign, amagn, rcond, triang, idist, iseed, a, lda)
SLATM4