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