420 SUBROUTINE zdrvgg( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
421 $ thrshn, nounit, a, lda, b, s, t, s2, t2, q,
422 $ ldq, z, alpha1, beta1, alpha2, beta2, vl, vr,
423 $ work, lwork, rwork, result, info )
431 INTEGER info, lda, ldq, lwork, nounit, nsizes, ntypes
432 DOUBLE PRECISION thresh, thrshn
439 INTEGER iseed( 4 ), nn( * )
440 DOUBLE PRECISION result( * ), rwork( * )
441 COMPLEX*16 a( lda, * ), alpha1( * ), alpha2( * ),
442 $ b( lda, * ), beta1( * ), beta2( * ),
443 $ q( ldq, * ), s( lda, * ), s2( lda, * ),
444 $ t( lda, * ), t2( lda, * ), vl( ldq, * ),
445 $ vr( ldq, * ), work( * ), z( ldq, * )
448 DOUBLE PRECISION zero, one
449 parameter( zero = 0.0d+0, one = 1.0d+0 )
450 COMPLEX*16 czero, cone
451 parameter( czero = ( 0.0d+0, 0.0d+0 ),
452 $ cone = ( 1.0d+0, 0.0d+0 ) )
454 parameter( maxtyp = 26 )
458 INTEGER i1, iadd, iinfo, in, j, jc, jr, jsize, jtype,
459 $ lwkopt, mtypes, n, n1, nb, nbz, nerrs, nmats,
460 $ nmax, ns, ntest, ntestt
461 DOUBLE PRECISION safmax, safmin, temp1, temp2, ulp, ulpinv
465 LOGICAL lasign( maxtyp ), lbsign( maxtyp )
466 INTEGER ioldsd( 4 ), kadd( 6 ), kamagn( maxtyp ),
467 $ katype( maxtyp ), kazero( maxtyp ),
468 $ kbmagn( maxtyp ), kbtype( maxtyp ),
469 $ kbzero( maxtyp ), kclass( maxtyp ),
470 $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
471 DOUBLE PRECISION dumma( 4 ), rmagn( 0: 3 )
484 INTRINSIC abs, dble, dconjg, dimag, max, min, sign
487 DOUBLE PRECISION abs1
490 abs1( x ) = abs( dble( x ) ) + abs( dimag( x ) )
493 DATA kclass / 15*1, 10*2, 1*3 /
494 DATA kz1 / 0, 1, 2, 1, 3, 3 /
495 DATA kz2 / 0, 0, 1, 2, 1, 1 /
496 DATA kadd / 0, 0, 0, 0, 3, 2 /
497 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
498 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
499 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
500 $ 1, 1, -4, 2, -4, 8*8, 0 /
501 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
503 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
505 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
507 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
509 DATA ktrian / 16*0, 10*1 /
510 DATA lasign / 6*.false., .true., .false., 2*.true.,
511 $ 2*.false., 3*.true., .false., .true.,
512 $ 3*.false., 5*.true., .false. /
513 DATA lbsign / 7*.false., .true., 2*.false.,
514 $ 2*.true., 2*.false., .true., .false., .true.,
526 nmax = max( nmax, nn( j ) )
534 nb = max( 1,
ilaenv( 1,
'ZGEQRF',
' ', nmax, nmax, -1, -1 ),
535 $
ilaenv( 1,
'ZUNMQR',
'LC', nmax, nmax, nmax, -1 ),
536 $
ilaenv( 1,
'ZUNGQR',
' ', nmax, nmax, nmax, -1 ) )
537 nbz =
ilaenv( 1,
'ZHGEQZ',
'SII', nmax, 1, nmax, 0 )
538 ns =
ilaenv( 4,
'ZHGEQZ',
'SII', nmax, 1, nmax, 0 )
540 lwkopt = max( 2*nmax, nmax*( nb+1 ), ( 2*i1+nmax+1 )*( i1+1 ) )
544 IF( nsizes.LT.0 )
THEN
546 ELSE IF( badnn )
THEN
548 ELSE IF( ntypes.LT.0 )
THEN
550 ELSE IF( thresh.LT.zero )
THEN
552 ELSE IF( lda.LE.1 .OR. lda.LT.nmax )
THEN
554 ELSE IF( ldq.LE.1 .OR. ldq.LT.nmax )
THEN
556 ELSE IF( lwkopt.GT.lwork )
THEN
561 CALL
xerbla(
'ZDRVGG', -info )
567 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
570 ulp =
dlamch(
'Precision' )
571 safmin =
dlamch(
'Safe minimum' )
572 safmin = safmin / ulp
573 safmax = one / safmin
574 CALL
dlabad( safmin, safmax )
588 DO 160 jsize = 1, nsizes
591 rmagn( 2 ) = safmax*ulp / dble( n1 )
592 rmagn( 3 ) = safmin*ulpinv*n1
594 IF( nsizes.NE.1 )
THEN
595 mtypes = min( maxtyp, ntypes )
597 mtypes = min( maxtyp+1, ntypes )
600 DO 150 jtype = 1, mtypes
601 IF( .NOT.dotype( jtype ) )
609 ioldsd( j ) = iseed( j )
639 IF( mtypes.GT.maxtyp )
642 IF( kclass( jtype ).LT.3 )
THEN
646 IF( abs( katype( jtype ) ).EQ.3 )
THEN
647 in = 2*( ( n-1 ) / 2 ) + 1
649 $ CALL
zlaset(
'Full', n, n, czero, czero, a, lda )
653 CALL
zlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
654 $ kz2( kazero( jtype ) ), lasign( jtype ),
655 $ rmagn( kamagn( jtype ) ), ulp,
656 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
658 iadd = kadd( kazero( jtype ) )
659 IF( iadd.GT.0 .AND. iadd.LE.n )
660 $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
664 IF( abs( kbtype( jtype ) ).EQ.3 )
THEN
665 in = 2*( ( n-1 ) / 2 ) + 1
667 $ CALL
zlaset(
'Full', n, n, czero, czero, b, lda )
671 CALL
zlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
672 $ kz2( kbzero( jtype ) ), lbsign( jtype ),
673 $ rmagn( kbmagn( jtype ) ), one,
674 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
676 iadd = kadd( kbzero( jtype ) )
677 IF( iadd.NE.0 .AND. iadd.LE.n )
678 $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
680 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 )
THEN
689 q( jr, jc ) =
zlarnd( 3, iseed )
690 z( jr, jc ) =
zlarnd( 3, iseed )
692 CALL
zlarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
694 work( 2*n+jc ) = sign( one, dble( q( jc, jc ) ) )
696 CALL
zlarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
698 work( 3*n+jc ) = sign( one, dble( z( jc, jc ) ) )
701 ctemp =
zlarnd( 3, iseed )
704 work( 3*n ) = ctemp / abs( ctemp )
705 ctemp =
zlarnd( 3, iseed )
708 work( 4*n ) = ctemp / abs( ctemp )
714 a( jr, jc ) = work( 2*n+jr )*
715 $ dconjg( work( 3*n+jc ) )*
717 b( jr, jc ) = work( 2*n+jr )*
718 $ dconjg( work( 3*n+jc ) )*
722 CALL
zunm2r(
'L',
'N', n, n, n-1, q, ldq, work, a,
723 $ lda, work( 2*n+1 ), iinfo )
726 CALL
zunm2r(
'R',
'C', n, n, n-1, z, ldq, work( n+1 ),
727 $ a, lda, work( 2*n+1 ), iinfo )
730 CALL
zunm2r(
'L',
'N', n, n, n-1, q, ldq, work, b,
731 $ lda, work( 2*n+1 ), iinfo )
734 CALL
zunm2r(
'R',
'C', n, n, n-1, z, ldq, work( n+1 ),
735 $ b, lda, work( 2*n+1 ), iinfo )
745 a( jr, jc ) = rmagn( kamagn( jtype ) )*
747 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
755 IF( iinfo.NE.0 )
THEN
756 WRITE( nounit, fmt = 9999 )
'Generator', iinfo, n, jtype,
766 CALL
zlacpy(
' ', n, n, a, lda, s, lda )
767 CALL
zlacpy(
' ', n, n, b, lda, t, lda )
771 CALL
zgegs(
'V',
'V', n, s, lda, t, lda, alpha1, beta1, q,
772 $ ldq, z, ldq, work, lwork, rwork, iinfo )
773 IF( iinfo.NE.0 )
THEN
774 WRITE( nounit, fmt = 9999 )
'ZGEGS', iinfo, n, jtype,
784 CALL
zget51( 1, n, a, lda, s, lda, q, ldq, z, ldq, work,
785 $ rwork, result( 1 ) )
786 CALL
zget51( 1, n, b, lda, t, lda, q, ldq, z, ldq, work,
787 $ rwork, result( 2 ) )
788 CALL
zget51( 3, n, b, lda, t, lda, q, ldq, q, ldq, work,
789 $ rwork, result( 3 ) )
790 CALL
zget51( 3, n, b, lda, t, lda, z, ldq, z, ldq, work,
791 $ rwork, result( 4 ) )
798 temp2 = ( abs1( alpha1( j )-s( j, j ) ) /
799 $ max( safmin, abs1( alpha1( j ) ), abs1( s( j,
800 $ j ) ) )+abs1( beta1( j )-t( j, j ) ) /
801 $ max( safmin, abs1( beta1( j ) ), abs1( t( j,
803 temp1 = max( temp1, temp2 )
811 CALL
zlacpy(
' ', n, n, a, lda, s2, lda )
812 CALL
zlacpy(
' ', n, n, b, lda, t2, lda )
816 CALL
zgegv(
'V',
'V', n, s2, lda, t2, lda, alpha2, beta2,
817 $ vl, ldq, vr, ldq, work, lwork, rwork, iinfo )
818 IF( iinfo.NE.0 )
THEN
819 WRITE( nounit, fmt = 9999 )
'ZGEGV', iinfo, n, jtype,
829 CALL
zget52( .true., n, a, lda, b, lda, vl, ldq, alpha2,
830 $ beta2, work, rwork, dumma( 1 ) )
831 result( 6 ) = dumma( 1 )
832 IF( dumma( 2 ).GT.thrshn )
THEN
833 WRITE( nounit, fmt = 9998 )
'Left',
'ZGEGV', dumma( 2 ),
837 CALL
zget52( .false., n, a, lda, b, lda, vr, ldq, alpha2,
838 $ beta2, work, rwork, dumma( 1 ) )
839 result( 7 ) = dumma( 1 )
840 IF( dumma( 2 ).GT.thresh )
THEN
841 WRITE( nounit, fmt = 9998 )
'Right',
'ZGEGV', dumma( 2 ),
849 ntestt = ntestt + ntest
854 IF( result( jr ).GE.thresh )
THEN
859 IF( nerrs.EQ.0 )
THEN
860 WRITE( nounit, fmt = 9997 )
'ZGG'
864 WRITE( nounit, fmt = 9996 )
865 WRITE( nounit, fmt = 9995 )
866 WRITE( nounit, fmt = 9994 )
'Unitary'
870 WRITE( nounit, fmt = 9993 )
'unitary',
'*',
871 $
'conjugate transpose', (
'*', j = 1, 5 )
875 IF( result( jr ).LT.10000.0d0 )
THEN
876 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
879 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
890 CALL
alasvm(
'ZGG', nounit, nerrs, ntestt, 0 )
893 9999 format(
' ZDRVGG: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
894 $ i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
')' )
896 9998 format(
' ZDRVGG: ', a,
' Eigenvectors from ', a,
' incorrectly ',
897 $
'normalized.', /
' Bits of error=', 0p, g10.3,
',', 9x,
898 $
'N=', i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ), i5,
901 9997 format( / 1x, a3,
902 $
' -- Complex Generalized eigenvalue problem driver' )
904 9996 format(
' Matrix types (see ZDRVGG for details): ' )
906 9995 format(
' Special Matrices:', 23x,
907 $
'(J''=transposed Jordan block)',
908 $ /
' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
909 $
'6=(diag(J'',I), diag(I,J''))', /
' Diagonal Matrices: ( ',
910 $
'D=diag(0,1,2,...) )', /
' 7=(D,I) 9=(large*D, small*I',
911 $
') 11=(large*I, small*D) 13=(large*D, large*I)', /
912 $
' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
913 $
' 14=(small*D, small*I)', /
' 15=(D, reversed D)' )
914 9994 format(
' Matrices Rotated by Random ', a,
' Matrices U, V:',
915 $ /
' 16=Transposed Jordan Blocks 19=geometric ',
916 $
'alpha, beta=0,1', /
' 17=arithm. alpha&beta ',
917 $
' 20=arithmetic alpha, beta=0,1', /
' 18=clustered ',
918 $
'alpha, beta=0,1 21=random alpha, beta=0,1',
919 $ /
' Large & Small Matrices:', /
' 22=(large, small) ',
920 $
'23=(small,large) 24=(small,small) 25=(large,large)',
921 $ /
' 26=random O(1) matrices.' )
923 9993 format( /
' Tests performed: (S is Schur, T is triangular, ',
924 $
'Q and Z are ', a,
',', / 20x,
925 $
'l and r are the appropriate left and right', / 19x,
926 $
'eigenvectors, resp., a is alpha, b is beta, and', / 19x, a,
927 $
' means ', a,
'.)', /
' 1 = | A - Q S Z', a,
928 $
' | / ( |A| n ulp ) 2 = | B - Q T Z', a,
929 $
' | / ( |B| n ulp )', /
' 3 = | I - QQ', a,
930 $
' | / ( n ulp ) 4 = | I - ZZ', a,
931 $
' | / ( n ulp )', /
932 $
' 5 = difference between (alpha,beta) and diagonals of',
933 $
' (S,T)', /
' 6 = max | ( b A - a B )', a,
934 $
' l | / const. 7 = max | ( b A - a B ) r | / const.',
936 9992 format(
' Matrix order=', i5,
', type=', i2,
', seed=',
937 $ 4( i4,
',' ),
' result ', i3,
' is', 0p, f8.2 )
938 9991 format(
' Matrix order=', i5,
', type=', i2,
', seed=',
939 $ 4( i4,
',' ),
' result ', i3,
' is', 1p, d10.3 )