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