 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE DCHKGT( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, 00002 $ A, AF, B, X, XACT, WORK, RWORK, IWORK, NOUT ) 00003 * 00004 * -- LAPACK test routine (version 3.1) -- 00005 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00006 * November 2006 00007 * 00008 * .. Scalar Arguments .. 00009 LOGICAL TSTERR 00010 INTEGER NN, NNS, NOUT 00011 DOUBLE PRECISION THRESH 00012 * .. 00013 * .. Array Arguments .. 00014 LOGICAL DOTYPE( * ) 00015 INTEGER IWORK( * ), NSVAL( * ), NVAL( * ) 00016 DOUBLE PRECISION A( * ), AF( * ), B( * ), RWORK( * ), WORK( * ), 00017 $ X( * ), XACT( * ) 00018 * .. 00019 * 00020 * Purpose 00021 * ======= 00022 * 00023 * DCHKGT tests DGTTRF, -TRS, -RFS, and -CON 00024 * 00025 * Arguments 00026 * ========= 00027 * 00028 * DOTYPE (input) LOGICAL array, dimension (NTYPES) 00029 * The matrix types to be used for testing. Matrices of type j 00030 * (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) = 00031 * .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used. 00032 * 00033 * NN (input) INTEGER 00034 * The number of values of N contained in the vector NVAL. 00035 * 00036 * NVAL (input) INTEGER array, dimension (NN) 00037 * The values of the matrix dimension N. 00038 * 00039 * NNS (input) INTEGER 00040 * The number of values of NRHS contained in the vector NSVAL. 00041 * 00042 * NSVAL (input) INTEGER array, dimension (NNS) 00043 * The values of the number of right hand sides NRHS. 00044 * 00045 * THRESH (input) DOUBLE PRECISION 00046 * The threshold value for the test ratios. A result is 00047 * included in the output file if RESULT >= THRESH. To have 00048 * every test ratio printed, use THRESH = 0. 00049 * 00050 * TSTERR (input) LOGICAL 00051 * Flag that indicates whether error exits are to be tested. 00052 * 00053 * A (workspace) DOUBLE PRECISION array, dimension (NMAX*4) 00054 * 00055 * AF (workspace) DOUBLE PRECISION array, dimension (NMAX*4) 00056 * 00057 * B (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00058 * where NSMAX is the largest entry in NSVAL. 00059 * 00060 * X (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00061 * 00062 * XACT (workspace) DOUBLE PRECISION array, dimension (NMAX*NSMAX) 00063 * 00064 * WORK (workspace) DOUBLE PRECISION array, dimension 00065 * (NMAX*max(3,NSMAX)) 00066 * 00067 * RWORK (workspace) DOUBLE PRECISION array, dimension 00068 * (max(NMAX,2*NSMAX)) 00069 * 00070 * IWORK (workspace) INTEGER array, dimension (2*NMAX) 00071 * 00072 * NOUT (input) INTEGER 00073 * The unit number for output. 00074 * 00075 * ===================================================================== 00076 * 00077 * .. Parameters .. 00078 DOUBLE PRECISION ONE, ZERO 00079 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 00080 INTEGER NTYPES 00081 PARAMETER ( NTYPES = 12 ) 00082 INTEGER NTESTS 00083 PARAMETER ( NTESTS = 7 ) 00084 * .. 00085 * .. Local Scalars .. 00086 LOGICAL TRFCON, ZEROT 00087 CHARACTER DIST, NORM, TRANS, TYPE 00088 CHARACTER*3 PATH 00089 INTEGER I, IMAT, IN, INFO, IRHS, ITRAN, IX, IZERO, J, 00090 $ K, KL, KOFF, KU, LDA, M, MODE, N, NERRS, NFAIL, 00091 $ NIMAT, NRHS, NRUN 00092 DOUBLE PRECISION AINVNM, ANORM, COND, RCOND, RCONDC, RCONDI, 00093 $ RCONDO 00094 * .. 00095 * .. Local Arrays .. 00096 CHARACTER TRANSS( 3 ) 00097 INTEGER ISEED( 4 ), ISEEDY( 4 ) 00098 DOUBLE PRECISION RESULT( NTESTS ), Z( 3 ) 00099 * .. 00100 * .. External Functions .. 00101 DOUBLE PRECISION DASUM, DGET06, DLANGT 00102 EXTERNAL DASUM, DGET06, DLANGT 00103 * .. 00104 * .. External Subroutines .. 00105 EXTERNAL ALAERH, ALAHD, ALASUM, DCOPY, DERRGE, DGET04, 00106 $ DGTCON, DGTRFS, DGTT01, DGTT02, DGTT05, DGTTRF, 00107 $ DGTTRS, DLACPY, DLAGTM, DLARNV, DLATB4, DLATMS, 00108 $ DSCAL 00109 * .. 00110 * .. Intrinsic Functions .. 00111 INTRINSIC MAX 00112 * .. 00113 * .. Scalars in Common .. 00114 LOGICAL LERR, OK 00115 CHARACTER*32 SRNAMT 00116 INTEGER INFOT, NUNIT 00117 * .. 00118 * .. Common blocks .. 00119 COMMON / INFOC / INFOT, NUNIT, OK, LERR 00120 COMMON / SRNAMC / SRNAMT 00121 * .. 00122 * .. Data statements .. 00123 DATA ISEEDY / 0, 0, 0, 1 / , TRANSS / 'N', 'T', 00124 $ 'C' / 00125 * .. 00126 * .. Executable Statements .. 00127 * 00128 PATH( 1: 1 ) = 'Double precision' 00129 PATH( 2: 3 ) = 'GT' 00130 NRUN = 0 00131 NFAIL = 0 00132 NERRS = 0 00133 DO 10 I = 1, 4 00134 ISEED( I ) = ISEEDY( I ) 00135 10 CONTINUE 00136 * 00137 * Test the error exits 00138 * 00139 IF( TSTERR ) 00140 $ CALL DERRGE( PATH, NOUT ) 00141 INFOT = 0 00142 * 00143 DO 110 IN = 1, NN 00144 * 00145 * Do for each value of N in NVAL. 00146 * 00147 N = NVAL( IN ) 00148 M = MAX( N-1, 0 ) 00149 LDA = MAX( 1, N ) 00150 NIMAT = NTYPES 00151 IF( N.LE.0 ) 00152 $ NIMAT = 1 00153 * 00154 DO 100 IMAT = 1, NIMAT 00155 * 00156 * Do the tests only if DOTYPE( IMAT ) is true. 00157 * 00158 IF( .NOT.DOTYPE( IMAT ) ) 00159 $ GO TO 100 00160 * 00161 * Set up parameters with DLATB4. 00162 * 00163 CALL DLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE, 00164 $ COND, DIST ) 00165 * 00166 ZEROT = IMAT.GE.8 .AND. IMAT.LE.10 00167 IF( IMAT.LE.6 ) THEN 00168 * 00169 * Types 1-6: generate matrices of known condition number. 00170 * 00171 KOFF = MAX( 2-KU, 3-MAX( 1, N ) ) 00172 SRNAMT = 'DLATMS' 00173 CALL DLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE, COND, 00174 $ ANORM, KL, KU, 'Z', AF( KOFF ), 3, WORK, 00175 $ INFO ) 00176 * 00177 * Check the error code from DLATMS. 00178 * 00179 IF( INFO.NE.0 ) THEN 00180 CALL ALAERH( PATH, 'DLATMS', INFO, 0, ' ', N, N, KL, 00181 $ KU, -1, IMAT, NFAIL, NERRS, NOUT ) 00182 GO TO 100 00183 END IF 00184 IZERO = 0 00185 * 00186 IF( N.GT.1 ) THEN 00187 CALL DCOPY( N-1, AF( 4 ), 3, A, 1 ) 00188 CALL DCOPY( N-1, AF( 3 ), 3, A( N+M+1 ), 1 ) 00189 END IF 00190 CALL DCOPY( N, AF( 2 ), 3, A( M+1 ), 1 ) 00191 ELSE 00192 * 00193 * Types 7-12: generate tridiagonal matrices with 00194 * unknown condition numbers. 00195 * 00196 IF( .NOT.ZEROT .OR. .NOT.DOTYPE( 7 ) ) THEN 00197 * 00198 * Generate a matrix with elements from [-1,1]. 00199 * 00200 CALL DLARNV( 2, ISEED, N+2*M, A ) 00201 IF( ANORM.NE.ONE ) 00202 $ CALL DSCAL( N+2*M, ANORM, A, 1 ) 00203 ELSE IF( IZERO.GT.0 ) THEN 00204 * 00205 * Reuse the last matrix by copying back the zeroed out 00206 * elements. 00207 * 00208 IF( IZERO.EQ.1 ) THEN 00209 A( N ) = Z( 2 ) 00210 IF( N.GT.1 ) 00211 $ A( 1 ) = Z( 3 ) 00212 ELSE IF( IZERO.EQ.N ) THEN 00213 A( 3*N-2 ) = Z( 1 ) 00214 A( 2*N-1 ) = Z( 2 ) 00215 ELSE 00216 A( 2*N-2+IZERO ) = Z( 1 ) 00217 A( N-1+IZERO ) = Z( 2 ) 00218 A( IZERO ) = Z( 3 ) 00219 END IF 00220 END IF 00221 * 00222 * If IMAT > 7, set one column of the matrix to 0. 00223 * 00224 IF( .NOT.ZEROT ) THEN 00225 IZERO = 0 00226 ELSE IF( IMAT.EQ.8 ) THEN 00227 IZERO = 1 00228 Z( 2 ) = A( N ) 00229 A( N ) = ZERO 00230 IF( N.GT.1 ) THEN 00231 Z( 3 ) = A( 1 ) 00232 A( 1 ) = ZERO 00233 END IF 00234 ELSE IF( IMAT.EQ.9 ) THEN 00235 IZERO = N 00236 Z( 1 ) = A( 3*N-2 ) 00237 Z( 2 ) = A( 2*N-1 ) 00238 A( 3*N-2 ) = ZERO 00239 A( 2*N-1 ) = ZERO 00240 ELSE 00241 IZERO = ( N+1 ) / 2 00242 DO 20 I = IZERO, N - 1 00243 A( 2*N-2+I ) = ZERO 00244 A( N-1+I ) = ZERO 00245 A( I ) = ZERO 00246 20 CONTINUE 00247 A( 3*N-2 ) = ZERO 00248 A( 2*N-1 ) = ZERO 00249 END IF 00250 END IF 00251 * 00252 *+ TEST 1 00253 * Factor A as L*U and compute the ratio 00254 * norm(L*U - A) / (n * norm(A) * EPS ) 00255 * 00256 CALL DCOPY( N+2*M, A, 1, AF, 1 ) 00257 SRNAMT = 'DGTTRF' 00258 CALL DGTTRF( N, AF, AF( M+1 ), AF( N+M+1 ), AF( N+2*M+1 ), 00259 $ IWORK, INFO ) 00260 * 00261 * Check error code from DGTTRF. 00262 * 00263 IF( INFO.NE.IZERO ) 00264 $ CALL ALAERH( PATH, 'DGTTRF', INFO, IZERO, ' ', N, N, 1, 00265 $ 1, -1, IMAT, NFAIL, NERRS, NOUT ) 00266 TRFCON = INFO.NE.0 00267 * 00268 CALL DGTT01( N, A, A( M+1 ), A( N+M+1 ), AF, AF( M+1 ), 00269 $ AF( N+M+1 ), AF( N+2*M+1 ), IWORK, WORK, LDA, 00270 $ RWORK, RESULT( 1 ) ) 00271 * 00272 * Print the test ratio if it is .GE. THRESH. 00273 * 00274 IF( RESULT( 1 ).GE.THRESH ) THEN 00275 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00276 $ CALL ALAHD( NOUT, PATH ) 00277 WRITE( NOUT, FMT = 9999 )N, IMAT, 1, RESULT( 1 ) 00278 NFAIL = NFAIL + 1 00279 END IF 00280 NRUN = NRUN + 1 00281 * 00282 DO 50 ITRAN = 1, 2 00283 TRANS = TRANSS( ITRAN ) 00284 IF( ITRAN.EQ.1 ) THEN 00285 NORM = 'O' 00286 ELSE 00287 NORM = 'I' 00288 END IF 00289 ANORM = DLANGT( NORM, N, A, A( M+1 ), A( N+M+1 ) ) 00290 * 00291 IF( .NOT.TRFCON ) THEN 00292 * 00293 * Use DGTTRS to solve for one column at a time of inv(A) 00294 * or inv(A^T), computing the maximum column sum as we 00295 * go. 00296 * 00297 AINVNM = ZERO 00298 DO 40 I = 1, N 00299 DO 30 J = 1, N 00300 X( J ) = ZERO 00301 30 CONTINUE 00302 X( I ) = ONE 00303 CALL DGTTRS( TRANS, N, 1, AF, AF( M+1 ), 00304 $ AF( N+M+1 ), AF( N+2*M+1 ), IWORK, X, 00305 $ LDA, INFO ) 00306 AINVNM = MAX( AINVNM, DASUM( N, X, 1 ) ) 00307 40 CONTINUE 00308 * 00309 * Compute RCONDC = 1 / (norm(A) * norm(inv(A)) 00310 * 00311 IF( ANORM.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN 00312 RCONDC = ONE 00313 ELSE 00314 RCONDC = ( ONE / ANORM ) / AINVNM 00315 END IF 00316 IF( ITRAN.EQ.1 ) THEN 00317 RCONDO = RCONDC 00318 ELSE 00319 RCONDI = RCONDC 00320 END IF 00321 ELSE 00322 RCONDC = ZERO 00323 END IF 00324 * 00325 *+ TEST 7 00326 * Estimate the reciprocal of the condition number of the 00327 * matrix. 00328 * 00329 SRNAMT = 'DGTCON' 00330 CALL DGTCON( NORM, N, AF, AF( M+1 ), AF( N+M+1 ), 00331 $ AF( N+2*M+1 ), IWORK, ANORM, RCOND, WORK, 00332 $ IWORK( N+1 ), INFO ) 00333 * 00334 * Check error code from DGTCON. 00335 * 00336 IF( INFO.NE.0 ) 00337 $ CALL ALAERH( PATH, 'DGTCON', INFO, 0, NORM, N, N, -1, 00338 $ -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00339 * 00340 RESULT( 7 ) = DGET06( RCOND, RCONDC ) 00341 * 00342 * Print the test ratio if it is .GE. THRESH. 00343 * 00344 IF( RESULT( 7 ).GE.THRESH ) THEN 00345 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00346 $ CALL ALAHD( NOUT, PATH ) 00347 WRITE( NOUT, FMT = 9997 )NORM, N, IMAT, 7, 00348 $ RESULT( 7 ) 00349 NFAIL = NFAIL + 1 00350 END IF 00351 NRUN = NRUN + 1 00352 50 CONTINUE 00353 * 00354 * Skip the remaining tests if the matrix is singular. 00355 * 00356 IF( TRFCON ) 00357 $ GO TO 100 00358 * 00359 DO 90 IRHS = 1, NNS 00360 NRHS = NSVAL( IRHS ) 00361 * 00362 * Generate NRHS random solution vectors. 00363 * 00364 IX = 1 00365 DO 60 J = 1, NRHS 00366 CALL DLARNV( 2, ISEED, N, XACT( IX ) ) 00367 IX = IX + LDA 00368 60 CONTINUE 00369 * 00370 DO 80 ITRAN = 1, 3 00371 TRANS = TRANSS( ITRAN ) 00372 IF( ITRAN.EQ.1 ) THEN 00373 RCONDC = RCONDO 00374 ELSE 00375 RCONDC = RCONDI 00376 END IF 00377 * 00378 * Set the right hand side. 00379 * 00380 CALL DLAGTM( TRANS, N, NRHS, ONE, A, A( M+1 ), 00381 $ A( N+M+1 ), XACT, LDA, ZERO, B, LDA ) 00382 * 00383 *+ TEST 2 00384 * Solve op(A) * X = B and compute the residual. 00385 * 00386 CALL DLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 00387 SRNAMT = 'DGTTRS' 00388 CALL DGTTRS( TRANS, N, NRHS, AF, AF( M+1 ), 00389 $ AF( N+M+1 ), AF( N+2*M+1 ), IWORK, X, 00390 $ LDA, INFO ) 00391 * 00392 * Check error code from DGTTRS. 00393 * 00394 IF( INFO.NE.0 ) 00395 $ CALL ALAERH( PATH, 'DGTTRS', INFO, 0, TRANS, N, N, 00396 $ -1, -1, NRHS, IMAT, NFAIL, NERRS, 00397 $ NOUT ) 00398 * 00399 CALL DLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA ) 00400 CALL DGTT02( TRANS, N, NRHS, A, A( M+1 ), A( N+M+1 ), 00401 $ X, LDA, WORK, LDA, RWORK, RESULT( 2 ) ) 00402 * 00403 *+ TEST 3 00404 * Check solution from generated exact solution. 00405 * 00406 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00407 $ RESULT( 3 ) ) 00408 * 00409 *+ TESTS 4, 5, and 6 00410 * Use iterative refinement to improve the solution. 00411 * 00412 SRNAMT = 'DGTRFS' 00413 CALL DGTRFS( TRANS, N, NRHS, A, A( M+1 ), A( N+M+1 ), 00414 $ AF, AF( M+1 ), AF( N+M+1 ), 00415 $ AF( N+2*M+1 ), IWORK, B, LDA, X, LDA, 00416 $ RWORK, RWORK( NRHS+1 ), WORK, 00417 $ IWORK( N+1 ), INFO ) 00418 * 00419 * Check error code from DGTRFS. 00420 * 00421 IF( INFO.NE.0 ) 00422 $ CALL ALAERH( PATH, 'DGTRFS', INFO, 0, TRANS, N, N, 00423 $ -1, -1, NRHS, IMAT, NFAIL, NERRS, 00424 $ NOUT ) 00425 * 00426 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00427 $ RESULT( 4 ) ) 00428 CALL DGTT05( TRANS, N, NRHS, A, A( M+1 ), A( N+M+1 ), 00429 $ B, LDA, X, LDA, XACT, LDA, RWORK, 00430 $ RWORK( NRHS+1 ), RESULT( 5 ) ) 00431 * 00432 * Print information about the tests that did not pass 00433 * the threshold. 00434 * 00435 DO 70 K = 2, 6 00436 IF( RESULT( K ).GE.THRESH ) THEN 00437 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00438 $ CALL ALAHD( NOUT, PATH ) 00439 WRITE( NOUT, FMT = 9998 )TRANS, N, NRHS, IMAT, 00440 $ K, RESULT( K ) 00441 NFAIL = NFAIL + 1 00442 END IF 00443 70 CONTINUE 00444 NRUN = NRUN + 5 00445 80 CONTINUE 00446 90 CONTINUE 00447 * 00448 100 CONTINUE 00449 110 CONTINUE 00450 * 00451 * Print a summary of the results. 00452 * 00453 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS ) 00454 * 00455 9999 FORMAT( 12X, 'N =', I5, ',', 10X, ' type ', I2, ', test(', I2, 00456 $ ') = ', G12.5 ) 00457 9998 FORMAT( ' TRANS=''', A1, ''', N =', I5, ', NRHS=', I3, ', type ', 00458 $ I2, ', test(', I2, ') = ', G12.5 ) 00459 9997 FORMAT( ' NORM =''', A1, ''', N =', I5, ',', 10X, ' type ', I2, 00460 $ ', test(', I2, ') = ', G12.5 ) 00461 RETURN 00462 * 00463 * End of DCHKGT 00464 * 00465 END
1.7.3 
