![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DCHKTP 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 DCHKTP( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, 00012 * NMAX, AP, AINVP, 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 AINVP( * ), AP( * ), B( * ), RWORK( * ), 00024 * $ WORK( * ), X( * ), XACT( * ) 00025 * .. 00026 * 00027 * 00028 *> \par Purpose: 00029 * ============= 00030 *> 00031 *> \verbatim 00032 *> 00033 *> DCHKTP tests DTPTRI, -TRS, -RFS, and -CON, and DLATPS 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 column 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 leading dimension of the work arrays. NMAX >= the 00089 *> maximumm value of N in NVAL. 00090 *> \endverbatim 00091 *> 00092 *> \param[out] AP 00093 *> \verbatim 00094 *> AP is DOUBLE PRECISION array, dimension 00095 *> (NMAX*(NMAX+1)/2) 00096 *> \endverbatim 00097 *> 00098 *> \param[out] AINVP 00099 *> \verbatim 00100 *> AINVP is DOUBLE PRECISION array, dimension 00101 *> (NMAX*(NMAX+1)/2) 00102 *> \endverbatim 00103 *> 00104 *> \param[out] B 00105 *> \verbatim 00106 *> B is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00107 *> where NSMAX is the largest entry in NSVAL. 00108 *> \endverbatim 00109 *> 00110 *> \param[out] X 00111 *> \verbatim 00112 *> X is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00113 *> \endverbatim 00114 *> 00115 *> \param[out] XACT 00116 *> \verbatim 00117 *> XACT is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00118 *> \endverbatim 00119 *> 00120 *> \param[out] WORK 00121 *> \verbatim 00122 *> WORK is DOUBLE PRECISION array, dimension 00123 *> (NMAX*max(3,NSMAX)) 00124 *> \endverbatim 00125 *> 00126 *> \param[out] IWORK 00127 *> \verbatim 00128 *> IWORK is INTEGER array, dimension (NMAX) 00129 *> \endverbatim 00130 *> 00131 *> \param[out] RWORK 00132 *> \verbatim 00133 *> RWORK is DOUBLE PRECISION array, dimension 00134 *> (max(NMAX,2*NSMAX)) 00135 *> \endverbatim 00136 *> 00137 *> \param[in] NOUT 00138 *> \verbatim 00139 *> NOUT is INTEGER 00140 *> The unit number for output. 00141 *> \endverbatim 00142 * 00143 * Authors: 00144 * ======== 00145 * 00146 *> \author Univ. of Tennessee 00147 *> \author Univ. of California Berkeley 00148 *> \author Univ. of Colorado Denver 00149 *> \author NAG Ltd. 00150 * 00151 *> \date November 2011 00152 * 00153 *> \ingroup double_lin 00154 * 00155 * ===================================================================== 00156 SUBROUTINE DCHKTP( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, 00157 $ NMAX, AP, AINVP, B, X, XACT, WORK, RWORK, 00158 $ IWORK, NOUT ) 00159 * 00160 * -- LAPACK test routine (version 3.4.0) -- 00161 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00162 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00163 * November 2011 00164 * 00165 * .. Scalar Arguments .. 00166 LOGICAL TSTERR 00167 INTEGER NMAX, NN, NNS, NOUT 00168 DOUBLE PRECISION THRESH 00169 * .. 00170 * .. Array Arguments .. 00171 LOGICAL DOTYPE( * ) 00172 INTEGER IWORK( * ), NSVAL( * ), NVAL( * ) 00173 DOUBLE PRECISION AINVP( * ), AP( * ), B( * ), RWORK( * ), 00174 $ WORK( * ), X( * ), XACT( * ) 00175 * .. 00176 * 00177 * ===================================================================== 00178 * 00179 * .. Parameters .. 00180 INTEGER NTYPE1, NTYPES 00181 PARAMETER ( NTYPE1 = 10, NTYPES = 18 ) 00182 INTEGER NTESTS 00183 PARAMETER ( NTESTS = 9 ) 00184 INTEGER NTRAN 00185 PARAMETER ( NTRAN = 3 ) 00186 DOUBLE PRECISION ONE, ZERO 00187 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00188 * .. 00189 * .. Local Scalars .. 00190 CHARACTER DIAG, NORM, TRANS, UPLO, XTYPE 00191 CHARACTER*3 PATH 00192 INTEGER I, IDIAG, IMAT, IN, INFO, IRHS, ITRAN, IUPLO, 00193 $ K, LAP, LDA, N, NERRS, NFAIL, NRHS, NRUN 00194 DOUBLE PRECISION AINVNM, ANORM, RCOND, RCONDC, RCONDI, RCONDO, 00195 $ SCALE 00196 * .. 00197 * .. Local Arrays .. 00198 CHARACTER TRANSS( NTRAN ), UPLOS( 2 ) 00199 INTEGER ISEED( 4 ), ISEEDY( 4 ) 00200 DOUBLE PRECISION RESULT( NTESTS ) 00201 * .. 00202 * .. External Functions .. 00203 LOGICAL LSAME 00204 DOUBLE PRECISION DLANTP 00205 EXTERNAL LSAME, DLANTP 00206 * .. 00207 * .. External Subroutines .. 00208 EXTERNAL ALAERH, ALAHD, ALASUM, DCOPY, DERRTR, DGET04, 00209 $ DLACPY, DLARHS, DLATPS, DLATTP, DTPCON, DTPRFS, 00210 $ DTPT01, DTPT02, DTPT03, DTPT05, DTPT06, DTPTRI, 00211 $ DTPTRS 00212 * .. 00213 * .. Scalars in Common .. 00214 LOGICAL LERR, OK 00215 CHARACTER*32 SRNAMT 00216 INTEGER INFOT, IOUNIT 00217 * .. 00218 * .. Common blocks .. 00219 COMMON / INFOC / INFOT, IOUNIT, OK, LERR 00220 COMMON / SRNAMC / SRNAMT 00221 * .. 00222 * .. Intrinsic Functions .. 00223 INTRINSIC MAX 00224 * .. 00225 * .. Data statements .. 00226 DATA ISEEDY / 1988, 1989, 1990, 1991 / 00227 DATA UPLOS / 'U', 'L' / , TRANSS / 'N', 'T', 'C' / 00228 * .. 00229 * .. Executable Statements .. 00230 * 00231 * Initialize constants and the random number seed. 00232 * 00233 PATH( 1: 1 ) = 'Double precision' 00234 PATH( 2: 3 ) = 'TP' 00235 NRUN = 0 00236 NFAIL = 0 00237 NERRS = 0 00238 DO 10 I = 1, 4 00239 ISEED( I ) = ISEEDY( I ) 00240 10 CONTINUE 00241 * 00242 * Test the error exits 00243 * 00244 IF( TSTERR ) 00245 $ CALL DERRTR( PATH, NOUT ) 00246 INFOT = 0 00247 * 00248 DO 110 IN = 1, NN 00249 * 00250 * Do for each value of N in NVAL 00251 * 00252 N = NVAL( IN ) 00253 LDA = MAX( 1, N ) 00254 LAP = LDA*( LDA+1 ) / 2 00255 XTYPE = 'N' 00256 * 00257 DO 70 IMAT = 1, NTYPE1 00258 * 00259 * Do the tests only if DOTYPE( IMAT ) is true. 00260 * 00261 IF( .NOT.DOTYPE( IMAT ) ) 00262 $ GO TO 70 00263 * 00264 DO 60 IUPLO = 1, 2 00265 * 00266 * Do first for UPLO = 'U', then for UPLO = 'L' 00267 * 00268 UPLO = UPLOS( IUPLO ) 00269 * 00270 * Call DLATTP to generate a triangular test matrix. 00271 * 00272 SRNAMT = 'DLATTP' 00273 CALL DLATTP( IMAT, UPLO, 'No transpose', DIAG, ISEED, N, 00274 $ AP, X, WORK, INFO ) 00275 * 00276 * Set IDIAG = 1 for non-unit matrices, 2 for unit. 00277 * 00278 IF( LSAME( DIAG, 'N' ) ) THEN 00279 IDIAG = 1 00280 ELSE 00281 IDIAG = 2 00282 END IF 00283 * 00284 *+ TEST 1 00285 * Form the inverse of A. 00286 * 00287 IF( N.GT.0 ) 00288 $ CALL DCOPY( LAP, AP, 1, AINVP, 1 ) 00289 SRNAMT = 'DTPTRI' 00290 CALL DTPTRI( UPLO, DIAG, N, AINVP, INFO ) 00291 * 00292 * Check error code from DTPTRI. 00293 * 00294 IF( INFO.NE.0 ) 00295 $ CALL ALAERH( PATH, 'DTPTRI', INFO, 0, UPLO // DIAG, N, 00296 $ N, -1, -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00297 * 00298 * Compute the infinity-norm condition number of A. 00299 * 00300 ANORM = DLANTP( 'I', UPLO, DIAG, N, AP, RWORK ) 00301 AINVNM = DLANTP( 'I', UPLO, DIAG, N, AINVP, RWORK ) 00302 IF( ANORM.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN 00303 RCONDI = ONE 00304 ELSE 00305 RCONDI = ( ONE / ANORM ) / AINVNM 00306 END IF 00307 * 00308 * Compute the residual for the triangular matrix times its 00309 * inverse. Also compute the 1-norm condition number of A. 00310 * 00311 CALL DTPT01( UPLO, DIAG, N, AP, AINVP, RCONDO, RWORK, 00312 $ RESULT( 1 ) ) 00313 * 00314 * Print the test ratio if it is .GE. THRESH. 00315 * 00316 IF( RESULT( 1 ).GE.THRESH ) THEN 00317 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00318 $ CALL ALAHD( NOUT, PATH ) 00319 WRITE( NOUT, FMT = 9999 )UPLO, DIAG, N, IMAT, 1, 00320 $ RESULT( 1 ) 00321 NFAIL = NFAIL + 1 00322 END IF 00323 NRUN = NRUN + 1 00324 * 00325 DO 40 IRHS = 1, NNS 00326 NRHS = NSVAL( IRHS ) 00327 XTYPE = 'N' 00328 * 00329 DO 30 ITRAN = 1, NTRAN 00330 * 00331 * Do for op(A) = A, A**T, or A**H. 00332 * 00333 TRANS = TRANSS( ITRAN ) 00334 IF( ITRAN.EQ.1 ) THEN 00335 NORM = 'O' 00336 RCONDC = RCONDO 00337 ELSE 00338 NORM = 'I' 00339 RCONDC = RCONDI 00340 END IF 00341 * 00342 *+ TEST 2 00343 * Solve and compute residual for op(A)*x = b. 00344 * 00345 SRNAMT = 'DLARHS' 00346 CALL DLARHS( PATH, XTYPE, UPLO, TRANS, N, N, 0, 00347 $ IDIAG, NRHS, AP, LAP, XACT, LDA, B, 00348 $ LDA, ISEED, INFO ) 00349 XTYPE = 'C' 00350 CALL DLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 00351 * 00352 SRNAMT = 'DTPTRS' 00353 CALL DTPTRS( UPLO, TRANS, DIAG, N, NRHS, AP, X, 00354 $ LDA, INFO ) 00355 * 00356 * Check error code from DTPTRS. 00357 * 00358 IF( INFO.NE.0 ) 00359 $ CALL ALAERH( PATH, 'DTPTRS', INFO, 0, 00360 $ UPLO // TRANS // DIAG, N, N, -1, 00361 $ -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00362 * 00363 CALL DTPT02( UPLO, TRANS, DIAG, N, NRHS, AP, X, 00364 $ LDA, B, LDA, WORK, RESULT( 2 ) ) 00365 * 00366 *+ TEST 3 00367 * Check solution from generated exact solution. 00368 * 00369 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00370 $ RESULT( 3 ) ) 00371 * 00372 *+ TESTS 4, 5, and 6 00373 * Use iterative refinement to improve the solution and 00374 * compute error bounds. 00375 * 00376 SRNAMT = 'DTPRFS' 00377 CALL DTPRFS( UPLO, TRANS, DIAG, N, NRHS, AP, B, 00378 $ LDA, X, LDA, RWORK, RWORK( NRHS+1 ), 00379 $ WORK, IWORK, INFO ) 00380 * 00381 * Check error code from DTPRFS. 00382 * 00383 IF( INFO.NE.0 ) 00384 $ CALL ALAERH( PATH, 'DTPRFS', INFO, 0, 00385 $ UPLO // TRANS // DIAG, N, N, -1, 00386 $ -1, NRHS, IMAT, NFAIL, NERRS, 00387 $ NOUT ) 00388 * 00389 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00390 $ RESULT( 4 ) ) 00391 CALL DTPT05( UPLO, TRANS, DIAG, N, NRHS, AP, B, 00392 $ LDA, X, LDA, XACT, LDA, RWORK, 00393 $ RWORK( NRHS+1 ), RESULT( 5 ) ) 00394 * 00395 * Print information about the tests that did not pass 00396 * the threshold. 00397 * 00398 DO 20 K = 2, 6 00399 IF( RESULT( K ).GE.THRESH ) THEN 00400 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00401 $ CALL ALAHD( NOUT, PATH ) 00402 WRITE( NOUT, FMT = 9998 )UPLO, TRANS, DIAG, 00403 $ N, NRHS, IMAT, K, RESULT( K ) 00404 NFAIL = NFAIL + 1 00405 END IF 00406 20 CONTINUE 00407 NRUN = NRUN + 5 00408 30 CONTINUE 00409 40 CONTINUE 00410 * 00411 *+ TEST 7 00412 * Get an estimate of RCOND = 1/CNDNUM. 00413 * 00414 DO 50 ITRAN = 1, 2 00415 IF( ITRAN.EQ.1 ) THEN 00416 NORM = 'O' 00417 RCONDC = RCONDO 00418 ELSE 00419 NORM = 'I' 00420 RCONDC = RCONDI 00421 END IF 00422 * 00423 SRNAMT = 'DTPCON' 00424 CALL DTPCON( NORM, UPLO, DIAG, N, AP, RCOND, WORK, 00425 $ IWORK, INFO ) 00426 * 00427 * Check error code from DTPCON. 00428 * 00429 IF( INFO.NE.0 ) 00430 $ CALL ALAERH( PATH, 'DTPCON', INFO, 0, 00431 $ NORM // UPLO // DIAG, N, N, -1, -1, 00432 $ -1, IMAT, NFAIL, NERRS, NOUT ) 00433 * 00434 CALL DTPT06( RCOND, RCONDC, UPLO, DIAG, N, AP, RWORK, 00435 $ RESULT( 7 ) ) 00436 * 00437 * Print the test ratio if it is .GE. THRESH. 00438 * 00439 IF( RESULT( 7 ).GE.THRESH ) THEN 00440 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00441 $ CALL ALAHD( NOUT, PATH ) 00442 WRITE( NOUT, FMT = 9997 ) 'DTPCON', NORM, UPLO, 00443 $ DIAG, N, IMAT, 7, RESULT( 7 ) 00444 NFAIL = NFAIL + 1 00445 END IF 00446 NRUN = NRUN + 1 00447 50 CONTINUE 00448 60 CONTINUE 00449 70 CONTINUE 00450 * 00451 * Use pathological test matrices to test DLATPS. 00452 * 00453 DO 100 IMAT = NTYPE1 + 1, NTYPES 00454 * 00455 * Do the tests only if DOTYPE( IMAT ) is true. 00456 * 00457 IF( .NOT.DOTYPE( IMAT ) ) 00458 $ GO TO 100 00459 * 00460 DO 90 IUPLO = 1, 2 00461 * 00462 * Do first for UPLO = 'U', then for UPLO = 'L' 00463 * 00464 UPLO = UPLOS( IUPLO ) 00465 DO 80 ITRAN = 1, NTRAN 00466 * 00467 * Do for op(A) = A, A**T, or A**H. 00468 * 00469 TRANS = TRANSS( ITRAN ) 00470 * 00471 * Call DLATTP to generate a triangular test matrix. 00472 * 00473 SRNAMT = 'DLATTP' 00474 CALL DLATTP( IMAT, UPLO, TRANS, DIAG, ISEED, N, AP, X, 00475 $ WORK, INFO ) 00476 * 00477 *+ TEST 8 00478 * Solve the system op(A)*x = b. 00479 * 00480 SRNAMT = 'DLATPS' 00481 CALL DCOPY( N, X, 1, B, 1 ) 00482 CALL DLATPS( UPLO, TRANS, DIAG, 'N', N, AP, B, SCALE, 00483 $ RWORK, INFO ) 00484 * 00485 * Check error code from DLATPS. 00486 * 00487 IF( INFO.NE.0 ) 00488 $ CALL ALAERH( PATH, 'DLATPS', INFO, 0, 00489 $ UPLO // TRANS // DIAG // 'N', N, N, 00490 $ -1, -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00491 * 00492 CALL DTPT03( UPLO, TRANS, DIAG, N, 1, AP, SCALE, 00493 $ RWORK, ONE, B, LDA, X, LDA, WORK, 00494 $ RESULT( 8 ) ) 00495 * 00496 *+ TEST 9 00497 * Solve op(A)*x = b again with NORMIN = 'Y'. 00498 * 00499 CALL DCOPY( N, X, 1, B( N+1 ), 1 ) 00500 CALL DLATPS( UPLO, TRANS, DIAG, 'Y', N, AP, B( N+1 ), 00501 $ SCALE, RWORK, INFO ) 00502 * 00503 * Check error code from DLATPS. 00504 * 00505 IF( INFO.NE.0 ) 00506 $ CALL ALAERH( PATH, 'DLATPS', INFO, 0, 00507 $ UPLO // TRANS // DIAG // 'Y', N, N, 00508 $ -1, -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00509 * 00510 CALL DTPT03( UPLO, TRANS, DIAG, N, 1, AP, SCALE, 00511 $ RWORK, ONE, B( N+1 ), LDA, X, LDA, WORK, 00512 $ RESULT( 9 ) ) 00513 * 00514 * Print information about the tests that did not pass 00515 * the threshold. 00516 * 00517 IF( RESULT( 8 ).GE.THRESH ) THEN 00518 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00519 $ CALL ALAHD( NOUT, PATH ) 00520 WRITE( NOUT, FMT = 9996 )'DLATPS', UPLO, TRANS, 00521 $ DIAG, 'N', N, IMAT, 8, RESULT( 8 ) 00522 NFAIL = NFAIL + 1 00523 END IF 00524 IF( RESULT( 9 ).GE.THRESH ) THEN 00525 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00526 $ CALL ALAHD( NOUT, PATH ) 00527 WRITE( NOUT, FMT = 9996 )'DLATPS', UPLO, TRANS, 00528 $ DIAG, 'Y', N, IMAT, 9, RESULT( 9 ) 00529 NFAIL = NFAIL + 1 00530 END IF 00531 NRUN = NRUN + 2 00532 80 CONTINUE 00533 90 CONTINUE 00534 100 CONTINUE 00535 110 CONTINUE 00536 * 00537 * Print a summary of the results. 00538 * 00539 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS ) 00540 * 00541 9999 FORMAT( ' UPLO=''', A1, ''', DIAG=''', A1, ''', N=', I5, 00542 $ ', type ', I2, ', test(', I2, ')= ', G12.5 ) 00543 9998 FORMAT( ' UPLO=''', A1, ''', TRANS=''', A1, ''', DIAG=''', A1, 00544 $ ''', N=', I5, ''', NRHS=', I5, ', type ', I2, ', test(', 00545 $ I2, ')= ', G12.5 ) 00546 9997 FORMAT( 1X, A, '( ''', A1, ''', ''', A1, ''', ''', A1, ''',', 00547 $ I5, ', ... ), type ', I2, ', test(', I2, ')=', G12.5 ) 00548 9996 FORMAT( 1X, A, '( ''', A1, ''', ''', A1, ''', ''', A1, ''', ''', 00549 $ A1, ''',', I5, ', ... ), type ', I2, ', test(', I2, ')=', 00550 $ G12.5 ) 00551 RETURN 00552 * 00553 * End of DCHKTP 00554 * 00555 END