 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE SCHKGE( DOTYPE, NM, MVAL, NN, NVAL, NNB, NBVAL, NNS, 00002 $ NSVAL, THRESH, TSTERR, NMAX, A, AFAC, AINV, B, 00003 $ X, XACT, WORK, RWORK, IWORK, NOUT ) 00004 * 00005 * -- LAPACK test routine (version 3.1.1) -- 00006 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00007 * January 2007 00008 * 00009 * .. Scalar Arguments .. 00010 LOGICAL TSTERR 00011 INTEGER NM, NMAX, NN, NNB, NNS, NOUT 00012 REAL THRESH 00013 * .. 00014 * .. Array Arguments .. 00015 LOGICAL DOTYPE( * ) 00016 INTEGER IWORK( * ), MVAL( * ), NBVAL( * ), NSVAL( * ), 00017 $ NVAL( * ) 00018 REAL A( * ), AFAC( * ), AINV( * ), B( * ), 00019 $ RWORK( * ), WORK( * ), X( * ), XACT( * ) 00020 * .. 00021 * 00022 * Purpose 00023 * ======= 00024 * 00025 * SCHKGE tests SGETRF, -TRI, -TRS, -RFS, and -CON. 00026 * 00027 * Arguments 00028 * ========= 00029 * 00030 * DOTYPE (input) LOGICAL array, dimension (NTYPES) 00031 * The matrix types to be used for testing. Matrices of type j 00032 * (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) = 00033 * .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used. 00034 * 00035 * NM (input) INTEGER 00036 * The number of values of M contained in the vector MVAL. 00037 * 00038 * MVAL (input) INTEGER array, dimension (NM) 00039 * The values of the matrix row dimension M. 00040 * 00041 * NN (input) INTEGER 00042 * The number of values of N contained in the vector NVAL. 00043 * 00044 * NVAL (input) INTEGER array, dimension (NN) 00045 * The values of the matrix column dimension N. 00046 * 00047 * NNB (input) INTEGER 00048 * The number of values of NB contained in the vector NBVAL. 00049 * 00050 * NBVAL (input) INTEGER array, dimension (NBVAL) 00051 * The values of the blocksize NB. 00052 * 00053 * NNS (input) INTEGER 00054 * The number of values of NRHS contained in the vector NSVAL. 00055 * 00056 * NSVAL (input) INTEGER array, dimension (NNS) 00057 * The values of the number of right hand sides NRHS. 00058 * 00059 * THRESH (input) REAL 00060 * The threshold value for the test ratios. A result is 00061 * included in the output file if RESULT >= THRESH. To have 00062 * every test ratio printed, use THRESH = 0. 00063 * 00064 * TSTERR (input) LOGICAL 00065 * Flag that indicates whether error exits are to be tested. 00066 * 00067 * NMAX (input) INTEGER 00068 * The maximum value permitted for M or N, used in dimensioning 00069 * the work arrays. 00070 * 00071 * A (workspace) REAL array, dimension (NMAX*NMAX) 00072 * 00073 * AFAC (workspace) REAL array, dimension (NMAX*NMAX) 00074 * 00075 * AINV (workspace) REAL array, dimension (NMAX*NMAX) 00076 * 00077 * B (workspace) REAL array, dimension (NMAX*NSMAX) 00078 * where NSMAX is the largest entry in NSVAL. 00079 * 00080 * X (workspace) REAL array, dimension (NMAX*NSMAX) 00081 * 00082 * XACT (workspace) REAL array, dimension (NMAX*NSMAX) 00083 * 00084 * WORK (workspace) REAL array, dimension 00085 * (NMAX*max(3,NSMAX)) 00086 * 00087 * RWORK (workspace) REAL array, dimension 00088 * (max(2*NMAX,2*NSMAX+NWORK)) 00089 * 00090 * IWORK (workspace) INTEGER array, dimension (2*NMAX) 00091 * 00092 * NOUT (input) INTEGER 00093 * The unit number for output. 00094 * 00095 * ===================================================================== 00096 * 00097 * .. Parameters .. 00098 REAL ONE, ZERO 00099 PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) 00100 INTEGER NTYPES 00101 PARAMETER ( NTYPES = 11 ) 00102 INTEGER NTESTS 00103 PARAMETER ( NTESTS = 8 ) 00104 INTEGER NTRAN 00105 PARAMETER ( NTRAN = 3 ) 00106 * .. 00107 * .. Local Scalars .. 00108 LOGICAL TRFCON, ZEROT 00109 CHARACTER DIST, NORM, TRANS, TYPE, XTYPE 00110 CHARACTER*3 PATH 00111 INTEGER I, IM, IMAT, IN, INB, INFO, IOFF, IRHS, ITRAN, 00112 $ IZERO, K, KL, KU, LDA, LWORK, M, MODE, N, NB, 00113 $ NERRS, NFAIL, NIMAT, NRHS, NRUN, NT 00114 REAL AINVNM, ANORM, ANORMI, ANORMO, CNDNUM, DUMMY, 00115 $ RCOND, RCONDC, RCONDI, RCONDO 00116 * .. 00117 * .. Local Arrays .. 00118 CHARACTER TRANSS( NTRAN ) 00119 INTEGER ISEED( 4 ), ISEEDY( 4 ) 00120 REAL RESULT( NTESTS ) 00121 * .. 00122 * .. External Functions .. 00123 REAL SGET06, SLANGE 00124 EXTERNAL SGET06, SLANGE 00125 * .. 00126 * .. External Subroutines .. 00127 EXTERNAL ALAERH, ALAHD, ALASUM, SERRGE, SGECON, SGERFS, 00128 $ SGET01, SGET02, SGET03, SGET04, SGET07, SGETRF, 00129 $ SGETRI, SGETRS, SLACPY, SLARHS, SLASET, SLATB4, 00130 $ SLATMS, XLAENV 00131 * .. 00132 * .. Intrinsic Functions .. 00133 INTRINSIC MAX, MIN 00134 * .. 00135 * .. Scalars in Common .. 00136 LOGICAL LERR, OK 00137 CHARACTER*32 SRNAMT 00138 INTEGER INFOT, NUNIT 00139 * .. 00140 * .. Common blocks .. 00141 COMMON / INFOC / INFOT, NUNIT, OK, LERR 00142 COMMON / SRNAMC / SRNAMT 00143 * .. 00144 * .. Data statements .. 00145 DATA ISEEDY / 1988, 1989, 1990, 1991 / , 00146 $ TRANSS / 'N', 'T', 'C' / 00147 * .. 00148 * .. Executable Statements .. 00149 * 00150 * Initialize constants and the random number seed. 00151 * 00152 PATH( 1: 1 ) = 'Single precision' 00153 PATH( 2: 3 ) = 'GE' 00154 NRUN = 0 00155 NFAIL = 0 00156 NERRS = 0 00157 DO 10 I = 1, 4 00158 ISEED( I ) = ISEEDY( I ) 00159 10 CONTINUE 00160 * 00161 * Test the error exits 00162 * 00163 CALL XLAENV( 1, 1 ) 00164 IF( TSTERR ) 00165 $ CALL SERRGE( PATH, NOUT ) 00166 INFOT = 0 00167 CALL XLAENV( 2, 2 ) 00168 * 00169 * Do for each value of M in MVAL 00170 * 00171 DO 120 IM = 1, NM 00172 M = MVAL( IM ) 00173 LDA = MAX( 1, M ) 00174 * 00175 * Do for each value of N in NVAL 00176 * 00177 DO 110 IN = 1, NN 00178 N = NVAL( IN ) 00179 XTYPE = 'N' 00180 NIMAT = NTYPES 00181 IF( M.LE.0 .OR. N.LE.0 ) 00182 $ NIMAT = 1 00183 * 00184 DO 100 IMAT = 1, NIMAT 00185 * 00186 * Do the tests only if DOTYPE( IMAT ) is true. 00187 * 00188 IF( .NOT.DOTYPE( IMAT ) ) 00189 $ GO TO 100 00190 * 00191 * Skip types 5, 6, or 7 if the matrix size is too small. 00192 * 00193 ZEROT = IMAT.GE.5 .AND. IMAT.LE.7 00194 IF( ZEROT .AND. N.LT.IMAT-4 ) 00195 $ GO TO 100 00196 * 00197 * Set up parameters with SLATB4 and generate a test matrix 00198 * with SLATMS. 00199 * 00200 CALL SLATB4( PATH, IMAT, M, N, TYPE, KL, KU, ANORM, MODE, 00201 $ CNDNUM, DIST ) 00202 * 00203 SRNAMT = 'SLATMS' 00204 CALL SLATMS( M, N, DIST, ISEED, TYPE, RWORK, MODE, 00205 $ CNDNUM, ANORM, KL, KU, 'No packing', A, LDA, 00206 $ WORK, INFO ) 00207 * 00208 * Check error code from SLATMS. 00209 * 00210 IF( INFO.NE.0 ) THEN 00211 CALL ALAERH( PATH, 'SLATMS', INFO, 0, ' ', M, N, -1, 00212 $ -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00213 GO TO 100 00214 END IF 00215 * 00216 * For types 5-7, zero one or more columns of the matrix to 00217 * test that INFO is returned correctly. 00218 * 00219 IF( ZEROT ) THEN 00220 IF( IMAT.EQ.5 ) THEN 00221 IZERO = 1 00222 ELSE IF( IMAT.EQ.6 ) THEN 00223 IZERO = MIN( M, N ) 00224 ELSE 00225 IZERO = MIN( M, N ) / 2 + 1 00226 END IF 00227 IOFF = ( IZERO-1 )*LDA 00228 IF( IMAT.LT.7 ) THEN 00229 DO 20 I = 1, M 00230 A( IOFF+I ) = ZERO 00231 20 CONTINUE 00232 ELSE 00233 CALL SLASET( 'Full', M, N-IZERO+1, ZERO, ZERO, 00234 $ A( IOFF+1 ), LDA ) 00235 END IF 00236 ELSE 00237 IZERO = 0 00238 END IF 00239 * 00240 * These lines, if used in place of the calls in the DO 60 00241 * loop, cause the code to bomb on a Sun SPARCstation. 00242 * 00243 * ANORMO = SLANGE( 'O', M, N, A, LDA, RWORK ) 00244 * ANORMI = SLANGE( 'I', M, N, A, LDA, RWORK ) 00245 * 00246 * Do for each blocksize in NBVAL 00247 * 00248 DO 90 INB = 1, NNB 00249 NB = NBVAL( INB ) 00250 CALL XLAENV( 1, NB ) 00251 * 00252 * Compute the LU factorization of the matrix. 00253 * 00254 CALL SLACPY( 'Full', M, N, A, LDA, AFAC, LDA ) 00255 SRNAMT = 'SGETRF' 00256 CALL SGETRF( M, N, AFAC, LDA, IWORK, INFO ) 00257 * 00258 * Check error code from SGETRF. 00259 * 00260 IF( INFO.NE.IZERO ) 00261 $ CALL ALAERH( PATH, 'SGETRF', INFO, IZERO, ' ', M, 00262 $ N, -1, -1, NB, IMAT, NFAIL, NERRS, 00263 $ NOUT ) 00264 TRFCON = .FALSE. 00265 * 00266 *+ TEST 1 00267 * Reconstruct matrix from factors and compute residual. 00268 * 00269 CALL SLACPY( 'Full', M, N, AFAC, LDA, AINV, LDA ) 00270 CALL SGET01( M, N, A, LDA, AINV, LDA, IWORK, RWORK, 00271 $ RESULT( 1 ) ) 00272 NT = 1 00273 * 00274 *+ TEST 2 00275 * Form the inverse if the factorization was successful 00276 * and compute the residual. 00277 * 00278 IF( M.EQ.N .AND. INFO.EQ.0 ) THEN 00279 CALL SLACPY( 'Full', N, N, AFAC, LDA, AINV, LDA ) 00280 SRNAMT = 'SGETRI' 00281 NRHS = NSVAL( 1 ) 00282 LWORK = NMAX*MAX( 3, NRHS ) 00283 CALL SGETRI( N, AINV, LDA, IWORK, WORK, LWORK, 00284 $ INFO ) 00285 * 00286 * Check error code from SGETRI. 00287 * 00288 IF( INFO.NE.0 ) 00289 $ CALL ALAERH( PATH, 'SGETRI', INFO, 0, ' ', N, N, 00290 $ -1, -1, NB, IMAT, NFAIL, NERRS, 00291 $ NOUT ) 00292 * 00293 * Compute the residual for the matrix times its 00294 * inverse. Also compute the 1-norm condition number 00295 * of A. 00296 * 00297 CALL SGET03( N, A, LDA, AINV, LDA, WORK, LDA, 00298 $ RWORK, RCONDO, RESULT( 2 ) ) 00299 ANORMO = SLANGE( 'O', M, N, A, LDA, RWORK ) 00300 * 00301 * Compute the infinity-norm condition number of A. 00302 * 00303 ANORMI = SLANGE( 'I', M, N, A, LDA, RWORK ) 00304 AINVNM = SLANGE( 'I', N, N, AINV, LDA, RWORK ) 00305 IF( ANORMI.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN 00306 RCONDI = ONE 00307 ELSE 00308 RCONDI = ( ONE / ANORMI ) / AINVNM 00309 END IF 00310 NT = 2 00311 ELSE 00312 * 00313 * Do only the condition estimate if INFO > 0. 00314 * 00315 TRFCON = .TRUE. 00316 ANORMO = SLANGE( 'O', M, N, A, LDA, RWORK ) 00317 ANORMI = SLANGE( 'I', M, N, A, LDA, RWORK ) 00318 RCONDO = ZERO 00319 RCONDI = ZERO 00320 END IF 00321 * 00322 * Print information about the tests so far that did not 00323 * pass the threshold. 00324 * 00325 DO 30 K = 1, NT 00326 IF( RESULT( K ).GE.THRESH ) THEN 00327 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00328 $ CALL ALAHD( NOUT, PATH ) 00329 WRITE( NOUT, FMT = 9999 )M, N, NB, IMAT, K, 00330 $ RESULT( K ) 00331 NFAIL = NFAIL + 1 00332 END IF 00333 30 CONTINUE 00334 NRUN = NRUN + NT 00335 * 00336 * Skip the remaining tests if this is not the first 00337 * block size or if M .ne. N. Skip the solve tests if 00338 * the matrix is singular. 00339 * 00340 IF( INB.GT.1 .OR. M.NE.N ) 00341 $ GO TO 90 00342 IF( TRFCON ) 00343 $ GO TO 70 00344 * 00345 DO 60 IRHS = 1, NNS 00346 NRHS = NSVAL( IRHS ) 00347 XTYPE = 'N' 00348 * 00349 DO 50 ITRAN = 1, NTRAN 00350 TRANS = TRANSS( ITRAN ) 00351 IF( ITRAN.EQ.1 ) THEN 00352 RCONDC = RCONDO 00353 ELSE 00354 RCONDC = RCONDI 00355 END IF 00356 * 00357 *+ TEST 3 00358 * Solve and compute residual for A * X = B. 00359 * 00360 SRNAMT = 'SLARHS' 00361 CALL SLARHS( PATH, XTYPE, ' ', TRANS, N, N, KL, 00362 $ KU, NRHS, A, LDA, XACT, LDA, B, 00363 $ LDA, ISEED, INFO ) 00364 XTYPE = 'C' 00365 * 00366 CALL SLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 00367 SRNAMT = 'SGETRS' 00368 CALL SGETRS( TRANS, N, NRHS, AFAC, LDA, IWORK, 00369 $ X, LDA, INFO ) 00370 * 00371 * Check error code from SGETRS. 00372 * 00373 IF( INFO.NE.0 ) 00374 $ CALL ALAERH( PATH, 'SGETRS', INFO, 0, TRANS, 00375 $ N, N, -1, -1, NRHS, IMAT, NFAIL, 00376 $ NERRS, NOUT ) 00377 * 00378 CALL SLACPY( 'Full', N, NRHS, B, LDA, WORK, 00379 $ LDA ) 00380 CALL SGET02( TRANS, N, N, NRHS, A, LDA, X, LDA, 00381 $ WORK, LDA, RWORK, RESULT( 3 ) ) 00382 * 00383 *+ TEST 4 00384 * Check solution from generated exact solution. 00385 * 00386 CALL SGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00387 $ RESULT( 4 ) ) 00388 * 00389 *+ TESTS 5, 6, and 7 00390 * Use iterative refinement to improve the 00391 * solution. 00392 * 00393 SRNAMT = 'SGERFS' 00394 CALL SGERFS( TRANS, N, NRHS, A, LDA, AFAC, LDA, 00395 $ IWORK, B, LDA, X, LDA, RWORK, 00396 $ RWORK( NRHS+1 ), WORK, 00397 $ IWORK( N+1 ), INFO ) 00398 * 00399 * Check error code from SGERFS. 00400 * 00401 IF( INFO.NE.0 ) 00402 $ CALL ALAERH( PATH, 'SGERFS', INFO, 0, TRANS, 00403 $ N, N, -1, -1, NRHS, IMAT, NFAIL, 00404 $ NERRS, NOUT ) 00405 * 00406 CALL SGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00407 $ RESULT( 5 ) ) 00408 CALL SGET07( TRANS, N, NRHS, A, LDA, B, LDA, X, 00409 $ LDA, XACT, LDA, RWORK, .TRUE., 00410 $ RWORK( NRHS+1 ), RESULT( 6 ) ) 00411 * 00412 * Print information about the tests that did not 00413 * pass the threshold. 00414 * 00415 DO 40 K = 3, 7 00416 IF( RESULT( K ).GE.THRESH ) THEN 00417 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00418 $ CALL ALAHD( NOUT, PATH ) 00419 WRITE( NOUT, FMT = 9998 )TRANS, N, NRHS, 00420 $ IMAT, K, RESULT( K ) 00421 NFAIL = NFAIL + 1 00422 END IF 00423 40 CONTINUE 00424 NRUN = NRUN + 5 00425 50 CONTINUE 00426 60 CONTINUE 00427 * 00428 *+ TEST 8 00429 * Get an estimate of RCOND = 1/CNDNUM. 00430 * 00431 70 CONTINUE 00432 DO 80 ITRAN = 1, 2 00433 IF( ITRAN.EQ.1 ) THEN 00434 ANORM = ANORMO 00435 RCONDC = RCONDO 00436 NORM = 'O' 00437 ELSE 00438 ANORM = ANORMI 00439 RCONDC = RCONDI 00440 NORM = 'I' 00441 END IF 00442 SRNAMT = 'SGECON' 00443 CALL SGECON( NORM, N, AFAC, LDA, ANORM, RCOND, 00444 $ WORK, IWORK( N+1 ), INFO ) 00445 * 00446 * Check error code from SGECON. 00447 * 00448 IF( INFO.NE.0 ) 00449 $ CALL ALAERH( PATH, 'SGECON', INFO, 0, NORM, N, 00450 $ N, -1, -1, -1, IMAT, NFAIL, NERRS, 00451 $ NOUT ) 00452 * 00453 * This line is needed on a Sun SPARCstation. 00454 * 00455 DUMMY = RCOND 00456 * 00457 RESULT( 8 ) = SGET06( RCOND, RCONDC ) 00458 * 00459 * Print information about the tests that did not pass 00460 * the threshold. 00461 * 00462 IF( RESULT( 8 ).GE.THRESH ) THEN 00463 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00464 $ CALL ALAHD( NOUT, PATH ) 00465 WRITE( NOUT, FMT = 9997 )NORM, N, IMAT, 8, 00466 $ RESULT( 8 ) 00467 NFAIL = NFAIL + 1 00468 END IF 00469 NRUN = NRUN + 1 00470 80 CONTINUE 00471 90 CONTINUE 00472 100 CONTINUE 00473 110 CONTINUE 00474 120 CONTINUE 00475 * 00476 * Print a summary of the results. 00477 * 00478 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS ) 00479 * 00480 9999 FORMAT( ' M = ', I5, ', N =', I5, ', NB =', I4, ', type ', I2, 00481 $ ', test(', I2, ') =', G12.5 ) 00482 9998 FORMAT( ' TRANS=''', A1, ''', N =', I5, ', NRHS=', I3, ', type ', 00483 $ I2, ', test(', I2, ') =', G12.5 ) 00484 9997 FORMAT( ' NORM =''', A1, ''', N =', I5, ',', 10X, ' type ', I2, 00485 $ ', test(', I2, ') =', G12.5 ) 00486 RETURN 00487 * 00488 * End of SCHKGE 00489 * 00490 END
1.7.3 
