LAPACK: TESTING/EIG/zget36.f Source File

Go to the documentation of this file.
00001       SUBROUTINE ZGET36( RMAX, LMAX, NINFO, KNT, NIN )
00002 *
00003 *  -- LAPACK test routine (version 3.1) --
00004 *     Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
00005 *     November 2006
00006 *
00007 *     .. Scalar Arguments ..
00008       INTEGER            KNT, LMAX, NIN, NINFO
00009       DOUBLE PRECISION   RMAX
00010 *     ..
00011 *
00012 *  Purpose
00013 *  =======
00014 *
00015 *  ZGET36 tests ZTREXC, a routine for reordering diagonal entries of a
00016 *  matrix in complex Schur form. Thus, ZLAEXC computes a unitary matrix
00017 *  Q such that
00018 *
00019 *     Q' * T1 * Q  = T2
00020 *
00021 *  and where one of the diagonal blocks of T1 (the one at row IFST) has
00022 *  been moved to position ILST.
00023 *
00024 *  The test code verifies that the residual Q'*T1*Q-T2 is small, that T2
00025 *  is in Schur form, and that the final position of the IFST block is
00026 *  ILST.
00027 *
00028 *  The test matrices are read from a file with logical unit number NIN.
00029 *
00030 *  Arguments
00031 *  ==========
00032 *
00033 *  RMAX    (output) DOUBLE PRECISION
00034 *          Value of the largest test ratio.
00035 *
00036 *  LMAX    (output) INTEGER
00037 *          Example number where largest test ratio achieved.
00038 *
00039 *  NINFO   (output) INTEGER
00040 *          Number of examples where INFO is nonzero.
00041 *
00042 *  KNT     (output) INTEGER
00043 *          Total number of examples tested.
00044 *
00045 *  NIN     (input) INTEGER
00046 *          Input logical unit number.
00047 *
00048 *  =====================================================================
00049 *
00050 *     .. Parameters ..
00051       DOUBLE PRECISION   ZERO, ONE
00052       PARAMETER          ( ZERO = 0.0D+0, ONE = 1.0D+0 )
00053       COMPLEX*16         CZERO, CONE
00054       PARAMETER          ( CZERO = ( 0.0D+0, 0.0D+0 ),
00055      $                   CONE = ( 1.0D+0, 0.0D+0 ) )
00056       INTEGER            LDT, LWORK
00057       PARAMETER          ( LDT = 10, LWORK = 2*LDT*LDT )
00058 *     ..
00059 *     .. Local Scalars ..
00060       INTEGER            I, IFST, ILST, INFO1, INFO2, J, N
00061       DOUBLE PRECISION   EPS, RES
00062       COMPLEX*16         CTEMP
00063 *     ..
00064 *     .. Local Arrays ..
00065       DOUBLE PRECISION   RESULT( 2 ), RWORK( LDT )
00066       COMPLEX*16         DIAG( LDT ), Q( LDT, LDT ), T1( LDT, LDT ),
00067      $                   T2( LDT, LDT ), TMP( LDT, LDT ), WORK( LWORK )
00068 *     ..
00069 *     .. External Functions ..
00070       DOUBLE PRECISION   DLAMCH
00071       EXTERNAL           DLAMCH
00072 *     ..
00073 *     .. External Subroutines ..
00074       EXTERNAL           ZCOPY, ZHST01, ZLACPY, ZLASET, ZTREXC
00075 *     ..
00076 *     .. Executable Statements ..
00077 *
00078       EPS = DLAMCH( 'P' )
00079       RMAX = ZERO
00080       LMAX = 0
00081       KNT = 0
00082       NINFO = 0
00083 *
00084 *     Read input data until N=0
00085 *
00086    10 CONTINUE
00087       READ( NIN, FMT = * )N, IFST, ILST
00088       IF( N.EQ.0 )
00089      $   RETURN
00090       KNT = KNT + 1
00091       DO 20 I = 1, N
00092          READ( NIN, FMT = * )( TMP( I, J ), J = 1, N )
00093    20 CONTINUE
00094       CALL ZLACPY( 'F', N, N, TMP, LDT, T1, LDT )
00095       CALL ZLACPY( 'F', N, N, TMP, LDT, T2, LDT )
00096       RES = ZERO
00097 *
00098 *     Test without accumulating Q
00099 *
00100       CALL ZLASET( 'Full', N, N, CZERO, CONE, Q, LDT )
00101       CALL ZTREXC( 'N', N, T1, LDT, Q, LDT, IFST, ILST, INFO1 )
00102       DO 40 I = 1, N
00103          DO 30 J = 1, N
00104             IF( I.EQ.J .AND. Q( I, J ).NE.CONE )
00105      $         RES = RES + ONE / EPS
00106             IF( I.NE.J .AND. Q( I, J ).NE.CZERO )
00107      $         RES = RES + ONE / EPS
00108    30    CONTINUE
00109    40 CONTINUE
00110 *
00111 *     Test with accumulating Q
00112 *
00113       CALL ZLASET( 'Full', N, N, CZERO, CONE, Q, LDT )
00114       CALL ZTREXC( 'V', N, T2, LDT, Q, LDT, IFST, ILST, INFO2 )
00115 *
00116 *     Compare T1 with T2
00117 *
00118       DO 60 I = 1, N
00119          DO 50 J = 1, N
00120             IF( T1( I, J ).NE.T2( I, J ) )
00121      $         RES = RES + ONE / EPS
00122    50    CONTINUE
00123    60 CONTINUE
00124       IF( INFO1.NE.0 .OR. INFO2.NE.0 )
00125      $   NINFO = NINFO + 1
00126       IF( INFO1.NE.INFO2 )
00127      $   RES = RES + ONE / EPS
00128 *
00129 *     Test for successful reordering of T2
00130 *
00131       CALL ZCOPY( N, TMP, LDT+1, DIAG, 1 )
00132       IF( IFST.LT.ILST ) THEN
00133          DO 70 I = IFST + 1, ILST
00134             CTEMP = DIAG( I )
00135             DIAG( I ) = DIAG( I-1 )
00136             DIAG( I-1 ) = CTEMP
00137    70    CONTINUE
00138       ELSE IF( IFST.GT.ILST ) THEN
00139          DO 80 I = IFST - 1, ILST, -1
00140             CTEMP = DIAG( I+1 )
00141             DIAG( I+1 ) = DIAG( I )
00142             DIAG( I ) = CTEMP
00143    80    CONTINUE
00144       END IF
00145       DO 90 I = 1, N
00146          IF( T2( I, I ).NE.DIAG( I ) )
00147      $      RES = RES + ONE / EPS
00148    90 CONTINUE
00149 *
00150 *     Test for small residual, and orthogonality of Q
00151 *
00152       CALL ZHST01( N, 1, N, TMP, LDT, T2, LDT, Q, LDT, WORK, LWORK,
00153      $             RWORK, RESULT )
00154       RES = RES + RESULT( 1 ) + RESULT( 2 )
00155 *
00156 *     Test for T2 being in Schur form
00157 *
00158       DO 110 J = 1, N - 1
00159          DO 100 I = J + 1, N
00160             IF( T2( I, J ).NE.CZERO )
00161      $         RES = RES + ONE / EPS
00162   100    CONTINUE
00163   110 CONTINUE
00164       IF( RES.GT.RMAX ) THEN
00165          RMAX = RES
00166          LMAX = KNT
00167       END IF
00168       GO TO 10
00169 *
00170 *     End of ZGET36
00171 *
00172       END