406 SUBROUTINE sdrgev( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
407 $ nounit, a, lda,
b, s, t, q, ldq, z, qe, ldqe,
408 $ alphar, alphai, beta, alphr1, alphi1, beta1,
409 $ work, lwork, result, info )
417 INTEGER info, lda, ldq, ldqe, lwork, nounit, nsizes,
423 INTEGER iseed( 4 ), nn( * )
424 REAL a( lda, * ), alphai( * ), alphi1( * ),
425 $ alphar( * ), alphr1( * ),
b( lda, * ),
426 $ beta( * ), beta1( * ), q( ldq, * ),
427 $ qe( ldqe, * ), result( * ), s( lda, * ),
428 $ t( lda, * ), work( * ), z( ldq, * )
435 parameter( zero = 0.0e+0, one = 1.0e+0 )
437 parameter( maxtyp = 26 )
441 INTEGER i, iadd, ierr, in,
j, jc, jr, jsize, jtype,
442 $ maxwrk, minwrk, mtypes, n, n1, nerrs, nmats,
444 REAL safmax, safmin, ulp, ulpinv
447 INTEGER iasign( maxtyp ), ibsign( maxtyp ),
448 $ ioldsd( 4 ), kadd( 6 ), kamagn( maxtyp ),
449 $ katype( maxtyp ), kazero( maxtyp ),
450 $ kbmagn( maxtyp ), kbtype( maxtyp ),
451 $ kbzero( maxtyp ), kclass( maxtyp ),
452 $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
465 INTRINSIC abs, max, min,
REAL, sign
468 DATA kclass / 15*1, 10*2, 1*3 /
469 DATA kz1 / 0, 1, 2, 1, 3, 3 /
470 DATA kz2 / 0, 0, 1, 2, 1, 1 /
471 DATA kadd / 0, 0, 0, 0, 3, 2 /
472 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
473 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
474 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
475 $ 1, 1, -4, 2, -4, 8*8, 0 /
476 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
478 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
480 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
482 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
484 DATA ktrian / 16*0, 10*1 /
485 DATA iasign / 6*0, 2, 0, 2*2, 2*0, 3*2, 0, 2, 3*0,
487 DATA ibsign / 7*0, 2, 2*0, 2*2, 2*0, 2, 0, 2, 9*0 /
498 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( ldq.LE.1 .OR. ldq.LT.nmax )
THEN
515 ELSE IF( ldqe.LE.1 .OR. ldqe.LT.nmax )
THEN
527 IF( info.EQ.0 .AND. lwork.GE.1 )
THEN
528 minwrk = max( 1, 8*nmax, nmax*( nmax+1 ) )
529 maxwrk = 7*nmax + nmax*
ilaenv( 1,
'SGEQRF',
' ', nmax, 1, nmax,
531 maxwrk = max( maxwrk, nmax*( nmax+1 ) )
535 IF( lwork.LT.minwrk )
539 CALL
xerbla(
'SDRGEV', -info )
545 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
548 safmin =
slamch(
'Safe minimum' )
550 safmin = safmin / ulp
551 safmax = one / safmin
552 CALL
slabad( safmin, safmax )
566 DO 220 jsize = 1, nsizes
569 rmagn( 2 ) = safmax*ulp /
REAL( n1 )
570 rmagn( 3 ) = safmin*ulpinv*n1
572 IF( nsizes.NE.1 )
THEN
573 mtypes = min( maxtyp, ntypes )
575 mtypes = min( maxtyp+1, ntypes )
578 DO 210 jtype = 1, mtypes
579 IF( .NOT.dotype( jtype ) )
586 ioldsd(
j ) = iseed(
j )
612 IF( mtypes.GT.maxtyp )
615 IF( kclass( jtype ).LT.3 )
THEN
619 IF( abs( katype( jtype ) ).EQ.3 )
THEN
620 in = 2*( ( n-1 ) / 2 ) + 1
622 $ CALL
slaset(
'Full', n, n, zero, zero, a, lda )
626 CALL
slatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
627 $ kz2( kazero( jtype ) ), iasign( jtype ),
628 $ rmagn( kamagn( jtype ) ), ulp,
629 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
631 iadd = kadd( kazero( jtype ) )
632 IF( iadd.GT.0 .AND. iadd.LE.n )
633 $ a( iadd, iadd ) = one
637 IF( abs( kbtype( jtype ) ).EQ.3 )
THEN
638 in = 2*( ( n-1 ) / 2 ) + 1
640 $ CALL
slaset(
'Full', n, n, zero, zero,
b, lda )
644 CALL
slatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
645 $ kz2( kbzero( jtype ) ), ibsign( jtype ),
646 $ rmagn( kbmagn( jtype ) ), one,
647 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
649 iadd = kadd( kbzero( jtype ) )
650 IF( iadd.NE.0 .AND. iadd.LE.n )
651 $
b( iadd, iadd ) = one
653 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 )
THEN
662 q( jr, jc ) =
slarnd( 3, iseed )
663 z( jr, jc ) =
slarnd( 3, iseed )
665 CALL
slarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
667 work( 2*n+jc ) = sign( one, q( jc, jc ) )
669 CALL
slarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
671 work( 3*n+jc ) = sign( one, z( jc, jc ) )
676 work( 3*n ) = sign( one,
slarnd( 2, iseed ) )
679 work( 4*n ) = sign( one,
slarnd( 2, iseed ) )
685 a( jr, jc ) = work( 2*n+jr )*work( 3*n+jc )*
687 b( jr, jc ) = work( 2*n+jr )*work( 3*n+jc )*
691 CALL
sorm2r(
'L',
'N', n, n, n-1, q, ldq, work, a,
692 $ lda, work( 2*n+1 ), ierr )
695 CALL
sorm2r(
'R',
'T', n, n, n-1, z, ldq, work( n+1 ),
696 $ a, lda, work( 2*n+1 ), ierr )
699 CALL
sorm2r(
'L',
'N', n, n, n-1, q, ldq, work,
b,
700 $ lda, work( 2*n+1 ), ierr )
703 CALL
sorm2r(
'R',
'T', n, n, n-1, z, ldq, work( n+1 ),
704 $
b, lda, work( 2*n+1 ), ierr )
714 a( jr, jc ) = rmagn( kamagn( jtype ) )*
716 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
725 WRITE( nounit, fmt = 9999 )
'Generator', ierr, n, jtype,
739 CALL
slacpy(
' ', n, n, a, lda, s, lda )
740 CALL
slacpy(
' ', n, n,
b, lda, t, lda )
741 CALL
sggev(
'V',
'V', n, s, lda, t, lda, alphar, alphai,
742 $ beta, q, ldq, z, ldq, work, lwork, ierr )
743 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
745 WRITE( nounit, fmt = 9999 )
'SGGEV1', ierr, n, jtype,
753 CALL
sget52( .true., n, a, lda,
b, lda, q, ldq, alphar,
754 $ alphai, beta, work, result( 1 ) )
755 IF( result( 2 ).GT.thresh )
THEN
756 WRITE( nounit, fmt = 9998 )
'Left',
'SGGEV1',
757 $ result( 2 ), n, jtype, ioldsd
762 CALL
sget52( .false., n, a, lda,
b, lda, z, ldq, alphar,
763 $ alphai, beta, work, result( 3 ) )
764 IF( result( 4 ).GT.thresh )
THEN
765 WRITE( nounit, fmt = 9998 )
'Right',
'SGGEV1',
766 $ result( 4 ), n, jtype, ioldsd
771 CALL
slacpy(
' ', n, n, a, lda, s, lda )
772 CALL
slacpy(
' ', n, n,
b, lda, t, lda )
773 CALL
sggev(
'N',
'N', n, s, lda, t, lda, alphr1, alphi1,
774 $ beta1, q, ldq, z, ldq, work, lwork, ierr )
775 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
777 WRITE( nounit, fmt = 9999 )
'SGGEV2', ierr, n, jtype,
784 IF( alphar(
j ).NE.alphr1(
j ) .OR. alphai(
j ).NE.
785 $ alphi1(
j ) .OR. beta(
j ).NE.beta1(
j ) )
786 $ result( 5 ) = ulpinv
792 CALL
slacpy(
' ', n, n, a, lda, s, lda )
793 CALL
slacpy(
' ', n, n,
b, lda, t, lda )
794 CALL
sggev(
'V',
'N', n, s, lda, t, lda, alphr1, alphi1,
795 $ beta1, qe, ldqe, z, ldq, work, lwork, ierr )
796 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
798 WRITE( nounit, fmt = 9999 )
'SGGEV3', ierr, n, jtype,
805 IF( alphar(
j ).NE.alphr1(
j ) .OR. alphai(
j ).NE.
806 $ alphi1(
j ) .OR. beta(
j ).NE.beta1(
j ) )
807 $ result( 6 ) = ulpinv
812 IF( q(
j, jc ).NE.qe(
j, jc ) )
813 $ result( 6 ) = ulpinv
820 CALL
slacpy(
' ', n, n, a, lda, s, lda )
821 CALL
slacpy(
' ', n, n,
b, lda, t, lda )
822 CALL
sggev(
'N',
'V', n, s, lda, t, lda, alphr1, alphi1,
823 $ beta1, q, ldq, qe, ldqe, work, lwork, ierr )
824 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
826 WRITE( nounit, fmt = 9999 )
'SGGEV4', ierr, n, jtype,
833 IF( alphar(
j ).NE.alphr1(
j ) .OR. alphai(
j ).NE.
834 $ alphi1(
j ) .OR. beta(
j ).NE.beta1(
j ) )
835 $ result( 7 ) = ulpinv
840 IF( z(
j, jc ).NE.qe(
j, jc ) )
841 $ result( 7 ) = ulpinv
854 IF( result( jr ).GE.thresh )
THEN
859 IF( nerrs.EQ.0 )
THEN
860 WRITE( nounit, fmt = 9997 )
'SGV'
864 WRITE( nounit, fmt = 9996 )
865 WRITE( nounit, fmt = 9995 )
866 WRITE( nounit, fmt = 9994 )
'Orthogonal'
870 WRITE( nounit, fmt = 9993 )
874 IF( result( jr ).LT.10000.0 )
THEN
875 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
878 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
889 CALL
alasvm(
'SGV', nounit, nerrs, ntestt, 0 )
895 9999
FORMAT(
' SDRGEV: ', a,
' returned INFO=', i6,
'.', / 3
x,
'N=',
896 $ i6,
', JTYPE=', i6,
', ISEED=(', 4( i4,
',' ), i5,
')' )
898 9998
FORMAT(
' SDRGEV: ', a,
' Eigenvectors from ', a,
' incorrectly ',
899 $
'normalized.', /
' Bits of error=', 0p, g10.3,
',', 3
x,
900 $
'N=', i4,
', JTYPE=', i3,
', ISEED=(', 4( i4,
',' ), i5,
903 9997
FORMAT( / 1
x, a3,
' -- Real Generalized eigenvalue problem driver'
906 9996
FORMAT(
' Matrix types (see SDRGEV for details): ' )
908 9995
FORMAT(
' Special Matrices:', 23
x,
909 $
'(J''=transposed Jordan block)',
910 $ /
' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
911 $
'6=(diag(J'',I), diag(I,J''))', /
' Diagonal Matrices: ( ',
912 $
'D=diag(0,1,2,...) )', /
' 7=(D,I) 9=(large*D, small*I',
913 $
') 11=(large*I, small*D) 13=(large*D, large*I)', /
914 $
' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
915 $
' 14=(small*D, small*I)', /
' 15=(D, reversed D)' )
916 9994
FORMAT(
' Matrices Rotated by Random ', a,
' Matrices U, V:',
917 $ /
' 16=Transposed Jordan Blocks 19=geometric ',
918 $
'alpha, beta=0,1', /
' 17=arithm. alpha&beta ',
919 $
' 20=arithmetic alpha, beta=0,1', /
' 18=clustered ',
920 $
'alpha, beta=0,1 21=random alpha, beta=0,1',
921 $ /
' Large & Small Matrices:', /
' 22=(large, small) ',
922 $
'23=(small,large) 24=(small,small) 25=(large,large)',
923 $ /
' 26=random O(1) matrices.' )
925 9993
FORMAT( /
' Tests performed: ',
926 $ /
' 1 = max | ( b A - a B )''*l | / const.,',
927 $ /
' 2 = | |VR(i)| - 1 | / ulp,',
928 $ /
' 3 = max | ( b A - a B )*r | / const.',
929 $ /
' 4 = | |VL(i)| - 1 | / ulp,',
930 $ /
' 5 = 0 if W same no matter if r or l computed,',
931 $ /
' 6 = 0 if l same no matter if l computed,',
932 $ /
' 7 = 0 if r same no matter if r computed,', / 1
x )
933 9992
FORMAT(
' Matrix order=', i5,
', type=', i2,
', seed=',
934 $ 4( i4,
',' ),
' result ', i2,
' is', 0p, f8.2 )
935 9991
FORMAT(
' Matrix order=', i5,
', type=', i2,
', seed=',
936 $ 4( i4,
',' ),
' result ', i2,
' is', 1p, e10.3 )
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 alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
REAL function slamch(CMACH)
SLAMCH
REAL function slarnd(IDIST, ISEED)
SLARND
subroutine sget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHAR, ALPHAI, BETA, WORK, RESULT)
SGET52
subroutine xerbla(SRNAME, INFO)
XERBLA
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real b(3) integer i
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
INTEGER function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
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 slarfg(N, ALPHA, X, INCX, TAU)
SLARFG generates an elementary reflector (Householder matrix).
set ue cd $ADTTMP cat<< EOF > tmp f Program LinearEquations Implicit none Real j
subroutine slatm4(ITYPE, N, NZ1, NZ2, ISIGN, AMAGN, RCOND, TRIANG, IDIST, ISEED, A, LDA)
SLATM4
subroutine slabad(SMALL, LARGE)
SLABAD
subroutine sdrgev(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)
SDRGEV
subroutine sggev(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, INFO)
SGGEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices ...