395 SUBROUTINE zdrgev( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
396 $ NOUNIT, A, LDA, B, S, T, Q, LDQ, Z, QE, LDQE,
397 $ ALPHA, BETA, ALPHA1, BETA1, WORK, LWORK, RWORK,
405 INTEGER INFO, LDA, LDQ, LDQE, LWORK, NOUNIT, NSIZES,
407 DOUBLE PRECISION THRESH
411 INTEGER ISEED( 4 ), NN( * )
412 DOUBLE PRECISION RESULT( * ), RWORK( * )
413 COMPLEX*16 A( LDA, * ), ALPHA( * ), ALPHA1( * ),
414 $ b( lda, * ), beta( * ), beta1( * ),
415 $ q( ldq, * ), qe( ldqe, * ), s( lda, * ),
416 $ t( lda, * ), work( * ), z( ldq, * )
422 DOUBLE PRECISION ZERO, ONE
423 PARAMETER ( ZERO = 0.0d+0, one = 1.0d+0 )
424 COMPLEX*16 CZERO, CONE
425 parameter( czero = ( 0.0d+0, 0.0d+0 ),
426 $ cone = ( 1.0d+0, 0.0d+0 ) )
428 parameter( maxtyp = 26 )
432 INTEGER I, IADD, IERR, IN, J, JC, JR, JSIZE, JTYPE,
433 $ MAXWRK, MINWRK, MTYPES, N, N1, NB, NERRS,
434 $ nmats, nmax, ntestt
435 DOUBLE PRECISION SAFMAX, SAFMIN, ULP, ULPINV
439 LOGICAL LASIGN( MAXTYP ), LBSIGN( MAXTYP )
440 INTEGER IOLDSD( 4 ), KADD( 6 ), KAMAGN( MAXTYP ),
441 $ KATYPE( MAXTYP ), KAZERO( MAXTYP ),
442 $ kbmagn( maxtyp ), kbtype( maxtyp ),
443 $ kbzero( maxtyp ), kclass( maxtyp ),
444 $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
445 DOUBLE PRECISION RMAGN( 0: 3 )
449 DOUBLE PRECISION DLAMCH
451 EXTERNAL ilaenv, dlamch, zlarnd
458 INTRINSIC abs, dble, dconjg, max, min, sign
461 DATA kclass / 15*1, 10*2, 1*3 /
462 DATA kz1 / 0, 1, 2, 1, 3, 3 /
463 DATA kz2 / 0, 0, 1, 2, 1, 1 /
464 DATA kadd / 0, 0, 0, 0, 3, 2 /
465 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
466 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
467 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
468 $ 1, 1, -4, 2, -4, 8*8, 0 /
469 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
471 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
473 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
475 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
477 DATA ktrian / 16*0, 10*1 /
478 DATA lasign / 6*.false., .true., .false., 2*.true.,
479 $ 2*.false., 3*.true., .false., .true.,
480 $ 3*.false., 5*.true., .false. /
481 DATA lbsign / 7*.false., .true., 2*.false.,
482 $ 2*.true., 2*.false., .true., .false., .true.,
494 nmax = max( nmax, nn( j ) )
499 IF( nsizes.LT.0 )
THEN
501 ELSE IF( badnn )
THEN
503 ELSE IF( ntypes.LT.0 )
THEN
505 ELSE IF( thresh.LT.zero )
THEN
507 ELSE IF( lda.LE.1 .OR. lda.LT.nmax )
THEN
509 ELSE IF( ldq.LE.1 .OR. ldq.LT.nmax )
THEN
511 ELSE IF( ldqe.LE.1 .OR. ldqe.LT.nmax )
THEN
523 IF( info.EQ.0 .AND. lwork.GE.1 )
THEN
524 minwrk = nmax*( nmax+1 )
525 nb = max( 1, ilaenv( 1,
'ZGEQRF',
' ', nmax, nmax, -1, -1 ),
526 $ ilaenv( 1,
'ZUNMQR',
'LC', nmax, nmax, nmax, -1 ),
527 $ ilaenv( 1,
'ZUNGQR',
' ', nmax, nmax, nmax, -1 ) )
528 maxwrk = max( 2*nmax, nmax*( nb+1 ), nmax*( nmax+1 ) )
532 IF( lwork.LT.minwrk )
536 CALL xerbla(
'ZDRGEV', -info )
542 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
545 ulp = dlamch(
'Precision' )
546 safmin = dlamch(
'Safe minimum' )
547 safmin = safmin / ulp
548 safmax = one / safmin
562 DO 220 jsize = 1, nsizes
565 rmagn( 2 ) = safmax*ulp / dble( 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 )
606 IF( mtypes.GT.maxtyp )
609 IF( kclass( jtype ).LT.3 )
THEN
613 IF( abs( katype( jtype ) ).EQ.3 )
THEN
614 in = 2*( ( n-1 ) / 2 ) + 1
616 $
CALL zlaset(
'Full', n, n, czero, czero, a, lda )
620 CALL zlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
621 $ kz2( kazero( jtype ) ), lasign( jtype ),
622 $ rmagn( kamagn( jtype ) ), ulp,
623 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
625 iadd = kadd( kazero( jtype ) )
626 IF( iadd.GT.0 .AND. iadd.LE.n )
627 $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
631 IF( abs( kbtype( jtype ) ).EQ.3 )
THEN
632 in = 2*( ( n-1 ) / 2 ) + 1
634 $
CALL zlaset(
'Full', n, n, czero, czero, b, lda )
638 CALL zlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
639 $ kz2( kbzero( jtype ) ), lbsign( jtype ),
640 $ rmagn( kbmagn( jtype ) ), one,
641 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
643 iadd = kadd( kbzero( jtype ) )
644 IF( iadd.NE.0 .AND. iadd.LE.n )
645 $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
647 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 )
THEN
656 q( jr, jc ) = zlarnd( 3, iseed )
657 z( jr, jc ) = zlarnd( 3, iseed )
659 CALL zlarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
661 work( 2*n+jc ) = sign( one, dble( q( jc, jc ) ) )
663 CALL zlarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
665 work( 3*n+jc ) = sign( one, dble( z( jc, jc ) ) )
668 ctemp = zlarnd( 3, iseed )
671 work( 3*n ) = ctemp / abs( ctemp )
672 ctemp = zlarnd( 3, iseed )
675 work( 4*n ) = ctemp / abs( ctemp )
681 a( jr, jc ) = work( 2*n+jr )*
682 $ dconjg( work( 3*n+jc ) )*
684 b( jr, jc ) = work( 2*n+jr )*
685 $ dconjg( work( 3*n+jc ) )*
689 CALL zunm2r(
'L',
'N', n, n, n-1, q, ldq, work, a,
690 $ lda, work( 2*n+1 ), ierr )
693 CALL zunm2r(
'R',
'C', n, n, n-1, z, ldq, work( n+1 ),
694 $ a, lda, work( 2*n+1 ), ierr )
697 CALL zunm2r(
'L',
'N', n, n, n-1, q, ldq, work, b,
698 $ lda, work( 2*n+1 ), ierr )
701 CALL zunm2r(
'R',
'C', n, n, n-1, z, ldq, work( n+1 ),
702 $ b, lda, work( 2*n+1 ), ierr )
712 a( jr, jc ) = rmagn( kamagn( jtype ) )*
714 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
723 WRITE( nounit, fmt = 9999 )
'Generator', ierr, n, jtype,
737 CALL zlacpy(
' ', n, n, a, lda, s, lda )
738 CALL zlacpy(
' ', n, n, b, lda, t, lda )
739 CALL zggev(
'V',
'V', n, s, lda, t, lda, alpha, beta, q,
740 $ ldq, z, ldq, work, lwork, rwork, ierr )
741 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
743 WRITE( nounit, fmt = 9999 )
'ZGGEV1', ierr, n, jtype,
751 CALL zget52( .true., n, a, lda, b, lda, q, ldq, alpha, beta,
752 $ work, rwork, result( 1 ) )
753 IF( result( 2 ).GT.thresh )
THEN
754 WRITE( nounit, fmt = 9998 )
'Left',
'ZGGEV1',
755 $ result( 2 ), n, jtype, ioldsd
760 CALL zget52( .false., n, a, lda, b, lda, z, ldq, alpha,
761 $ beta, work, rwork, result( 3 ) )
762 IF( result( 4 ).GT.thresh )
THEN
763 WRITE( nounit, fmt = 9998 )
'Right',
'ZGGEV1',
764 $ result( 4 ), n, jtype, ioldsd
769 CALL zlacpy(
' ', n, n, a, lda, s, lda )
770 CALL zlacpy(
' ', n, n, b, lda, t, lda )
771 CALL zggev(
'N',
'N', n, s, lda, t, lda, alpha1, beta1, q,
772 $ ldq, z, ldq, work, lwork, rwork, ierr )
773 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
775 WRITE( nounit, fmt = 9999 )
'ZGGEV2', ierr, n, jtype,
782 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
783 $ beta1( j ) )result( 5 ) = ulpinv
789 CALL zlacpy(
' ', n, n, a, lda, s, lda )
790 CALL zlacpy(
' ', n, n, b, lda, t, lda )
791 CALL zggev(
'V',
'N', n, s, lda, t, lda, alpha1, beta1, qe,
792 $ ldqe, z, ldq, work, lwork, rwork, ierr )
793 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
795 WRITE( nounit, fmt = 9999 )
'ZGGEV3', ierr, n, jtype,
802 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
803 $ beta1( j ) )result( 6 ) = ulpinv
808 IF( q( j, jc ).NE.qe( j, jc ) )
809 $ result( 6 ) = ulpinv
816 CALL zlacpy(
' ', n, n, a, lda, s, lda )
817 CALL zlacpy(
' ', n, n, b, lda, t, lda )
818 CALL zggev(
'N',
'V', n, s, lda, t, lda, alpha1, beta1, q,
819 $ ldq, qe, ldqe, work, lwork, rwork, ierr )
820 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
822 WRITE( nounit, fmt = 9999 )
'ZGGEV4', ierr, n, jtype,
829 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
830 $ beta1( j ) )result( 7 ) = ulpinv
835 IF( z( j, jc ).NE.qe( j, jc ) )
836 $ result( 7 ) = ulpinv
849 IF( result( jr ).GE.thresh )
THEN
854 IF( nerrs.EQ.0 )
THEN
855 WRITE( nounit, fmt = 9997 )
'ZGV'
859 WRITE( nounit, fmt = 9996 )
860 WRITE( nounit, fmt = 9995 )
861 WRITE( nounit, fmt = 9994 )
'Orthogonal'
865 WRITE( nounit, fmt = 9993 )
869 IF( result( jr ).LT.10000.0d0 )
THEN
870 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
873 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
884 CALL alasvm(
'ZGV', nounit, nerrs, ntestt, 0 )
890 9999
FORMAT(
' ZDRGEV: ', a,
' returned INFO=', i6,
'.', / 3x,
'N=',
891 $ i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
893 9998
FORMAT(
' ZDRGEV: ', a,
' Eigenvectors from ', a,
' incorrectly ',
894 $
'normalized.', /
' Bits of error=', 0p, g10.3,
',', 3x,
895 $
'N=', i4,
', JTYPE=', i3,
', ISEED=(', 3( i4,
',' ), i5,
898 9997
FORMAT( / 1x, a3,
' -- Complex Generalized eigenvalue problem ',
901 9996
FORMAT(
' Matrix types (see ZDRGEV for details): ' )
903 9995
FORMAT(
' Special Matrices:', 23x,
904 $
'(J''=transposed Jordan block)',
905 $ /
' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
906 $
'6=(diag(J'',I), diag(I,J''))', /
' Diagonal Matrices: ( ',
907 $
'D=diag(0,1,2,...) )', /
' 7=(D,I) 9=(large*D, small*I',
908 $
') 11=(large*I, small*D) 13=(large*D, large*I)', /
909 $
' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
910 $
' 14=(small*D, small*I)', /
' 15=(D, reversed D)' )
911 9994
FORMAT(
' Matrices Rotated by Random ', a,
' Matrices U, V:',
912 $ /
' 16=Transposed Jordan Blocks 19=geometric ',
913 $
'alpha, beta=0,1', /
' 17=arithm. alpha&beta ',
914 $
' 20=arithmetic alpha, beta=0,1', /
' 18=clustered ',
915 $
'alpha, beta=0,1 21=random alpha, beta=0,1',
916 $ /
' Large & Small Matrices:', /
' 22=(large, small) ',
917 $
'23=(small,large) 24=(small,small) 25=(large,large)',
918 $ /
' 26=random O(1) matrices.' )
920 9993
FORMAT( /
' Tests performed: ',
921 $ /
' 1 = max | ( b A - a B )''*l | / const.,',
922 $ /
' 2 = | |VR(i)| - 1 | / ulp,',
923 $ /
' 3 = max | ( b A - a B )*r | / const.',
924 $ /
' 4 = | |VL(i)| - 1 | / ulp,',
925 $ /
' 5 = 0 if W same no matter if r or l computed,',
926 $ /
' 6 = 0 if l same no matter if l computed,',
927 $ /
' 7 = 0 if r same no matter if r computed,', / 1x )
928 9992
FORMAT(
' Matrix order=', i5,
', type=', i2,
', seed=',
929 $ 4( i4,
',' ),
' result ', i2,
' is', 0p, f8.2 )
930 9991
FORMAT(
' Matrix order=', i5,
', type=', i2,
', seed=',
931 $ 4( i4,
',' ),
' result ', i2,
' is', 1p, d10.3 )
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
subroutine xerbla(srname, info)
subroutine zggev(jobvl, jobvr, n, a, lda, b, ldb, alpha, beta, vl, ldvl, vr, ldvr, work, lwork, rwork, info)
ZGGEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
subroutine zlarfg(n, alpha, x, incx, tau)
ZLARFG generates an elementary reflector (Householder matrix).
subroutine zlaset(uplo, m, n, alpha, beta, a, lda)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
subroutine zunm2r(side, trans, m, n, k, a, lda, tau, c, ldc, work, info)
ZUNM2R multiplies a general matrix by the unitary matrix from a QR factorization determined by cgeqrf...
subroutine zdrgev(nsizes, nn, ntypes, dotype, iseed, thresh, nounit, a, lda, b, s, t, q, ldq, z, qe, ldqe, alpha, beta, alpha1, beta1, work, lwork, rwork, result, info)
ZDRGEV
subroutine zget52(left, n, a, lda, b, ldb, e, lde, alpha, beta, work, rwork, result)
ZGET52
subroutine zlatm4(itype, n, nz1, nz2, rsign, amagn, rcond, triang, idist, iseed, a, lda)
ZLATM4