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