1*> \brief \b DCHKPT 2* 3* =========== DOCUMENTATION =========== 4* 5* Online html documentation available at 6* http://www.netlib.org/lapack/explore-html/ 7* 8* Definition: 9* =========== 10* 11* SUBROUTINE DCHKPT( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, 12* A, D, E, B, X, XACT, WORK, RWORK, NOUT ) 13* 14* .. Scalar Arguments .. 15* LOGICAL TSTERR 16* INTEGER NN, NNS, NOUT 17* DOUBLE PRECISION THRESH 18* .. 19* .. Array Arguments .. 20* LOGICAL DOTYPE( * ) 21* INTEGER NSVAL( * ), NVAL( * ) 22* DOUBLE PRECISION A( * ), B( * ), D( * ), E( * ), RWORK( * ), 23* $ WORK( * ), X( * ), XACT( * ) 24* .. 25* 26* 27*> \par Purpose: 28* ============= 29*> 30*> \verbatim 31*> 32*> DCHKPT tests DPTTRF, -TRS, -RFS, and -CON 33*> \endverbatim 34* 35* Arguments: 36* ========== 37* 38*> \param[in] DOTYPE 39*> \verbatim 40*> DOTYPE is LOGICAL array, dimension (NTYPES) 41*> The matrix types to be used for testing. Matrices of type j 42*> (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) = 43*> .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used. 44*> \endverbatim 45*> 46*> \param[in] NN 47*> \verbatim 48*> NN is INTEGER 49*> The number of values of N contained in the vector NVAL. 50*> \endverbatim 51*> 52*> \param[in] NVAL 53*> \verbatim 54*> NVAL is INTEGER array, dimension (NN) 55*> The values of the matrix dimension N. 56*> \endverbatim 57*> 58*> \param[in] NNS 59*> \verbatim 60*> NNS is INTEGER 61*> The number of values of NRHS contained in the vector NSVAL. 62*> \endverbatim 63*> 64*> \param[in] NSVAL 65*> \verbatim 66*> NSVAL is INTEGER array, dimension (NNS) 67*> The values of the number of right hand sides NRHS. 68*> \endverbatim 69*> 70*> \param[in] THRESH 71*> \verbatim 72*> THRESH is DOUBLE PRECISION 73*> The threshold value for the test ratios. A result is 74*> included in the output file if RESULT >= THRESH. To have 75*> every test ratio printed, use THRESH = 0. 76*> \endverbatim 77*> 78*> \param[in] TSTERR 79*> \verbatim 80*> TSTERR is LOGICAL 81*> Flag that indicates whether error exits are to be tested. 82*> \endverbatim 83*> 84*> \param[out] A 85*> \verbatim 86*> A is DOUBLE PRECISION array, dimension (NMAX*2) 87*> \endverbatim 88*> 89*> \param[out] D 90*> \verbatim 91*> D is DOUBLE PRECISION array, dimension (NMAX*2) 92*> \endverbatim 93*> 94*> \param[out] E 95*> \verbatim 96*> E is DOUBLE PRECISION array, dimension (NMAX*2) 97*> \endverbatim 98*> 99*> \param[out] B 100*> \verbatim 101*> B is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 102*> where NSMAX is the largest entry in NSVAL. 103*> \endverbatim 104*> 105*> \param[out] X 106*> \verbatim 107*> X is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 108*> \endverbatim 109*> 110*> \param[out] XACT 111*> \verbatim 112*> XACT is DOUBLE PRECISION array, dimension (NMAX*NSMAX) 113*> \endverbatim 114*> 115*> \param[out] WORK 116*> \verbatim 117*> WORK is DOUBLE PRECISION array, dimension 118*> (NMAX*max(3,NSMAX)) 119*> \endverbatim 120*> 121*> \param[out] RWORK 122*> \verbatim 123*> RWORK is DOUBLE PRECISION array, dimension 124*> (max(NMAX,2*NSMAX)) 125*> \endverbatim 126*> 127*> \param[in] NOUT 128*> \verbatim 129*> NOUT is INTEGER 130*> The unit number for output. 131*> \endverbatim 132* 133* Authors: 134* ======== 135* 136*> \author Univ. of Tennessee 137*> \author Univ. of California Berkeley 138*> \author Univ. of Colorado Denver 139*> \author NAG Ltd. 140* 141*> \ingroup double_lin 142* 143* ===================================================================== 144 SUBROUTINE DCHKPT( DOTYPE, NN, NVAL, NNS, NSVAL, THRESH, TSTERR, 145 $ A, D, E, B, X, XACT, WORK, RWORK, NOUT ) 146* 147* -- LAPACK test routine -- 148* -- LAPACK is a software package provided by Univ. of Tennessee, -- 149* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 150* 151* .. Scalar Arguments .. 152 LOGICAL TSTERR 153 INTEGER NN, NNS, NOUT 154 DOUBLE PRECISION THRESH 155* .. 156* .. Array Arguments .. 157 LOGICAL DOTYPE( * ) 158 INTEGER NSVAL( * ), NVAL( * ) 159 DOUBLE PRECISION A( * ), B( * ), D( * ), E( * ), RWORK( * ), 160 $ WORK( * ), X( * ), XACT( * ) 161* .. 162* 163* ===================================================================== 164* 165* .. Parameters .. 166 DOUBLE PRECISION ONE, ZERO 167 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 ) 168 INTEGER NTYPES 169 PARAMETER ( NTYPES = 12 ) 170 INTEGER NTESTS 171 PARAMETER ( NTESTS = 7 ) 172* .. 173* .. Local Scalars .. 174 LOGICAL ZEROT 175 CHARACTER DIST, TYPE 176 CHARACTER*3 PATH 177 INTEGER I, IA, IMAT, IN, INFO, IRHS, IX, IZERO, J, K, 178 $ KL, KU, LDA, MODE, N, NERRS, NFAIL, NIMAT, 179 $ NRHS, NRUN 180 DOUBLE PRECISION AINVNM, ANORM, COND, DMAX, RCOND, RCONDC 181* .. 182* .. Local Arrays .. 183 INTEGER ISEED( 4 ), ISEEDY( 4 ) 184 DOUBLE PRECISION RESULT( NTESTS ), Z( 3 ) 185* .. 186* .. External Functions .. 187 INTEGER IDAMAX 188 DOUBLE PRECISION DASUM, DGET06, DLANST 189 EXTERNAL IDAMAX, DASUM, DGET06, DLANST 190* .. 191* .. External Subroutines .. 192 EXTERNAL ALAERH, ALAHD, ALASUM, DCOPY, DERRGT, DGET04, 193 $ DLACPY, DLAPTM, DLARNV, DLATB4, DLATMS, DPTCON, 194 $ DPTRFS, DPTT01, DPTT02, DPTT05, DPTTRF, DPTTRS, 195 $ DSCAL 196* .. 197* .. Intrinsic Functions .. 198 INTRINSIC ABS, MAX 199* .. 200* .. Scalars in Common .. 201 LOGICAL LERR, OK 202 CHARACTER*32 SRNAMT 203 INTEGER INFOT, NUNIT 204* .. 205* .. Common blocks .. 206 COMMON / INFOC / INFOT, NUNIT, OK, LERR 207 COMMON / SRNAMC / SRNAMT 208* .. 209* .. Data statements .. 210 DATA ISEEDY / 0, 0, 0, 1 / 211* .. 212* .. Executable Statements .. 213* 214 PATH( 1: 1 ) = 'Double precision' 215 PATH( 2: 3 ) = 'PT' 216 NRUN = 0 217 NFAIL = 0 218 NERRS = 0 219 DO 10 I = 1, 4 220 ISEED( I ) = ISEEDY( I ) 221 10 CONTINUE 222* 223* Test the error exits 224* 225 IF( TSTERR ) 226 $ CALL DERRGT( PATH, NOUT ) 227 INFOT = 0 228* 229 DO 110 IN = 1, NN 230* 231* Do for each value of N in NVAL. 232* 233 N = NVAL( IN ) 234 LDA = MAX( 1, N ) 235 NIMAT = NTYPES 236 IF( N.LE.0 ) 237 $ NIMAT = 1 238* 239 DO 100 IMAT = 1, NIMAT 240* 241* Do the tests only if DOTYPE( IMAT ) is true. 242* 243 IF( N.GT.0 .AND. .NOT.DOTYPE( IMAT ) ) 244 $ GO TO 100 245* 246* Set up parameters with DLATB4. 247* 248 CALL DLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE, 249 $ COND, DIST ) 250* 251 ZEROT = IMAT.GE.8 .AND. IMAT.LE.10 252 IF( IMAT.LE.6 ) THEN 253* 254* Type 1-6: generate a symmetric tridiagonal matrix of 255* known condition number in lower triangular band storage. 256* 257 SRNAMT = 'DLATMS' 258 CALL DLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE, COND, 259 $ ANORM, KL, KU, 'B', A, 2, WORK, INFO ) 260* 261* Check the error code from DLATMS. 262* 263 IF( INFO.NE.0 ) THEN 264 CALL ALAERH( PATH, 'DLATMS', INFO, 0, ' ', N, N, KL, 265 $ KU, -1, IMAT, NFAIL, NERRS, NOUT ) 266 GO TO 100 267 END IF 268 IZERO = 0 269* 270* Copy the matrix to D and E. 271* 272 IA = 1 273 DO 20 I = 1, N - 1 274 D( I ) = A( IA ) 275 E( I ) = A( IA+1 ) 276 IA = IA + 2 277 20 CONTINUE 278 IF( N.GT.0 ) 279 $ D( N ) = A( IA ) 280 ELSE 281* 282* Type 7-12: generate a diagonally dominant matrix with 283* unknown condition number in the vectors D and E. 284* 285 IF( .NOT.ZEROT .OR. .NOT.DOTYPE( 7 ) ) THEN 286* 287* Let D and E have values from [-1,1]. 288* 289 CALL DLARNV( 2, ISEED, N, D ) 290 CALL DLARNV( 2, ISEED, N-1, E ) 291* 292* Make the tridiagonal matrix diagonally dominant. 293* 294 IF( N.EQ.1 ) THEN 295 D( 1 ) = ABS( D( 1 ) ) 296 ELSE 297 D( 1 ) = ABS( D( 1 ) ) + ABS( E( 1 ) ) 298 D( N ) = ABS( D( N ) ) + ABS( E( N-1 ) ) 299 DO 30 I = 2, N - 1 300 D( I ) = ABS( D( I ) ) + ABS( E( I ) ) + 301 $ ABS( E( I-1 ) ) 302 30 CONTINUE 303 END IF 304* 305* Scale D and E so the maximum element is ANORM. 306* 307 IX = IDAMAX( N, D, 1 ) 308 DMAX = D( IX ) 309 CALL DSCAL( N, ANORM / DMAX, D, 1 ) 310 CALL DSCAL( N-1, ANORM / DMAX, E, 1 ) 311* 312 ELSE IF( IZERO.GT.0 ) THEN 313* 314* Reuse the last matrix by copying back the zeroed out 315* elements. 316* 317 IF( IZERO.EQ.1 ) THEN 318 D( 1 ) = Z( 2 ) 319 IF( N.GT.1 ) 320 $ E( 1 ) = Z( 3 ) 321 ELSE IF( IZERO.EQ.N ) THEN 322 E( N-1 ) = Z( 1 ) 323 D( N ) = Z( 2 ) 324 ELSE 325 E( IZERO-1 ) = Z( 1 ) 326 D( IZERO ) = Z( 2 ) 327 E( IZERO ) = Z( 3 ) 328 END IF 329 END IF 330* 331* For types 8-10, set one row and column of the matrix to 332* zero. 333* 334 IZERO = 0 335 IF( IMAT.EQ.8 ) THEN 336 IZERO = 1 337 Z( 2 ) = D( 1 ) 338 D( 1 ) = ZERO 339 IF( N.GT.1 ) THEN 340 Z( 3 ) = E( 1 ) 341 E( 1 ) = ZERO 342 END IF 343 ELSE IF( IMAT.EQ.9 ) THEN 344 IZERO = N 345 IF( N.GT.1 ) THEN 346 Z( 1 ) = E( N-1 ) 347 E( N-1 ) = ZERO 348 END IF 349 Z( 2 ) = D( N ) 350 D( N ) = ZERO 351 ELSE IF( IMAT.EQ.10 ) THEN 352 IZERO = ( N+1 ) / 2 353 IF( IZERO.GT.1 ) THEN 354 Z( 1 ) = E( IZERO-1 ) 355 E( IZERO-1 ) = ZERO 356 Z( 3 ) = E( IZERO ) 357 E( IZERO ) = ZERO 358 END IF 359 Z( 2 ) = D( IZERO ) 360 D( IZERO ) = ZERO 361 END IF 362 END IF 363* 364 CALL DCOPY( N, D, 1, D( N+1 ), 1 ) 365 IF( N.GT.1 ) 366 $ CALL DCOPY( N-1, E, 1, E( N+1 ), 1 ) 367* 368*+ TEST 1 369* Factor A as L*D*L' and compute the ratio 370* norm(L*D*L' - A) / (n * norm(A) * EPS ) 371* 372 CALL DPTTRF( N, D( N+1 ), E( N+1 ), INFO ) 373* 374* Check error code from DPTTRF. 375* 376 IF( INFO.NE.IZERO ) THEN 377 CALL ALAERH( PATH, 'DPTTRF', INFO, IZERO, ' ', N, N, -1, 378 $ -1, -1, IMAT, NFAIL, NERRS, NOUT ) 379 GO TO 100 380 END IF 381* 382 IF( INFO.GT.0 ) THEN 383 RCONDC = ZERO 384 GO TO 90 385 END IF 386* 387 CALL DPTT01( N, D, E, D( N+1 ), E( N+1 ), WORK, 388 $ RESULT( 1 ) ) 389* 390* Print the test ratio if greater than or equal to THRESH. 391* 392 IF( RESULT( 1 ).GE.THRESH ) THEN 393 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 394 $ CALL ALAHD( NOUT, PATH ) 395 WRITE( NOUT, FMT = 9999 )N, IMAT, 1, RESULT( 1 ) 396 NFAIL = NFAIL + 1 397 END IF 398 NRUN = NRUN + 1 399* 400* Compute RCONDC = 1 / (norm(A) * norm(inv(A)) 401* 402* Compute norm(A). 403* 404 ANORM = DLANST( '1', N, D, E ) 405* 406* Use DPTTRS to solve for one column at a time of inv(A), 407* computing the maximum column sum as we go. 408* 409 AINVNM = ZERO 410 DO 50 I = 1, N 411 DO 40 J = 1, N 412 X( J ) = ZERO 413 40 CONTINUE 414 X( I ) = ONE 415 CALL DPTTRS( N, 1, D( N+1 ), E( N+1 ), X, LDA, INFO ) 416 AINVNM = MAX( AINVNM, DASUM( N, X, 1 ) ) 417 50 CONTINUE 418 RCONDC = ONE / MAX( ONE, ANORM*AINVNM ) 419* 420 DO 80 IRHS = 1, NNS 421 NRHS = NSVAL( IRHS ) 422* 423* Generate NRHS random solution vectors. 424* 425 IX = 1 426 DO 60 J = 1, NRHS 427 CALL DLARNV( 2, ISEED, N, XACT( IX ) ) 428 IX = IX + LDA 429 60 CONTINUE 430* 431* Set the right hand side. 432* 433 CALL DLAPTM( N, NRHS, ONE, D, E, XACT, LDA, ZERO, B, 434 $ LDA ) 435* 436*+ TEST 2 437* Solve A*x = b and compute the residual. 438* 439 CALL DLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 440 CALL DPTTRS( N, NRHS, D( N+1 ), E( N+1 ), X, LDA, INFO ) 441* 442* Check error code from DPTTRS. 443* 444 IF( INFO.NE.0 ) 445 $ CALL ALAERH( PATH, 'DPTTRS', INFO, 0, ' ', N, N, -1, 446 $ -1, NRHS, IMAT, NFAIL, NERRS, NOUT ) 447* 448 CALL DLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA ) 449 CALL DPTT02( N, NRHS, D, E, X, LDA, WORK, LDA, 450 $ RESULT( 2 ) ) 451* 452*+ TEST 3 453* Check solution from generated exact solution. 454* 455 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 456 $ RESULT( 3 ) ) 457* 458*+ TESTS 4, 5, and 6 459* Use iterative refinement to improve the solution. 460* 461 SRNAMT = 'DPTRFS' 462 CALL DPTRFS( N, NRHS, D, E, D( N+1 ), E( N+1 ), B, LDA, 463 $ X, LDA, RWORK, RWORK( NRHS+1 ), WORK, INFO ) 464* 465* Check error code from DPTRFS. 466* 467 IF( INFO.NE.0 ) 468 $ CALL ALAERH( PATH, 'DPTRFS', INFO, 0, ' ', N, N, -1, 469 $ -1, NRHS, IMAT, NFAIL, NERRS, NOUT ) 470* 471 CALL DGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 472 $ RESULT( 4 ) ) 473 CALL DPTT05( N, NRHS, D, E, B, LDA, X, LDA, XACT, LDA, 474 $ RWORK, RWORK( NRHS+1 ), RESULT( 5 ) ) 475* 476* Print information about the tests that did not pass the 477* threshold. 478* 479 DO 70 K = 2, 6 480 IF( RESULT( K ).GE.THRESH ) THEN 481 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 482 $ CALL ALAHD( NOUT, PATH ) 483 WRITE( NOUT, FMT = 9998 )N, NRHS, IMAT, K, 484 $ RESULT( K ) 485 NFAIL = NFAIL + 1 486 END IF 487 70 CONTINUE 488 NRUN = NRUN + 5 489 80 CONTINUE 490* 491*+ TEST 7 492* Estimate the reciprocal of the condition number of the 493* matrix. 494* 495 90 CONTINUE 496 SRNAMT = 'DPTCON' 497 CALL DPTCON( N, D( N+1 ), E( N+1 ), ANORM, RCOND, RWORK, 498 $ INFO ) 499* 500* Check error code from DPTCON. 501* 502 IF( INFO.NE.0 ) 503 $ CALL ALAERH( PATH, 'DPTCON', INFO, 0, ' ', N, N, -1, -1, 504 $ -1, IMAT, NFAIL, NERRS, NOUT ) 505* 506 RESULT( 7 ) = DGET06( RCOND, RCONDC ) 507* 508* Print the test ratio if greater than or equal to THRESH. 509* 510 IF( RESULT( 7 ).GE.THRESH ) THEN 511 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 512 $ CALL ALAHD( NOUT, PATH ) 513 WRITE( NOUT, FMT = 9999 )N, IMAT, 7, RESULT( 7 ) 514 NFAIL = NFAIL + 1 515 END IF 516 NRUN = NRUN + 1 517 100 CONTINUE 518 110 CONTINUE 519* 520* Print a summary of the results. 521* 522 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS ) 523* 524 9999 FORMAT( ' N =', I5, ', type ', I2, ', test ', I2, ', ratio = ', 525 $ G12.5 ) 526 9998 FORMAT( ' N =', I5, ', NRHS=', I3, ', type ', I2, ', test(', I2, 527 $ ') = ', G12.5 ) 528 RETURN 529* 530* End of DCHKPT 531* 532 END 533