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