 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE DCHKPP( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, 00002 $ NMAX, A, AFAC, AINV, B, X, XACT, WORK, RWORK, 00003 $ IWORK, NOUT ) 00004 * 00005 * -- LAPACK test routine (version 3.1) -- 00006 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00007 * November 2006 00008 * 00009 * .. Scalar Arguments .. 00010 LOGICAL TSTERR 00011 INTEGER NMAX, NN, NNS, NOUT 00012 DOUBLE PRECISION THRESH 00013 * .. 00014 * .. Array Arguments .. 00015 LOGICAL DOTYPE( * ) 00016 INTEGER IWORK( * ), NSVAL( * ), NVAL( * ) 00017 DOUBLE PRECISION A( * ), AFAC( * ), AINV( * ), B( * ), 00018 $ RWORK( * ), WORK( * ), X( * ), XACT( * ) 00019 * .. 00020 * 00021 * Purpose 00022 * ======= 00023 * 00024 * DCHKPP tests DPPTRF, -TRI, -TRS, -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 * NNS (input) INTEGER 00041 * The number of values of NRHS contained in the vector NSVAL. 00042 * 00043 * NSVAL (input) INTEGER array, dimension (NNS) 00044 * The values of the number of right hand sides NRHS. 00045 * 00046 * THRESH (input) DOUBLE PRECISION 00047 * The threshold value for the test ratios. A result is 00048 * included in the output file if RESULT >= THRESH. To have 00049 * every test ratio printed, use THRESH = 0. 00050 * 00051 * TSTERR (input) LOGICAL 00052 * Flag that indicates whether error exits are to be tested. 00053 * 00054 * NMAX (input) INTEGER 00055 * The maximum value permitted for N, used in dimensioning the 00056 * work arrays. 00057 * 00058 * A (workspace) DOUBLE PRECISION array, dimension 00059 * (NMAX*(NMAX+1)/2) 00060 * 00061 * AFAC (workspace) DOUBLE PRECISION array, dimension 00062 * (NMAX*(NMAX+1)/2) 00063 * 00064 * AINV (workspace) DOUBLE PRECISION array, dimension 00065 * (NMAX*(NMAX+1)/2) 00066 * 00067 * B (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00068 * where NSMAX is the largest entry in NSVAL. 00069 * 00070 * X (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00071 * 00072 * XACT (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00073 * 00074 * WORK (workspace) DOUBLE PRECISION array, dimension 00075 * (NMAX*max(3,NSMAX)) 00076 * 00077 * RWORK (workspace) DOUBLE PRECISION array, dimension 00078 * (max(NMAX,2*NSMAX)) 00079 * 00080 * IWORK (workspace) INTEGER array, dimension (NMAX) 00081 * 00082 * NOUT (input) INTEGER 00083 * The unit number for output. 00084 * 00085 * ===================================================================== 00086 * 00087 * .. Parameters .. 00088 DOUBLE PRECISION ZERO 00089 PARAMETER ( ZERO = 0.0D+0 ) 00090 INTEGER NTYPES 00091 PARAMETER ( NTYPES = 9 ) 00092 INTEGER NTESTS 00093 PARAMETER ( NTESTS = 8 ) 00094 * .. 00095 * .. Local Scalars .. 00096 LOGICAL ZEROT 00097 CHARACTER DIST, PACKIT, TYPE, UPLO, XTYPE 00098 CHARACTER*3 PATH 00099 INTEGER I, IMAT, IN, INFO, IOFF, IRHS, IUPLO, IZERO, K, 00100 $ KL, KU, LDA, MODE, N, NERRS, NFAIL, NIMAT, NPP, 00101 $ NRHS, NRUN 00102 DOUBLE PRECISION ANORM, CNDNUM, RCOND, RCONDC 00103 * .. 00104 * .. Local Arrays .. 00105 CHARACTER PACKS( 2 ), UPLOS( 2 ) 00106 INTEGER ISEED( 4 ), ISEEDY( 4 ) 00107 DOUBLE PRECISION RESULT( NTESTS ) 00108 * .. 00109 * .. External Functions .. 00110 DOUBLE PRECISION DGET06, DLANSP 00111 EXTERNAL DGET06, DLANSP 00112 * .. 00113 * .. External Subroutines .. 00114 EXTERNAL ALAERH, ALAHD, ALASUM, DCOPY, DERRPO, DGET04, 00115 $ DLACPY, DLARHS, DLATB4, DLATMS, DPPCON, DPPRFS, 00116 $ DPPT01, DPPT02, DPPT03, DPPT05, DPPTRF, DPPTRI, 00117 $ DPPTRS 00118 * .. 00119 * .. Scalars in Common .. 00120 LOGICAL LERR, OK 00121 CHARACTER*32 SRNAMT 00122 INTEGER INFOT, NUNIT 00123 * .. 00124 * .. Common blocks .. 00125 COMMON / INFOC / INFOT, NUNIT, OK, LERR 00126 COMMON / SRNAMC / SRNAMT 00127 * .. 00128 * .. Intrinsic Functions .. 00129 INTRINSIC MAX 00130 * .. 00131 * .. Data statements .. 00132 DATA ISEEDY / 1988, 1989, 1990, 1991 / 00133 DATA UPLOS / 'U', 'L' / , PACKS / 'C', 'R' / 00134 * .. 00135 * .. Executable Statements .. 00136 * 00137 * Initialize constants and the random number seed. 00138 * 00139 PATH( 1: 1 ) = 'Double precision' 00140 PATH( 2: 3 ) = 'PP' 00141 NRUN = 0 00142 NFAIL = 0 00143 NERRS = 0 00144 DO 10 I = 1, 4 00145 ISEED( I ) = ISEEDY( I ) 00146 10 CONTINUE 00147 * 00148 * Test the error exits 00149 * 00150 IF( TSTERR ) 00151 $ CALL DERRPO( PATH, NOUT ) 00152 INFOT = 0 00153 * 00154 * Do for each value of N in NVAL 00155 * 00156 DO 110 IN = 1, NN 00157 N = NVAL( IN ) 00158 LDA = MAX( N, 1 ) 00159 XTYPE = 'N' 00160 NIMAT = NTYPES 00161 IF( N.LE.0 ) 00162 $ NIMAT = 1 00163 * 00164 DO 100 IMAT = 1, NIMAT 00165 * 00166 * Do the tests only if DOTYPE( IMAT ) is true. 00167 * 00168 IF( .NOT.DOTYPE( IMAT ) ) 00169 $ GO TO 100 00170 * 00171 * Skip types 3, 4, or 5 if the matrix size is too small. 00172 * 00173 ZEROT = IMAT.GE.3 .AND. IMAT.LE.5 00174 IF( ZEROT .AND. N.LT.IMAT-2 ) 00175 $ GO TO 100 00176 * 00177 * Do first for UPLO = 'U', then for UPLO = 'L' 00178 * 00179 DO 90 IUPLO = 1, 2 00180 UPLO = UPLOS( IUPLO ) 00181 PACKIT = PACKS( IUPLO ) 00182 * 00183 * Set up parameters with DLATB4 and generate a test matrix 00184 * with DLATMS. 00185 * 00186 CALL DLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE, 00187 $ CNDNUM, DIST ) 00188 * 00189 SRNAMT = 'DLATMS' 00190 CALL DLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE, 00191 $ CNDNUM, ANORM, KL, KU, PACKIT, A, LDA, WORK, 00192 $ INFO ) 00193 * 00194 * Check error code from DLATMS. 00195 * 00196 IF( INFO.NE.0 ) THEN 00197 CALL ALAERH( PATH, 'DLATMS', INFO, 0, UPLO, N, N, -1, 00198 $ -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00199 GO TO 90 00200 END IF 00201 * 00202 * For types 3-5, zero one row and column of the matrix to 00203 * test that INFO is returned correctly. 00204 * 00205 IF( ZEROT ) THEN 00206 IF( IMAT.EQ.3 ) THEN 00207 IZERO = 1 00208 ELSE IF( IMAT.EQ.4 ) THEN 00209 IZERO = N 00210 ELSE 00211 IZERO = N / 2 + 1 00212 END IF 00213 * 00214 * Set row and column IZERO of A to 0. 00215 * 00216 IF( IUPLO.EQ.1 ) THEN 00217 IOFF = ( IZERO-1 )*IZERO / 2 00218 DO 20 I = 1, IZERO - 1 00219 A( IOFF+I ) = ZERO 00220 20 CONTINUE 00221 IOFF = IOFF + IZERO 00222 DO 30 I = IZERO, N 00223 A( IOFF ) = ZERO 00224 IOFF = IOFF + I 00225 30 CONTINUE 00226 ELSE 00227 IOFF = IZERO 00228 DO 40 I = 1, IZERO - 1 00229 A( IOFF ) = ZERO 00230 IOFF = IOFF + N - I 00231 40 CONTINUE 00232 IOFF = IOFF - IZERO 00233 DO 50 I = IZERO, N 00234 A( IOFF+I ) = ZERO 00235 50 CONTINUE 00236 END IF 00237 ELSE 00238 IZERO = 0 00239 END IF 00240 * 00241 * Compute the L*L' or U'*U factorization of the matrix. 00242 * 00243 NPP = N*( N+1 ) / 2 00244 CALL DCOPY( NPP, A, 1, AFAC, 1 ) 00245 SRNAMT = 'DPPTRF' 00246 CALL DPPTRF( UPLO, N, AFAC, INFO ) 00247 * 00248 * Check error code from DPPTRF. 00249 * 00250 IF( INFO.NE.IZERO ) THEN 00251 CALL ALAERH( PATH, 'DPPTRF', INFO, IZERO, UPLO, N, N, 00252 $ -1, -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00253 GO TO 90 00254 END IF 00255 * 00256 * Skip the tests if INFO is not 0. 00257 * 00258 IF( INFO.NE.0 ) 00259 $ GO TO 90 00260 * 00261 *+ TEST 1 00262 * Reconstruct matrix from factors and compute residual. 00263 * 00264 CALL DCOPY( NPP, AFAC, 1, AINV, 1 ) 00265 CALL DPPT01( UPLO, N, A, AINV, RWORK, RESULT( 1 ) ) 00266 * 00267 *+ TEST 2 00268 * Form the inverse and compute the residual. 00269 * 00270 CALL DCOPY( NPP, AFAC, 1, AINV, 1 ) 00271 SRNAMT = 'DPPTRI' 00272 CALL DPPTRI( UPLO, N, AINV, INFO ) 00273 * 00274 * Check error code from DPPTRI. 00275 * 00276 IF( INFO.NE.0 ) 00277 $ CALL ALAERH( PATH, 'DPPTRI', INFO, 0, UPLO, N, N, -1, 00278 $ -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00279 * 00280 CALL DPPT03( UPLO, N, A, AINV, WORK, LDA, RWORK, RCONDC, 00281 $ RESULT( 2 ) ) 00282 * 00283 * Print information about the tests that did not pass 00284 * the threshold. 00285 * 00286 DO 60 K = 1, 2 00287 IF( RESULT( K ).GE.THRESH ) THEN 00288 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00289 $ CALL ALAHD( NOUT, PATH ) 00290 WRITE( NOUT, FMT = 9999 )UPLO, N, IMAT, K, 00291 $ RESULT( K ) 00292 NFAIL = NFAIL + 1 00293 END IF 00294 60 CONTINUE 00295 NRUN = NRUN + 2 00296 * 00297 DO 80 IRHS = 1, NNS 00298 NRHS = NSVAL( IRHS ) 00299 * 00300 *+ TEST 3 00301 * Solve and compute residual for A * X = B. 00302 * 00303 SRNAMT = 'DLARHS' 00304 CALL DLARHS( PATH, XTYPE, UPLO, ' ', N, N, KL, KU, 00305 $ NRHS, A, LDA, XACT, LDA, B, LDA, ISEED, 00306 $ INFO ) 00307 CALL DLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 00308 * 00309 SRNAMT = 'DPPTRS' 00310 CALL DPPTRS( UPLO, N, NRHS, AFAC, X, LDA, INFO ) 00311 * 00312 * Check error code from DPPTRS. 00313 * 00314 IF( INFO.NE.0 ) 00315 $ CALL ALAERH( PATH, 'DPPTRS', INFO, 0, UPLO, N, N, 00316 $ -1, -1, NRHS, IMAT, NFAIL, NERRS, 00317 $ NOUT ) 00318 * 00319 CALL DLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA ) 00320 CALL DPPT02( UPLO, N, NRHS, A, X, LDA, WORK, LDA, 00321 $ RWORK, RESULT( 3 ) ) 00322 * 00323 *+ TEST 4 00324 * Check solution from generated exact solution. 00325 * 00326 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00327 $ RESULT( 4 ) ) 00328 * 00329 *+ TESTS 5, 6, and 7 00330 * Use iterative refinement to improve the solution. 00331 * 00332 SRNAMT = 'DPPRFS' 00333 CALL DPPRFS( UPLO, N, NRHS, A, AFAC, B, LDA, X, LDA, 00334 $ RWORK, RWORK( NRHS+1 ), WORK, IWORK, 00335 $ INFO ) 00336 * 00337 * Check error code from DPPRFS. 00338 * 00339 IF( INFO.NE.0 ) 00340 $ CALL ALAERH( PATH, 'DPPRFS', INFO, 0, UPLO, N, N, 00341 $ -1, -1, NRHS, IMAT, NFAIL, NERRS, 00342 $ NOUT ) 00343 * 00344 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00345 $ RESULT( 5 ) ) 00346 CALL DPPT05( UPLO, N, NRHS, A, B, LDA, X, LDA, XACT, 00347 $ LDA, RWORK, RWORK( NRHS+1 ), 00348 $ RESULT( 6 ) ) 00349 * 00350 * Print information about the tests that did not pass 00351 * the threshold. 00352 * 00353 DO 70 K = 3, 7 00354 IF( RESULT( K ).GE.THRESH ) THEN 00355 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00356 $ CALL ALAHD( NOUT, PATH ) 00357 WRITE( NOUT, FMT = 9998 )UPLO, N, NRHS, IMAT, 00358 $ K, RESULT( K ) 00359 NFAIL = NFAIL + 1 00360 END IF 00361 70 CONTINUE 00362 NRUN = NRUN + 5 00363 80 CONTINUE 00364 * 00365 *+ TEST 8 00366 * Get an estimate of RCOND = 1/CNDNUM. 00367 * 00368 ANORM = DLANSP( '1', UPLO, N, A, RWORK ) 00369 SRNAMT = 'DPPCON' 00370 CALL DPPCON( UPLO, N, AFAC, ANORM, RCOND, WORK, IWORK, 00371 $ INFO ) 00372 * 00373 * Check error code from DPPCON. 00374 * 00375 IF( INFO.NE.0 ) 00376 $ CALL ALAERH( PATH, 'DPPCON', INFO, 0, UPLO, N, N, -1, 00377 $ -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00378 * 00379 RESULT( 8 ) = DGET06( RCOND, RCONDC ) 00380 * 00381 * Print the test ratio if greater than or equal to THRESH. 00382 * 00383 IF( RESULT( 8 ).GE.THRESH ) THEN 00384 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00385 $ CALL ALAHD( NOUT, PATH ) 00386 WRITE( NOUT, FMT = 9999 )UPLO, N, IMAT, 8, 00387 $ RESULT( 8 ) 00388 NFAIL = NFAIL + 1 00389 END IF 00390 NRUN = NRUN + 1 00391 90 CONTINUE 00392 100 CONTINUE 00393 110 CONTINUE 00394 * 00395 * Print a summary of the results. 00396 * 00397 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS ) 00398 * 00399 9999 FORMAT( ' UPLO = ''', A1, ''', N =', I5, ', type ', I2, ', test ', 00400 $ I2, ', ratio =', G12.5 ) 00401 9998 FORMAT( ' UPLO = ''', A1, ''', N =', I5, ', NRHS=', I3, ', type ', 00402 $ I2, ', test(', I2, ') =', G12.5 ) 00403 RETURN 00404 * 00405 * End of DCHKPP 00406 * 00407 END
1.7.3 
