 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE CDRVSP( DOTYPE, NN, NVAL, NRHS, THRESH, TSTERR, NMAX, 00002 $ A, AFAC, AINV, B, X, XACT, WORK, RWORK, IWORK, 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, NOUT, NRHS 00012 REAL THRESH 00013 * .. 00014 * .. Array Arguments .. 00015 LOGICAL DOTYPE( * ) 00016 INTEGER IWORK( * ), NVAL( * ) 00017 REAL RWORK( * ) 00018 COMPLEX A( * ), AFAC( * ), AINV( * ), B( * ), 00019 $ WORK( * ), X( * ), XACT( * ) 00020 * .. 00021 * 00022 * Purpose 00023 * ======= 00024 * 00025 * CDRVSP tests the driver routines CSPSV and -SVX. 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 * NRHS (input) INTEGER 00042 * The number of right hand side vectors to be generated for 00043 * each linear system. 00044 * 00045 * THRESH (input) REAL 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 * NMAX (input) INTEGER 00054 * The maximum value permitted for N, used in dimensioning the 00055 * work arrays. 00056 * 00057 * A (workspace) COMPLEX array, dimension 00058 * (NMAX*(NMAX+1)/2) 00059 * 00060 * AFAC (workspace) COMPLEX array, dimension 00061 * (NMAX*(NMAX+1)/2) 00062 * 00063 * AINV (workspace) COMPLEX array, dimension 00064 * (NMAX*(NMAX+1)/2) 00065 * 00066 * B (workspace) COMPLEX array, dimension (NMAX*NRHS) 00067 * 00068 * X (workspace) COMPLEX array, dimension (NMAX*NRHS) 00069 * 00070 * XACT (workspace) COMPLEX array, dimension (NMAX*NRHS) 00071 * 00072 * WORK (workspace) COMPLEX array, dimension 00073 * (NMAX*max(2,NRHS)) 00074 * 00075 * RWORK (workspace) REAL array, dimension (NMAX+2*NRHS) 00076 * 00077 * IWORK (workspace) INTEGER array, dimension (NMAX) 00078 * 00079 * NOUT (input) INTEGER 00080 * The unit number for output. 00081 * 00082 * ===================================================================== 00083 * 00084 * .. Parameters .. 00085 REAL ONE, ZERO 00086 PARAMETER ( ONE = 1.0E+0, ZERO = 0.0E+0 ) 00087 INTEGER NTYPES, NTESTS 00088 PARAMETER ( NTYPES = 11, NTESTS = 6 ) 00089 INTEGER NFACT 00090 PARAMETER ( NFACT = 2 ) 00091 * .. 00092 * .. Local Scalars .. 00093 LOGICAL ZEROT 00094 CHARACTER DIST, FACT, PACKIT, TYPE, UPLO, XTYPE 00095 CHARACTER*3 PATH 00096 INTEGER I, I1, I2, IFACT, IMAT, IN, INFO, IOFF, IUPLO, 00097 $ IZERO, J, K, K1, KL, KU, LDA, MODE, N, NB, 00098 $ NBMIN, NERRS, NFAIL, NIMAT, NPP, NRUN, NT 00099 REAL AINVNM, ANORM, CNDNUM, RCOND, RCONDC 00100 * .. 00101 * .. Local Arrays .. 00102 CHARACTER FACTS( NFACT ) 00103 INTEGER ISEED( 4 ), ISEEDY( 4 ) 00104 REAL RESULT( NTESTS ) 00105 * .. 00106 * .. External Functions .. 00107 REAL CLANSP, SGET06 00108 EXTERNAL CLANSP, SGET06 00109 * .. 00110 * .. External Subroutines .. 00111 EXTERNAL ALADHD, ALAERH, ALASVM, CCOPY, CERRVX, CGET04, 00112 $ CLACPY, CLARHS, CLASET, CLATB4, CLATMS, CLATSP, 00113 $ CPPT05, CSPSV, CSPSVX, CSPT01, CSPT02, CSPTRF, 00114 $ CSPTRI, XLAENV 00115 * .. 00116 * .. Scalars in Common .. 00117 LOGICAL LERR, OK 00118 CHARACTER*32 SRNAMT 00119 INTEGER INFOT, NUNIT 00120 * .. 00121 * .. Common blocks .. 00122 COMMON / INFOC / INFOT, NUNIT, OK, LERR 00123 COMMON / SRNAMC / SRNAMT 00124 * .. 00125 * .. Intrinsic Functions .. 00126 INTRINSIC CMPLX, MAX, MIN 00127 * .. 00128 * .. Data statements .. 00129 DATA ISEEDY / 1988, 1989, 1990, 1991 / 00130 DATA FACTS / 'F', 'N' / 00131 * .. 00132 * .. Executable Statements .. 00133 * 00134 * Initialize constants and the random number seed. 00135 * 00136 PATH( 1: 1 ) = 'Complex precision' 00137 PATH( 2: 3 ) = 'SP' 00138 NRUN = 0 00139 NFAIL = 0 00140 NERRS = 0 00141 DO 10 I = 1, 4 00142 ISEED( I ) = ISEEDY( I ) 00143 10 CONTINUE 00144 * 00145 * Test the error exits 00146 * 00147 IF( TSTERR ) 00148 $ CALL CERRVX( PATH, NOUT ) 00149 INFOT = 0 00150 * 00151 * Set the block size and minimum block size for testing. 00152 * 00153 NB = 1 00154 NBMIN = 2 00155 CALL XLAENV( 1, NB ) 00156 CALL XLAENV( 2, NBMIN ) 00157 * 00158 * Do for each value of N in NVAL 00159 * 00160 DO 180 IN = 1, NN 00161 N = NVAL( IN ) 00162 LDA = MAX( N, 1 ) 00163 NPP = N*( N+1 ) / 2 00164 XTYPE = 'N' 00165 NIMAT = NTYPES 00166 IF( N.LE.0 ) 00167 $ NIMAT = 1 00168 * 00169 DO 170 IMAT = 1, NIMAT 00170 * 00171 * Do the tests only if DOTYPE( IMAT ) is true. 00172 * 00173 IF( .NOT.DOTYPE( IMAT ) ) 00174 $ GO TO 170 00175 * 00176 * Skip types 3, 4, 5, or 6 if the matrix size is too small. 00177 * 00178 ZEROT = IMAT.GE.3 .AND. IMAT.LE.6 00179 IF( ZEROT .AND. N.LT.IMAT-2 ) 00180 $ GO TO 170 00181 * 00182 * Do first for UPLO = 'U', then for UPLO = 'L' 00183 * 00184 DO 160 IUPLO = 1, 2 00185 IF( IUPLO.EQ.1 ) THEN 00186 UPLO = 'U' 00187 PACKIT = 'C' 00188 ELSE 00189 UPLO = 'L' 00190 PACKIT = 'R' 00191 END IF 00192 * 00193 IF( IMAT.NE.NTYPES ) THEN 00194 * 00195 * Set up parameters with CLATB4 and generate a test 00196 * matrix with CLATMS. 00197 * 00198 CALL CLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, 00199 $ MODE, CNDNUM, DIST ) 00200 * 00201 SRNAMT = 'CLATMS' 00202 CALL CLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE, 00203 $ CNDNUM, ANORM, KL, KU, PACKIT, A, LDA, 00204 $ WORK, INFO ) 00205 * 00206 * Check error code from CLATMS. 00207 * 00208 IF( INFO.NE.0 ) THEN 00209 CALL ALAERH( PATH, 'CLATMS', INFO, 0, UPLO, N, N, 00210 $ -1, -1, -1, IMAT, NFAIL, NERRS, NOUT ) 00211 GO TO 160 00212 END IF 00213 * 00214 * For types 3-6, zero one or more rows and columns of 00215 * the matrix to test that INFO is returned correctly. 00216 * 00217 IF( ZEROT ) THEN 00218 IF( IMAT.EQ.3 ) THEN 00219 IZERO = 1 00220 ELSE IF( IMAT.EQ.4 ) THEN 00221 IZERO = N 00222 ELSE 00223 IZERO = N / 2 + 1 00224 END IF 00225 * 00226 IF( IMAT.LT.6 ) THEN 00227 * 00228 * Set row and column IZERO to zero. 00229 * 00230 IF( IUPLO.EQ.1 ) THEN 00231 IOFF = ( IZERO-1 )*IZERO / 2 00232 DO 20 I = 1, IZERO - 1 00233 A( IOFF+I ) = ZERO 00234 20 CONTINUE 00235 IOFF = IOFF + IZERO 00236 DO 30 I = IZERO, N 00237 A( IOFF ) = ZERO 00238 IOFF = IOFF + I 00239 30 CONTINUE 00240 ELSE 00241 IOFF = IZERO 00242 DO 40 I = 1, IZERO - 1 00243 A( IOFF ) = ZERO 00244 IOFF = IOFF + N - I 00245 40 CONTINUE 00246 IOFF = IOFF - IZERO 00247 DO 50 I = IZERO, N 00248 A( IOFF+I ) = ZERO 00249 50 CONTINUE 00250 END IF 00251 ELSE 00252 IF( IUPLO.EQ.1 ) THEN 00253 * 00254 * Set the first IZERO rows and columns to zero. 00255 * 00256 IOFF = 0 00257 DO 70 J = 1, N 00258 I2 = MIN( J, IZERO ) 00259 DO 60 I = 1, I2 00260 A( IOFF+I ) = ZERO 00261 60 CONTINUE 00262 IOFF = IOFF + J 00263 70 CONTINUE 00264 ELSE 00265 * 00266 * Set the last IZERO rows and columns to zero. 00267 * 00268 IOFF = 0 00269 DO 90 J = 1, N 00270 I1 = MAX( J, IZERO ) 00271 DO 80 I = I1, N 00272 A( IOFF+I ) = ZERO 00273 80 CONTINUE 00274 IOFF = IOFF + N - J 00275 90 CONTINUE 00276 END IF 00277 END IF 00278 ELSE 00279 IZERO = 0 00280 END IF 00281 ELSE 00282 * 00283 * Use a special block diagonal matrix to test alternate 00284 * code for the 2-by-2 blocks. 00285 * 00286 CALL CLATSP( UPLO, N, A, ISEED ) 00287 END IF 00288 * 00289 DO 150 IFACT = 1, NFACT 00290 * 00291 * Do first for FACT = 'F', then for other values. 00292 * 00293 FACT = FACTS( IFACT ) 00294 * 00295 * Compute the condition number for comparison with 00296 * the value returned by CSPSVX. 00297 * 00298 IF( ZEROT ) THEN 00299 IF( IFACT.EQ.1 ) 00300 $ GO TO 150 00301 RCONDC = ZERO 00302 * 00303 ELSE IF( IFACT.EQ.1 ) THEN 00304 * 00305 * Compute the 1-norm of A. 00306 * 00307 ANORM = CLANSP( '1', UPLO, N, A, RWORK ) 00308 * 00309 * Factor the matrix A. 00310 * 00311 CALL CCOPY( NPP, A, 1, AFAC, 1 ) 00312 CALL CSPTRF( UPLO, N, AFAC, IWORK, INFO ) 00313 * 00314 * Compute inv(A) and take its norm. 00315 * 00316 CALL CCOPY( NPP, AFAC, 1, AINV, 1 ) 00317 CALL CSPTRI( UPLO, N, AINV, IWORK, WORK, INFO ) 00318 AINVNM = CLANSP( '1', UPLO, N, AINV, RWORK ) 00319 * 00320 * Compute the 1-norm condition number of A. 00321 * 00322 IF( ANORM.LE.ZERO .OR. AINVNM.LE.ZERO ) THEN 00323 RCONDC = ONE 00324 ELSE 00325 RCONDC = ( ONE / ANORM ) / AINVNM 00326 END IF 00327 END IF 00328 * 00329 * Form an exact solution and set the right hand side. 00330 * 00331 SRNAMT = 'CLARHS' 00332 CALL CLARHS( PATH, XTYPE, UPLO, ' ', N, N, KL, KU, 00333 $ NRHS, A, LDA, XACT, LDA, B, LDA, ISEED, 00334 $ INFO ) 00335 XTYPE = 'C' 00336 * 00337 * --- Test CSPSV --- 00338 * 00339 IF( IFACT.EQ.2 ) THEN 00340 CALL CCOPY( NPP, A, 1, AFAC, 1 ) 00341 CALL CLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 00342 * 00343 * Factor the matrix and solve the system using CSPSV. 00344 * 00345 SRNAMT = 'CSPSV ' 00346 CALL CSPSV( UPLO, N, NRHS, AFAC, IWORK, X, LDA, 00347 $ INFO ) 00348 * 00349 * Adjust the expected value of INFO to account for 00350 * pivoting. 00351 * 00352 K = IZERO 00353 IF( K.GT.0 ) THEN 00354 100 CONTINUE 00355 IF( IWORK( K ).LT.0 ) THEN 00356 IF( IWORK( K ).NE.-K ) THEN 00357 K = -IWORK( K ) 00358 GO TO 100 00359 END IF 00360 ELSE IF( IWORK( K ).NE.K ) THEN 00361 K = IWORK( K ) 00362 GO TO 100 00363 END IF 00364 END IF 00365 * 00366 * Check error code from CSPSV . 00367 * 00368 IF( INFO.NE.K ) THEN 00369 CALL ALAERH( PATH, 'CSPSV ', INFO, K, UPLO, N, 00370 $ N, -1, -1, NRHS, IMAT, NFAIL, 00371 $ NERRS, NOUT ) 00372 GO TO 120 00373 ELSE IF( INFO.NE.0 ) THEN 00374 GO TO 120 00375 END IF 00376 * 00377 * Reconstruct matrix from factors and compute 00378 * residual. 00379 * 00380 CALL CSPT01( UPLO, N, A, AFAC, IWORK, AINV, LDA, 00381 $ RWORK, RESULT( 1 ) ) 00382 * 00383 * Compute residual of the computed solution. 00384 * 00385 CALL CLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA ) 00386 CALL CSPT02( UPLO, N, NRHS, A, X, LDA, WORK, LDA, 00387 $ RWORK, RESULT( 2 ) ) 00388 * 00389 * Check solution from generated exact solution. 00390 * 00391 CALL CGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00392 $ RESULT( 3 ) ) 00393 NT = 3 00394 * 00395 * Print information about the tests that did not pass 00396 * the threshold. 00397 * 00398 DO 110 K = 1, NT 00399 IF( RESULT( K ).GE.THRESH ) THEN 00400 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00401 $ CALL ALADHD( NOUT, PATH ) 00402 WRITE( NOUT, FMT = 9999 )'CSPSV ', UPLO, N, 00403 $ IMAT, K, RESULT( K ) 00404 NFAIL = NFAIL + 1 00405 END IF 00406 110 CONTINUE 00407 NRUN = NRUN + NT 00408 120 CONTINUE 00409 END IF 00410 * 00411 * --- Test CSPSVX --- 00412 * 00413 IF( IFACT.EQ.2 .AND. NPP.GT.0 ) 00414 $ CALL CLASET( 'Full', NPP, 1, CMPLX( ZERO ), 00415 $ CMPLX( ZERO ), AFAC, NPP ) 00416 CALL CLASET( 'Full', N, NRHS, CMPLX( ZERO ), 00417 $ CMPLX( ZERO ), X, LDA ) 00418 * 00419 * Solve the system and compute the condition number and 00420 * error bounds using CSPSVX. 00421 * 00422 SRNAMT = 'CSPSVX' 00423 CALL CSPSVX( FACT, UPLO, N, NRHS, A, AFAC, IWORK, B, 00424 $ LDA, X, LDA, RCOND, RWORK, 00425 $ RWORK( NRHS+1 ), WORK, RWORK( 2*NRHS+1 ), 00426 $ INFO ) 00427 * 00428 * Adjust the expected value of INFO to account for 00429 * pivoting. 00430 * 00431 K = IZERO 00432 IF( K.GT.0 ) THEN 00433 130 CONTINUE 00434 IF( IWORK( K ).LT.0 ) THEN 00435 IF( IWORK( K ).NE.-K ) THEN 00436 K = -IWORK( K ) 00437 GO TO 130 00438 END IF 00439 ELSE IF( IWORK( K ).NE.K ) THEN 00440 K = IWORK( K ) 00441 GO TO 130 00442 END IF 00443 END IF 00444 * 00445 * Check the error code from CSPSVX. 00446 * 00447 IF( INFO.NE.K ) THEN 00448 CALL ALAERH( PATH, 'CSPSVX', INFO, K, FACT // UPLO, 00449 $ N, N, -1, -1, NRHS, IMAT, NFAIL, 00450 $ NERRS, NOUT ) 00451 GO TO 150 00452 END IF 00453 * 00454 IF( INFO.EQ.0 ) THEN 00455 IF( IFACT.GE.2 ) THEN 00456 * 00457 * Reconstruct matrix from factors and compute 00458 * residual. 00459 * 00460 CALL CSPT01( UPLO, N, A, AFAC, IWORK, AINV, LDA, 00461 $ RWORK( 2*NRHS+1 ), RESULT( 1 ) ) 00462 K1 = 1 00463 ELSE 00464 K1 = 2 00465 END IF 00466 * 00467 * Compute residual of the computed solution. 00468 * 00469 CALL CLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA ) 00470 CALL CSPT02( UPLO, N, NRHS, A, X, LDA, WORK, LDA, 00471 $ RWORK( 2*NRHS+1 ), RESULT( 2 ) ) 00472 * 00473 * Check solution from generated exact solution. 00474 * 00475 CALL CGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 00476 $ RESULT( 3 ) ) 00477 * 00478 * Check the error bounds from iterative refinement. 00479 * 00480 CALL CPPT05( UPLO, N, NRHS, A, B, LDA, X, LDA, 00481 $ XACT, LDA, RWORK, RWORK( NRHS+1 ), 00482 $ RESULT( 4 ) ) 00483 ELSE 00484 K1 = 6 00485 END IF 00486 * 00487 * Compare RCOND from CSPSVX with the computed value 00488 * in RCONDC. 00489 * 00490 RESULT( 6 ) = SGET06( RCOND, RCONDC ) 00491 * 00492 * Print information about the tests that did not pass 00493 * the threshold. 00494 * 00495 DO 140 K = K1, 6 00496 IF( RESULT( K ).GE.THRESH ) THEN 00497 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 00498 $ CALL ALADHD( NOUT, PATH ) 00499 WRITE( NOUT, FMT = 9998 )'CSPSVX', FACT, UPLO, 00500 $ N, IMAT, K, RESULT( K ) 00501 NFAIL = NFAIL + 1 00502 END IF 00503 140 CONTINUE 00504 NRUN = NRUN + 7 - K1 00505 * 00506 150 CONTINUE 00507 * 00508 160 CONTINUE 00509 170 CONTINUE 00510 180 CONTINUE 00511 * 00512 * Print a summary of the results. 00513 * 00514 CALL ALASVM( PATH, NOUT, NFAIL, NRUN, NERRS ) 00515 * 00516 9999 FORMAT( 1X, A, ', UPLO=''', A1, ''', N =', I5, ', type ', I2, 00517 $ ', test ', I2, ', ratio =', G12.5 ) 00518 9998 FORMAT( 1X, A, ', FACT=''', A1, ''', UPLO=''', A1, ''', N =', I5, 00519 $ ', type ', I2, ', test ', I2, ', ratio =', G12.5 ) 00520 RETURN 00521 * 00522 * End of CDRVSP 00523 * 00524 END
1.7.3 
