1*> \brief \b SCHKST 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 SCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, 12* NOUNIT, A, LDA, AP, SD, SE, D1, D2, D3, D4, D5, 13* WA1, WA2, WA3, WR, U, LDU, V, VP, TAU, Z, WORK, 14* LWORK, IWORK, LIWORK, RESULT, INFO ) 15* 16* .. Scalar Arguments .. 17* INTEGER INFO, LDA, LDU, LIWORK, LWORK, NOUNIT, NSIZES, 18* $ NTYPES 19* REAL THRESH 20* .. 21* .. Array Arguments .. 22* LOGICAL DOTYPE( * ) 23* INTEGER ISEED( 4 ), IWORK( * ), NN( * ) 24* REAL A( LDA, * ), AP( * ), D1( * ), D2( * ), 25* $ D3( * ), D4( * ), D5( * ), RESULT( * ), 26* $ SD( * ), SE( * ), TAU( * ), U( LDU, * ), 27* $ V( LDU, * ), VP( * ), WA1( * ), WA2( * ), 28* $ WA3( * ), WORK( * ), WR( * ), Z( LDU, * ) 29* .. 30* 31* 32*> \par Purpose: 33* ============= 34*> 35*> \verbatim 36*> 37*> SCHKST checks the symmetric eigenvalue problem routines. 38*> 39*> SSYTRD factors A as U S U' , where ' means transpose, 40*> S is symmetric tridiagonal, and U is orthogonal. 41*> SSYTRD can use either just the lower or just the upper triangle 42*> of A; SCHKST checks both cases. 43*> U is represented as a product of Householder 44*> transformations, whose vectors are stored in the first 45*> n-1 columns of V, and whose scale factors are in TAU. 46*> 47*> SSPTRD does the same as SSYTRD, except that A and V are stored 48*> in "packed" format. 49*> 50*> SORGTR constructs the matrix U from the contents of V and TAU. 51*> 52*> SOPGTR constructs the matrix U from the contents of VP and TAU. 53*> 54*> SSTEQR factors S as Z D1 Z' , where Z is the orthogonal 55*> matrix of eigenvectors and D1 is a diagonal matrix with 56*> the eigenvalues on the diagonal. D2 is the matrix of 57*> eigenvalues computed when Z is not computed. 58*> 59*> SSTERF computes D3, the matrix of eigenvalues, by the 60*> PWK method, which does not yield eigenvectors. 61*> 62*> SPTEQR factors S as Z4 D4 Z4' , for a 63*> symmetric positive definite tridiagonal matrix. 64*> D5 is the matrix of eigenvalues computed when Z is not 65*> computed. 66*> 67*> SSTEBZ computes selected eigenvalues. WA1, WA2, and 68*> WA3 will denote eigenvalues computed to high 69*> absolute accuracy, with different range options. 70*> WR will denote eigenvalues computed to high relative 71*> accuracy. 72*> 73*> SSTEIN computes Y, the eigenvectors of S, given the 74*> eigenvalues. 75*> 76*> SSTEDC factors S as Z D1 Z' , where Z is the orthogonal 77*> matrix of eigenvectors and D1 is a diagonal matrix with 78*> the eigenvalues on the diagonal ('I' option). It may also 79*> update an input orthogonal matrix, usually the output 80*> from SSYTRD/SORGTR or SSPTRD/SOPGTR ('V' option). It may 81*> also just compute eigenvalues ('N' option). 82*> 83*> SSTEMR factors S as Z D1 Z' , where Z is the orthogonal 84*> matrix of eigenvectors and D1 is a diagonal matrix with 85*> the eigenvalues on the diagonal ('I' option). SSTEMR 86*> uses the Relatively Robust Representation whenever possible. 87*> 88*> When SCHKST is called, a number of matrix "sizes" ("n's") and a 89*> number of matrix "types" are specified. For each size ("n") 90*> and each type of matrix, one matrix will be generated and used 91*> to test the symmetric eigenroutines. For each matrix, a number 92*> of tests will be performed: 93*> 94*> (1) | A - V S V' | / ( |A| n ulp ) SSYTRD( UPLO='U', ... ) 95*> 96*> (2) | I - UV' | / ( n ulp ) SORGTR( UPLO='U', ... ) 97*> 98*> (3) | A - V S V' | / ( |A| n ulp ) SSYTRD( UPLO='L', ... ) 99*> 100*> (4) | I - UV' | / ( n ulp ) SORGTR( UPLO='L', ... ) 101*> 102*> (5-8) Same as 1-4, but for SSPTRD and SOPGTR. 103*> 104*> (9) | S - Z D Z' | / ( |S| n ulp ) SSTEQR('V',...) 105*> 106*> (10) | I - ZZ' | / ( n ulp ) SSTEQR('V',...) 107*> 108*> (11) | D1 - D2 | / ( |D1| ulp ) SSTEQR('N',...) 109*> 110*> (12) | D1 - D3 | / ( |D1| ulp ) SSTERF 111*> 112*> (13) 0 if the true eigenvalues (computed by sturm count) 113*> of S are within THRESH of 114*> those in D1. 2*THRESH if they are not. (Tested using 115*> SSTECH) 116*> 117*> For S positive definite, 118*> 119*> (14) | S - Z4 D4 Z4' | / ( |S| n ulp ) SPTEQR('V',...) 120*> 121*> (15) | I - Z4 Z4' | / ( n ulp ) SPTEQR('V',...) 122*> 123*> (16) | D4 - D5 | / ( 100 |D4| ulp ) SPTEQR('N',...) 124*> 125*> When S is also diagonally dominant by the factor gamma < 1, 126*> 127*> (17) max | D4(i) - WR(i) | / ( |D4(i)| omega ) , 128*> i 129*> omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4 130*> SSTEBZ( 'A', 'E', ...) 131*> 132*> (18) | WA1 - D3 | / ( |D3| ulp ) SSTEBZ( 'A', 'E', ...) 133*> 134*> (19) ( max { min | WA2(i)-WA3(j) | } + 135*> i j 136*> max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) 137*> i j 138*> SSTEBZ( 'I', 'E', ...) 139*> 140*> (20) | S - Y WA1 Y' | / ( |S| n ulp ) SSTEBZ, SSTEIN 141*> 142*> (21) | I - Y Y' | / ( n ulp ) SSTEBZ, SSTEIN 143*> 144*> (22) | S - Z D Z' | / ( |S| n ulp ) SSTEDC('I') 145*> 146*> (23) | I - ZZ' | / ( n ulp ) SSTEDC('I') 147*> 148*> (24) | S - Z D Z' | / ( |S| n ulp ) SSTEDC('V') 149*> 150*> (25) | I - ZZ' | / ( n ulp ) SSTEDC('V') 151*> 152*> (26) | D1 - D2 | / ( |D1| ulp ) SSTEDC('V') and 153*> SSTEDC('N') 154*> 155*> Test 27 is disabled at the moment because SSTEMR does not 156*> guarantee high relatvie accuracy. 157*> 158*> (27) max | D6(i) - WR(i) | / ( |D6(i)| omega ) , 159*> i 160*> omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4 161*> SSTEMR('V', 'A') 162*> 163*> (28) max | D6(i) - WR(i) | / ( |D6(i)| omega ) , 164*> i 165*> omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4 166*> SSTEMR('V', 'I') 167*> 168*> Tests 29 through 34 are disable at present because SSTEMR 169*> does not handle partial spectrum requests. 170*> 171*> (29) | S - Z D Z' | / ( |S| n ulp ) SSTEMR('V', 'I') 172*> 173*> (30) | I - ZZ' | / ( n ulp ) SSTEMR('V', 'I') 174*> 175*> (31) ( max { min | WA2(i)-WA3(j) | } + 176*> i j 177*> max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) 178*> i j 179*> SSTEMR('N', 'I') vs. SSTEMR('V', 'I') 180*> 181*> (32) | S - Z D Z' | / ( |S| n ulp ) SSTEMR('V', 'V') 182*> 183*> (33) | I - ZZ' | / ( n ulp ) SSTEMR('V', 'V') 184*> 185*> (34) ( max { min | WA2(i)-WA3(j) | } + 186*> i j 187*> max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) 188*> i j 189*> SSTEMR('N', 'V') vs. SSTEMR('V', 'V') 190*> 191*> (35) | S - Z D Z' | / ( |S| n ulp ) SSTEMR('V', 'A') 192*> 193*> (36) | I - ZZ' | / ( n ulp ) SSTEMR('V', 'A') 194*> 195*> (37) ( max { min | WA2(i)-WA3(j) | } + 196*> i j 197*> max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) 198*> i j 199*> SSTEMR('N', 'A') vs. SSTEMR('V', 'A') 200*> 201*> The "sizes" are specified by an array NN(1:NSIZES); the value of 202*> each element NN(j) specifies one size. 203*> The "types" are specified by a logical array DOTYPE( 1:NTYPES ); 204*> if DOTYPE(j) is .TRUE., then matrix type "j" will be generated. 205*> Currently, the list of possible types is: 206*> 207*> (1) The zero matrix. 208*> (2) The identity matrix. 209*> 210*> (3) A diagonal matrix with evenly spaced entries 211*> 1, ..., ULP and random signs. 212*> (ULP = (first number larger than 1) - 1 ) 213*> (4) A diagonal matrix with geometrically spaced entries 214*> 1, ..., ULP and random signs. 215*> (5) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP 216*> and random signs. 217*> 218*> (6) Same as (4), but multiplied by SQRT( overflow threshold ) 219*> (7) Same as (4), but multiplied by SQRT( underflow threshold ) 220*> 221*> (8) A matrix of the form U' D U, where U is orthogonal and 222*> D has evenly spaced entries 1, ..., ULP with random signs 223*> on the diagonal. 224*> 225*> (9) A matrix of the form U' D U, where U is orthogonal and 226*> D has geometrically spaced entries 1, ..., ULP with random 227*> signs on the diagonal. 228*> 229*> (10) A matrix of the form U' D U, where U is orthogonal and 230*> D has "clustered" entries 1, ULP,..., ULP with random 231*> signs on the diagonal. 232*> 233*> (11) Same as (8), but multiplied by SQRT( overflow threshold ) 234*> (12) Same as (8), but multiplied by SQRT( underflow threshold ) 235*> 236*> (13) Symmetric matrix with random entries chosen from (-1,1). 237*> (14) Same as (13), but multiplied by SQRT( overflow threshold ) 238*> (15) Same as (13), but multiplied by SQRT( underflow threshold ) 239*> (16) Same as (8), but diagonal elements are all positive. 240*> (17) Same as (9), but diagonal elements are all positive. 241*> (18) Same as (10), but diagonal elements are all positive. 242*> (19) Same as (16), but multiplied by SQRT( overflow threshold ) 243*> (20) Same as (16), but multiplied by SQRT( underflow threshold ) 244*> (21) A diagonally dominant tridiagonal matrix with geometrically 245*> spaced diagonal entries 1, ..., ULP. 246*> \endverbatim 247* 248* Arguments: 249* ========== 250* 251*> \param[in] NSIZES 252*> \verbatim 253*> NSIZES is INTEGER 254*> The number of sizes of matrices to use. If it is zero, 255*> SCHKST does nothing. It must be at least zero. 256*> \endverbatim 257*> 258*> \param[in] NN 259*> \verbatim 260*> NN is INTEGER array, dimension (NSIZES) 261*> An array containing the sizes to be used for the matrices. 262*> Zero values will be skipped. The values must be at least 263*> zero. 264*> \endverbatim 265*> 266*> \param[in] NTYPES 267*> \verbatim 268*> NTYPES is INTEGER 269*> The number of elements in DOTYPE. If it is zero, SCHKST 270*> does nothing. It must be at least zero. If it is MAXTYP+1 271*> and NSIZES is 1, then an additional type, MAXTYP+1 is 272*> defined, which is to use whatever matrix is in A. This 273*> is only useful if DOTYPE(1:MAXTYP) is .FALSE. and 274*> DOTYPE(MAXTYP+1) is .TRUE. . 275*> \endverbatim 276*> 277*> \param[in] DOTYPE 278*> \verbatim 279*> DOTYPE is LOGICAL array, dimension (NTYPES) 280*> If DOTYPE(j) is .TRUE., then for each size in NN a 281*> matrix of that size and of type j will be generated. 282*> If NTYPES is smaller than the maximum number of types 283*> defined (PARAMETER MAXTYP), then types NTYPES+1 through 284*> MAXTYP will not be generated. If NTYPES is larger 285*> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES) 286*> will be ignored. 287*> \endverbatim 288*> 289*> \param[in,out] ISEED 290*> \verbatim 291*> ISEED is INTEGER array, dimension (4) 292*> On entry ISEED specifies the seed of the random number 293*> generator. The array elements should be between 0 and 4095; 294*> if not they will be reduced mod 4096. Also, ISEED(4) must 295*> be odd. The random number generator uses a linear 296*> congruential sequence limited to small integers, and so 297*> should produce machine independent random numbers. The 298*> values of ISEED are changed on exit, and can be used in the 299*> next call to SCHKST to continue the same random number 300*> sequence. 301*> \endverbatim 302*> 303*> \param[in] THRESH 304*> \verbatim 305*> THRESH is REAL 306*> A test will count as "failed" if the "error", computed as 307*> described above, exceeds THRESH. Note that the error 308*> is scaled to be O(1), so THRESH should be a reasonably 309*> small multiple of 1, e.g., 10 or 100. In particular, 310*> it should not depend on the precision (single vs. double) 311*> or the size of the matrix. It must be at least zero. 312*> \endverbatim 313*> 314*> \param[in] NOUNIT 315*> \verbatim 316*> NOUNIT is INTEGER 317*> The FORTRAN unit number for printing out error messages 318*> (e.g., if a routine returns IINFO not equal to 0.) 319*> \endverbatim 320*> 321*> \param[in,out] A 322*> \verbatim 323*> A is REAL array of 324*> dimension ( LDA , max(NN) ) 325*> Used to hold the matrix whose eigenvalues are to be 326*> computed. On exit, A contains the last matrix actually 327*> used. 328*> \endverbatim 329*> 330*> \param[in] LDA 331*> \verbatim 332*> LDA is INTEGER 333*> The leading dimension of A. It must be at 334*> least 1 and at least max( NN ). 335*> \endverbatim 336*> 337*> \param[out] AP 338*> \verbatim 339*> AP is REAL array of 340*> dimension( max(NN)*max(NN+1)/2 ) 341*> The matrix A stored in packed format. 342*> \endverbatim 343*> 344*> \param[out] SD 345*> \verbatim 346*> SD is REAL array of 347*> dimension( max(NN) ) 348*> The diagonal of the tridiagonal matrix computed by SSYTRD. 349*> On exit, SD and SE contain the tridiagonal form of the 350*> matrix in A. 351*> \endverbatim 352*> 353*> \param[out] SE 354*> \verbatim 355*> SE is REAL array of 356*> dimension( max(NN) ) 357*> The off-diagonal of the tridiagonal matrix computed by 358*> SSYTRD. On exit, SD and SE contain the tridiagonal form of 359*> the matrix in A. 360*> \endverbatim 361*> 362*> \param[out] D1 363*> \verbatim 364*> D1 is REAL array of 365*> dimension( max(NN) ) 366*> The eigenvalues of A, as computed by SSTEQR simlutaneously 367*> with Z. On exit, the eigenvalues in D1 correspond with the 368*> matrix in A. 369*> \endverbatim 370*> 371*> \param[out] D2 372*> \verbatim 373*> D2 is REAL array of 374*> dimension( max(NN) ) 375*> The eigenvalues of A, as computed by SSTEQR if Z is not 376*> computed. On exit, the eigenvalues in D2 correspond with 377*> the matrix in A. 378*> \endverbatim 379*> 380*> \param[out] D3 381*> \verbatim 382*> D3 is REAL array of 383*> dimension( max(NN) ) 384*> The eigenvalues of A, as computed by SSTERF. On exit, the 385*> eigenvalues in D3 correspond with the matrix in A. 386*> \endverbatim 387*> 388*> \param[out] D4 389*> \verbatim 390*> D4 is REAL array of 391*> dimension( max(NN) ) 392*> The eigenvalues of A, as computed by SPTEQR(V). 393*> ZPTEQR factors S as Z4 D4 Z4* 394*> On exit, the eigenvalues in D4 correspond with the matrix in A. 395*> \endverbatim 396*> 397*> \param[out] D5 398*> \verbatim 399*> D5 is REAL array of 400*> dimension( max(NN) ) 401*> The eigenvalues of A, as computed by SPTEQR(N) 402*> when Z is not computed. On exit, the 403*> eigenvalues in D4 correspond with the matrix in A. 404*> \endverbatim 405*> 406*> \param[out] WA1 407*> \verbatim 408*> WA1 is REAL array of 409*> dimension( max(NN) ) 410*> All eigenvalues of A, computed to high 411*> absolute accuracy, with different range options. 412*> as computed by SSTEBZ. 413*> \endverbatim 414*> 415*> \param[out] WA2 416*> \verbatim 417*> WA2 is REAL array of 418*> dimension( max(NN) ) 419*> Selected eigenvalues of A, computed to high 420*> absolute accuracy, with different range options. 421*> as computed by SSTEBZ. 422*> Choose random values for IL and IU, and ask for the 423*> IL-th through IU-th eigenvalues. 424*> \endverbatim 425*> 426*> \param[out] WA3 427*> \verbatim 428*> WA3 is REAL array of 429*> dimension( max(NN) ) 430*> Selected eigenvalues of A, computed to high 431*> absolute accuracy, with different range options. 432*> as computed by SSTEBZ. 433*> Determine the values VL and VU of the IL-th and IU-th 434*> eigenvalues and ask for all eigenvalues in this range. 435*> \endverbatim 436*> 437*> \param[out] WR 438*> \verbatim 439*> WR is REAL array of 440*> dimension( max(NN) ) 441*> All eigenvalues of A, computed to high 442*> absolute accuracy, with different options. 443*> as computed by SSTEBZ. 444*> \endverbatim 445*> 446*> \param[out] U 447*> \verbatim 448*> U is REAL array of 449*> dimension( LDU, max(NN) ). 450*> The orthogonal matrix computed by SSYTRD + SORGTR. 451*> \endverbatim 452*> 453*> \param[in] LDU 454*> \verbatim 455*> LDU is INTEGER 456*> The leading dimension of U, Z, and V. It must be at least 1 457*> and at least max( NN ). 458*> \endverbatim 459*> 460*> \param[out] V 461*> \verbatim 462*> V is REAL array of 463*> dimension( LDU, max(NN) ). 464*> The Housholder vectors computed by SSYTRD in reducing A to 465*> tridiagonal form. The vectors computed with UPLO='U' are 466*> in the upper triangle, and the vectors computed with UPLO='L' 467*> are in the lower triangle. (As described in SSYTRD, the 468*> sub- and superdiagonal are not set to 1, although the 469*> true Householder vector has a 1 in that position. The 470*> routines that use V, such as SORGTR, set those entries to 471*> 1 before using them, and then restore them later.) 472*> \endverbatim 473*> 474*> \param[out] VP 475*> \verbatim 476*> VP is REAL array of 477*> dimension( max(NN)*max(NN+1)/2 ) 478*> The matrix V stored in packed format. 479*> \endverbatim 480*> 481*> \param[out] TAU 482*> \verbatim 483*> TAU is REAL array of 484*> dimension( max(NN) ) 485*> The Householder factors computed by SSYTRD in reducing A 486*> to tridiagonal form. 487*> \endverbatim 488*> 489*> \param[out] Z 490*> \verbatim 491*> Z is REAL array of 492*> dimension( LDU, max(NN) ). 493*> The orthogonal matrix of eigenvectors computed by SSTEQR, 494*> SPTEQR, and SSTEIN. 495*> \endverbatim 496*> 497*> \param[out] WORK 498*> \verbatim 499*> WORK is REAL array of 500*> dimension( LWORK ) 501*> \endverbatim 502*> 503*> \param[in] LWORK 504*> \verbatim 505*> LWORK is INTEGER 506*> The number of entries in WORK. This must be at least 507*> 1 + 4 * Nmax + 2 * Nmax * lg Nmax + 3 * Nmax**2 508*> where Nmax = max( NN(j), 2 ) and lg = log base 2. 509*> \endverbatim 510*> 511*> \param[out] IWORK 512*> \verbatim 513*> IWORK is INTEGER array, 514*> Workspace. 515*> \endverbatim 516*> 517*> \param[out] LIWORK 518*> \verbatim 519*> LIWORK is INTEGER 520*> The number of entries in IWORK. This must be at least 521*> 6 + 6*Nmax + 5 * Nmax * lg Nmax 522*> where Nmax = max( NN(j), 2 ) and lg = log base 2. 523*> \endverbatim 524*> 525*> \param[out] RESULT 526*> \verbatim 527*> RESULT is REAL array, dimension (26) 528*> The values computed by the tests described above. 529*> The values are currently limited to 1/ulp, to avoid 530*> overflow. 531*> \endverbatim 532*> 533*> \param[out] INFO 534*> \verbatim 535*> INFO is INTEGER 536*> If 0, then everything ran OK. 537*> -1: NSIZES < 0 538*> -2: Some NN(j) < 0 539*> -3: NTYPES < 0 540*> -5: THRESH < 0 541*> -9: LDA < 1 or LDA < NMAX, where NMAX is max( NN(j) ). 542*> -23: LDU < 1 or LDU < NMAX. 543*> -29: LWORK too small. 544*> If SLATMR, SLATMS, SSYTRD, SORGTR, SSTEQR, SSTERF, 545*> or SORMC2 returns an error code, the 546*> absolute value of it is returned. 547*> 548*>----------------------------------------------------------------------- 549*> 550*> Some Local Variables and Parameters: 551*> ---- ----- --------- --- ---------- 552*> ZERO, ONE Real 0 and 1. 553*> MAXTYP The number of types defined. 554*> NTEST The number of tests performed, or which can 555*> be performed so far, for the current matrix. 556*> NTESTT The total number of tests performed so far. 557*> NBLOCK Blocksize as returned by ENVIR. 558*> NMAX Largest value in NN. 559*> NMATS The number of matrices generated so far. 560*> NERRS The number of tests which have exceeded THRESH 561*> so far. 562*> COND, IMODE Values to be passed to the matrix generators. 563*> ANORM Norm of A; passed to matrix generators. 564*> 565*> OVFL, UNFL Overflow and underflow thresholds. 566*> ULP, ULPINV Finest relative precision and its inverse. 567*> RTOVFL, RTUNFL Square roots of the previous 2 values. 568*> The following four arrays decode JTYPE: 569*> KTYPE(j) The general type (1-10) for type "j". 570*> KMODE(j) The MODE value to be passed to the matrix 571*> generator for type "j". 572*> KMAGN(j) The order of magnitude ( O(1), 573*> O(overflow^(1/2) ), O(underflow^(1/2) ) 574*> \endverbatim 575* 576* Authors: 577* ======== 578* 579*> \author Univ. of Tennessee 580*> \author Univ. of California Berkeley 581*> \author Univ. of Colorado Denver 582*> \author NAG Ltd. 583* 584*> \ingroup single_eig 585* 586* ===================================================================== 587 SUBROUTINE SCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, 588 $ NOUNIT, A, LDA, AP, SD, SE, D1, D2, D3, D4, D5, 589 $ WA1, WA2, WA3, WR, U, LDU, V, VP, TAU, Z, WORK, 590 $ LWORK, IWORK, LIWORK, RESULT, INFO ) 591* 592* -- LAPACK test routine -- 593* -- LAPACK is a software package provided by Univ. of Tennessee, -- 594* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 595* 596* .. Scalar Arguments .. 597 INTEGER INFO, LDA, LDU, LIWORK, LWORK, NOUNIT, NSIZES, 598 $ NTYPES 599 REAL THRESH 600* .. 601* .. Array Arguments .. 602 LOGICAL DOTYPE( * ) 603 INTEGER ISEED( 4 ), IWORK( * ), NN( * ) 604 REAL A( LDA, * ), AP( * ), D1( * ), D2( * ), 605 $ D3( * ), D4( * ), D5( * ), RESULT( * ), 606 $ SD( * ), SE( * ), TAU( * ), U( LDU, * ), 607 $ V( LDU, * ), VP( * ), WA1( * ), WA2( * ), 608 $ WA3( * ), WORK( * ), WR( * ), Z( LDU, * ) 609* .. 610* 611* ===================================================================== 612* 613* .. Parameters .. 614 REAL ZERO, ONE, TWO, EIGHT, TEN, HUN 615 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0, 616 $ EIGHT = 8.0E0, TEN = 10.0E0, HUN = 100.0E0 ) 617 REAL HALF 618 PARAMETER ( HALF = ONE / TWO ) 619 INTEGER MAXTYP 620 PARAMETER ( MAXTYP = 21 ) 621 LOGICAL SRANGE 622 PARAMETER ( SRANGE = .FALSE. ) 623 LOGICAL SREL 624 PARAMETER ( SREL = .FALSE. ) 625* .. 626* .. Local Scalars .. 627 LOGICAL BADNN, TRYRAC 628 INTEGER I, IINFO, IL, IMODE, ITEMP, ITYPE, IU, J, JC, 629 $ JR, JSIZE, JTYPE, LGN, LIWEDC, LOG2UI, LWEDC, 630 $ M, M2, M3, MTYPES, N, NAP, NBLOCK, NERRS, 631 $ NMATS, NMAX, NSPLIT, NTEST, NTESTT 632 REAL ABSTOL, ANINV, ANORM, COND, OVFL, RTOVFL, 633 $ RTUNFL, TEMP1, TEMP2, TEMP3, TEMP4, ULP, 634 $ ULPINV, UNFL, VL, VU 635* .. 636* .. Local Arrays .. 637 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ), 638 $ KMAGN( MAXTYP ), KMODE( MAXTYP ), 639 $ KTYPE( MAXTYP ) 640 REAL DUMMA( 1 ) 641* .. 642* .. External Functions .. 643 INTEGER ILAENV 644 REAL SLAMCH, SLARND, SSXT1 645 EXTERNAL ILAENV, SLAMCH, SLARND, SSXT1 646* .. 647* .. External Subroutines .. 648 EXTERNAL SCOPY, SLABAD, SLACPY, SLASET, SLASUM, SLATMR, 649 $ SLATMS, SOPGTR, SORGTR, SPTEQR, SSPT21, SSPTRD, 650 $ SSTEBZ, SSTECH, SSTEDC, SSTEMR, SSTEIN, SSTEQR, 651 $ SSTERF, SSTT21, SSTT22, SSYT21, SSYTRD, XERBLA 652* .. 653* .. Intrinsic Functions .. 654 INTRINSIC ABS, INT, LOG, MAX, MIN, REAL, SQRT 655* .. 656* .. Data statements .. 657 DATA KTYPE / 1, 2, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 8, 658 $ 8, 8, 9, 9, 9, 9, 9, 10 / 659 DATA KMAGN / 1, 1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1, 660 $ 2, 3, 1, 1, 1, 2, 3, 1 / 661 DATA KMODE / 0, 0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0, 662 $ 0, 0, 4, 3, 1, 4, 4, 3 / 663* .. 664* .. Executable Statements .. 665* 666* Keep ftnchek happy 667 IDUMMA( 1 ) = 1 668* 669* Check for errors 670* 671 NTESTT = 0 672 INFO = 0 673* 674* Important constants 675* 676 BADNN = .FALSE. 677 TRYRAC = .TRUE. 678 NMAX = 1 679 DO 10 J = 1, NSIZES 680 NMAX = MAX( NMAX, NN( J ) ) 681 IF( NN( J ).LT.0 ) 682 $ BADNN = .TRUE. 683 10 CONTINUE 684* 685 NBLOCK = ILAENV( 1, 'SSYTRD', 'L', NMAX, -1, -1, -1 ) 686 NBLOCK = MIN( NMAX, MAX( 1, NBLOCK ) ) 687* 688* Check for errors 689* 690 IF( NSIZES.LT.0 ) THEN 691 INFO = -1 692 ELSE IF( BADNN ) THEN 693 INFO = -2 694 ELSE IF( NTYPES.LT.0 ) THEN 695 INFO = -3 696 ELSE IF( LDA.LT.NMAX ) THEN 697 INFO = -9 698 ELSE IF( LDU.LT.NMAX ) THEN 699 INFO = -23 700 ELSE IF( 2*MAX( 2, NMAX )**2.GT.LWORK ) THEN 701 INFO = -29 702 END IF 703* 704 IF( INFO.NE.0 ) THEN 705 CALL XERBLA( 'SCHKST', -INFO ) 706 RETURN 707 END IF 708* 709* Quick return if possible 710* 711 IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 ) 712 $ RETURN 713* 714* More Important constants 715* 716 UNFL = SLAMCH( 'Safe minimum' ) 717 OVFL = ONE / UNFL 718 CALL SLABAD( UNFL, OVFL ) 719 ULP = SLAMCH( 'Epsilon' )*SLAMCH( 'Base' ) 720 ULPINV = ONE / ULP 721 LOG2UI = INT( LOG( ULPINV ) / LOG( TWO ) ) 722 RTUNFL = SQRT( UNFL ) 723 RTOVFL = SQRT( OVFL ) 724* 725* Loop over sizes, types 726* 727 DO 20 I = 1, 4 728 ISEED2( I ) = ISEED( I ) 729 20 CONTINUE 730 NERRS = 0 731 NMATS = 0 732* 733 DO 310 JSIZE = 1, NSIZES 734 N = NN( JSIZE ) 735 IF( N.GT.0 ) THEN 736 LGN = INT( LOG( REAL( N ) ) / LOG( TWO ) ) 737 IF( 2**LGN.LT.N ) 738 $ LGN = LGN + 1 739 IF( 2**LGN.LT.N ) 740 $ LGN = LGN + 1 741 LWEDC = 1 + 4*N + 2*N*LGN + 4*N**2 742 LIWEDC = 6 + 6*N + 5*N*LGN 743 ELSE 744 LWEDC = 8 745 LIWEDC = 12 746 END IF 747 NAP = ( N*( N+1 ) ) / 2 748 ANINV = ONE / REAL( MAX( 1, N ) ) 749* 750 IF( NSIZES.NE.1 ) THEN 751 MTYPES = MIN( MAXTYP, NTYPES ) 752 ELSE 753 MTYPES = MIN( MAXTYP+1, NTYPES ) 754 END IF 755* 756 DO 300 JTYPE = 1, MTYPES 757 IF( .NOT.DOTYPE( JTYPE ) ) 758 $ GO TO 300 759 NMATS = NMATS + 1 760 NTEST = 0 761* 762 DO 30 J = 1, 4 763 IOLDSD( J ) = ISEED( J ) 764 30 CONTINUE 765* 766* Compute "A" 767* 768* Control parameters: 769* 770* KMAGN KMODE KTYPE 771* =1 O(1) clustered 1 zero 772* =2 large clustered 2 identity 773* =3 small exponential (none) 774* =4 arithmetic diagonal, (w/ eigenvalues) 775* =5 random log symmetric, w/ eigenvalues 776* =6 random (none) 777* =7 random diagonal 778* =8 random symmetric 779* =9 positive definite 780* =10 diagonally dominant tridiagonal 781* 782 IF( MTYPES.GT.MAXTYP ) 783 $ GO TO 100 784* 785 ITYPE = KTYPE( JTYPE ) 786 IMODE = KMODE( JTYPE ) 787* 788* Compute norm 789* 790 GO TO ( 40, 50, 60 )KMAGN( JTYPE ) 791* 792 40 CONTINUE 793 ANORM = ONE 794 GO TO 70 795* 796 50 CONTINUE 797 ANORM = ( RTOVFL*ULP )*ANINV 798 GO TO 70 799* 800 60 CONTINUE 801 ANORM = RTUNFL*N*ULPINV 802 GO TO 70 803* 804 70 CONTINUE 805* 806 CALL SLASET( 'Full', LDA, N, ZERO, ZERO, A, LDA ) 807 IINFO = 0 808 IF( JTYPE.LE.15 ) THEN 809 COND = ULPINV 810 ELSE 811 COND = ULPINV*ANINV / TEN 812 END IF 813* 814* Special Matrices -- Identity & Jordan block 815* 816* Zero 817* 818 IF( ITYPE.EQ.1 ) THEN 819 IINFO = 0 820* 821 ELSE IF( ITYPE.EQ.2 ) THEN 822* 823* Identity 824* 825 DO 80 JC = 1, N 826 A( JC, JC ) = ANORM 827 80 CONTINUE 828* 829 ELSE IF( ITYPE.EQ.4 ) THEN 830* 831* Diagonal Matrix, [Eigen]values Specified 832* 833 CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND, 834 $ ANORM, 0, 0, 'N', A, LDA, WORK( N+1 ), 835 $ IINFO ) 836* 837* 838 ELSE IF( ITYPE.EQ.5 ) THEN 839* 840* Symmetric, eigenvalues specified 841* 842 CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND, 843 $ ANORM, N, N, 'N', A, LDA, WORK( N+1 ), 844 $ IINFO ) 845* 846 ELSE IF( ITYPE.EQ.7 ) THEN 847* 848* Diagonal, random eigenvalues 849* 850 CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE, 851 $ 'T', 'N', WORK( N+1 ), 1, ONE, 852 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0, 853 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO ) 854* 855 ELSE IF( ITYPE.EQ.8 ) THEN 856* 857* Symmetric, random eigenvalues 858* 859 CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE, 860 $ 'T', 'N', WORK( N+1 ), 1, ONE, 861 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N, 862 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO ) 863* 864 ELSE IF( ITYPE.EQ.9 ) THEN 865* 866* Positive definite, eigenvalues specified. 867* 868 CALL SLATMS( N, N, 'S', ISEED, 'P', WORK, IMODE, COND, 869 $ ANORM, N, N, 'N', A, LDA, WORK( N+1 ), 870 $ IINFO ) 871* 872 ELSE IF( ITYPE.EQ.10 ) THEN 873* 874* Positive definite tridiagonal, eigenvalues specified. 875* 876 CALL SLATMS( N, N, 'S', ISEED, 'P', WORK, IMODE, COND, 877 $ ANORM, 1, 1, 'N', A, LDA, WORK( N+1 ), 878 $ IINFO ) 879 DO 90 I = 2, N 880 TEMP1 = ABS( A( I-1, I ) ) / 881 $ SQRT( ABS( A( I-1, I-1 )*A( I, I ) ) ) 882 IF( TEMP1.GT.HALF ) THEN 883 A( I-1, I ) = HALF*SQRT( ABS( A( I-1, I-1 )*A( I, 884 $ I ) ) ) 885 A( I, I-1 ) = A( I-1, I ) 886 END IF 887 90 CONTINUE 888* 889 ELSE 890* 891 IINFO = 1 892 END IF 893* 894 IF( IINFO.NE.0 ) THEN 895 WRITE( NOUNIT, FMT = 9999 )'Generator', IINFO, N, JTYPE, 896 $ IOLDSD 897 INFO = ABS( IINFO ) 898 RETURN 899 END IF 900* 901 100 CONTINUE 902* 903* Call SSYTRD and SORGTR to compute S and U from 904* upper triangle. 905* 906 CALL SLACPY( 'U', N, N, A, LDA, V, LDU ) 907* 908 NTEST = 1 909 CALL SSYTRD( 'U', N, V, LDU, SD, SE, TAU, WORK, LWORK, 910 $ IINFO ) 911* 912 IF( IINFO.NE.0 ) THEN 913 WRITE( NOUNIT, FMT = 9999 )'SSYTRD(U)', IINFO, N, JTYPE, 914 $ IOLDSD 915 INFO = ABS( IINFO ) 916 IF( IINFO.LT.0 ) THEN 917 RETURN 918 ELSE 919 RESULT( 1 ) = ULPINV 920 GO TO 280 921 END IF 922 END IF 923* 924 CALL SLACPY( 'U', N, N, V, LDU, U, LDU ) 925* 926 NTEST = 2 927 CALL SORGTR( 'U', N, U, LDU, TAU, WORK, LWORK, IINFO ) 928 IF( IINFO.NE.0 ) THEN 929 WRITE( NOUNIT, FMT = 9999 )'SORGTR(U)', IINFO, N, JTYPE, 930 $ IOLDSD 931 INFO = ABS( IINFO ) 932 IF( IINFO.LT.0 ) THEN 933 RETURN 934 ELSE 935 RESULT( 2 ) = ULPINV 936 GO TO 280 937 END IF 938 END IF 939* 940* Do tests 1 and 2 941* 942 CALL SSYT21( 2, 'Upper', N, 1, A, LDA, SD, SE, U, LDU, V, 943 $ LDU, TAU, WORK, RESULT( 1 ) ) 944 CALL SSYT21( 3, 'Upper', N, 1, A, LDA, SD, SE, U, LDU, V, 945 $ LDU, TAU, WORK, RESULT( 2 ) ) 946* 947* Call SSYTRD and SORGTR to compute S and U from 948* lower triangle, do tests. 949* 950 CALL SLACPY( 'L', N, N, A, LDA, V, LDU ) 951* 952 NTEST = 3 953 CALL SSYTRD( 'L', N, V, LDU, SD, SE, TAU, WORK, LWORK, 954 $ IINFO ) 955* 956 IF( IINFO.NE.0 ) THEN 957 WRITE( NOUNIT, FMT = 9999 )'SSYTRD(L)', IINFO, N, JTYPE, 958 $ IOLDSD 959 INFO = ABS( IINFO ) 960 IF( IINFO.LT.0 ) THEN 961 RETURN 962 ELSE 963 RESULT( 3 ) = ULPINV 964 GO TO 280 965 END IF 966 END IF 967* 968 CALL SLACPY( 'L', N, N, V, LDU, U, LDU ) 969* 970 NTEST = 4 971 CALL SORGTR( 'L', N, U, LDU, TAU, WORK, LWORK, IINFO ) 972 IF( IINFO.NE.0 ) THEN 973 WRITE( NOUNIT, FMT = 9999 )'SORGTR(L)', IINFO, N, JTYPE, 974 $ IOLDSD 975 INFO = ABS( IINFO ) 976 IF( IINFO.LT.0 ) THEN 977 RETURN 978 ELSE 979 RESULT( 4 ) = ULPINV 980 GO TO 280 981 END IF 982 END IF 983* 984 CALL SSYT21( 2, 'Lower', N, 1, A, LDA, SD, SE, U, LDU, V, 985 $ LDU, TAU, WORK, RESULT( 3 ) ) 986 CALL SSYT21( 3, 'Lower', N, 1, A, LDA, SD, SE, U, LDU, V, 987 $ LDU, TAU, WORK, RESULT( 4 ) ) 988* 989* Store the upper triangle of A in AP 990* 991 I = 0 992 DO 120 JC = 1, N 993 DO 110 JR = 1, JC 994 I = I + 1 995 AP( I ) = A( JR, JC ) 996 110 CONTINUE 997 120 CONTINUE 998* 999* Call SSPTRD and SOPGTR to compute S and U from AP 1000* 1001 CALL SCOPY( NAP, AP, 1, VP, 1 ) 1002* 1003 NTEST = 5 1004 CALL SSPTRD( 'U', N, VP, SD, SE, TAU, IINFO ) 1005* 1006 IF( IINFO.NE.0 ) THEN 1007 WRITE( NOUNIT, FMT = 9999 )'SSPTRD(U)', IINFO, N, JTYPE, 1008 $ IOLDSD 1009 INFO = ABS( IINFO ) 1010 IF( IINFO.LT.0 ) THEN 1011 RETURN 1012 ELSE 1013 RESULT( 5 ) = ULPINV 1014 GO TO 280 1015 END IF 1016 END IF 1017* 1018 NTEST = 6 1019 CALL SOPGTR( 'U', N, VP, TAU, U, LDU, WORK, IINFO ) 1020 IF( IINFO.NE.0 ) THEN 1021 WRITE( NOUNIT, FMT = 9999 )'SOPGTR(U)', IINFO, N, JTYPE, 1022 $ IOLDSD 1023 INFO = ABS( IINFO ) 1024 IF( IINFO.LT.0 ) THEN 1025 RETURN 1026 ELSE 1027 RESULT( 6 ) = ULPINV 1028 GO TO 280 1029 END IF 1030 END IF 1031* 1032* Do tests 5 and 6 1033* 1034 CALL SSPT21( 2, 'Upper', N, 1, AP, SD, SE, U, LDU, VP, TAU, 1035 $ WORK, RESULT( 5 ) ) 1036 CALL SSPT21( 3, 'Upper', N, 1, AP, SD, SE, U, LDU, VP, TAU, 1037 $ WORK, RESULT( 6 ) ) 1038* 1039* Store the lower triangle of A in AP 1040* 1041 I = 0 1042 DO 140 JC = 1, N 1043 DO 130 JR = JC, N 1044 I = I + 1 1045 AP( I ) = A( JR, JC ) 1046 130 CONTINUE 1047 140 CONTINUE 1048* 1049* Call SSPTRD and SOPGTR to compute S and U from AP 1050* 1051 CALL SCOPY( NAP, AP, 1, VP, 1 ) 1052* 1053 NTEST = 7 1054 CALL SSPTRD( 'L', N, VP, SD, SE, TAU, IINFO ) 1055* 1056 IF( IINFO.NE.0 ) THEN 1057 WRITE( NOUNIT, FMT = 9999 )'SSPTRD(L)', IINFO, N, JTYPE, 1058 $ IOLDSD 1059 INFO = ABS( IINFO ) 1060 IF( IINFO.LT.0 ) THEN 1061 RETURN 1062 ELSE 1063 RESULT( 7 ) = ULPINV 1064 GO TO 280 1065 END IF 1066 END IF 1067* 1068 NTEST = 8 1069 CALL SOPGTR( 'L', N, VP, TAU, U, LDU, WORK, IINFO ) 1070 IF( IINFO.NE.0 ) THEN 1071 WRITE( NOUNIT, FMT = 9999 )'SOPGTR(L)', IINFO, N, JTYPE, 1072 $ IOLDSD 1073 INFO = ABS( IINFO ) 1074 IF( IINFO.LT.0 ) THEN 1075 RETURN 1076 ELSE 1077 RESULT( 8 ) = ULPINV 1078 GO TO 280 1079 END IF 1080 END IF 1081* 1082 CALL SSPT21( 2, 'Lower', N, 1, AP, SD, SE, U, LDU, VP, TAU, 1083 $ WORK, RESULT( 7 ) ) 1084 CALL SSPT21( 3, 'Lower', N, 1, AP, SD, SE, U, LDU, VP, TAU, 1085 $ WORK, RESULT( 8 ) ) 1086* 1087* Call SSTEQR to compute D1, D2, and Z, do tests. 1088* 1089* Compute D1 and Z 1090* 1091 CALL SCOPY( N, SD, 1, D1, 1 ) 1092 IF( N.GT.0 ) 1093 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1094 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1095* 1096 NTEST = 9 1097 CALL SSTEQR( 'V', N, D1, WORK, Z, LDU, WORK( N+1 ), IINFO ) 1098 IF( IINFO.NE.0 ) THEN 1099 WRITE( NOUNIT, FMT = 9999 )'SSTEQR(V)', IINFO, N, JTYPE, 1100 $ IOLDSD 1101 INFO = ABS( IINFO ) 1102 IF( IINFO.LT.0 ) THEN 1103 RETURN 1104 ELSE 1105 RESULT( 9 ) = ULPINV 1106 GO TO 280 1107 END IF 1108 END IF 1109* 1110* Compute D2 1111* 1112 CALL SCOPY( N, SD, 1, D2, 1 ) 1113 IF( N.GT.0 ) 1114 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1115* 1116 NTEST = 11 1117 CALL SSTEQR( 'N', N, D2, WORK, WORK( N+1 ), LDU, 1118 $ WORK( N+1 ), IINFO ) 1119 IF( IINFO.NE.0 ) THEN 1120 WRITE( NOUNIT, FMT = 9999 )'SSTEQR(N)', IINFO, N, JTYPE, 1121 $ IOLDSD 1122 INFO = ABS( IINFO ) 1123 IF( IINFO.LT.0 ) THEN 1124 RETURN 1125 ELSE 1126 RESULT( 11 ) = ULPINV 1127 GO TO 280 1128 END IF 1129 END IF 1130* 1131* Compute D3 (using PWK method) 1132* 1133 CALL SCOPY( N, SD, 1, D3, 1 ) 1134 IF( N.GT.0 ) 1135 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1136* 1137 NTEST = 12 1138 CALL SSTERF( N, D3, WORK, IINFO ) 1139 IF( IINFO.NE.0 ) THEN 1140 WRITE( NOUNIT, FMT = 9999 )'SSTERF', IINFO, N, JTYPE, 1141 $ IOLDSD 1142 INFO = ABS( IINFO ) 1143 IF( IINFO.LT.0 ) THEN 1144 RETURN 1145 ELSE 1146 RESULT( 12 ) = ULPINV 1147 GO TO 280 1148 END IF 1149 END IF 1150* 1151* Do Tests 9 and 10 1152* 1153 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1154 $ RESULT( 9 ) ) 1155* 1156* Do Tests 11 and 12 1157* 1158 TEMP1 = ZERO 1159 TEMP2 = ZERO 1160 TEMP3 = ZERO 1161 TEMP4 = ZERO 1162* 1163 DO 150 J = 1, N 1164 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 1165 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1166 TEMP3 = MAX( TEMP3, ABS( D1( J ) ), ABS( D3( J ) ) ) 1167 TEMP4 = MAX( TEMP4, ABS( D1( J )-D3( J ) ) ) 1168 150 CONTINUE 1169* 1170 RESULT( 11 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 1171 RESULT( 12 ) = TEMP4 / MAX( UNFL, ULP*MAX( TEMP3, TEMP4 ) ) 1172* 1173* Do Test 13 -- Sturm Sequence Test of Eigenvalues 1174* Go up by factors of two until it succeeds 1175* 1176 NTEST = 13 1177 TEMP1 = THRESH*( HALF-ULP ) 1178* 1179 DO 160 J = 0, LOG2UI 1180 CALL SSTECH( N, SD, SE, D1, TEMP1, WORK, IINFO ) 1181 IF( IINFO.EQ.0 ) 1182 $ GO TO 170 1183 TEMP1 = TEMP1*TWO 1184 160 CONTINUE 1185* 1186 170 CONTINUE 1187 RESULT( 13 ) = TEMP1 1188* 1189* For positive definite matrices ( JTYPE.GT.15 ) call SPTEQR 1190* and do tests 14, 15, and 16 . 1191* 1192 IF( JTYPE.GT.15 ) THEN 1193* 1194* Compute D4 and Z4 1195* 1196 CALL SCOPY( N, SD, 1, D4, 1 ) 1197 IF( N.GT.0 ) 1198 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1199 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1200* 1201 NTEST = 14 1202 CALL SPTEQR( 'V', N, D4, WORK, Z, LDU, WORK( N+1 ), 1203 $ IINFO ) 1204 IF( IINFO.NE.0 ) THEN 1205 WRITE( NOUNIT, FMT = 9999 )'SPTEQR(V)', IINFO, N, 1206 $ JTYPE, IOLDSD 1207 INFO = ABS( IINFO ) 1208 IF( IINFO.LT.0 ) THEN 1209 RETURN 1210 ELSE 1211 RESULT( 14 ) = ULPINV 1212 GO TO 280 1213 END IF 1214 END IF 1215* 1216* Do Tests 14 and 15 1217* 1218 CALL SSTT21( N, 0, SD, SE, D4, DUMMA, Z, LDU, WORK, 1219 $ RESULT( 14 ) ) 1220* 1221* Compute D5 1222* 1223 CALL SCOPY( N, SD, 1, D5, 1 ) 1224 IF( N.GT.0 ) 1225 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1226* 1227 NTEST = 16 1228 CALL SPTEQR( 'N', N, D5, WORK, Z, LDU, WORK( N+1 ), 1229 $ IINFO ) 1230 IF( IINFO.NE.0 ) THEN 1231 WRITE( NOUNIT, FMT = 9999 )'SPTEQR(N)', IINFO, N, 1232 $ JTYPE, IOLDSD 1233 INFO = ABS( IINFO ) 1234 IF( IINFO.LT.0 ) THEN 1235 RETURN 1236 ELSE 1237 RESULT( 16 ) = ULPINV 1238 GO TO 280 1239 END IF 1240 END IF 1241* 1242* Do Test 16 1243* 1244 TEMP1 = ZERO 1245 TEMP2 = ZERO 1246 DO 180 J = 1, N 1247 TEMP1 = MAX( TEMP1, ABS( D4( J ) ), ABS( D5( J ) ) ) 1248 TEMP2 = MAX( TEMP2, ABS( D4( J )-D5( J ) ) ) 1249 180 CONTINUE 1250* 1251 RESULT( 16 ) = TEMP2 / MAX( UNFL, 1252 $ HUN*ULP*MAX( TEMP1, TEMP2 ) ) 1253 ELSE 1254 RESULT( 14 ) = ZERO 1255 RESULT( 15 ) = ZERO 1256 RESULT( 16 ) = ZERO 1257 END IF 1258* 1259* Call SSTEBZ with different options and do tests 17-18. 1260* 1261* If S is positive definite and diagonally dominant, 1262* ask for all eigenvalues with high relative accuracy. 1263* 1264 VL = ZERO 1265 VU = ZERO 1266 IL = 0 1267 IU = 0 1268 IF( JTYPE.EQ.21 ) THEN 1269 NTEST = 17 1270 ABSTOL = UNFL + UNFL 1271 CALL SSTEBZ( 'A', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 1272 $ M, NSPLIT, WR, IWORK( 1 ), IWORK( N+1 ), 1273 $ WORK, IWORK( 2*N+1 ), IINFO ) 1274 IF( IINFO.NE.0 ) THEN 1275 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A,rel)', IINFO, N, 1276 $ JTYPE, IOLDSD 1277 INFO = ABS( IINFO ) 1278 IF( IINFO.LT.0 ) THEN 1279 RETURN 1280 ELSE 1281 RESULT( 17 ) = ULPINV 1282 GO TO 280 1283 END IF 1284 END IF 1285* 1286* Do test 17 1287* 1288 TEMP2 = TWO*( TWO*N-ONE )*ULP*( ONE+EIGHT*HALF**2 ) / 1289 $ ( ONE-HALF )**4 1290* 1291 TEMP1 = ZERO 1292 DO 190 J = 1, N 1293 TEMP1 = MAX( TEMP1, ABS( D4( J )-WR( N-J+1 ) ) / 1294 $ ( ABSTOL+ABS( D4( J ) ) ) ) 1295 190 CONTINUE 1296* 1297 RESULT( 17 ) = TEMP1 / TEMP2 1298 ELSE 1299 RESULT( 17 ) = ZERO 1300 END IF 1301* 1302* Now ask for all eigenvalues with high absolute accuracy. 1303* 1304 NTEST = 18 1305 ABSTOL = UNFL + UNFL 1306 CALL SSTEBZ( 'A', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, M, 1307 $ NSPLIT, WA1, IWORK( 1 ), IWORK( N+1 ), WORK, 1308 $ IWORK( 2*N+1 ), IINFO ) 1309 IF( IINFO.NE.0 ) THEN 1310 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A)', IINFO, N, JTYPE, 1311 $ IOLDSD 1312 INFO = ABS( IINFO ) 1313 IF( IINFO.LT.0 ) THEN 1314 RETURN 1315 ELSE 1316 RESULT( 18 ) = ULPINV 1317 GO TO 280 1318 END IF 1319 END IF 1320* 1321* Do test 18 1322* 1323 TEMP1 = ZERO 1324 TEMP2 = ZERO 1325 DO 200 J = 1, N 1326 TEMP1 = MAX( TEMP1, ABS( D3( J ) ), ABS( WA1( J ) ) ) 1327 TEMP2 = MAX( TEMP2, ABS( D3( J )-WA1( J ) ) ) 1328 200 CONTINUE 1329* 1330 RESULT( 18 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 1331* 1332* Choose random values for IL and IU, and ask for the 1333* IL-th through IU-th eigenvalues. 1334* 1335 NTEST = 19 1336 IF( N.LE.1 ) THEN 1337 IL = 1 1338 IU = N 1339 ELSE 1340 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1341 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1342 IF( IU.LT.IL ) THEN 1343 ITEMP = IU 1344 IU = IL 1345 IL = ITEMP 1346 END IF 1347 END IF 1348* 1349 CALL SSTEBZ( 'I', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 1350 $ M2, NSPLIT, WA2, IWORK( 1 ), IWORK( N+1 ), 1351 $ WORK, IWORK( 2*N+1 ), IINFO ) 1352 IF( IINFO.NE.0 ) THEN 1353 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(I)', IINFO, N, JTYPE, 1354 $ IOLDSD 1355 INFO = ABS( IINFO ) 1356 IF( IINFO.LT.0 ) THEN 1357 RETURN 1358 ELSE 1359 RESULT( 19 ) = ULPINV 1360 GO TO 280 1361 END IF 1362 END IF 1363* 1364* Determine the values VL and VU of the IL-th and IU-th 1365* eigenvalues and ask for all eigenvalues in this range. 1366* 1367 IF( N.GT.0 ) THEN 1368 IF( IL.NE.1 ) THEN 1369 VL = WA1( IL ) - MAX( HALF*( WA1( IL )-WA1( IL-1 ) ), 1370 $ ULP*ANORM, TWO*RTUNFL ) 1371 ELSE 1372 VL = WA1( 1 ) - MAX( HALF*( WA1( N )-WA1( 1 ) ), 1373 $ ULP*ANORM, TWO*RTUNFL ) 1374 END IF 1375 IF( IU.NE.N ) THEN 1376 VU = WA1( IU ) + MAX( HALF*( WA1( IU+1 )-WA1( IU ) ), 1377 $ ULP*ANORM, TWO*RTUNFL ) 1378 ELSE 1379 VU = WA1( N ) + MAX( HALF*( WA1( N )-WA1( 1 ) ), 1380 $ ULP*ANORM, TWO*RTUNFL ) 1381 END IF 1382 ELSE 1383 VL = ZERO 1384 VU = ONE 1385 END IF 1386* 1387 CALL SSTEBZ( 'V', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 1388 $ M3, NSPLIT, WA3, IWORK( 1 ), IWORK( N+1 ), 1389 $ WORK, IWORK( 2*N+1 ), IINFO ) 1390 IF( IINFO.NE.0 ) THEN 1391 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(V)', IINFO, N, JTYPE, 1392 $ IOLDSD 1393 INFO = ABS( IINFO ) 1394 IF( IINFO.LT.0 ) THEN 1395 RETURN 1396 ELSE 1397 RESULT( 19 ) = ULPINV 1398 GO TO 280 1399 END IF 1400 END IF 1401* 1402 IF( M3.EQ.0 .AND. N.NE.0 ) THEN 1403 RESULT( 19 ) = ULPINV 1404 GO TO 280 1405 END IF 1406* 1407* Do test 19 1408* 1409 TEMP1 = SSXT1( 1, WA2, M2, WA3, M3, ABSTOL, ULP, UNFL ) 1410 TEMP2 = SSXT1( 1, WA3, M3, WA2, M2, ABSTOL, ULP, UNFL ) 1411 IF( N.GT.0 ) THEN 1412 TEMP3 = MAX( ABS( WA1( N ) ), ABS( WA1( 1 ) ) ) 1413 ELSE 1414 TEMP3 = ZERO 1415 END IF 1416* 1417 RESULT( 19 ) = ( TEMP1+TEMP2 ) / MAX( UNFL, TEMP3*ULP ) 1418* 1419* Call SSTEIN to compute eigenvectors corresponding to 1420* eigenvalues in WA1. (First call SSTEBZ again, to make sure 1421* it returns these eigenvalues in the correct order.) 1422* 1423 NTEST = 21 1424 CALL SSTEBZ( 'A', 'B', N, VL, VU, IL, IU, ABSTOL, SD, SE, M, 1425 $ NSPLIT, WA1, IWORK( 1 ), IWORK( N+1 ), WORK, 1426 $ IWORK( 2*N+1 ), IINFO ) 1427 IF( IINFO.NE.0 ) THEN 1428 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A,B)', IINFO, N, 1429 $ JTYPE, IOLDSD 1430 INFO = ABS( IINFO ) 1431 IF( IINFO.LT.0 ) THEN 1432 RETURN 1433 ELSE 1434 RESULT( 20 ) = ULPINV 1435 RESULT( 21 ) = ULPINV 1436 GO TO 280 1437 END IF 1438 END IF 1439* 1440 CALL SSTEIN( N, SD, SE, M, WA1, IWORK( 1 ), IWORK( N+1 ), Z, 1441 $ LDU, WORK, IWORK( 2*N+1 ), IWORK( 3*N+1 ), 1442 $ IINFO ) 1443 IF( IINFO.NE.0 ) THEN 1444 WRITE( NOUNIT, FMT = 9999 )'SSTEIN', IINFO, N, JTYPE, 1445 $ IOLDSD 1446 INFO = ABS( IINFO ) 1447 IF( IINFO.LT.0 ) THEN 1448 RETURN 1449 ELSE 1450 RESULT( 20 ) = ULPINV 1451 RESULT( 21 ) = ULPINV 1452 GO TO 280 1453 END IF 1454 END IF 1455* 1456* Do tests 20 and 21 1457* 1458 CALL SSTT21( N, 0, SD, SE, WA1, DUMMA, Z, LDU, WORK, 1459 $ RESULT( 20 ) ) 1460* 1461* Call SSTEDC(I) to compute D1 and Z, do tests. 1462* 1463* Compute D1 and Z 1464* 1465 CALL SCOPY( N, SD, 1, D1, 1 ) 1466 IF( N.GT.0 ) 1467 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1468 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1469* 1470 NTEST = 22 1471 CALL SSTEDC( 'I', N, D1, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 1472 $ IWORK, LIWEDC, IINFO ) 1473 IF( IINFO.NE.0 ) THEN 1474 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(I)', IINFO, N, JTYPE, 1475 $ IOLDSD 1476 INFO = ABS( IINFO ) 1477 IF( IINFO.LT.0 ) THEN 1478 RETURN 1479 ELSE 1480 RESULT( 22 ) = ULPINV 1481 GO TO 280 1482 END IF 1483 END IF 1484* 1485* Do Tests 22 and 23 1486* 1487 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1488 $ RESULT( 22 ) ) 1489* 1490* Call SSTEDC(V) to compute D1 and Z, do tests. 1491* 1492* Compute D1 and Z 1493* 1494 CALL SCOPY( N, SD, 1, D1, 1 ) 1495 IF( N.GT.0 ) 1496 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1497 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1498* 1499 NTEST = 24 1500 CALL SSTEDC( 'V', N, D1, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 1501 $ IWORK, LIWEDC, IINFO ) 1502 IF( IINFO.NE.0 ) THEN 1503 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(V)', IINFO, N, JTYPE, 1504 $ IOLDSD 1505 INFO = ABS( IINFO ) 1506 IF( IINFO.LT.0 ) THEN 1507 RETURN 1508 ELSE 1509 RESULT( 24 ) = ULPINV 1510 GO TO 280 1511 END IF 1512 END IF 1513* 1514* Do Tests 24 and 25 1515* 1516 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1517 $ RESULT( 24 ) ) 1518* 1519* Call SSTEDC(N) to compute D2, do tests. 1520* 1521* Compute D2 1522* 1523 CALL SCOPY( N, SD, 1, D2, 1 ) 1524 IF( N.GT.0 ) 1525 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1526 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1527* 1528 NTEST = 26 1529 CALL SSTEDC( 'N', N, D2, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 1530 $ IWORK, LIWEDC, IINFO ) 1531 IF( IINFO.NE.0 ) THEN 1532 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(N)', IINFO, N, JTYPE, 1533 $ IOLDSD 1534 INFO = ABS( IINFO ) 1535 IF( IINFO.LT.0 ) THEN 1536 RETURN 1537 ELSE 1538 RESULT( 26 ) = ULPINV 1539 GO TO 280 1540 END IF 1541 END IF 1542* 1543* Do Test 26 1544* 1545 TEMP1 = ZERO 1546 TEMP2 = ZERO 1547* 1548 DO 210 J = 1, N 1549 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 1550 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1551 210 CONTINUE 1552* 1553 RESULT( 26 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 1554* 1555* Only test SSTEMR if IEEE compliant 1556* 1557 IF( ILAENV( 10, 'SSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 .AND. 1558 $ ILAENV( 11, 'SSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 ) THEN 1559* 1560* Call SSTEMR, do test 27 (relative eigenvalue accuracy) 1561* 1562* If S is positive definite and diagonally dominant, 1563* ask for all eigenvalues with high relative accuracy. 1564* 1565 VL = ZERO 1566 VU = ZERO 1567 IL = 0 1568 IU = 0 1569 IF( JTYPE.EQ.21 .AND. SREL ) THEN 1570 NTEST = 27 1571 ABSTOL = UNFL + UNFL 1572 CALL SSTEMR( 'V', 'A', N, SD, SE, VL, VU, IL, IU, 1573 $ M, WR, Z, LDU, N, IWORK( 1 ), TRYRAC, 1574 $ WORK, LWORK, IWORK( 2*N+1 ), LWORK-2*N, 1575 $ IINFO ) 1576 IF( IINFO.NE.0 ) THEN 1577 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,A,rel)', 1578 $ IINFO, N, JTYPE, IOLDSD 1579 INFO = ABS( IINFO ) 1580 IF( IINFO.LT.0 ) THEN 1581 RETURN 1582 ELSE 1583 RESULT( 27 ) = ULPINV 1584 GO TO 270 1585 END IF 1586 END IF 1587* 1588* Do test 27 1589* 1590 TEMP2 = TWO*( TWO*N-ONE )*ULP*( ONE+EIGHT*HALF**2 ) / 1591 $ ( ONE-HALF )**4 1592* 1593 TEMP1 = ZERO 1594 DO 220 J = 1, N 1595 TEMP1 = MAX( TEMP1, ABS( D4( J )-WR( N-J+1 ) ) / 1596 $ ( ABSTOL+ABS( D4( J ) ) ) ) 1597 220 CONTINUE 1598* 1599 RESULT( 27 ) = TEMP1 / TEMP2 1600* 1601 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1602 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1603 IF( IU.LT.IL ) THEN 1604 ITEMP = IU 1605 IU = IL 1606 IL = ITEMP 1607 END IF 1608* 1609 IF( SRANGE ) THEN 1610 NTEST = 28 1611 ABSTOL = UNFL + UNFL 1612 CALL SSTEMR( 'V', 'I', N, SD, SE, VL, VU, IL, IU, 1613 $ M, WR, Z, LDU, N, IWORK( 1 ), TRYRAC, 1614 $ WORK, LWORK, IWORK( 2*N+1 ), 1615 $ LWORK-2*N, IINFO ) 1616* 1617 IF( IINFO.NE.0 ) THEN 1618 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,I,rel)', 1619 $ IINFO, N, JTYPE, IOLDSD 1620 INFO = ABS( IINFO ) 1621 IF( IINFO.LT.0 ) THEN 1622 RETURN 1623 ELSE 1624 RESULT( 28 ) = ULPINV 1625 GO TO 270 1626 END IF 1627 END IF 1628* 1629* 1630* Do test 28 1631* 1632 TEMP2 = TWO*( TWO*N-ONE )*ULP* 1633 $ ( ONE+EIGHT*HALF**2 ) / ( ONE-HALF )**4 1634* 1635 TEMP1 = ZERO 1636 DO 230 J = IL, IU 1637 TEMP1 = MAX( TEMP1, ABS( WR( J-IL+1 )-D4( N-J+ 1638 $ 1 ) ) / ( ABSTOL+ABS( WR( J-IL+1 ) ) ) ) 1639 230 CONTINUE 1640* 1641 RESULT( 28 ) = TEMP1 / TEMP2 1642 ELSE 1643 RESULT( 28 ) = ZERO 1644 END IF 1645 ELSE 1646 RESULT( 27 ) = ZERO 1647 RESULT( 28 ) = ZERO 1648 END IF 1649* 1650* Call SSTEMR(V,I) to compute D1 and Z, do tests. 1651* 1652* Compute D1 and Z 1653* 1654 CALL SCOPY( N, SD, 1, D5, 1 ) 1655 IF( N.GT.0 ) 1656 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1657 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1658* 1659 IF( SRANGE ) THEN 1660 NTEST = 29 1661 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1662 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1663 IF( IU.LT.IL ) THEN 1664 ITEMP = IU 1665 IU = IL 1666 IL = ITEMP 1667 END IF 1668 CALL SSTEMR( 'V', 'I', N, D5, WORK, VL, VU, IL, IU, 1669 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 1670 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1671 $ LIWORK-2*N, IINFO ) 1672 IF( IINFO.NE.0 ) THEN 1673 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,I)', IINFO, 1674 $ N, JTYPE, IOLDSD 1675 INFO = ABS( IINFO ) 1676 IF( IINFO.LT.0 ) THEN 1677 RETURN 1678 ELSE 1679 RESULT( 29 ) = ULPINV 1680 GO TO 280 1681 END IF 1682 END IF 1683* 1684* Do Tests 29 and 30 1685* 1686 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1687 $ M, RESULT( 29 ) ) 1688* 1689* Call SSTEMR to compute D2, do tests. 1690* 1691* Compute D2 1692* 1693 CALL SCOPY( N, SD, 1, D5, 1 ) 1694 IF( N.GT.0 ) 1695 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1696* 1697 NTEST = 31 1698 CALL SSTEMR( 'N', 'I', N, D5, WORK, VL, VU, IL, IU, 1699 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 1700 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1701 $ LIWORK-2*N, IINFO ) 1702 IF( IINFO.NE.0 ) THEN 1703 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,I)', IINFO, 1704 $ N, JTYPE, IOLDSD 1705 INFO = ABS( IINFO ) 1706 IF( IINFO.LT.0 ) THEN 1707 RETURN 1708 ELSE 1709 RESULT( 31 ) = ULPINV 1710 GO TO 280 1711 END IF 1712 END IF 1713* 1714* Do Test 31 1715* 1716 TEMP1 = ZERO 1717 TEMP2 = ZERO 1718* 1719 DO 240 J = 1, IU - IL + 1 1720 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), 1721 $ ABS( D2( J ) ) ) 1722 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1723 240 CONTINUE 1724* 1725 RESULT( 31 ) = TEMP2 / MAX( UNFL, 1726 $ ULP*MAX( TEMP1, TEMP2 ) ) 1727* 1728* 1729* Call SSTEMR(V,V) to compute D1 and Z, do tests. 1730* 1731* Compute D1 and Z 1732* 1733 CALL SCOPY( N, SD, 1, D5, 1 ) 1734 IF( N.GT.0 ) 1735 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1736 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1737* 1738 NTEST = 32 1739* 1740 IF( N.GT.0 ) THEN 1741 IF( IL.NE.1 ) THEN 1742 VL = D2( IL ) - MAX( HALF* 1743 $ ( D2( IL )-D2( IL-1 ) ), ULP*ANORM, 1744 $ TWO*RTUNFL ) 1745 ELSE 1746 VL = D2( 1 ) - MAX( HALF*( D2( N )-D2( 1 ) ), 1747 $ ULP*ANORM, TWO*RTUNFL ) 1748 END IF 1749 IF( IU.NE.N ) THEN 1750 VU = D2( IU ) + MAX( HALF* 1751 $ ( D2( IU+1 )-D2( IU ) ), ULP*ANORM, 1752 $ TWO*RTUNFL ) 1753 ELSE 1754 VU = D2( N ) + MAX( HALF*( D2( N )-D2( 1 ) ), 1755 $ ULP*ANORM, TWO*RTUNFL ) 1756 END IF 1757 ELSE 1758 VL = ZERO 1759 VU = ONE 1760 END IF 1761* 1762 CALL SSTEMR( 'V', 'V', N, D5, WORK, VL, VU, IL, IU, 1763 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 1764 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1765 $ LIWORK-2*N, IINFO ) 1766 IF( IINFO.NE.0 ) THEN 1767 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,V)', IINFO, 1768 $ N, JTYPE, IOLDSD 1769 INFO = ABS( IINFO ) 1770 IF( IINFO.LT.0 ) THEN 1771 RETURN 1772 ELSE 1773 RESULT( 32 ) = ULPINV 1774 GO TO 280 1775 END IF 1776 END IF 1777* 1778* Do Tests 32 and 33 1779* 1780 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1781 $ M, RESULT( 32 ) ) 1782* 1783* Call SSTEMR to compute D2, do tests. 1784* 1785* Compute D2 1786* 1787 CALL SCOPY( N, SD, 1, D5, 1 ) 1788 IF( N.GT.0 ) 1789 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1790* 1791 NTEST = 34 1792 CALL SSTEMR( 'N', 'V', N, D5, WORK, VL, VU, IL, IU, 1793 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 1794 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1795 $ LIWORK-2*N, IINFO ) 1796 IF( IINFO.NE.0 ) THEN 1797 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,V)', IINFO, 1798 $ N, JTYPE, IOLDSD 1799 INFO = ABS( IINFO ) 1800 IF( IINFO.LT.0 ) THEN 1801 RETURN 1802 ELSE 1803 RESULT( 34 ) = ULPINV 1804 GO TO 280 1805 END IF 1806 END IF 1807* 1808* Do Test 34 1809* 1810 TEMP1 = ZERO 1811 TEMP2 = ZERO 1812* 1813 DO 250 J = 1, IU - IL + 1 1814 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), 1815 $ ABS( D2( J ) ) ) 1816 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1817 250 CONTINUE 1818* 1819 RESULT( 34 ) = TEMP2 / MAX( UNFL, 1820 $ ULP*MAX( TEMP1, TEMP2 ) ) 1821 ELSE 1822 RESULT( 29 ) = ZERO 1823 RESULT( 30 ) = ZERO 1824 RESULT( 31 ) = ZERO 1825 RESULT( 32 ) = ZERO 1826 RESULT( 33 ) = ZERO 1827 RESULT( 34 ) = ZERO 1828 END IF 1829* 1830* 1831* Call SSTEMR(V,A) to compute D1 and Z, do tests. 1832* 1833* Compute D1 and Z 1834* 1835 CALL SCOPY( N, SD, 1, D5, 1 ) 1836 IF( N.GT.0 ) 1837 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1838* 1839 NTEST = 35 1840* 1841 CALL SSTEMR( 'V', 'A', N, D5, WORK, VL, VU, IL, IU, 1842 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 1843 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1844 $ LIWORK-2*N, IINFO ) 1845 IF( IINFO.NE.0 ) THEN 1846 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,A)', IINFO, N, 1847 $ JTYPE, IOLDSD 1848 INFO = ABS( IINFO ) 1849 IF( IINFO.LT.0 ) THEN 1850 RETURN 1851 ELSE 1852 RESULT( 35 ) = ULPINV 1853 GO TO 280 1854 END IF 1855 END IF 1856* 1857* Do Tests 35 and 36 1858* 1859 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, M, 1860 $ RESULT( 35 ) ) 1861* 1862* Call SSTEMR to compute D2, do tests. 1863* 1864* Compute D2 1865* 1866 CALL SCOPY( N, SD, 1, D5, 1 ) 1867 IF( N.GT.0 ) 1868 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1869* 1870 NTEST = 37 1871 CALL SSTEMR( 'N', 'A', N, D5, WORK, VL, VU, IL, IU, 1872 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 1873 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1874 $ LIWORK-2*N, IINFO ) 1875 IF( IINFO.NE.0 ) THEN 1876 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,A)', IINFO, N, 1877 $ JTYPE, IOLDSD 1878 INFO = ABS( IINFO ) 1879 IF( IINFO.LT.0 ) THEN 1880 RETURN 1881 ELSE 1882 RESULT( 37 ) = ULPINV 1883 GO TO 280 1884 END IF 1885 END IF 1886* 1887* Do Test 34 1888* 1889 TEMP1 = ZERO 1890 TEMP2 = ZERO 1891* 1892 DO 260 J = 1, N 1893 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 1894 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1895 260 CONTINUE 1896* 1897 RESULT( 37 ) = TEMP2 / MAX( UNFL, 1898 $ ULP*MAX( TEMP1, TEMP2 ) ) 1899 END IF 1900 270 CONTINUE 1901 280 CONTINUE 1902 NTESTT = NTESTT + NTEST 1903* 1904* End of Loop -- Check for RESULT(j) > THRESH 1905* 1906* 1907* Print out tests which fail. 1908* 1909 DO 290 JR = 1, NTEST 1910 IF( RESULT( JR ).GE.THRESH ) THEN 1911* 1912* If this is the first test to fail, 1913* print a header to the data file. 1914* 1915 IF( NERRS.EQ.0 ) THEN 1916 WRITE( NOUNIT, FMT = 9998 )'SST' 1917 WRITE( NOUNIT, FMT = 9997 ) 1918 WRITE( NOUNIT, FMT = 9996 ) 1919 WRITE( NOUNIT, FMT = 9995 )'Symmetric' 1920 WRITE( NOUNIT, FMT = 9994 ) 1921* 1922* Tests performed 1923* 1924 WRITE( NOUNIT, FMT = 9988 ) 1925 END IF 1926 NERRS = NERRS + 1 1927 WRITE( NOUNIT, FMT = 9990 )N, IOLDSD, JTYPE, JR, 1928 $ RESULT( JR ) 1929 END IF 1930 290 CONTINUE 1931 300 CONTINUE 1932 310 CONTINUE 1933* 1934* Summary 1935* 1936 CALL SLASUM( 'SST', NOUNIT, NERRS, NTESTT ) 1937 RETURN 1938* 1939 9999 FORMAT( ' SCHKST: ', A, ' returned INFO=', I6, '.', / 9X, 'N=', 1940 $ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' ) 1941* 1942 9998 FORMAT( / 1X, A3, ' -- Real Symmetric eigenvalue problem' ) 1943 9997 FORMAT( ' Matrix types (see SCHKST for details): ' ) 1944* 1945 9996 FORMAT( / ' Special Matrices:', 1946 $ / ' 1=Zero matrix. ', 1947 $ ' 5=Diagonal: clustered entries.', 1948 $ / ' 2=Identity matrix. ', 1949 $ ' 6=Diagonal: large, evenly spaced.', 1950 $ / ' 3=Diagonal: evenly spaced entries. ', 1951 $ ' 7=Diagonal: small, evenly spaced.', 1952 $ / ' 4=Diagonal: geometr. spaced entries.' ) 1953 9995 FORMAT( ' Dense ', A, ' Matrices:', 1954 $ / ' 8=Evenly spaced eigenvals. ', 1955 $ ' 12=Small, evenly spaced eigenvals.', 1956 $ / ' 9=Geometrically spaced eigenvals. ', 1957 $ ' 13=Matrix with random O(1) entries.', 1958 $ / ' 10=Clustered eigenvalues. ', 1959 $ ' 14=Matrix with large random entries.', 1960 $ / ' 11=Large, evenly spaced eigenvals. ', 1961 $ ' 15=Matrix with small random entries.' ) 1962 9994 FORMAT( ' 16=Positive definite, evenly spaced eigenvalues', 1963 $ / ' 17=Positive definite, geometrically spaced eigenvlaues', 1964 $ / ' 18=Positive definite, clustered eigenvalues', 1965 $ / ' 19=Positive definite, small evenly spaced eigenvalues', 1966 $ / ' 20=Positive definite, large evenly spaced eigenvalues', 1967 $ / ' 21=Diagonally dominant tridiagonal, geometrically', 1968 $ ' spaced eigenvalues' ) 1969* 1970 9990 FORMAT( ' N=', I5, ', seed=', 4( I4, ',' ), ' type ', I2, 1971 $ ', test(', I2, ')=', G10.3 ) 1972* 1973 9988 FORMAT( / 'Test performed: see SCHKST for details.', / ) 1974* End of SCHKST 1975* 1976 END 1977