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