380 SUBROUTINE zdrges3( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
381 $ nounit, a, lda, b, s, t, q, ldq, z, alpha,
382 $ beta, work, lwork, rwork, result, bwork,
391 INTEGER INFO, LDA, LDQ, LWORK, NOUNIT, NSIZES, NTYPES
392 DOUBLE PRECISION THRESH
395 LOGICAL BWORK( * ), DOTYPE( * )
396 INTEGER ISEED( 4 ), NN( * )
397 DOUBLE PRECISION RESULT( 13 ), RWORK( * )
398 COMPLEX*16 A( lda, * ), ALPHA( * ), B( lda, * ),
399 $ beta( * ), q( ldq, * ), s( lda, * ),
400 $ t( lda, * ), work( * ), z( ldq, * )
406 DOUBLE PRECISION ZERO, ONE
407 parameter ( zero = 0.0d+0, one = 1.0d+0 )
408 COMPLEX*16 CZERO, CONE
409 parameter ( czero = ( 0.0d+0, 0.0d+0 ),
410 $ cone = ( 1.0d+0, 0.0d+0 ) )
412 parameter ( maxtyp = 26 )
415 LOGICAL BADNN, ILABAD
417 INTEGER I, IADD, IINFO, IN, ISORT, J, JC, JR, JSIZE,
418 $ jtype, knteig, maxwrk, minwrk, mtypes, n, n1,
419 $ nb, nerrs, nmats, nmax, ntest, ntestt, rsub,
421 DOUBLE PRECISION SAFMAX, SAFMIN, TEMP1, TEMP2, ULP, ULPINV
425 LOGICAL LASIGN( maxtyp ), LBSIGN( maxtyp )
426 INTEGER IOLDSD( 4 ), KADD( 6 ), KAMAGN( maxtyp ),
427 $ katype( maxtyp ), kazero( maxtyp ),
428 $ kbmagn( maxtyp ), kbtype( maxtyp ),
429 $ kbzero( maxtyp ), kclass( maxtyp ),
430 $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
431 DOUBLE PRECISION RMAGN( 0: 3 )
436 DOUBLE PRECISION DLAMCH
438 EXTERNAL zlctes, ilaenv, dlamch, zlarnd
445 INTRINSIC abs, dble, dconjg, dimag, max, min, sign
448 DOUBLE PRECISION ABS1
451 abs1( x ) = abs( dble( x ) ) + abs( dimag( x ) )
454 DATA kclass / 15*1, 10*2, 1*3 /
455 DATA kz1 / 0, 1, 2, 1, 3, 3 /
456 DATA kz2 / 0, 0, 1, 2, 1, 1 /
457 DATA kadd / 0, 0, 0, 0, 3, 2 /
458 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
459 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
460 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
461 $ 1, 1, -4, 2, -4, 8*8, 0 /
462 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
464 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
466 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
468 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
470 DATA ktrian / 16*0, 10*1 /
471 DATA lasign / 6*.false., .true., .false., 2*.true.,
472 $ 2*.false., 3*.true., .false., .true.,
473 $ 3*.false., 5*.true., .false. /
474 DATA lbsign / 7*.false., .true., 2*.false.,
475 $ 2*.true., 2*.false., .true., .false., .true.,
487 nmax = max( nmax, nn( j ) )
492 IF( nsizes.LT.0 )
THEN
494 ELSE IF( badnn )
THEN
496 ELSE IF( ntypes.LT.0 )
THEN
498 ELSE IF( thresh.LT.zero )
THEN
500 ELSE IF( lda.LE.1 .OR. lda.LT.nmax )
THEN
502 ELSE IF( ldq.LE.1 .OR. ldq.LT.nmax )
THEN
514 IF( info.EQ.0 .AND. lwork.GE.1 )
THEN
516 nb = max( 1, ilaenv( 1,
'ZGEQRF',
' ', nmax, nmax, -1, -1 ),
517 $ ilaenv( 1,
'ZUNMQR',
'LC', nmax, nmax, nmax, -1 ),
518 $ ilaenv( 1,
'ZUNGQR',
' ', nmax, nmax, nmax, -1 ) )
519 maxwrk = max( nmax+nmax*nb, 3*nmax*nmax )
523 IF( lwork.LT.minwrk )
527 CALL xerbla(
'ZDRGES3', -info )
533 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
536 ulp = dlamch(
'Precision' )
537 safmin = dlamch(
'Safe minimum' )
538 safmin = safmin / ulp
539 safmax = one / safmin
540 CALL dlabad( safmin, safmax )
554 DO 190 jsize = 1, nsizes
557 rmagn( 2 ) = safmax*ulp / dble( n1 )
558 rmagn( 3 ) = safmin*ulpinv*dble( n1 )
560 IF( nsizes.NE.1 )
THEN
561 mtypes = min( maxtyp, ntypes )
563 mtypes = min( maxtyp+1, ntypes )
568 DO 180 jtype = 1, mtypes
569 IF( .NOT.dotype( jtype ) )
577 ioldsd( j ) = iseed( j )
607 IF( mtypes.GT.maxtyp )
610 IF( kclass( jtype ).LT.3 )
THEN
614 IF( abs( katype( jtype ) ).EQ.3 )
THEN
615 in = 2*( ( n-1 ) / 2 ) + 1
617 $
CALL zlaset(
'Full', n, n, czero, czero, a, lda )
621 CALL zlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
622 $ kz2( kazero( jtype ) ), lasign( jtype ),
623 $ rmagn( kamagn( jtype ) ), ulp,
624 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
626 iadd = kadd( kazero( jtype ) )
627 IF( iadd.GT.0 .AND. iadd.LE.n )
628 $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
632 IF( abs( kbtype( jtype ) ).EQ.3 )
THEN
633 in = 2*( ( n-1 ) / 2 ) + 1
635 $
CALL zlaset(
'Full', n, n, czero, czero, b, lda )
639 CALL zlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
640 $ kz2( kbzero( jtype ) ), lbsign( jtype ),
641 $ rmagn( kbmagn( jtype ) ), one,
642 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
644 iadd = kadd( kbzero( jtype ) )
645 IF( iadd.NE.0 .AND. iadd.LE.n )
646 $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
648 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 )
THEN
657 q( jr, jc ) = zlarnd( 3, iseed )
658 z( jr, jc ) = zlarnd( 3, iseed )
660 CALL zlarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
662 work( 2*n+jc ) = sign( one, dble( q( jc, jc ) ) )
664 CALL zlarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
666 work( 3*n+jc ) = sign( one, dble( z( jc, jc ) ) )
669 ctemp = zlarnd( 3, iseed )
672 work( 3*n ) = ctemp / abs( ctemp )
673 ctemp = zlarnd( 3, iseed )
676 work( 4*n ) = ctemp / abs( ctemp )
682 a( jr, jc ) = work( 2*n+jr )*
683 $ dconjg( work( 3*n+jc ) )*
685 b( jr, jc ) = work( 2*n+jr )*
686 $ dconjg( work( 3*n+jc ) )*
690 CALL zunm2r(
'L',
'N', n, n, n-1, q, ldq, work, a,
691 $ lda, work( 2*n+1 ), iinfo )
694 CALL zunm2r(
'R',
'C', n, n, n-1, z, ldq, work( n+1 ),
695 $ a, lda, work( 2*n+1 ), iinfo )
698 CALL zunm2r(
'L',
'N', n, n, n-1, q, ldq, work, b,
699 $ lda, work( 2*n+1 ), iinfo )
702 CALL zunm2r(
'R',
'C', n, n, n-1, z, ldq, work( n+1 ),
703 $ b, lda, work( 2*n+1 ), iinfo )
713 a( jr, jc ) = rmagn( kamagn( jtype ) )*
715 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
723 IF( iinfo.NE.0 )
THEN
724 WRITE( nounit, fmt = 9999 )
'Generator', iinfo, n, jtype,
739 IF( isort.EQ.0 )
THEN
749 CALL zlacpy(
'Full', n, n, a, lda, s, lda )
750 CALL zlacpy(
'Full', n, n, b, lda, t, lda )
751 ntest = 1 + rsub + isort
752 result( 1+rsub+isort ) = ulpinv
753 CALL zgges3(
'V',
'V', sort, zlctes, n, s, lda, t, lda,
754 $ sdim, alpha, beta, q, ldq, z, ldq, work,
755 $ lwork, rwork, bwork, iinfo )
756 IF( iinfo.NE.0 .AND. iinfo.NE.n+2 )
THEN
757 result( 1+rsub+isort ) = ulpinv
758 WRITE( nounit, fmt = 9999 )
'ZGGES3', iinfo, n, jtype,
768 IF( isort.EQ.0 )
THEN
769 CALL zget51( 1, n, a, lda, s, lda, q, ldq, z, ldq,
770 $ work, rwork, result( 1 ) )
771 CALL zget51( 1, n, b, lda, t, lda, q, ldq, z, ldq,
772 $ work, rwork, result( 2 ) )
774 CALL zget54( n, a, lda, b, lda, s, lda, t, lda, q,
775 $ ldq, z, ldq, work, result( 2+rsub ) )
778 CALL zget51( 3, n, b, lda, t, lda, q, ldq, q, ldq, work,
779 $ rwork, result( 3+rsub ) )
780 CALL zget51( 3, n, b, lda, t, lda, z, ldq, z, ldq, work,
781 $ rwork, result( 4+rsub ) )
792 temp2 = ( abs1( alpha( j )-s( j, j ) ) /
793 $ max( safmin, abs1( alpha( j ) ), abs1( s( j,
794 $ j ) ) )+abs1( beta( j )-t( j, j ) ) /
795 $ max( safmin, abs1( beta( j ) ), abs1( t( j,
799 IF( s( j+1, j ).NE.zero )
THEN
801 result( 5+rsub ) = ulpinv
805 IF( s( j, j-1 ).NE.zero )
THEN
807 result( 5+rsub ) = ulpinv
810 temp1 = max( temp1, temp2 )
812 WRITE( nounit, fmt = 9998 )j, n, jtype, ioldsd
815 result( 6+rsub ) = temp1
817 IF( isort.GE.1 )
THEN
825 IF( zlctes( alpha( i ), beta( i ) ) )
826 $ knteig = knteig + 1
829 $ result( 13 ) = ulpinv
838 ntestt = ntestt + ntest
843 IF( result( jr ).GE.thresh )
THEN
848 IF( nerrs.EQ.0 )
THEN
849 WRITE( nounit, fmt = 9997 )
'ZGS'
853 WRITE( nounit, fmt = 9996 )
854 WRITE( nounit, fmt = 9995 )
855 WRITE( nounit, fmt = 9994 )
'Unitary'
859 WRITE( nounit, fmt = 9993 )
'unitary',
'''',
860 $
'transpose', (
'''', j = 1, 8 )
864 IF( result( jr ).LT.10000.0d0 )
THEN
865 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
868 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
879 CALL alasvm(
'ZGS', nounit, nerrs, ntestt, 0 )
885 9999
FORMAT(
' ZDRGES3: ', a,
' returned INFO=', i6,
'.', / 9x,
'N=',
886 $ i6,
', JTYPE=', i6,
', ISEED=(', 4( i4,
',' ), i5,
')' )
888 9998
FORMAT(
' ZDRGES3: S not in Schur form at eigenvalue ', i6,
'.',
889 $ / 9x,
'N=', i6,
', JTYPE=', i6,
', ISEED=(', 3( i5,
',' ),
892 9997
FORMAT( / 1x, a3,
' -- Complex Generalized Schur from problem ',
895 9996
FORMAT(
' Matrix types (see ZDRGES3 for details): ' )
897 9995
FORMAT(
' Special Matrices:', 23x,
898 $
'(J''=transposed Jordan block)',
899 $ /
' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
900 $
'6=(diag(J'',I), diag(I,J''))', /
' Diagonal Matrices: ( ',
901 $
'D=diag(0,1,2,...) )', /
' 7=(D,I) 9=(large*D, small*I',
902 $
') 11=(large*I, small*D) 13=(large*D, large*I)', /
903 $
' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
904 $
' 14=(small*D, small*I)', /
' 15=(D, reversed D)' )
905 9994
FORMAT(
' Matrices Rotated by Random ', a,
' Matrices U, V:',
906 $ /
' 16=Transposed Jordan Blocks 19=geometric ',
907 $
'alpha, beta=0,1', /
' 17=arithm. alpha&beta ',
908 $
' 20=arithmetic alpha, beta=0,1', /
' 18=clustered ',
909 $
'alpha, beta=0,1 21=random alpha, beta=0,1',
910 $ /
' Large & Small Matrices:', /
' 22=(large, small) ',
911 $
'23=(small,large) 24=(small,small) 25=(large,large)',
912 $ /
' 26=random O(1) matrices.' )
914 9993
FORMAT( /
' Tests performed: (S is Schur, T is triangular, ',
915 $
'Q and Z are ', a,
',', / 19x,
916 $
'l and r are the appropriate left and right', / 19x,
917 $
'eigenvectors, resp., a is alpha, b is beta, and', / 19x, a,
918 $
' means ', a,
'.)', /
' Without ordering: ',
919 $ /
' 1 = | A - Q S Z', a,
920 $
' | / ( |A| n ulp ) 2 = | B - Q T Z', a,
921 $
' | / ( |B| n ulp )', /
' 3 = | I - QQ', a,
922 $
' | / ( n ulp ) 4 = | I - ZZ', a,
923 $
' | / ( n ulp )', /
' 5 = A is in Schur form S',
924 $ /
' 6 = difference between (alpha,beta)',
925 $
' and diagonals of (S,T)', /
' With ordering: ',
926 $ /
' 7 = | (A,B) - Q (S,T) Z', a,
' | / ( |(A,B)| n ulp )',
927 $ /
' 8 = | I - QQ', a,
928 $
' | / ( n ulp ) 9 = | I - ZZ', a,
929 $
' | / ( n ulp )', /
' 10 = A is in Schur form S',
930 $ /
' 11 = difference between (alpha,beta) and diagonals',
931 $
' of (S,T)', /
' 12 = SDIM is the correct number of ',
932 $
'selected eigenvalues', / )
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, 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 zget54(N, A, LDA, B, LDB, S, LDS, T, LDT, U, LDU, V, LDV, WORK, RESULT)
ZGET54
subroutine zgges3(JOBVSL, JOBVSR, SORT, SELCTG, N, A, LDA, B, LDB, SDIM, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK, RWORK, BWORK, INFO)
ZGGES3 computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE...
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 zlarfg(N, ALPHA, X, INCX, TAU)
ZLARFG generates an elementary reflector (Householder matrix).
subroutine zget51(ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK, RWORK, RESULT)
ZGET51
subroutine zdrges3(NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, NOUNIT, A, LDA, B, S, T, Q, LDQ, Z, ALPHA, BETA, WORK, LWORK, RWORK, RESULT, BWORK, INFO)
ZDRGES3