|
LAPACK 3.3.0
|
|
LAPACK 3.3.0
|
00001 SUBROUTINE SGET36( 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 00009 REAL RMAX 00010 * .. 00011 * .. Array Arguments .. 00012 INTEGER NINFO( 3 ) 00013 * .. 00014 * 00015 * Purpose 00016 * ======= 00017 * 00018 * SGET36 tests STREXC, a routine for moving blocks (either 1 by 1 or 00019 * 2 by 2) on the diagonal of a matrix in real Schur form. Thus, SLAEXC 00020 * computes an orthogonal matrix Q such that 00021 * 00022 * Q' * T1 * Q = T2 00023 * 00024 * and where one of the diagonal blocks of T1 (the one at row IFST) has 00025 * been moved to position ILST. 00026 * 00027 * The test code verifies that the residual Q'*T1*Q-T2 is small, that T2 00028 * is in Schur form, and that the final position of the IFST block is 00029 * ILST (within +-1). 00030 * 00031 * The test matrices are read from a file with logical unit number NIN. 00032 * 00033 * Arguments 00034 * ========== 00035 * 00036 * RMAX (output) REAL 00037 * Value of the largest test ratio. 00038 * 00039 * LMAX (output) INTEGER 00040 * Example number where largest test ratio achieved. 00041 * 00042 * NINFO (output) INTEGER array, dimension (3) 00043 * NINFO(J) is the number of examples where INFO=J. 00044 * 00045 * KNT (output) INTEGER 00046 * Total number of examples tested. 00047 * 00048 * NIN (input) INTEGER 00049 * Input logical unit number. 00050 * 00051 * ===================================================================== 00052 * 00053 * .. Parameters .. 00054 REAL ZERO, ONE 00055 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0 ) 00056 INTEGER LDT, LWORK 00057 PARAMETER ( LDT = 10, LWORK = 2*LDT*LDT ) 00058 * .. 00059 * .. Local Scalars .. 00060 INTEGER I, IFST, IFST1, IFST2, IFSTSV, ILST, ILST1, 00061 $ ILST2, ILSTSV, INFO1, INFO2, J, LOC, N 00062 REAL EPS, RES 00063 * .. 00064 * .. Local Arrays .. 00065 REAL Q( LDT, LDT ), RESULT( 2 ), T1( LDT, LDT ), 00066 $ T2( LDT, LDT ), TMP( LDT, LDT ), WORK( LWORK ) 00067 * .. 00068 * .. External Functions .. 00069 REAL SLAMCH 00070 EXTERNAL SLAMCH 00071 * .. 00072 * .. External Subroutines .. 00073 EXTERNAL SHST01, SLACPY, SLASET, STREXC 00074 * .. 00075 * .. Intrinsic Functions .. 00076 INTRINSIC ABS, SIGN 00077 * .. 00078 * .. Executable Statements .. 00079 * 00080 EPS = SLAMCH( 'P' ) 00081 RMAX = ZERO 00082 LMAX = 0 00083 KNT = 0 00084 NINFO( 1 ) = 0 00085 NINFO( 2 ) = 0 00086 NINFO( 3 ) = 0 00087 * 00088 * Read input data until N=0 00089 * 00090 10 CONTINUE 00091 READ( NIN, FMT = * )N, IFST, ILST 00092 IF( N.EQ.0 ) 00093 $ RETURN 00094 KNT = KNT + 1 00095 DO 20 I = 1, N 00096 READ( NIN, FMT = * )( TMP( I, J ), J = 1, N ) 00097 20 CONTINUE 00098 CALL SLACPY( 'F', N, N, TMP, LDT, T1, LDT ) 00099 CALL SLACPY( 'F', N, N, TMP, LDT, T2, LDT ) 00100 IFSTSV = IFST 00101 ILSTSV = ILST 00102 IFST1 = IFST 00103 ILST1 = ILST 00104 IFST2 = IFST 00105 ILST2 = ILST 00106 RES = ZERO 00107 * 00108 * Test without accumulating Q 00109 * 00110 CALL SLASET( 'Full', N, N, ZERO, ONE, Q, LDT ) 00111 CALL STREXC( 'N', N, T1, LDT, Q, LDT, IFST1, ILST1, WORK, INFO1 ) 00112 DO 40 I = 1, N 00113 DO 30 J = 1, N 00114 IF( I.EQ.J .AND. Q( I, J ).NE.ONE ) 00115 $ RES = RES + ONE / EPS 00116 IF( I.NE.J .AND. Q( I, J ).NE.ZERO ) 00117 $ RES = RES + ONE / EPS 00118 30 CONTINUE 00119 40 CONTINUE 00120 * 00121 * Test with accumulating Q 00122 * 00123 CALL SLASET( 'Full', N, N, ZERO, ONE, Q, LDT ) 00124 CALL STREXC( 'V', N, T2, LDT, Q, LDT, IFST2, ILST2, WORK, INFO2 ) 00125 * 00126 * Compare T1 with T2 00127 * 00128 DO 60 I = 1, N 00129 DO 50 J = 1, N 00130 IF( T1( I, J ).NE.T2( I, J ) ) 00131 $ RES = RES + ONE / EPS 00132 50 CONTINUE 00133 60 CONTINUE 00134 IF( IFST1.NE.IFST2 ) 00135 $ RES = RES + ONE / EPS 00136 IF( ILST1.NE.ILST2 ) 00137 $ RES = RES + ONE / EPS 00138 IF( INFO1.NE.INFO2 ) 00139 $ RES = RES + ONE / EPS 00140 * 00141 * Test for successful reordering of T2 00142 * 00143 IF( INFO2.NE.0 ) THEN 00144 NINFO( INFO2 ) = NINFO( INFO2 ) + 1 00145 ELSE 00146 IF( ABS( IFST2-IFSTSV ).GT.1 ) 00147 $ RES = RES + ONE / EPS 00148 IF( ABS( ILST2-ILSTSV ).GT.1 ) 00149 $ RES = RES + ONE / EPS 00150 END IF 00151 * 00152 * Test for small residual, and orthogonality of Q 00153 * 00154 CALL SHST01( N, 1, N, TMP, LDT, T2, LDT, Q, LDT, WORK, LWORK, 00155 $ RESULT ) 00156 RES = RES + RESULT( 1 ) + RESULT( 2 ) 00157 * 00158 * Test for T2 being in Schur form 00159 * 00160 LOC = 1 00161 70 CONTINUE 00162 IF( T2( LOC+1, LOC ).NE.ZERO ) THEN 00163 * 00164 * 2 by 2 block 00165 * 00166 IF( T2( LOC, LOC+1 ).EQ.ZERO .OR. T2( LOC, LOC ).NE. 00167 $ T2( LOC+1, LOC+1 ) .OR. SIGN( ONE, T2( LOC, LOC+1 ) ).EQ. 00168 $ SIGN( ONE, T2( LOC+1, LOC ) ) )RES = RES + ONE / EPS 00169 DO 80 I = LOC + 2, N 00170 IF( T2( I, LOC ).NE.ZERO ) 00171 $ RES = RES + ONE / RES 00172 IF( T2( I, LOC+1 ).NE.ZERO ) 00173 $ RES = RES + ONE / RES 00174 80 CONTINUE 00175 LOC = LOC + 2 00176 ELSE 00177 * 00178 * 1 by 1 block 00179 * 00180 DO 90 I = LOC + 1, N 00181 IF( T2( I, LOC ).NE.ZERO ) 00182 $ RES = RES + ONE / RES 00183 90 CONTINUE 00184 LOC = LOC + 1 00185 END IF 00186 IF( LOC.LT.N ) 00187 $ GO TO 70 00188 IF( RES.GT.RMAX ) THEN 00189 RMAX = RES 00190 LMAX = KNT 00191 END IF 00192 GO TO 10 00193 * 00194 * End of SGET36 00195 * 00196 END
1.7.3 
