LAPACK 3.3.0
|
00001 SUBROUTINE CDRVRF3( NOUT, NN, NVAL, THRESH, A, LDA, ARF, B1, B2, 00002 + S_WORK_CLANGE, C_WORK_CGEQRF, TAU ) 00003 * 00004 * -- LAPACK test routine (version 3.2.0) -- 00005 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00006 * November 2008 00007 * 00008 * .. Scalar Arguments .. 00009 INTEGER LDA, NN, NOUT 00010 REAL THRESH 00011 * .. 00012 * .. Array Arguments .. 00013 INTEGER NVAL( NN ) 00014 REAL S_WORK_CLANGE( * ) 00015 COMPLEX A( LDA, * ), ARF( * ), B1( LDA, * ), 00016 + B2( LDA, * ) 00017 COMPLEX C_WORK_CGEQRF( * ), TAU( * ) 00018 * .. 00019 * 00020 * Purpose 00021 * ======= 00022 * 00023 * CDRVRF3 tests the LAPACK RFP routines: 00024 * CTFSM 00025 * 00026 * Arguments 00027 * ========= 00028 * 00029 * NOUT (input) INTEGER 00030 * The unit number for output. 00031 * 00032 * NN (input) INTEGER 00033 * The number of values of N contained in the vector NVAL. 00034 * 00035 * NVAL (input) INTEGER array, dimension (NN) 00036 * The values of the matrix dimension N. 00037 * 00038 * THRESH (input) DOUBLE PRECISION 00039 * The threshold value for the test ratios. A result is 00040 * included in the output file if RESULT >= THRESH. To have 00041 * every test ratio printed, use THRESH = 0. 00042 * 00043 * A (workspace) COMPLEX*16 array, dimension (LDA,NMAX) 00044 * 00045 * LDA (input) INTEGER 00046 * The leading dimension of the array A. LDA >= max(1,NMAX). 00047 * 00048 * ARF (workspace) COMPLEX array, dimension ((NMAX*(NMAX+1))/2). 00049 * 00050 * B1 (workspace) COMPLEX array, dimension (LDA,NMAX) 00051 * 00052 * B2 (workspace) COMPLEX array, dimension (LDA,NMAX) 00053 * 00054 * S_WORK_CLANGE (workspace) REAL array, dimension (NMAX) 00055 * 00056 * C_WORK_CGEQRF (workspace) COMPLEX array, dimension (NMAX) 00057 * 00058 * TAU (workspace) COMPLEX array, dimension (NMAX) 00059 * 00060 * ===================================================================== 00061 * .. 00062 * .. Parameters .. 00063 COMPLEX ZERO, ONE 00064 PARAMETER ( ZERO = ( 0.0E+0, 0.0E+0 ) , 00065 + ONE = ( 1.0E+0, 0.0E+0 ) ) 00066 INTEGER NTESTS 00067 PARAMETER ( NTESTS = 1 ) 00068 * .. 00069 * .. Local Scalars .. 00070 CHARACTER UPLO, CFORM, DIAG, TRANS, SIDE 00071 INTEGER I, IFORM, IIM, IIN, INFO, IUPLO, J, M, N, NA, 00072 + NFAIL, NRUN, ISIDE, IDIAG, IALPHA, ITRANS 00073 COMPLEX ALPHA 00074 REAL EPS 00075 * .. 00076 * .. Local Arrays .. 00077 CHARACTER UPLOS( 2 ), FORMS( 2 ), TRANSS( 2 ), 00078 + DIAGS( 2 ), SIDES( 2 ) 00079 INTEGER ISEED( 4 ), ISEEDY( 4 ) 00080 REAL RESULT( NTESTS ) 00081 * .. 00082 * .. External Functions .. 00083 REAL SLAMCH, CLANGE 00084 COMPLEX CLARND 00085 EXTERNAL SLAMCH, CLARND, CLANGE 00086 * .. 00087 * .. External Subroutines .. 00088 EXTERNAL CTRTTF, CGEQRF, CGEQLF, CTFSM, CTRSM 00089 * .. 00090 * .. Intrinsic Functions .. 00091 INTRINSIC MAX, SQRT 00092 * .. 00093 * .. Scalars in Common .. 00094 CHARACTER*32 SRNAMT 00095 * .. 00096 * .. Common blocks .. 00097 COMMON / SRNAMC / SRNAMT 00098 * .. 00099 * .. Data statements .. 00100 DATA ISEEDY / 1988, 1989, 1990, 1991 / 00101 DATA UPLOS / 'U', 'L' / 00102 DATA FORMS / 'N', 'C' / 00103 DATA SIDES / 'L', 'R' / 00104 DATA TRANSS / 'N', 'C' / 00105 DATA DIAGS / 'N', 'U' / 00106 * .. 00107 * .. Executable Statements .. 00108 * 00109 * Initialize constants and the random number seed. 00110 * 00111 NRUN = 0 00112 NFAIL = 0 00113 INFO = 0 00114 DO 10 I = 1, 4 00115 ISEED( I ) = ISEEDY( I ) 00116 10 CONTINUE 00117 EPS = SLAMCH( 'Precision' ) 00118 * 00119 DO 170 IIM = 1, NN 00120 * 00121 M = NVAL( IIM ) 00122 * 00123 DO 160 IIN = 1, NN 00124 * 00125 N = NVAL( IIN ) 00126 * 00127 DO 150 IFORM = 1, 2 00128 * 00129 CFORM = FORMS( IFORM ) 00130 * 00131 DO 140 IUPLO = 1, 2 00132 * 00133 UPLO = UPLOS( IUPLO ) 00134 * 00135 DO 130 ISIDE = 1, 2 00136 * 00137 SIDE = SIDES( ISIDE ) 00138 * 00139 DO 120 ITRANS = 1, 2 00140 * 00141 TRANS = TRANSS( ITRANS ) 00142 * 00143 DO 110 IDIAG = 1, 2 00144 * 00145 DIAG = DIAGS( IDIAG ) 00146 * 00147 DO 100 IALPHA = 1, 3 00148 * 00149 IF ( IALPHA.EQ. 1) THEN 00150 ALPHA = ZERO 00151 ELSE IF ( IALPHA.EQ. 1) THEN 00152 ALPHA = ONE 00153 ELSE 00154 ALPHA = CLARND( 4, ISEED ) 00155 END IF 00156 * 00157 * All the parameters are set: 00158 * CFORM, SIDE, UPLO, TRANS, DIAG, M, N, 00159 * and ALPHA 00160 * READY TO TEST! 00161 * 00162 NRUN = NRUN + 1 00163 * 00164 IF ( ISIDE.EQ.1 ) THEN 00165 * 00166 * The case ISIDE.EQ.1 is when SIDE.EQ.'L' 00167 * -> A is M-by-M ( B is M-by-N ) 00168 * 00169 NA = M 00170 * 00171 ELSE 00172 * 00173 * The case ISIDE.EQ.2 is when SIDE.EQ.'R' 00174 * -> A is N-by-N ( B is M-by-N ) 00175 * 00176 NA = N 00177 * 00178 END IF 00179 * 00180 * Generate A our NA--by--NA triangular 00181 * matrix. 00182 * Our test is based on forward error so we 00183 * do want A to be well conditionned! To get 00184 * a well-conditionned triangular matrix, we 00185 * take the R factor of the QR/LQ factorization 00186 * of a random matrix. 00187 * 00188 DO J = 1, NA 00189 DO I = 1, NA 00190 A( I, J) = CLARND( 4, ISEED ) 00191 END DO 00192 END DO 00193 * 00194 IF ( IUPLO.EQ.1 ) THEN 00195 * 00196 * The case IUPLO.EQ.1 is when SIDE.EQ.'U' 00197 * -> QR factorization. 00198 * 00199 SRNAMT = 'CGEQRF' 00200 CALL CGEQRF( NA, NA, A, LDA, TAU, 00201 + C_WORK_CGEQRF, LDA, 00202 + INFO ) 00203 ELSE 00204 * 00205 * The case IUPLO.EQ.2 is when SIDE.EQ.'L' 00206 * -> QL factorization. 00207 * 00208 SRNAMT = 'CGELQF' 00209 CALL CGELQF( NA, NA, A, LDA, TAU, 00210 + C_WORK_CGEQRF, LDA, 00211 + INFO ) 00212 END IF 00213 * 00214 * After the QR factorization, the diagonal 00215 * of A is made of real numbers, we multiply 00216 * by a random complex number of absolute 00217 * value 1.0E+00. 00218 * 00219 DO J = 1, NA 00220 A( J, J) = A(J,J) * CLARND( 5, ISEED ) 00221 END DO 00222 * 00223 * Store a copy of A in RFP format (in ARF). 00224 * 00225 SRNAMT = 'CTRTTF' 00226 CALL CTRTTF( CFORM, UPLO, NA, A, LDA, ARF, 00227 + INFO ) 00228 * 00229 * Generate B1 our M--by--N right-hand side 00230 * and store a copy in B2. 00231 * 00232 DO J = 1, N 00233 DO I = 1, M 00234 B1( I, J) = CLARND( 4, ISEED ) 00235 B2( I, J) = B1( I, J) 00236 END DO 00237 END DO 00238 * 00239 * Solve op( A ) X = B or X op( A ) = B 00240 * with CTRSM 00241 * 00242 SRNAMT = 'CTRSM' 00243 CALL CTRSM( SIDE, UPLO, TRANS, DIAG, M, N, 00244 + ALPHA, A, LDA, B1, LDA ) 00245 * 00246 * Solve op( A ) X = B or X op( A ) = B 00247 * with CTFSM 00248 * 00249 SRNAMT = 'CTFSM' 00250 CALL CTFSM( CFORM, SIDE, UPLO, TRANS, 00251 + DIAG, M, N, ALPHA, ARF, B2, 00252 + LDA ) 00253 * 00254 * Check that the result agrees. 00255 * 00256 DO J = 1, N 00257 DO I = 1, M 00258 B1( I, J) = B2( I, J ) - B1( I, J ) 00259 END DO 00260 END DO 00261 * 00262 RESULT(1) = CLANGE( 'I', M, N, B1, LDA, 00263 + S_WORK_CLANGE ) 00264 * 00265 RESULT(1) = RESULT(1) / SQRT( EPS ) 00266 + / MAX ( MAX( M, N), 1 ) 00267 * 00268 IF( RESULT(1).GE.THRESH ) THEN 00269 IF( NFAIL.EQ.0 ) THEN 00270 WRITE( NOUT, * ) 00271 WRITE( NOUT, FMT = 9999 ) 00272 END IF 00273 WRITE( NOUT, FMT = 9997 ) 'CTFSM', 00274 + CFORM, SIDE, UPLO, TRANS, DIAG, M, 00275 + N, RESULT(1) 00276 NFAIL = NFAIL + 1 00277 END IF 00278 * 00279 100 CONTINUE 00280 110 CONTINUE 00281 120 CONTINUE 00282 130 CONTINUE 00283 140 CONTINUE 00284 150 CONTINUE 00285 160 CONTINUE 00286 170 CONTINUE 00287 * 00288 * Print a summary of the results. 00289 * 00290 IF ( NFAIL.EQ.0 ) THEN 00291 WRITE( NOUT, FMT = 9996 ) 'CTFSM', NRUN 00292 ELSE 00293 WRITE( NOUT, FMT = 9995 ) 'CTFSM', NFAIL, NRUN 00294 END IF 00295 * 00296 9999 FORMAT( 1X, 00297 ' *** Error(s) or Failure(s) while testing CTFSM + ***') 00298 9997 FORMAT( 1X, ' Failure in ',A5,', CFORM=''',A1,''',', 00299 + ' SIDE=''',A1,''',',' UPLO=''',A1,''',',' TRANS=''',A1,''',', 00300 + ' DIAG=''',A1,''',',' M=',I3,', N =', I3,', test=',G12.5) 00301 9996 FORMAT( 1X, 'All tests for ',A5,' auxiliary routine passed the ', 00302 + 'threshold (',I5,' tests run)') 00303 9995 FORMAT( 1X, A6, ' auxiliary routine:',I5,' out of ',I5, 00304 + ' tests failed to pass the threshold') 00305 * 00306 RETURN 00307 * 00308 * End of CDRVRF3 00309 * 00310 END