397 SUBROUTINE zdrgev( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
398 $ nounit, a, lda, b, s, t, q, ldq, z, qe, ldqe,
399 $ alpha, beta, alpha1, beta1, work, lwork, rwork,
408 INTEGER INFO, LDA, LDQ, LDQE, LWORK, NOUNIT, NSIZES,
410 DOUBLE PRECISION THRESH
414 INTEGER ISEED( 4 ), NN( * )
415 DOUBLE PRECISION RESULT( * ), RWORK( * )
416 COMPLEX*16 A( lda, * ), ALPHA( * ), ALPHA1( * ),
417 $ b( lda, * ), beta( * ), beta1( * ),
418 $ q( ldq, * ), qe( ldqe, * ), s( lda, * ),
419 $ t( lda, * ), work( * ), z( ldq, * )
425 DOUBLE PRECISION ZERO, ONE
426 parameter ( zero = 0.0d+0, one = 1.0d+0 )
427 COMPLEX*16 CZERO, CONE
428 parameter ( czero = ( 0.0d+0, 0.0d+0 ),
429 $ cone = ( 1.0d+0, 0.0d+0 ) )
431 parameter ( maxtyp = 26 )
435 INTEGER I, IADD, IERR, IN, J, JC, JR, JSIZE, JTYPE,
436 $ maxwrk, minwrk, mtypes, n, n1, nb, nerrs,
437 $ nmats, nmax, ntestt
438 DOUBLE PRECISION SAFMAX, SAFMIN, ULP, ULPINV
442 LOGICAL LASIGN( maxtyp ), LBSIGN( maxtyp )
443 INTEGER IOLDSD( 4 ), KADD( 6 ), KAMAGN( maxtyp ),
444 $ katype( maxtyp ), kazero( maxtyp ),
445 $ kbmagn( maxtyp ), kbtype( maxtyp ),
446 $ kbzero( maxtyp ), kclass( maxtyp ),
447 $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
448 DOUBLE PRECISION RMAGN( 0: 3 )
452 DOUBLE PRECISION DLAMCH
454 EXTERNAL ilaenv, dlamch, zlarnd
461 INTRINSIC abs, dble, dconjg, max, min, sign
464 DATA kclass / 15*1, 10*2, 1*3 /
465 DATA kz1 / 0, 1, 2, 1, 3, 3 /
466 DATA kz2 / 0, 0, 1, 2, 1, 1 /
467 DATA kadd / 0, 0, 0, 0, 3, 2 /
468 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
469 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
470 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
471 $ 1, 1, -4, 2, -4, 8*8, 0 /
472 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
474 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
476 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
478 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
480 DATA ktrian / 16*0, 10*1 /
481 DATA lasign / 6*.false., .true., .false., 2*.true.,
482 $ 2*.false., 3*.true., .false., .true.,
483 $ 3*.false., 5*.true., .false. /
484 DATA lbsign / 7*.false., .true., 2*.false.,
485 $ 2*.true., 2*.false., .true., .false., .true.,
497 nmax = max( nmax, nn( j ) )
502 IF( nsizes.LT.0 )
THEN
504 ELSE IF( badnn )
THEN
506 ELSE IF( ntypes.LT.0 )
THEN
508 ELSE IF( thresh.LT.zero )
THEN
510 ELSE IF( lda.LE.1 .OR. lda.LT.nmax )
THEN
512 ELSE IF( ldq.LE.1 .OR. ldq.LT.nmax )
THEN
514 ELSE IF( ldqe.LE.1 .OR. ldqe.LT.nmax )
THEN
526 IF( info.EQ.0 .AND. lwork.GE.1 )
THEN
527 minwrk = nmax*( nmax+1 )
528 nb = max( 1, ilaenv( 1,
'ZGEQRF',
' ', nmax, nmax, -1, -1 ),
529 $ ilaenv( 1,
'ZUNMQR',
'LC', nmax, nmax, nmax, -1 ),
530 $ ilaenv( 1,
'ZUNGQR',
' ', nmax, nmax, nmax, -1 ) )
531 maxwrk = max( 2*nmax, nmax*( nb+1 ), nmax*( nmax+1 ) )
535 IF( lwork.LT.minwrk )
539 CALL xerbla(
'ZDRGEV', -info )
545 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
548 ulp = dlamch(
'Precision' )
549 safmin = dlamch(
'Safe minimum' )
550 safmin = safmin / ulp
551 safmax = one / safmin
552 CALL dlabad( safmin, safmax )
566 DO 220 jsize = 1, nsizes
569 rmagn( 2 ) = safmax*ulp / dble( 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 )
610 IF( mtypes.GT.maxtyp )
613 IF( kclass( jtype ).LT.3 )
THEN
617 IF( abs( katype( jtype ) ).EQ.3 )
THEN
618 in = 2*( ( n-1 ) / 2 ) + 1
620 $
CALL zlaset(
'Full', n, n, czero, czero, a, lda )
624 CALL zlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
625 $ kz2( kazero( jtype ) ), lasign( jtype ),
626 $ rmagn( kamagn( jtype ) ), ulp,
627 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
629 iadd = kadd( kazero( jtype ) )
630 IF( iadd.GT.0 .AND. iadd.LE.n )
631 $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
635 IF( abs( kbtype( jtype ) ).EQ.3 )
THEN
636 in = 2*( ( n-1 ) / 2 ) + 1
638 $
CALL zlaset(
'Full', n, n, czero, czero, b, lda )
642 CALL zlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
643 $ kz2( kbzero( jtype ) ), lbsign( jtype ),
644 $ rmagn( kbmagn( jtype ) ), one,
645 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
647 iadd = kadd( kbzero( jtype ) )
648 IF( iadd.NE.0 .AND. iadd.LE.n )
649 $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
651 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 )
THEN
660 q( jr, jc ) = zlarnd( 3, iseed )
661 z( jr, jc ) = zlarnd( 3, iseed )
663 CALL zlarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
665 work( 2*n+jc ) = sign( one, dble( q( jc, jc ) ) )
667 CALL zlarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
669 work( 3*n+jc ) = sign( one, dble( z( jc, jc ) ) )
672 ctemp = zlarnd( 3, iseed )
675 work( 3*n ) = ctemp / abs( ctemp )
676 ctemp = zlarnd( 3, iseed )
679 work( 4*n ) = ctemp / abs( ctemp )
685 a( jr, jc ) = work( 2*n+jr )*
686 $ dconjg( work( 3*n+jc ) )*
688 b( jr, jc ) = work( 2*n+jr )*
689 $ dconjg( work( 3*n+jc ) )*
693 CALL zunm2r(
'L',
'N', n, n, n-1, q, ldq, work, a,
694 $ lda, work( 2*n+1 ), ierr )
697 CALL zunm2r(
'R',
'C', n, n, n-1, z, ldq, work( n+1 ),
698 $ a, lda, work( 2*n+1 ), ierr )
701 CALL zunm2r(
'L',
'N', n, n, n-1, q, ldq, work, b,
702 $ lda, work( 2*n+1 ), ierr )
705 CALL zunm2r(
'R',
'C', n, n, n-1, z, ldq, work( n+1 ),
706 $ b, lda, work( 2*n+1 ), ierr )
716 a( jr, jc ) = rmagn( kamagn( jtype ) )*
718 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
727 WRITE( nounit, fmt = 9999 )
'Generator', ierr, n, jtype,
741 CALL zlacpy(
' ', n, n, a, lda, s, lda )
742 CALL zlacpy(
' ', n, n, b, lda, t, lda )
743 CALL zggev(
'V',
'V', n, s, lda, t, lda, alpha, beta, q,
744 $ ldq, z, ldq, work, lwork, rwork, ierr )
745 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
747 WRITE( nounit, fmt = 9999 )
'ZGGEV1', ierr, n, jtype,
755 CALL zget52( .true., n, a, lda, b, lda, q, ldq, alpha, beta,
756 $ work, rwork, result( 1 ) )
757 IF( result( 2 ).GT.thresh )
THEN
758 WRITE( nounit, fmt = 9998 )
'Left',
'ZGGEV1',
759 $ result( 2 ), n, jtype, ioldsd
764 CALL zget52( .false., n, a, lda, b, lda, z, ldq, alpha,
765 $ beta, work, rwork, result( 3 ) )
766 IF( result( 4 ).GT.thresh )
THEN
767 WRITE( nounit, fmt = 9998 )
'Right',
'ZGGEV1',
768 $ result( 4 ), n, jtype, ioldsd
773 CALL zlacpy(
' ', n, n, a, lda, s, lda )
774 CALL zlacpy(
' ', n, n, b, lda, t, lda )
775 CALL zggev(
'N',
'N', n, s, lda, t, lda, alpha1, beta1, q,
776 $ ldq, z, ldq, work, lwork, rwork, ierr )
777 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
779 WRITE( nounit, fmt = 9999 )
'ZGGEV2', ierr, n, jtype,
786 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
787 $ beta1( j ) )result( 5 ) = ulpinv
793 CALL zlacpy(
' ', n, n, a, lda, s, lda )
794 CALL zlacpy(
' ', n, n, b, lda, t, lda )
795 CALL zggev(
'V',
'N', n, s, lda, t, lda, alpha1, beta1, qe,
796 $ ldqe, z, ldq, work, lwork, rwork, ierr )
797 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
799 WRITE( nounit, fmt = 9999 )
'ZGGEV3', ierr, n, jtype,
806 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
807 $ beta1( j ) )result( 6 ) = ulpinv
812 IF( q( j, jc ).NE.qe( j, jc ) )
813 $ result( 6 ) = ulpinv
820 CALL zlacpy(
' ', n, n, a, lda, s, lda )
821 CALL zlacpy(
' ', n, n, b, lda, t, lda )
822 CALL zggev(
'N',
'V', n, s, lda, t, lda, alpha1, beta1, q,
823 $ ldq, qe, ldqe, work, lwork, rwork, ierr )
824 IF( ierr.NE.0 .AND. ierr.NE.n+1 )
THEN
826 WRITE( nounit, fmt = 9999 )
'ZGGEV4', ierr, n, jtype,
833 IF( alpha( j ).NE.alpha1( j ) .OR. beta( j ).NE.
834 $ beta1( j ) )result( 7 ) = ulpinv
839 IF( z( j, jc ).NE.qe( j, jc ) )
840 $ result( 7 ) = ulpinv
853 IF( result( jr ).GE.thresh )
THEN
858 IF( nerrs.EQ.0 )
THEN
859 WRITE( nounit, fmt = 9997 )
'ZGV'
863 WRITE( nounit, fmt = 9996 )
864 WRITE( nounit, fmt = 9995 )
865 WRITE( nounit, fmt = 9994 )
'Orthogonal'
869 WRITE( nounit, fmt = 9993 )
873 IF( result( jr ).LT.10000.0d0 )
THEN
874 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
877 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
888 CALL alasvm(
'ZGV', nounit, nerrs, ntestt, 0 )
894 9999
FORMAT(
' ZDRGEV: ', a,
' returned INFO=', i6,
'.', / 3x,
'N=',
895 $ i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
897 9998
FORMAT(
' ZDRGEV: ', a,
' Eigenvectors from ', a,
' incorrectly ',
898 $
'normalized.', /
' Bits of error=', 0p, g10.3,
',', 3x,
899 $
'N=', i4,
', JTYPE=', i3,
', ISEED=(', 3( i4,
',' ), i5,
902 9997
FORMAT( / 1x, a3,
' -- Complex Generalized eigenvalue problem ',
905 9996
FORMAT(
' Matrix types (see ZDRGEV for details): ' )
907 9995
FORMAT(
' Special Matrices:', 23x,
908 $
'(J''=transposed Jordan block)',
909 $ /
' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
910 $
'6=(diag(J'',I), diag(I,J''))', /
' Diagonal Matrices: ( ',
911 $
'D=diag(0,1,2,...) )', /
' 7=(D,I) 9=(large*D, small*I',
912 $
') 11=(large*I, small*D) 13=(large*D, large*I)', /
913 $
' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
914 $
' 14=(small*D, small*I)', /
' 15=(D, reversed D)' )
915 9994
FORMAT(
' Matrices Rotated by Random ', a,
' Matrices U, V:',
916 $ /
' 16=Transposed Jordan Blocks 19=geometric ',
917 $
'alpha, beta=0,1', /
' 17=arithm. alpha&beta ',
918 $
' 20=arithmetic alpha, beta=0,1', /
' 18=clustered ',
919 $
'alpha, beta=0,1 21=random alpha, beta=0,1',
920 $ /
' Large & Small Matrices:', /
' 22=(large, small) ',
921 $
'23=(small,large) 24=(small,small) 25=(large,large)',
922 $ /
' 26=random O(1) matrices.' )
924 9993
FORMAT( /
' Tests performed: ',
925 $ /
' 1 = max | ( b A - a B )''*l | / const.,',
926 $ /
' 2 = | |VR(i)| - 1 | / ulp,',
927 $ /
' 3 = max | ( b A - a B )*r | / const.',
928 $ /
' 4 = | |VL(i)| - 1 | / ulp,',
929 $ /
' 5 = 0 if W same no matter if r or l computed,',
930 $ /
' 6 = 0 if l same no matter if l computed,',
931 $ /
' 7 = 0 if r same no matter if r computed,', / 1x )
932 9992
FORMAT(
' Matrix order=', i5,
', type=', i2,
', seed=',
933 $ 4( i4,
',' ),
' result ', i2,
' is', 0p, f8.2 )
934 9991
FORMAT(
' Matrix order=', i5,
', type=', i2,
', seed=',
935 $ 4( i4,
',' ),
' result ', i2,
' is', 1p, d10.3 )
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 alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
subroutine zget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHA, BETA, WORK, RWORK, RESULT)
ZGET52
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
subroutine xerbla(SRNAME, INFO)
XERBLA
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 zlatm4(ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND, TRIANG, IDIST, ISEED, A, LDA)
ZLATM4
subroutine dlabad(SMALL, LARGE)
DLABAD
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 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 zlarfg(N, ALPHA, X, INCX, TAU)
ZLARFG generates an elementary reflector (Householder matrix).