LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE CGET35( 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 REAL RMAX 00010 * .. 00011 * 00012 * Purpose 00013 * ======= 00014 * 00015 * CGET35 tests CTRSYL, a routine for solving the Sylvester matrix 00016 * equation 00017 * 00018 * op(A)*X + ISGN*X*op(B) = scale*C, 00019 * 00020 * A and B are assumed to be in Schur canonical form, op() represents an 00021 * optional transpose, and ISGN can be -1 or +1. Scale is an output 00022 * less than or equal to 1, chosen to avoid overflow in X. 00023 * 00024 * The test code verifies that the following residual is order 1: 00025 * 00026 * norm(op(A)*X + ISGN*X*op(B) - scale*C) / 00027 * (EPS*max(norm(A),norm(B))*norm(X)) 00028 * 00029 * Arguments 00030 * ========== 00031 * 00032 * RMAX (output) REAL 00033 * Value of the largest test ratio. 00034 * 00035 * LMAX (output) INTEGER 00036 * Example number where largest test ratio achieved. 00037 * 00038 * NINFO (output) INTEGER 00039 * Number of examples where INFO is nonzero. 00040 * 00041 * KNT (output) INTEGER 00042 * Total number of examples tested. 00043 * 00044 * NIN (input) INTEGER 00045 * Input logical unit number. 00046 * 00047 * ===================================================================== 00048 * 00049 * .. Parameters .. 00050 INTEGER LDT 00051 PARAMETER ( LDT = 10 ) 00052 REAL ZERO, ONE, TWO 00053 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0 ) 00054 REAL LARGE 00055 PARAMETER ( LARGE = 1.0E6 ) 00056 COMPLEX CONE 00057 PARAMETER ( CONE = 1.0E0 ) 00058 * .. 00059 * .. Local Scalars .. 00060 CHARACTER TRANA, TRANB 00061 INTEGER I, IMLA, IMLAD, IMLB, IMLC, INFO, ISGN, ITRANA, 00062 $ ITRANB, J, M, N 00063 REAL BIGNUM, EPS, RES, RES1, SCALE, SMLNUM, TNRM, 00064 $ XNRM 00065 COMPLEX RMUL 00066 * .. 00067 * .. Local Arrays .. 00068 REAL DUM( 1 ), VM1( 3 ), VM2( 3 ) 00069 COMPLEX A( LDT, LDT ), ATMP( LDT, LDT ), B( LDT, LDT ), 00070 $ BTMP( LDT, LDT ), C( LDT, LDT ), 00071 $ CSAV( LDT, LDT ), CTMP( LDT, LDT ) 00072 * .. 00073 * .. External Functions .. 00074 REAL CLANGE, SLAMCH 00075 EXTERNAL CLANGE, SLAMCH 00076 * .. 00077 * .. External Subroutines .. 00078 EXTERNAL CGEMM, CTRSYL 00079 * .. 00080 * .. Intrinsic Functions .. 00081 INTRINSIC ABS, MAX, REAL, SQRT 00082 * .. 00083 * .. Executable Statements .. 00084 * 00085 * Get machine parameters 00086 * 00087 EPS = SLAMCH( 'P' ) 00088 SMLNUM = SLAMCH( 'S' ) / EPS 00089 BIGNUM = ONE / SMLNUM 00090 CALL SLABAD( SMLNUM, BIGNUM ) 00091 * 00092 * Set up test case parameters 00093 * 00094 VM1( 1 ) = SQRT( SMLNUM ) 00095 VM1( 2 ) = ONE 00096 VM1( 3 ) = LARGE 00097 VM2( 1 ) = ONE 00098 VM2( 2 ) = ONE + TWO*EPS 00099 VM2( 3 ) = TWO 00100 * 00101 KNT = 0 00102 NINFO = 0 00103 LMAX = 0 00104 RMAX = ZERO 00105 * 00106 * Begin test loop 00107 * 00108 10 CONTINUE 00109 READ( NIN, FMT = * )M, N 00110 IF( N.EQ.0 ) 00111 $ RETURN 00112 DO 20 I = 1, M 00113 READ( NIN, FMT = * )( ATMP( I, J ), J = 1, M ) 00114 20 CONTINUE 00115 DO 30 I = 1, N 00116 READ( NIN, FMT = * )( BTMP( I, J ), J = 1, N ) 00117 30 CONTINUE 00118 DO 40 I = 1, M 00119 READ( NIN, FMT = * )( CTMP( I, J ), J = 1, N ) 00120 40 CONTINUE 00121 DO 170 IMLA = 1, 3 00122 DO 160 IMLAD = 1, 3 00123 DO 150 IMLB = 1, 3 00124 DO 140 IMLC = 1, 3 00125 DO 130 ITRANA = 1, 2 00126 DO 120 ITRANB = 1, 2 00127 DO 110 ISGN = -1, 1, 2 00128 IF( ITRANA.EQ.1 ) 00129 $ TRANA = 'N' 00130 IF( ITRANA.EQ.2 ) 00131 $ TRANA = 'C' 00132 IF( ITRANB.EQ.1 ) 00133 $ TRANB = 'N' 00134 IF( ITRANB.EQ.2 ) 00135 $ TRANB = 'C' 00136 TNRM = ZERO 00137 DO 60 I = 1, M 00138 DO 50 J = 1, M 00139 A( I, J ) = ATMP( I, J )*VM1( IMLA ) 00140 TNRM = MAX( TNRM, ABS( A( I, J ) ) ) 00141 50 CONTINUE 00142 A( I, I ) = A( I, I )*VM2( IMLAD ) 00143 TNRM = MAX( TNRM, ABS( A( I, I ) ) ) 00144 60 CONTINUE 00145 DO 80 I = 1, N 00146 DO 70 J = 1, N 00147 B( I, J ) = BTMP( I, J )*VM1( IMLB ) 00148 TNRM = MAX( TNRM, ABS( B( I, J ) ) ) 00149 70 CONTINUE 00150 80 CONTINUE 00151 IF( TNRM.EQ.ZERO ) 00152 $ TNRM = ONE 00153 DO 100 I = 1, M 00154 DO 90 J = 1, N 00155 C( I, J ) = CTMP( I, J )*VM1( IMLC ) 00156 CSAV( I, J ) = C( I, J ) 00157 90 CONTINUE 00158 100 CONTINUE 00159 KNT = KNT + 1 00160 CALL CTRSYL( TRANA, TRANB, ISGN, M, N, A, 00161 $ LDT, B, LDT, C, LDT, SCALE, 00162 $ INFO ) 00163 IF( INFO.NE.0 ) 00164 $ NINFO = NINFO + 1 00165 XNRM = CLANGE( 'M', M, N, C, LDT, DUM ) 00166 RMUL = CONE 00167 IF( XNRM.GT.ONE .AND. TNRM.GT.ONE ) THEN 00168 IF( XNRM.GT.BIGNUM / TNRM ) THEN 00169 RMUL = MAX( XNRM, TNRM ) 00170 RMUL = CONE / RMUL 00171 END IF 00172 END IF 00173 CALL CGEMM( TRANA, 'N', M, N, M, RMUL, A, 00174 $ LDT, C, LDT, -SCALE*RMUL, CSAV, 00175 $ LDT ) 00176 CALL CGEMM( 'N', TRANB, M, N, N, 00177 $ REAL( ISGN )*RMUL, C, LDT, B, 00178 $ LDT, CONE, CSAV, LDT ) 00179 RES1 = CLANGE( 'M', M, N, CSAV, LDT, DUM ) 00180 RES = RES1 / MAX( SMLNUM, SMLNUM*XNRM, 00181 $ ( ( ABS( RMUL )*TNRM )*EPS )*XNRM ) 00182 IF( RES.GT.RMAX ) THEN 00183 LMAX = KNT 00184 RMAX = RES 00185 END IF 00186 110 CONTINUE 00187 120 CONTINUE 00188 130 CONTINUE 00189 140 CONTINUE 00190 150 CONTINUE 00191 160 CONTINUE 00192 170 CONTINUE 00193 GO TO 10 00194 * 00195 * End of CGET35 00196 * 00197 END