1*> \brief \b SCHKSY 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 SCHKSY( DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, 12* THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, 13* XACT, WORK, RWORK, IWORK, NOUT ) 14* 15* .. Scalar Arguments .. 16* LOGICAL TSTERR 17* INTEGER NMAX, NN, NNB, NNS, NOUT 18* REAL THRESH 19* .. 20* .. Array Arguments .. 21* LOGICAL DOTYPE( * ) 22* INTEGER IWORK( * ), NBVAL( * ), NSVAL( * ), NVAL( * ) 23* REAL A( * ), AFAC( * ), AINV( * ), B( * ), 24* $ RWORK( * ), WORK( * ), X( * ), XACT( * ) 25* .. 26* 27* 28*> \par Purpose: 29* ============= 30*> 31*> \verbatim 32*> 33*> SCHKSY tests SSYTRF, -TRI2, -TRS, -TRS2, -RFS, and -CON. 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] NNB 60*> \verbatim 61*> NNB is INTEGER 62*> The number of values of NB contained in the vector NBVAL. 63*> \endverbatim 64*> 65*> \param[in] NBVAL 66*> \verbatim 67*> NBVAL is INTEGER array, dimension (NBVAL) 68*> The values of the blocksize NB. 69*> \endverbatim 70*> 71*> \param[in] NNS 72*> \verbatim 73*> NNS is INTEGER 74*> The number of values of NRHS contained in the vector NSVAL. 75*> \endverbatim 76*> 77*> \param[in] NSVAL 78*> \verbatim 79*> NSVAL is INTEGER array, dimension (NNS) 80*> The values of the number of right hand sides NRHS. 81*> \endverbatim 82*> 83*> \param[in] THRESH 84*> \verbatim 85*> THRESH is REAL 86*> The threshold value for the test ratios. A result is 87*> included in the output file if RESULT >= THRESH. To have 88*> every test ratio printed, use THRESH = 0. 89*> \endverbatim 90*> 91*> \param[in] TSTERR 92*> \verbatim 93*> TSTERR is LOGICAL 94*> Flag that indicates whether error exits are to be tested. 95*> \endverbatim 96*> 97*> \param[in] NMAX 98*> \verbatim 99*> NMAX is INTEGER 100*> The maximum value permitted for N, used in dimensioning the 101*> work arrays. 102*> \endverbatim 103*> 104*> \param[out] A 105*> \verbatim 106*> A is REAL array, dimension (NMAX*NMAX) 107*> \endverbatim 108*> 109*> \param[out] AFAC 110*> \verbatim 111*> AFAC is REAL array, dimension (NMAX*NMAX) 112*> \endverbatim 113*> 114*> \param[out] AINV 115*> \verbatim 116*> AINV is REAL array, dimension (NMAX*NMAX) 117*> \endverbatim 118*> 119*> \param[out] B 120*> \verbatim 121*> B is REAL array, dimension (NMAX*NSMAX) 122*> where NSMAX is the largest entry in NSVAL. 123*> \endverbatim 124*> 125*> \param[out] X 126*> \verbatim 127*> X is REAL array, dimension (NMAX*NSMAX) 128*> \endverbatim 129*> 130*> \param[out] XACT 131*> \verbatim 132*> XACT is REAL array, dimension (NMAX*NSMAX) 133*> \endverbatim 134*> 135*> \param[out] WORK 136*> \verbatim 137*> WORK is REAL array, dimension 138*> (NMAX*max(3,NSMAX)) 139*> \endverbatim 140*> 141*> \param[out] RWORK 142*> \verbatim 143*> RWORK is REAL array, dimension 144*> (max(NMAX,2*NSMAX)) 145*> \endverbatim 146*> 147*> \param[out] IWORK 148*> \verbatim 149*> IWORK is INTEGER array, dimension (2*NMAX) 150*> \endverbatim 151*> 152*> \param[in] NOUT 153*> \verbatim 154*> NOUT is INTEGER 155*> The unit number for output. 156*> \endverbatim 157* 158* Authors: 159* ======== 160* 161*> \author Univ. of Tennessee 162*> \author Univ. of California Berkeley 163*> \author Univ. of Colorado Denver 164*> \author NAG Ltd. 165* 166*> \date April 2012 167* 168*> \ingroup single_lin 169* 170* ===================================================================== 171 SUBROUTINE SCHKSY( DOTYPE, NN, NVAL, NNB, NBVAL, NNS, NSVAL, 172 $ THRESH, TSTERR, NMAX, A, AFAC, AINV, B, X, 173 $ XACT, WORK, RWORK, IWORK, NOUT ) 174* 175* -- LAPACK test routine (version 3.4.1) -- 176* -- LAPACK is a software package provided by Univ. of Tennessee, -- 177* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 178* April 2012 179* 180* .. Scalar Arguments .. 181 LOGICAL TSTERR 182 INTEGER NMAX, NN, NNB, NNS, NOUT 183 REAL THRESH 184* .. 185* .. Array Arguments .. 186 LOGICAL DOTYPE( * ) 187 INTEGER IWORK( * ), NBVAL( * ), NSVAL( * ), NVAL( * ) 188 REAL A( * ), AFAC( * ), AINV( * ), B( * ), 189 $ RWORK( * ), WORK( * ), X( * ), XACT( * ) 190* .. 191* 192* ===================================================================== 193* 194* .. Parameters .. 195 REAL ZERO 196 PARAMETER ( ZERO = 0.0E+0 ) 197 INTEGER NTYPES 198 PARAMETER ( NTYPES = 10 ) 199 INTEGER NTESTS 200 PARAMETER ( NTESTS = 9 ) 201* .. 202* .. Local Scalars .. 203 LOGICAL TRFCON, ZEROT 204 CHARACTER DIST, TYPE, UPLO, XTYPE 205 CHARACTER*3 PATH 206 INTEGER I, I1, I2, IMAT, IN, INB, INFO, IOFF, IRHS, 207 $ IUPLO, IZERO, J, K, KL, KU, LDA, LWORK, MODE, 208 $ N, NB, NERRS, NFAIL, NIMAT, NRHS, NRUN, NT 209 REAL ANORM, CNDNUM, RCOND, RCONDC 210* .. 211* .. Local Arrays .. 212 CHARACTER UPLOS( 2 ) 213 INTEGER ISEED( 4 ), ISEEDY( 4 ) 214 REAL RESULT( NTESTS ) 215* .. 216* .. External Functions .. 217 REAL SGET06, SLANSY 218 EXTERNAL SGET06, SLANSY 219* .. 220* .. External Subroutines .. 221 EXTERNAL ALAERH, ALAHD, ALASUM, SERRSY, SGET04, SLACPY, 222 $ SLARHS, SLATB4, SLATMS, SPOT02, SPOT03, SPOT05, 223 $ SSYCON, SSYCONV, SSYRFS, SSYT01, SSYTRF, 224 $ SSYTRI2, SSYTRS, SSYTRS2, XLAENV 225* .. 226* .. Intrinsic Functions .. 227 INTRINSIC MAX, MIN 228* .. 229* .. Scalars in Common .. 230 LOGICAL LERR, OK 231 CHARACTER*32 SRNAMT 232 INTEGER INFOT, NUNIT 233* .. 234* .. Common blocks .. 235 COMMON / INFOC / INFOT, NUNIT, OK, LERR 236 COMMON / SRNAMC / SRNAMT 237* .. 238* .. Data statements .. 239 DATA ISEEDY / 1988, 1989, 1990, 1991 / 240 DATA UPLOS / 'U', 'L' / 241* .. 242* .. Executable Statements .. 243* 244* Initialize constants and the random number seed. 245* 246 PATH( 1: 1 ) = 'Single precision' 247 PATH( 2: 3 ) = 'SY' 248 NRUN = 0 249 NFAIL = 0 250 NERRS = 0 251 DO 10 I = 1, 4 252 ISEED( I ) = ISEEDY( I ) 253 10 CONTINUE 254* 255* Test the error exits 256* 257 IF( TSTERR ) 258 $ CALL SERRSY( PATH, NOUT ) 259 INFOT = 0 260* 261* Set the minimum block size for which the block routine should 262* be used, which will be later returned by ILAENV 263* 264 CALL XLAENV( 2, 2 ) 265* 266* Do for each value of N in NVAL 267* 268 DO 180 IN = 1, NN 269 N = NVAL( IN ) 270 LDA = MAX( N, 1 ) 271 XTYPE = 'N' 272 NIMAT = NTYPES 273 IF( N.LE.0 ) 274 $ NIMAT = 1 275* 276 IZERO = 0 277* 278* Do for each value of matrix type IMAT 279* 280 DO 170 IMAT = 1, NIMAT 281* 282* Do the tests only if DOTYPE( IMAT ) is true. 283* 284 IF( .NOT.DOTYPE( IMAT ) ) 285 $ GO TO 170 286* 287* Skip types 3, 4, 5, or 6 if the matrix size is too small. 288* 289 ZEROT = IMAT.GE.3 .AND. IMAT.LE.6 290 IF( ZEROT .AND. N.LT.IMAT-2 ) 291 $ GO TO 170 292* 293* Do first for UPLO = 'U', then for UPLO = 'L' 294* 295 DO 160 IUPLO = 1, 2 296 UPLO = UPLOS( IUPLO ) 297* 298* Begin generate the test matrix A. 299* 300* Set up parameters with SLATB4 for the matrix generator 301* based on the type of matrix to be generated. 302* 303 CALL SLATB4( PATH, IMAT, N, N, TYPE, KL, KU, ANORM, MODE, 304 $ CNDNUM, DIST ) 305* 306* Generate a matrix with SLATMS. 307* 308 SRNAMT = 'SLATMS' 309 CALL SLATMS( N, N, DIST, ISEED, TYPE, RWORK, MODE, 310 $ CNDNUM, ANORM, KL, KU, UPLO, A, LDA, WORK, 311 $ INFO ) 312* 313* Check error code from SLATMS and handle error. 314* 315 IF( INFO.NE.0 ) THEN 316 CALL ALAERH( PATH, 'SLATMS', INFO, 0, UPLO, N, N, -1, 317 $ -1, -1, IMAT, NFAIL, NERRS, NOUT ) 318 GO TO 160 319 END IF 320* 321* For matrix types 3-6, zero one or more rows and 322* columns of the matrix to test that INFO is returned 323* correctly. 324* 325 IF( ZEROT ) THEN 326 IF( IMAT.EQ.3 ) THEN 327 IZERO = 1 328 ELSE IF( IMAT.EQ.4 ) THEN 329 IZERO = N 330 ELSE 331 IZERO = N / 2 + 1 332 END IF 333* 334 IF( IMAT.LT.6 ) THEN 335* 336* Set row and column IZERO to zero. 337* 338 IF( IUPLO.EQ.1 ) THEN 339 IOFF = ( IZERO-1 )*LDA 340 DO 20 I = 1, IZERO - 1 341 A( IOFF+I ) = ZERO 342 20 CONTINUE 343 IOFF = IOFF + IZERO 344 DO 30 I = IZERO, N 345 A( IOFF ) = ZERO 346 IOFF = IOFF + LDA 347 30 CONTINUE 348 ELSE 349 IOFF = IZERO 350 DO 40 I = 1, IZERO - 1 351 A( IOFF ) = ZERO 352 IOFF = IOFF + LDA 353 40 CONTINUE 354 IOFF = IOFF - IZERO 355 DO 50 I = IZERO, N 356 A( IOFF+I ) = ZERO 357 50 CONTINUE 358 END IF 359 ELSE 360 IOFF = 0 361 IF( IUPLO.EQ.1 ) THEN 362* 363* Set the first IZERO rows and columns to zero. 364* 365 DO 70 J = 1, N 366 I2 = MIN( J, IZERO ) 367 DO 60 I = 1, I2 368 A( IOFF+I ) = ZERO 369 60 CONTINUE 370 IOFF = IOFF + LDA 371 70 CONTINUE 372 ELSE 373* 374* Set the last IZERO rows and columns to zero. 375* 376 DO 90 J = 1, N 377 I1 = MAX( J, IZERO ) 378 DO 80 I = I1, N 379 A( IOFF+I ) = ZERO 380 80 CONTINUE 381 IOFF = IOFF + LDA 382 90 CONTINUE 383 END IF 384 END IF 385 ELSE 386 IZERO = 0 387 END IF 388* 389* End generate the test matrix A. 390* 391* Do for each value of NB in NBVAL 392* 393 DO 150 INB = 1, NNB 394* 395* Set the optimal blocksize, which will be later 396* returned by ILAENV. 397* 398 NB = NBVAL( INB ) 399 CALL XLAENV( 1, NB ) 400* 401* Copy the test matrix A into matrix AFAC which 402* will be factorized in place. This is needed to 403* preserve the test matrix A for subsequent tests. 404* 405 CALL SLACPY( UPLO, N, N, A, LDA, AFAC, LDA ) 406* 407* Compute the L*D*L**T or U*D*U**T factorization of the 408* matrix. IWORK stores details of the interchanges and 409* the block structure of D. AINV is a work array for 410* block factorization, LWORK is the length of AINV. 411* 412 LWORK = MAX( 2, NB )*LDA 413 SRNAMT = 'SSYTRF' 414 CALL SSYTRF( UPLO, N, AFAC, LDA, IWORK, AINV, LWORK, 415 $ INFO ) 416* 417* Adjust the expected value of INFO to account for 418* pivoting. 419* 420 K = IZERO 421 IF( K.GT.0 ) THEN 422 100 CONTINUE 423 IF( IWORK( K ).LT.0 ) THEN 424 IF( IWORK( K ).NE.-K ) THEN 425 K = -IWORK( K ) 426 GO TO 100 427 END IF 428 ELSE IF( IWORK( K ).NE.K ) THEN 429 K = IWORK( K ) 430 GO TO 100 431 END IF 432 END IF 433* 434* Check error code from SSYTRF and handle error. 435* 436 IF( INFO.NE.K ) 437 $ CALL ALAERH( PATH, 'SSYTRF', INFO, K, UPLO, N, N, 438 $ -1, -1, NB, IMAT, NFAIL, NERRS, NOUT ) 439* 440* Set the condition estimate flag if the INFO is not 0. 441* 442 IF( INFO.NE.0 ) THEN 443 TRFCON = .TRUE. 444 ELSE 445 TRFCON = .FALSE. 446 END IF 447* 448*+ TEST 1 449* Reconstruct matrix from factors and compute residual. 450* 451 CALL SSYT01( UPLO, N, A, LDA, AFAC, LDA, IWORK, AINV, 452 $ LDA, RWORK, RESULT( 1 ) ) 453 NT = 1 454* 455*+ TEST 2 456* Form the inverse and compute the residual, 457* if the factorization was competed without INFO > 0 458* (i.e. there is no zero rows and columns). 459* Do it only for the first block size. 460* 461 IF( INB.EQ.1 .AND. .NOT.TRFCON ) THEN 462 CALL SLACPY( UPLO, N, N, AFAC, LDA, AINV, LDA ) 463 SRNAMT = 'SSYTRI2' 464 LWORK = (N+NB+1)*(NB+3) 465 CALL SSYTRI2( UPLO, N, AINV, LDA, IWORK, WORK, 466 $ LWORK, INFO ) 467* 468* Check error code from SSYTRI2 and handle error. 469* 470 IF( INFO.NE.0 ) 471 $ CALL ALAERH( PATH, 'SSYTRI2', INFO, -1, UPLO, N, 472 $ N, -1, -1, -1, IMAT, NFAIL, NERRS, 473 $ NOUT ) 474* 475* Compute the residual for a symmetric matrix times 476* its inverse. 477* 478 CALL SPOT03( UPLO, N, A, LDA, AINV, LDA, WORK, LDA, 479 $ RWORK, RCONDC, RESULT( 2 ) ) 480 NT = 2 481 END IF 482* 483* Print information about the tests that did not pass 484* the threshold. 485* 486 DO 110 K = 1, NT 487 IF( RESULT( K ).GE.THRESH ) THEN 488 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 489 $ CALL ALAHD( NOUT, PATH ) 490 WRITE( NOUT, FMT = 9999 )UPLO, N, NB, IMAT, K, 491 $ RESULT( K ) 492 NFAIL = NFAIL + 1 493 END IF 494 110 CONTINUE 495 NRUN = NRUN + NT 496* 497* Skip the other tests if this is not the first block 498* size. 499* 500 IF( INB.GT.1 ) 501 $ GO TO 150 502* 503* Do only the condition estimate if INFO is not 0. 504* 505 IF( TRFCON ) THEN 506 RCONDC = ZERO 507 GO TO 140 508 END IF 509* 510 DO 130 IRHS = 1, NNS 511 NRHS = NSVAL( IRHS ) 512* 513*+ TEST 3 (Using DSYTRS) 514* Solve and compute residual for A * X = B. 515* 516* Choose a set of NRHS random solution vectors 517* stored in XACT and set up the right hand side B 518* 519 SRNAMT = 'SLARHS' 520 CALL SLARHS( PATH, XTYPE, UPLO, ' ', N, N, KL, KU, 521 $ NRHS, A, LDA, XACT, LDA, B, LDA, 522 $ ISEED, INFO ) 523 CALL SLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 524* 525 SRNAMT = 'SSYTRS' 526 CALL SSYTRS( UPLO, N, NRHS, AFAC, LDA, IWORK, X, 527 $ LDA, INFO ) 528* 529* Check error code from SSYTRS and handle error. 530* 531 IF( INFO.NE.0 ) 532 $ CALL ALAERH( PATH, 'SSYTRS', INFO, 0, UPLO, N, 533 $ N, -1, -1, NRHS, IMAT, NFAIL, 534 $ NERRS, NOUT ) 535* 536 CALL SLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA ) 537* 538* Compute the residual for the solution 539* 540 CALL SPOT02( UPLO, N, NRHS, A, LDA, X, LDA, WORK, 541 $ LDA, RWORK, RESULT( 3 ) ) 542* 543*+ TEST 4 (Using DSYTRS2) 544* Solve and compute residual for A * X = B. 545* 546* Choose a set of NRHS random solution vectors 547* stored in XACT and set up the right hand side B 548* 549 SRNAMT = 'SLARHS' 550 CALL SLARHS( PATH, XTYPE, UPLO, ' ', N, N, KL, KU, 551 $ NRHS, A, LDA, XACT, LDA, B, LDA, 552 $ ISEED, INFO ) 553 CALL SLACPY( 'Full', N, NRHS, B, LDA, X, LDA ) 554* 555 SRNAMT = 'DSYTRS2' 556 CALL SSYTRS2( UPLO, N, NRHS, AFAC, LDA, IWORK, X, 557 $ LDA, WORK, INFO ) 558* 559* Check error code from SSYTRS2 and handle error. 560* 561 IF( INFO.NE.0 ) 562 $ CALL ALAERH( PATH, 'SSYTRS2', INFO, 0, UPLO, N, 563 $ N, -1, -1, NRHS, IMAT, NFAIL, 564 $ NERRS, NOUT ) 565* 566 CALL SLACPY( 'Full', N, NRHS, B, LDA, WORK, LDA ) 567* 568* Compute the residual for the solution 569* 570 CALL SPOT02( UPLO, N, NRHS, A, LDA, X, LDA, WORK, 571 $ LDA, RWORK, RESULT( 4 ) ) 572* 573*+ TEST 5 574* Check solution from generated exact solution. 575* 576 CALL SGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 577 $ RESULT( 5 ) ) 578* 579*+ TESTS 6, 7, and 8 580* Use iterative refinement to improve the solution. 581* 582 SRNAMT = 'SSYRFS' 583 CALL SSYRFS( UPLO, N, NRHS, A, LDA, AFAC, LDA, 584 $ IWORK, B, LDA, X, LDA, RWORK, 585 $ RWORK( NRHS+1 ), WORK, IWORK( N+1 ), 586 $ INFO ) 587* 588* Check error code from SSYRFS. 589* 590 IF( INFO.NE.0 ) 591 $ CALL ALAERH( PATH, 'SSYRFS', INFO, 0, UPLO, N, 592 $ N, -1, -1, NRHS, IMAT, NFAIL, 593 $ NERRS, NOUT ) 594* 595 CALL SGET04( N, NRHS, X, LDA, XACT, LDA, RCONDC, 596 $ RESULT( 6 ) ) 597 CALL SPOT05( UPLO, N, NRHS, A, LDA, B, LDA, X, LDA, 598 $ XACT, LDA, RWORK, RWORK( NRHS+1 ), 599 $ RESULT( 7 ) ) 600* 601* Print information about the tests that did not pass 602* the threshold. 603* 604 DO 120 K = 3, 8 605 IF( RESULT( K ).GE.THRESH ) THEN 606 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 607 $ CALL ALAHD( NOUT, PATH ) 608 WRITE( NOUT, FMT = 9998 )UPLO, N, NRHS, 609 $ IMAT, K, RESULT( K ) 610 NFAIL = NFAIL + 1 611 END IF 612 120 CONTINUE 613 NRUN = NRUN + 6 614 130 CONTINUE 615* 616*+ TEST 9 617* Get an estimate of RCOND = 1/CNDNUM. 618* 619 140 CONTINUE 620 ANORM = SLANSY( '1', UPLO, N, A, LDA, RWORK ) 621 SRNAMT = 'SSYCON' 622 CALL SSYCON( UPLO, N, AFAC, LDA, IWORK, ANORM, RCOND, 623 $ WORK, IWORK( N+1 ), INFO ) 624* 625* Check error code from SSYCON and handle error. 626* 627 IF( INFO.NE.0 ) 628 $ CALL ALAERH( PATH, 'SSYCON', INFO, 0, UPLO, N, N, 629 $ -1, -1, -1, IMAT, NFAIL, NERRS, NOUT ) 630* 631* Compute the test ratio to compare to values of RCOND 632* 633 RESULT( 9 ) = SGET06( RCOND, RCONDC ) 634* 635* Print information about the tests that did not pass 636* the threshold. 637* 638 IF( RESULT( 9 ).GE.THRESH ) THEN 639 IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 ) 640 $ CALL ALAHD( NOUT, PATH ) 641 WRITE( NOUT, FMT = 9997 )UPLO, N, IMAT, 9, 642 $ RESULT( 9 ) 643 NFAIL = NFAIL + 1 644 END IF 645 NRUN = NRUN + 1 646 150 CONTINUE 647* 648 160 CONTINUE 649 170 CONTINUE 650 180 CONTINUE 651* 652* Print a summary of the results. 653* 654 CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS ) 655* 656 9999 FORMAT( ' UPLO = ''', A1, ''', N =', I5, ', NB =', I4, ', type ', 657 $ I2, ', test ', I2, ', ratio =', G12.5 ) 658 9998 FORMAT( ' UPLO = ''', A1, ''', N =', I5, ', NRHS=', I3, ', type ', 659 $ I2, ', test(', I2, ') =', G12.5 ) 660 9997 FORMAT( ' UPLO = ''', A1, ''', N =', I5, ',', 10X, ' type ', I2, 661 $ ', test(', I2, ') =', G12.5 ) 662 RETURN 663* 664* End of SCHKSY 665* 666 END 667