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*> \date December 2016 585* 586*> \ingroup single_eig 587* 588* ===================================================================== 589 SUBROUTINE SCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, 590 $ NOUNIT, A, LDA, AP, SD, SE, D1, D2, D3, D4, D5, 591 $ WA1, WA2, WA3, WR, U, LDU, V, VP, TAU, Z, WORK, 592 $ LWORK, IWORK, LIWORK, RESULT, INFO ) 593* 594* -- LAPACK test routine (version 3.7.0) -- 595* -- LAPACK is a software package provided by Univ. of Tennessee, -- 596* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 597* December 2016 598* 599* .. Scalar Arguments .. 600 INTEGER INFO, LDA, LDU, LIWORK, LWORK, NOUNIT, NSIZES, 601 $ NTYPES 602 REAL THRESH 603* .. 604* .. Array Arguments .. 605 LOGICAL DOTYPE( * ) 606 INTEGER ISEED( 4 ), IWORK( * ), NN( * ) 607 REAL A( LDA, * ), AP( * ), D1( * ), D2( * ), 608 $ D3( * ), D4( * ), D5( * ), RESULT( * ), 609 $ SD( * ), SE( * ), TAU( * ), U( LDU, * ), 610 $ V( LDU, * ), VP( * ), WA1( * ), WA2( * ), 611 $ WA3( * ), WORK( * ), WR( * ), Z( LDU, * ) 612* .. 613* 614* ===================================================================== 615* 616* .. Parameters .. 617 REAL ZERO, ONE, TWO, EIGHT, TEN, HUN 618 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0, 619 $ EIGHT = 8.0E0, TEN = 10.0E0, HUN = 100.0E0 ) 620 REAL HALF 621 PARAMETER ( HALF = ONE / TWO ) 622 INTEGER MAXTYP 623 PARAMETER ( MAXTYP = 21 ) 624 LOGICAL SRANGE 625 PARAMETER ( SRANGE = .FALSE. ) 626 LOGICAL SREL 627 PARAMETER ( SREL = .FALSE. ) 628* .. 629* .. Local Scalars .. 630 LOGICAL BADNN, TRYRAC 631 INTEGER I, IINFO, IL, IMODE, ITEMP, ITYPE, IU, J, JC, 632 $ JR, JSIZE, JTYPE, LGN, LIWEDC, LOG2UI, LWEDC, 633 $ M, M2, M3, MTYPES, N, NAP, NBLOCK, NERRS, 634 $ NMATS, NMAX, NSPLIT, NTEST, NTESTT 635 REAL ABSTOL, ANINV, ANORM, COND, OVFL, RTOVFL, 636 $ RTUNFL, TEMP1, TEMP2, TEMP3, TEMP4, ULP, 637 $ ULPINV, UNFL, VL, VU 638* .. 639* .. Local Arrays .. 640 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ), 641 $ KMAGN( MAXTYP ), KMODE( MAXTYP ), 642 $ KTYPE( MAXTYP ) 643 REAL DUMMA( 1 ) 644* .. 645* .. External Functions .. 646 INTEGER ILAENV 647 REAL SLAMCH, SLARND, SSXT1 648 EXTERNAL ILAENV, SLAMCH, SLARND, SSXT1 649* .. 650* .. External Subroutines .. 651 EXTERNAL SCOPY, SLABAD, SLACPY, SLASET, SLASUM, SLATMR, 652 $ SLATMS, SOPGTR, SORGTR, SPTEQR, SSPT21, SSPTRD, 653 $ SSTEBZ, SSTECH, SSTEDC, SSTEMR, SSTEIN, SSTEQR, 654 $ SSTERF, SSTT21, SSTT22, SSYT21, SSYTRD, XERBLA 655* .. 656* .. Intrinsic Functions .. 657 INTRINSIC ABS, INT, LOG, MAX, MIN, REAL, SQRT 658* .. 659* .. Data statements .. 660 DATA KTYPE / 1, 2, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 8, 661 $ 8, 8, 9, 9, 9, 9, 9, 10 / 662 DATA KMAGN / 1, 1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1, 663 $ 2, 3, 1, 1, 1, 2, 3, 1 / 664 DATA KMODE / 0, 0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0, 665 $ 0, 0, 4, 3, 1, 4, 4, 3 / 666* .. 667* .. Executable Statements .. 668* 669* Keep ftnchek happy 670 IDUMMA( 1 ) = 1 671* 672* Check for errors 673* 674 NTESTT = 0 675 INFO = 0 676* 677* Important constants 678* 679 BADNN = .FALSE. 680 TRYRAC = .TRUE. 681 NMAX = 1 682 DO 10 J = 1, NSIZES 683 NMAX = MAX( NMAX, NN( J ) ) 684 IF( NN( J ).LT.0 ) 685 $ BADNN = .TRUE. 686 10 CONTINUE 687* 688 NBLOCK = ILAENV( 1, 'SSYTRD', 'L', NMAX, -1, -1, -1 ) 689 NBLOCK = MIN( NMAX, MAX( 1, NBLOCK ) ) 690* 691* Check for errors 692* 693 IF( NSIZES.LT.0 ) THEN 694 INFO = -1 695 ELSE IF( BADNN ) THEN 696 INFO = -2 697 ELSE IF( NTYPES.LT.0 ) THEN 698 INFO = -3 699 ELSE IF( LDA.LT.NMAX ) THEN 700 INFO = -9 701 ELSE IF( LDU.LT.NMAX ) THEN 702 INFO = -23 703 ELSE IF( 2*MAX( 2, NMAX )**2.GT.LWORK ) THEN 704 INFO = -29 705 END IF 706* 707 IF( INFO.NE.0 ) THEN 708 CALL XERBLA( 'SCHKST', -INFO ) 709 RETURN 710 END IF 711* 712* Quick return if possible 713* 714 IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 ) 715 $ RETURN 716* 717* More Important constants 718* 719 UNFL = SLAMCH( 'Safe minimum' ) 720 OVFL = ONE / UNFL 721 CALL SLABAD( UNFL, OVFL ) 722 ULP = SLAMCH( 'Epsilon' )*SLAMCH( 'Base' ) 723 ULPINV = ONE / ULP 724 LOG2UI = INT( LOG( ULPINV ) / LOG( TWO ) ) 725 RTUNFL = SQRT( UNFL ) 726 RTOVFL = SQRT( OVFL ) 727* 728* Loop over sizes, types 729* 730 DO 20 I = 1, 4 731 ISEED2( I ) = ISEED( I ) 732 20 CONTINUE 733 NERRS = 0 734 NMATS = 0 735* 736 DO 310 JSIZE = 1, NSIZES 737 N = NN( JSIZE ) 738 IF( N.GT.0 ) THEN 739 LGN = INT( LOG( REAL( N ) ) / LOG( TWO ) ) 740 IF( 2**LGN.LT.N ) 741 $ LGN = LGN + 1 742 IF( 2**LGN.LT.N ) 743 $ LGN = LGN + 1 744 LWEDC = 1 + 4*N + 2*N*LGN + 4*N**2 745 LIWEDC = 6 + 6*N + 5*N*LGN 746 ELSE 747 LWEDC = 8 748 LIWEDC = 12 749 END IF 750 NAP = ( N*( N+1 ) ) / 2 751 ANINV = ONE / REAL( MAX( 1, N ) ) 752* 753 IF( NSIZES.NE.1 ) THEN 754 MTYPES = MIN( MAXTYP, NTYPES ) 755 ELSE 756 MTYPES = MIN( MAXTYP+1, NTYPES ) 757 END IF 758* 759 DO 300 JTYPE = 1, MTYPES 760 IF( .NOT.DOTYPE( JTYPE ) ) 761 $ GO TO 300 762 NMATS = NMATS + 1 763 NTEST = 0 764* 765 DO 30 J = 1, 4 766 IOLDSD( J ) = ISEED( J ) 767 30 CONTINUE 768* 769* Compute "A" 770* 771* Control parameters: 772* 773* KMAGN KMODE KTYPE 774* =1 O(1) clustered 1 zero 775* =2 large clustered 2 identity 776* =3 small exponential (none) 777* =4 arithmetic diagonal, (w/ eigenvalues) 778* =5 random log symmetric, w/ eigenvalues 779* =6 random (none) 780* =7 random diagonal 781* =8 random symmetric 782* =9 positive definite 783* =10 diagonally dominant tridiagonal 784* 785 IF( MTYPES.GT.MAXTYP ) 786 $ GO TO 100 787* 788 ITYPE = KTYPE( JTYPE ) 789 IMODE = KMODE( JTYPE ) 790* 791* Compute norm 792* 793 GO TO ( 40, 50, 60 )KMAGN( JTYPE ) 794* 795 40 CONTINUE 796 ANORM = ONE 797 GO TO 70 798* 799 50 CONTINUE 800 ANORM = ( RTOVFL*ULP )*ANINV 801 GO TO 70 802* 803 60 CONTINUE 804 ANORM = RTUNFL*N*ULPINV 805 GO TO 70 806* 807 70 CONTINUE 808* 809 CALL SLASET( 'Full', LDA, N, ZERO, ZERO, A, LDA ) 810 IINFO = 0 811 IF( JTYPE.LE.15 ) THEN 812 COND = ULPINV 813 ELSE 814 COND = ULPINV*ANINV / TEN 815 END IF 816* 817* Special Matrices -- Identity & Jordan block 818* 819* Zero 820* 821 IF( ITYPE.EQ.1 ) THEN 822 IINFO = 0 823* 824 ELSE IF( ITYPE.EQ.2 ) THEN 825* 826* Identity 827* 828 DO 80 JC = 1, N 829 A( JC, JC ) = ANORM 830 80 CONTINUE 831* 832 ELSE IF( ITYPE.EQ.4 ) THEN 833* 834* Diagonal Matrix, [Eigen]values Specified 835* 836 CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND, 837 $ ANORM, 0, 0, 'N', A, LDA, WORK( N+1 ), 838 $ IINFO ) 839* 840* 841 ELSE IF( ITYPE.EQ.5 ) THEN 842* 843* Symmetric, eigenvalues specified 844* 845 CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND, 846 $ ANORM, N, N, 'N', A, LDA, WORK( N+1 ), 847 $ IINFO ) 848* 849 ELSE IF( ITYPE.EQ.7 ) THEN 850* 851* Diagonal, random eigenvalues 852* 853 CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE, 854 $ 'T', 'N', WORK( N+1 ), 1, ONE, 855 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0, 856 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO ) 857* 858 ELSE IF( ITYPE.EQ.8 ) THEN 859* 860* Symmetric, random eigenvalues 861* 862 CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE, 863 $ 'T', 'N', WORK( N+1 ), 1, ONE, 864 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N, 865 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO ) 866* 867 ELSE IF( ITYPE.EQ.9 ) THEN 868* 869* Positive definite, eigenvalues specified. 870* 871 CALL SLATMS( N, N, 'S', ISEED, 'P', WORK, IMODE, COND, 872 $ ANORM, N, N, 'N', A, LDA, WORK( N+1 ), 873 $ IINFO ) 874* 875 ELSE IF( ITYPE.EQ.10 ) THEN 876* 877* Positive definite tridiagonal, eigenvalues specified. 878* 879 CALL SLATMS( N, N, 'S', ISEED, 'P', WORK, IMODE, COND, 880 $ ANORM, 1, 1, 'N', A, LDA, WORK( N+1 ), 881 $ IINFO ) 882 DO 90 I = 2, N 883 TEMP1 = ABS( A( I-1, I ) ) / 884 $ SQRT( ABS( A( I-1, I-1 )*A( I, I ) ) ) 885 IF( TEMP1.GT.HALF ) THEN 886 A( I-1, I ) = HALF*SQRT( ABS( A( I-1, I-1 )*A( I, 887 $ I ) ) ) 888 A( I, I-1 ) = A( I-1, I ) 889 END IF 890 90 CONTINUE 891* 892 ELSE 893* 894 IINFO = 1 895 END IF 896* 897 IF( IINFO.NE.0 ) THEN 898 WRITE( NOUNIT, FMT = 9999 )'Generator', IINFO, N, JTYPE, 899 $ IOLDSD 900 INFO = ABS( IINFO ) 901 RETURN 902 END IF 903* 904 100 CONTINUE 905* 906* Call SSYTRD and SORGTR to compute S and U from 907* upper triangle. 908* 909 CALL SLACPY( 'U', N, N, A, LDA, V, LDU ) 910* 911 NTEST = 1 912 CALL SSYTRD( 'U', N, V, LDU, SD, SE, TAU, WORK, LWORK, 913 $ IINFO ) 914* 915 IF( IINFO.NE.0 ) THEN 916 WRITE( NOUNIT, FMT = 9999 )'SSYTRD(U)', IINFO, N, JTYPE, 917 $ IOLDSD 918 INFO = ABS( IINFO ) 919 IF( IINFO.LT.0 ) THEN 920 RETURN 921 ELSE 922 RESULT( 1 ) = ULPINV 923 GO TO 280 924 END IF 925 END IF 926* 927 CALL SLACPY( 'U', N, N, V, LDU, U, LDU ) 928* 929 NTEST = 2 930 CALL SORGTR( 'U', N, U, LDU, TAU, WORK, LWORK, IINFO ) 931 IF( IINFO.NE.0 ) THEN 932 WRITE( NOUNIT, FMT = 9999 )'SORGTR(U)', IINFO, N, JTYPE, 933 $ IOLDSD 934 INFO = ABS( IINFO ) 935 IF( IINFO.LT.0 ) THEN 936 RETURN 937 ELSE 938 RESULT( 2 ) = ULPINV 939 GO TO 280 940 END IF 941 END IF 942* 943* Do tests 1 and 2 944* 945 CALL SSYT21( 2, 'Upper', N, 1, A, LDA, SD, SE, U, LDU, V, 946 $ LDU, TAU, WORK, RESULT( 1 ) ) 947 CALL SSYT21( 3, 'Upper', N, 1, A, LDA, SD, SE, U, LDU, V, 948 $ LDU, TAU, WORK, RESULT( 2 ) ) 949* 950* Call SSYTRD and SORGTR to compute S and U from 951* lower triangle, do tests. 952* 953 CALL SLACPY( 'L', N, N, A, LDA, V, LDU ) 954* 955 NTEST = 3 956 CALL SSYTRD( 'L', N, V, LDU, SD, SE, TAU, WORK, LWORK, 957 $ IINFO ) 958* 959 IF( IINFO.NE.0 ) THEN 960 WRITE( NOUNIT, FMT = 9999 )'SSYTRD(L)', IINFO, N, JTYPE, 961 $ IOLDSD 962 INFO = ABS( IINFO ) 963 IF( IINFO.LT.0 ) THEN 964 RETURN 965 ELSE 966 RESULT( 3 ) = ULPINV 967 GO TO 280 968 END IF 969 END IF 970* 971 CALL SLACPY( 'L', N, N, V, LDU, U, LDU ) 972* 973 NTEST = 4 974 CALL SORGTR( 'L', N, U, LDU, TAU, WORK, LWORK, IINFO ) 975 IF( IINFO.NE.0 ) THEN 976 WRITE( NOUNIT, FMT = 9999 )'SORGTR(L)', IINFO, N, JTYPE, 977 $ IOLDSD 978 INFO = ABS( IINFO ) 979 IF( IINFO.LT.0 ) THEN 980 RETURN 981 ELSE 982 RESULT( 4 ) = ULPINV 983 GO TO 280 984 END IF 985 END IF 986* 987 CALL SSYT21( 2, 'Lower', N, 1, A, LDA, SD, SE, U, LDU, V, 988 $ LDU, TAU, WORK, RESULT( 3 ) ) 989 CALL SSYT21( 3, 'Lower', N, 1, A, LDA, SD, SE, U, LDU, V, 990 $ LDU, TAU, WORK, RESULT( 4 ) ) 991* 992* Store the upper triangle of A in AP 993* 994 I = 0 995 DO 120 JC = 1, N 996 DO 110 JR = 1, JC 997 I = I + 1 998 AP( I ) = A( JR, JC ) 999 110 CONTINUE 1000 120 CONTINUE 1001* 1002* Call SSPTRD and SOPGTR to compute S and U from AP 1003* 1004 CALL SCOPY( NAP, AP, 1, VP, 1 ) 1005* 1006 NTEST = 5 1007 CALL SSPTRD( 'U', N, VP, SD, SE, TAU, IINFO ) 1008* 1009 IF( IINFO.NE.0 ) THEN 1010 WRITE( NOUNIT, FMT = 9999 )'SSPTRD(U)', IINFO, N, JTYPE, 1011 $ IOLDSD 1012 INFO = ABS( IINFO ) 1013 IF( IINFO.LT.0 ) THEN 1014 RETURN 1015 ELSE 1016 RESULT( 5 ) = ULPINV 1017 GO TO 280 1018 END IF 1019 END IF 1020* 1021 NTEST = 6 1022 CALL SOPGTR( 'U', N, VP, TAU, U, LDU, WORK, IINFO ) 1023 IF( IINFO.NE.0 ) THEN 1024 WRITE( NOUNIT, FMT = 9999 )'SOPGTR(U)', IINFO, N, JTYPE, 1025 $ IOLDSD 1026 INFO = ABS( IINFO ) 1027 IF( IINFO.LT.0 ) THEN 1028 RETURN 1029 ELSE 1030 RESULT( 6 ) = ULPINV 1031 GO TO 280 1032 END IF 1033 END IF 1034* 1035* Do tests 5 and 6 1036* 1037 CALL SSPT21( 2, 'Upper', N, 1, AP, SD, SE, U, LDU, VP, TAU, 1038 $ WORK, RESULT( 5 ) ) 1039 CALL SSPT21( 3, 'Upper', N, 1, AP, SD, SE, U, LDU, VP, TAU, 1040 $ WORK, RESULT( 6 ) ) 1041* 1042* Store the lower triangle of A in AP 1043* 1044 I = 0 1045 DO 140 JC = 1, N 1046 DO 130 JR = JC, N 1047 I = I + 1 1048 AP( I ) = A( JR, JC ) 1049 130 CONTINUE 1050 140 CONTINUE 1051* 1052* Call SSPTRD and SOPGTR to compute S and U from AP 1053* 1054 CALL SCOPY( NAP, AP, 1, VP, 1 ) 1055* 1056 NTEST = 7 1057 CALL SSPTRD( 'L', N, VP, SD, SE, TAU, IINFO ) 1058* 1059 IF( IINFO.NE.0 ) THEN 1060 WRITE( NOUNIT, FMT = 9999 )'SSPTRD(L)', IINFO, N, JTYPE, 1061 $ IOLDSD 1062 INFO = ABS( IINFO ) 1063 IF( IINFO.LT.0 ) THEN 1064 RETURN 1065 ELSE 1066 RESULT( 7 ) = ULPINV 1067 GO TO 280 1068 END IF 1069 END IF 1070* 1071 NTEST = 8 1072 CALL SOPGTR( 'L', N, VP, TAU, U, LDU, WORK, IINFO ) 1073 IF( IINFO.NE.0 ) THEN 1074 WRITE( NOUNIT, FMT = 9999 )'SOPGTR(L)', IINFO, N, JTYPE, 1075 $ IOLDSD 1076 INFO = ABS( IINFO ) 1077 IF( IINFO.LT.0 ) THEN 1078 RETURN 1079 ELSE 1080 RESULT( 8 ) = ULPINV 1081 GO TO 280 1082 END IF 1083 END IF 1084* 1085 CALL SSPT21( 2, 'Lower', N, 1, AP, SD, SE, U, LDU, VP, TAU, 1086 $ WORK, RESULT( 7 ) ) 1087 CALL SSPT21( 3, 'Lower', N, 1, AP, SD, SE, U, LDU, VP, TAU, 1088 $ WORK, RESULT( 8 ) ) 1089* 1090* Call SSTEQR to compute D1, D2, and Z, do tests. 1091* 1092* Compute D1 and Z 1093* 1094 CALL SCOPY( N, SD, 1, D1, 1 ) 1095 IF( N.GT.0 ) 1096 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1097 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1098* 1099 NTEST = 9 1100 CALL SSTEQR( 'V', N, D1, WORK, Z, LDU, WORK( N+1 ), IINFO ) 1101 IF( IINFO.NE.0 ) THEN 1102 WRITE( NOUNIT, FMT = 9999 )'SSTEQR(V)', IINFO, N, JTYPE, 1103 $ IOLDSD 1104 INFO = ABS( IINFO ) 1105 IF( IINFO.LT.0 ) THEN 1106 RETURN 1107 ELSE 1108 RESULT( 9 ) = ULPINV 1109 GO TO 280 1110 END IF 1111 END IF 1112* 1113* Compute D2 1114* 1115 CALL SCOPY( N, SD, 1, D2, 1 ) 1116 IF( N.GT.0 ) 1117 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1118* 1119 NTEST = 11 1120 CALL SSTEQR( 'N', N, D2, WORK, WORK( N+1 ), LDU, 1121 $ WORK( N+1 ), IINFO ) 1122 IF( IINFO.NE.0 ) THEN 1123 WRITE( NOUNIT, FMT = 9999 )'SSTEQR(N)', IINFO, N, JTYPE, 1124 $ IOLDSD 1125 INFO = ABS( IINFO ) 1126 IF( IINFO.LT.0 ) THEN 1127 RETURN 1128 ELSE 1129 RESULT( 11 ) = ULPINV 1130 GO TO 280 1131 END IF 1132 END IF 1133* 1134* Compute D3 (using PWK method) 1135* 1136 CALL SCOPY( N, SD, 1, D3, 1 ) 1137 IF( N.GT.0 ) 1138 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1139* 1140 NTEST = 12 1141 CALL SSTERF( N, D3, WORK, IINFO ) 1142 IF( IINFO.NE.0 ) THEN 1143 WRITE( NOUNIT, FMT = 9999 )'SSTERF', IINFO, N, JTYPE, 1144 $ IOLDSD 1145 INFO = ABS( IINFO ) 1146 IF( IINFO.LT.0 ) THEN 1147 RETURN 1148 ELSE 1149 RESULT( 12 ) = ULPINV 1150 GO TO 280 1151 END IF 1152 END IF 1153* 1154* Do Tests 9 and 10 1155* 1156 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1157 $ RESULT( 9 ) ) 1158* 1159* Do Tests 11 and 12 1160* 1161 TEMP1 = ZERO 1162 TEMP2 = ZERO 1163 TEMP3 = ZERO 1164 TEMP4 = ZERO 1165* 1166 DO 150 J = 1, N 1167 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 1168 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1169 TEMP3 = MAX( TEMP3, ABS( D1( J ) ), ABS( D3( J ) ) ) 1170 TEMP4 = MAX( TEMP4, ABS( D1( J )-D3( J ) ) ) 1171 150 CONTINUE 1172* 1173 RESULT( 11 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 1174 RESULT( 12 ) = TEMP4 / MAX( UNFL, ULP*MAX( TEMP3, TEMP4 ) ) 1175* 1176* Do Test 13 -- Sturm Sequence Test of Eigenvalues 1177* Go up by factors of two until it succeeds 1178* 1179 NTEST = 13 1180 TEMP1 = THRESH*( HALF-ULP ) 1181* 1182 DO 160 J = 0, LOG2UI 1183 CALL SSTECH( N, SD, SE, D1, TEMP1, WORK, IINFO ) 1184 IF( IINFO.EQ.0 ) 1185 $ GO TO 170 1186 TEMP1 = TEMP1*TWO 1187 160 CONTINUE 1188* 1189 170 CONTINUE 1190 RESULT( 13 ) = TEMP1 1191* 1192* For positive definite matrices ( JTYPE.GT.15 ) call SPTEQR 1193* and do tests 14, 15, and 16 . 1194* 1195 IF( JTYPE.GT.15 ) THEN 1196* 1197* Compute D4 and Z4 1198* 1199 CALL SCOPY( N, SD, 1, D4, 1 ) 1200 IF( N.GT.0 ) 1201 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1202 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1203* 1204 NTEST = 14 1205 CALL SPTEQR( 'V', N, D4, WORK, Z, LDU, WORK( N+1 ), 1206 $ IINFO ) 1207 IF( IINFO.NE.0 ) THEN 1208 WRITE( NOUNIT, FMT = 9999 )'SPTEQR(V)', IINFO, N, 1209 $ JTYPE, IOLDSD 1210 INFO = ABS( IINFO ) 1211 IF( IINFO.LT.0 ) THEN 1212 RETURN 1213 ELSE 1214 RESULT( 14 ) = ULPINV 1215 GO TO 280 1216 END IF 1217 END IF 1218* 1219* Do Tests 14 and 15 1220* 1221 CALL SSTT21( N, 0, SD, SE, D4, DUMMA, Z, LDU, WORK, 1222 $ RESULT( 14 ) ) 1223* 1224* Compute D5 1225* 1226 CALL SCOPY( N, SD, 1, D5, 1 ) 1227 IF( N.GT.0 ) 1228 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1229* 1230 NTEST = 16 1231 CALL SPTEQR( 'N', N, D5, WORK, Z, LDU, WORK( N+1 ), 1232 $ IINFO ) 1233 IF( IINFO.NE.0 ) THEN 1234 WRITE( NOUNIT, FMT = 9999 )'SPTEQR(N)', IINFO, N, 1235 $ JTYPE, IOLDSD 1236 INFO = ABS( IINFO ) 1237 IF( IINFO.LT.0 ) THEN 1238 RETURN 1239 ELSE 1240 RESULT( 16 ) = ULPINV 1241 GO TO 280 1242 END IF 1243 END IF 1244* 1245* Do Test 16 1246* 1247 TEMP1 = ZERO 1248 TEMP2 = ZERO 1249 DO 180 J = 1, N 1250 TEMP1 = MAX( TEMP1, ABS( D4( J ) ), ABS( D5( J ) ) ) 1251 TEMP2 = MAX( TEMP2, ABS( D4( J )-D5( J ) ) ) 1252 180 CONTINUE 1253* 1254 RESULT( 16 ) = TEMP2 / MAX( UNFL, 1255 $ HUN*ULP*MAX( TEMP1, TEMP2 ) ) 1256 ELSE 1257 RESULT( 14 ) = ZERO 1258 RESULT( 15 ) = ZERO 1259 RESULT( 16 ) = ZERO 1260 END IF 1261* 1262* Call SSTEBZ with different options and do tests 17-18. 1263* 1264* If S is positive definite and diagonally dominant, 1265* ask for all eigenvalues with high relative accuracy. 1266* 1267 VL = ZERO 1268 VU = ZERO 1269 IL = 0 1270 IU = 0 1271 IF( JTYPE.EQ.21 ) THEN 1272 NTEST = 17 1273 ABSTOL = UNFL + UNFL 1274 CALL SSTEBZ( 'A', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 1275 $ M, NSPLIT, WR, IWORK( 1 ), IWORK( N+1 ), 1276 $ WORK, IWORK( 2*N+1 ), IINFO ) 1277 IF( IINFO.NE.0 ) THEN 1278 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A,rel)', IINFO, N, 1279 $ JTYPE, IOLDSD 1280 INFO = ABS( IINFO ) 1281 IF( IINFO.LT.0 ) THEN 1282 RETURN 1283 ELSE 1284 RESULT( 17 ) = ULPINV 1285 GO TO 280 1286 END IF 1287 END IF 1288* 1289* Do test 17 1290* 1291 TEMP2 = TWO*( TWO*N-ONE )*ULP*( ONE+EIGHT*HALF**2 ) / 1292 $ ( ONE-HALF )**4 1293* 1294 TEMP1 = ZERO 1295 DO 190 J = 1, N 1296 TEMP1 = MAX( TEMP1, ABS( D4( J )-WR( N-J+1 ) ) / 1297 $ ( ABSTOL+ABS( D4( J ) ) ) ) 1298 190 CONTINUE 1299* 1300 RESULT( 17 ) = TEMP1 / TEMP2 1301 ELSE 1302 RESULT( 17 ) = ZERO 1303 END IF 1304* 1305* Now ask for all eigenvalues with high absolute accuracy. 1306* 1307 NTEST = 18 1308 ABSTOL = UNFL + UNFL 1309 CALL SSTEBZ( 'A', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, M, 1310 $ NSPLIT, WA1, IWORK( 1 ), IWORK( N+1 ), WORK, 1311 $ IWORK( 2*N+1 ), IINFO ) 1312 IF( IINFO.NE.0 ) THEN 1313 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A)', IINFO, N, JTYPE, 1314 $ IOLDSD 1315 INFO = ABS( IINFO ) 1316 IF( IINFO.LT.0 ) THEN 1317 RETURN 1318 ELSE 1319 RESULT( 18 ) = ULPINV 1320 GO TO 280 1321 END IF 1322 END IF 1323* 1324* Do test 18 1325* 1326 TEMP1 = ZERO 1327 TEMP2 = ZERO 1328 DO 200 J = 1, N 1329 TEMP1 = MAX( TEMP1, ABS( D3( J ) ), ABS( WA1( J ) ) ) 1330 TEMP2 = MAX( TEMP2, ABS( D3( J )-WA1( J ) ) ) 1331 200 CONTINUE 1332* 1333 RESULT( 18 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 1334* 1335* Choose random values for IL and IU, and ask for the 1336* IL-th through IU-th eigenvalues. 1337* 1338 NTEST = 19 1339 IF( N.LE.1 ) THEN 1340 IL = 1 1341 IU = N 1342 ELSE 1343 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1344 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1345 IF( IU.LT.IL ) THEN 1346 ITEMP = IU 1347 IU = IL 1348 IL = ITEMP 1349 END IF 1350 END IF 1351* 1352 CALL SSTEBZ( 'I', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 1353 $ M2, NSPLIT, WA2, IWORK( 1 ), IWORK( N+1 ), 1354 $ WORK, IWORK( 2*N+1 ), IINFO ) 1355 IF( IINFO.NE.0 ) THEN 1356 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(I)', IINFO, N, JTYPE, 1357 $ IOLDSD 1358 INFO = ABS( IINFO ) 1359 IF( IINFO.LT.0 ) THEN 1360 RETURN 1361 ELSE 1362 RESULT( 19 ) = ULPINV 1363 GO TO 280 1364 END IF 1365 END IF 1366* 1367* Determine the values VL and VU of the IL-th and IU-th 1368* eigenvalues and ask for all eigenvalues in this range. 1369* 1370 IF( N.GT.0 ) THEN 1371 IF( IL.NE.1 ) THEN 1372 VL = WA1( IL ) - MAX( HALF*( WA1( IL )-WA1( IL-1 ) ), 1373 $ ULP*ANORM, TWO*RTUNFL ) 1374 ELSE 1375 VL = WA1( 1 ) - MAX( HALF*( WA1( N )-WA1( 1 ) ), 1376 $ ULP*ANORM, TWO*RTUNFL ) 1377 END IF 1378 IF( IU.NE.N ) THEN 1379 VU = WA1( IU ) + MAX( HALF*( WA1( IU+1 )-WA1( IU ) ), 1380 $ ULP*ANORM, TWO*RTUNFL ) 1381 ELSE 1382 VU = WA1( N ) + MAX( HALF*( WA1( N )-WA1( 1 ) ), 1383 $ ULP*ANORM, TWO*RTUNFL ) 1384 END IF 1385 ELSE 1386 VL = ZERO 1387 VU = ONE 1388 END IF 1389* 1390 CALL SSTEBZ( 'V', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 1391 $ M3, NSPLIT, WA3, IWORK( 1 ), IWORK( N+1 ), 1392 $ WORK, IWORK( 2*N+1 ), IINFO ) 1393 IF( IINFO.NE.0 ) THEN 1394 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(V)', IINFO, N, JTYPE, 1395 $ IOLDSD 1396 INFO = ABS( IINFO ) 1397 IF( IINFO.LT.0 ) THEN 1398 RETURN 1399 ELSE 1400 RESULT( 19 ) = ULPINV 1401 GO TO 280 1402 END IF 1403 END IF 1404* 1405 IF( M3.EQ.0 .AND. N.NE.0 ) THEN 1406 RESULT( 19 ) = ULPINV 1407 GO TO 280 1408 END IF 1409* 1410* Do test 19 1411* 1412 TEMP1 = SSXT1( 1, WA2, M2, WA3, M3, ABSTOL, ULP, UNFL ) 1413 TEMP2 = SSXT1( 1, WA3, M3, WA2, M2, ABSTOL, ULP, UNFL ) 1414 IF( N.GT.0 ) THEN 1415 TEMP3 = MAX( ABS( WA1( N ) ), ABS( WA1( 1 ) ) ) 1416 ELSE 1417 TEMP3 = ZERO 1418 END IF 1419* 1420 RESULT( 19 ) = ( TEMP1+TEMP2 ) / MAX( UNFL, TEMP3*ULP ) 1421* 1422* Call SSTEIN to compute eigenvectors corresponding to 1423* eigenvalues in WA1. (First call SSTEBZ again, to make sure 1424* it returns these eigenvalues in the correct order.) 1425* 1426 NTEST = 21 1427 CALL SSTEBZ( 'A', 'B', N, VL, VU, IL, IU, ABSTOL, SD, SE, M, 1428 $ NSPLIT, WA1, IWORK( 1 ), IWORK( N+1 ), WORK, 1429 $ IWORK( 2*N+1 ), IINFO ) 1430 IF( IINFO.NE.0 ) THEN 1431 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A,B)', IINFO, N, 1432 $ JTYPE, IOLDSD 1433 INFO = ABS( IINFO ) 1434 IF( IINFO.LT.0 ) THEN 1435 RETURN 1436 ELSE 1437 RESULT( 20 ) = ULPINV 1438 RESULT( 21 ) = ULPINV 1439 GO TO 280 1440 END IF 1441 END IF 1442* 1443 CALL SSTEIN( N, SD, SE, M, WA1, IWORK( 1 ), IWORK( N+1 ), Z, 1444 $ LDU, WORK, IWORK( 2*N+1 ), IWORK( 3*N+1 ), 1445 $ IINFO ) 1446 IF( IINFO.NE.0 ) THEN 1447 WRITE( NOUNIT, FMT = 9999 )'SSTEIN', IINFO, N, JTYPE, 1448 $ IOLDSD 1449 INFO = ABS( IINFO ) 1450 IF( IINFO.LT.0 ) THEN 1451 RETURN 1452 ELSE 1453 RESULT( 20 ) = ULPINV 1454 RESULT( 21 ) = ULPINV 1455 GO TO 280 1456 END IF 1457 END IF 1458* 1459* Do tests 20 and 21 1460* 1461 CALL SSTT21( N, 0, SD, SE, WA1, DUMMA, Z, LDU, WORK, 1462 $ RESULT( 20 ) ) 1463* 1464* Call SSTEDC(I) to compute D1 and Z, do tests. 1465* 1466* Compute D1 and Z 1467* 1468 CALL SCOPY( N, SD, 1, D1, 1 ) 1469 IF( N.GT.0 ) 1470 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1471 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1472* 1473 NTEST = 22 1474 CALL SSTEDC( 'I', N, D1, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 1475 $ IWORK, LIWEDC, IINFO ) 1476 IF( IINFO.NE.0 ) THEN 1477 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(I)', IINFO, N, JTYPE, 1478 $ IOLDSD 1479 INFO = ABS( IINFO ) 1480 IF( IINFO.LT.0 ) THEN 1481 RETURN 1482 ELSE 1483 RESULT( 22 ) = ULPINV 1484 GO TO 280 1485 END IF 1486 END IF 1487* 1488* Do Tests 22 and 23 1489* 1490 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1491 $ RESULT( 22 ) ) 1492* 1493* Call SSTEDC(V) to compute D1 and Z, do tests. 1494* 1495* Compute D1 and Z 1496* 1497 CALL SCOPY( N, SD, 1, D1, 1 ) 1498 IF( N.GT.0 ) 1499 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1500 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1501* 1502 NTEST = 24 1503 CALL SSTEDC( 'V', N, D1, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 1504 $ IWORK, LIWEDC, IINFO ) 1505 IF( IINFO.NE.0 ) THEN 1506 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(V)', IINFO, N, JTYPE, 1507 $ IOLDSD 1508 INFO = ABS( IINFO ) 1509 IF( IINFO.LT.0 ) THEN 1510 RETURN 1511 ELSE 1512 RESULT( 24 ) = ULPINV 1513 GO TO 280 1514 END IF 1515 END IF 1516* 1517* Do Tests 24 and 25 1518* 1519 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1520 $ RESULT( 24 ) ) 1521* 1522* Call SSTEDC(N) to compute D2, do tests. 1523* 1524* Compute D2 1525* 1526 CALL SCOPY( N, SD, 1, D2, 1 ) 1527 IF( N.GT.0 ) 1528 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1529 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1530* 1531 NTEST = 26 1532 CALL SSTEDC( 'N', N, D2, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 1533 $ IWORK, LIWEDC, IINFO ) 1534 IF( IINFO.NE.0 ) THEN 1535 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(N)', IINFO, N, JTYPE, 1536 $ IOLDSD 1537 INFO = ABS( IINFO ) 1538 IF( IINFO.LT.0 ) THEN 1539 RETURN 1540 ELSE 1541 RESULT( 26 ) = ULPINV 1542 GO TO 280 1543 END IF 1544 END IF 1545* 1546* Do Test 26 1547* 1548 TEMP1 = ZERO 1549 TEMP2 = ZERO 1550* 1551 DO 210 J = 1, N 1552 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 1553 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1554 210 CONTINUE 1555* 1556 RESULT( 26 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 1557* 1558* Only test SSTEMR if IEEE compliant 1559* 1560 IF( ILAENV( 10, 'SSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 .AND. 1561 $ ILAENV( 11, 'SSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 ) THEN 1562* 1563* Call SSTEMR, do test 27 (relative eigenvalue accuracy) 1564* 1565* If S is positive definite and diagonally dominant, 1566* ask for all eigenvalues with high relative accuracy. 1567* 1568 VL = ZERO 1569 VU = ZERO 1570 IL = 0 1571 IU = 0 1572 IF( JTYPE.EQ.21 .AND. SREL ) THEN 1573 NTEST = 27 1574 ABSTOL = UNFL + UNFL 1575 CALL SSTEMR( 'V', 'A', N, SD, SE, VL, VU, IL, IU, 1576 $ M, WR, Z, LDU, N, IWORK( 1 ), TRYRAC, 1577 $ WORK, LWORK, IWORK( 2*N+1 ), LWORK-2*N, 1578 $ IINFO ) 1579 IF( IINFO.NE.0 ) THEN 1580 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,A,rel)', 1581 $ IINFO, N, JTYPE, IOLDSD 1582 INFO = ABS( IINFO ) 1583 IF( IINFO.LT.0 ) THEN 1584 RETURN 1585 ELSE 1586 RESULT( 27 ) = ULPINV 1587 GO TO 270 1588 END IF 1589 END IF 1590* 1591* Do test 27 1592* 1593 TEMP2 = TWO*( TWO*N-ONE )*ULP*( ONE+EIGHT*HALF**2 ) / 1594 $ ( ONE-HALF )**4 1595* 1596 TEMP1 = ZERO 1597 DO 220 J = 1, N 1598 TEMP1 = MAX( TEMP1, ABS( D4( J )-WR( N-J+1 ) ) / 1599 $ ( ABSTOL+ABS( D4( J ) ) ) ) 1600 220 CONTINUE 1601* 1602 RESULT( 27 ) = TEMP1 / TEMP2 1603* 1604 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1605 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1606 IF( IU.LT.IL ) THEN 1607 ITEMP = IU 1608 IU = IL 1609 IL = ITEMP 1610 END IF 1611* 1612 IF( SRANGE ) THEN 1613 NTEST = 28 1614 ABSTOL = UNFL + UNFL 1615 CALL SSTEMR( 'V', 'I', N, SD, SE, VL, VU, IL, IU, 1616 $ M, WR, Z, LDU, N, IWORK( 1 ), TRYRAC, 1617 $ WORK, LWORK, IWORK( 2*N+1 ), 1618 $ LWORK-2*N, IINFO ) 1619* 1620 IF( IINFO.NE.0 ) THEN 1621 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,I,rel)', 1622 $ IINFO, N, JTYPE, IOLDSD 1623 INFO = ABS( IINFO ) 1624 IF( IINFO.LT.0 ) THEN 1625 RETURN 1626 ELSE 1627 RESULT( 28 ) = ULPINV 1628 GO TO 270 1629 END IF 1630 END IF 1631* 1632* 1633* Do test 28 1634* 1635 TEMP2 = TWO*( TWO*N-ONE )*ULP* 1636 $ ( ONE+EIGHT*HALF**2 ) / ( ONE-HALF )**4 1637* 1638 TEMP1 = ZERO 1639 DO 230 J = IL, IU 1640 TEMP1 = MAX( TEMP1, ABS( WR( J-IL+1 )-D4( N-J+ 1641 $ 1 ) ) / ( ABSTOL+ABS( WR( J-IL+1 ) ) ) ) 1642 230 CONTINUE 1643* 1644 RESULT( 28 ) = TEMP1 / TEMP2 1645 ELSE 1646 RESULT( 28 ) = ZERO 1647 END IF 1648 ELSE 1649 RESULT( 27 ) = ZERO 1650 RESULT( 28 ) = ZERO 1651 END IF 1652* 1653* Call SSTEMR(V,I) to compute D1 and Z, do tests. 1654* 1655* Compute D1 and Z 1656* 1657 CALL SCOPY( N, SD, 1, D5, 1 ) 1658 IF( N.GT.0 ) 1659 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1660 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1661* 1662 IF( SRANGE ) THEN 1663 NTEST = 29 1664 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1665 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 1666 IF( IU.LT.IL ) THEN 1667 ITEMP = IU 1668 IU = IL 1669 IL = ITEMP 1670 END IF 1671 CALL SSTEMR( 'V', 'I', N, D5, WORK, VL, VU, IL, IU, 1672 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 1673 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1674 $ LIWORK-2*N, IINFO ) 1675 IF( IINFO.NE.0 ) THEN 1676 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,I)', IINFO, 1677 $ N, JTYPE, IOLDSD 1678 INFO = ABS( IINFO ) 1679 IF( IINFO.LT.0 ) THEN 1680 RETURN 1681 ELSE 1682 RESULT( 29 ) = ULPINV 1683 GO TO 280 1684 END IF 1685 END IF 1686* 1687* Do Tests 29 and 30 1688* 1689 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1690 $ M, RESULT( 29 ) ) 1691* 1692* Call SSTEMR to compute D2, do tests. 1693* 1694* Compute D2 1695* 1696 CALL SCOPY( N, SD, 1, D5, 1 ) 1697 IF( N.GT.0 ) 1698 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1699* 1700 NTEST = 31 1701 CALL SSTEMR( 'N', 'I', N, D5, WORK, VL, VU, IL, IU, 1702 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 1703 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1704 $ LIWORK-2*N, IINFO ) 1705 IF( IINFO.NE.0 ) THEN 1706 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,I)', IINFO, 1707 $ N, JTYPE, IOLDSD 1708 INFO = ABS( IINFO ) 1709 IF( IINFO.LT.0 ) THEN 1710 RETURN 1711 ELSE 1712 RESULT( 31 ) = ULPINV 1713 GO TO 280 1714 END IF 1715 END IF 1716* 1717* Do Test 31 1718* 1719 TEMP1 = ZERO 1720 TEMP2 = ZERO 1721* 1722 DO 240 J = 1, IU - IL + 1 1723 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), 1724 $ ABS( D2( J ) ) ) 1725 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1726 240 CONTINUE 1727* 1728 RESULT( 31 ) = TEMP2 / MAX( UNFL, 1729 $ ULP*MAX( TEMP1, TEMP2 ) ) 1730* 1731* 1732* Call SSTEMR(V,V) to compute D1 and Z, do tests. 1733* 1734* Compute D1 and Z 1735* 1736 CALL SCOPY( N, SD, 1, D5, 1 ) 1737 IF( N.GT.0 ) 1738 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1739 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 1740* 1741 NTEST = 32 1742* 1743 IF( N.GT.0 ) THEN 1744 IF( IL.NE.1 ) THEN 1745 VL = D2( IL ) - MAX( HALF* 1746 $ ( D2( IL )-D2( IL-1 ) ), ULP*ANORM, 1747 $ TWO*RTUNFL ) 1748 ELSE 1749 VL = D2( 1 ) - MAX( HALF*( D2( N )-D2( 1 ) ), 1750 $ ULP*ANORM, TWO*RTUNFL ) 1751 END IF 1752 IF( IU.NE.N ) THEN 1753 VU = D2( IU ) + MAX( HALF* 1754 $ ( D2( IU+1 )-D2( IU ) ), ULP*ANORM, 1755 $ TWO*RTUNFL ) 1756 ELSE 1757 VU = D2( N ) + MAX( HALF*( D2( N )-D2( 1 ) ), 1758 $ ULP*ANORM, TWO*RTUNFL ) 1759 END IF 1760 ELSE 1761 VL = ZERO 1762 VU = ONE 1763 END IF 1764* 1765 CALL SSTEMR( 'V', 'V', N, D5, WORK, VL, VU, IL, IU, 1766 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 1767 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1768 $ LIWORK-2*N, IINFO ) 1769 IF( IINFO.NE.0 ) THEN 1770 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,V)', IINFO, 1771 $ N, JTYPE, IOLDSD 1772 INFO = ABS( IINFO ) 1773 IF( IINFO.LT.0 ) THEN 1774 RETURN 1775 ELSE 1776 RESULT( 32 ) = ULPINV 1777 GO TO 280 1778 END IF 1779 END IF 1780* 1781* Do Tests 32 and 33 1782* 1783 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 1784 $ M, RESULT( 32 ) ) 1785* 1786* Call SSTEMR to compute D2, do tests. 1787* 1788* Compute D2 1789* 1790 CALL SCOPY( N, SD, 1, D5, 1 ) 1791 IF( N.GT.0 ) 1792 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1793* 1794 NTEST = 34 1795 CALL SSTEMR( 'N', 'V', N, D5, WORK, VL, VU, IL, IU, 1796 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 1797 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1798 $ LIWORK-2*N, IINFO ) 1799 IF( IINFO.NE.0 ) THEN 1800 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,V)', IINFO, 1801 $ N, JTYPE, IOLDSD 1802 INFO = ABS( IINFO ) 1803 IF( IINFO.LT.0 ) THEN 1804 RETURN 1805 ELSE 1806 RESULT( 34 ) = ULPINV 1807 GO TO 280 1808 END IF 1809 END IF 1810* 1811* Do Test 34 1812* 1813 TEMP1 = ZERO 1814 TEMP2 = ZERO 1815* 1816 DO 250 J = 1, IU - IL + 1 1817 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), 1818 $ ABS( D2( J ) ) ) 1819 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1820 250 CONTINUE 1821* 1822 RESULT( 34 ) = TEMP2 / MAX( UNFL, 1823 $ ULP*MAX( TEMP1, TEMP2 ) ) 1824 ELSE 1825 RESULT( 29 ) = ZERO 1826 RESULT( 30 ) = ZERO 1827 RESULT( 31 ) = ZERO 1828 RESULT( 32 ) = ZERO 1829 RESULT( 33 ) = ZERO 1830 RESULT( 34 ) = ZERO 1831 END IF 1832* 1833* 1834* Call SSTEMR(V,A) to compute D1 and Z, do tests. 1835* 1836* Compute D1 and Z 1837* 1838 CALL SCOPY( N, SD, 1, D5, 1 ) 1839 IF( N.GT.0 ) 1840 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1841* 1842 NTEST = 35 1843* 1844 CALL SSTEMR( 'V', 'A', N, D5, WORK, VL, VU, IL, IU, 1845 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 1846 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1847 $ LIWORK-2*N, IINFO ) 1848 IF( IINFO.NE.0 ) THEN 1849 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,A)', IINFO, N, 1850 $ JTYPE, IOLDSD 1851 INFO = ABS( IINFO ) 1852 IF( IINFO.LT.0 ) THEN 1853 RETURN 1854 ELSE 1855 RESULT( 35 ) = ULPINV 1856 GO TO 280 1857 END IF 1858 END IF 1859* 1860* Do Tests 35 and 36 1861* 1862 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, M, 1863 $ RESULT( 35 ) ) 1864* 1865* Call SSTEMR to compute D2, do tests. 1866* 1867* Compute D2 1868* 1869 CALL SCOPY( N, SD, 1, D5, 1 ) 1870 IF( N.GT.0 ) 1871 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 1872* 1873 NTEST = 37 1874 CALL SSTEMR( 'N', 'A', N, D5, WORK, VL, VU, IL, IU, 1875 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 1876 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 1877 $ LIWORK-2*N, IINFO ) 1878 IF( IINFO.NE.0 ) THEN 1879 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,A)', IINFO, N, 1880 $ JTYPE, IOLDSD 1881 INFO = ABS( IINFO ) 1882 IF( IINFO.LT.0 ) THEN 1883 RETURN 1884 ELSE 1885 RESULT( 37 ) = ULPINV 1886 GO TO 280 1887 END IF 1888 END IF 1889* 1890* Do Test 34 1891* 1892 TEMP1 = ZERO 1893 TEMP2 = ZERO 1894* 1895 DO 260 J = 1, N 1896 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 1897 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 1898 260 CONTINUE 1899* 1900 RESULT( 37 ) = TEMP2 / MAX( UNFL, 1901 $ ULP*MAX( TEMP1, TEMP2 ) ) 1902 END IF 1903 270 CONTINUE 1904 280 CONTINUE 1905 NTESTT = NTESTT + NTEST 1906* 1907* End of Loop -- Check for RESULT(j) > THRESH 1908* 1909* 1910* Print out tests which fail. 1911* 1912 DO 290 JR = 1, NTEST 1913 IF( RESULT( JR ).GE.THRESH ) THEN 1914* 1915* If this is the first test to fail, 1916* print a header to the data file. 1917* 1918 IF( NERRS.EQ.0 ) THEN 1919 WRITE( NOUNIT, FMT = 9998 )'SST' 1920 WRITE( NOUNIT, FMT = 9997 ) 1921 WRITE( NOUNIT, FMT = 9996 ) 1922 WRITE( NOUNIT, FMT = 9995 )'Symmetric' 1923 WRITE( NOUNIT, FMT = 9994 ) 1924* 1925* Tests performed 1926* 1927 WRITE( NOUNIT, FMT = 9988 ) 1928 END IF 1929 NERRS = NERRS + 1 1930 WRITE( NOUNIT, FMT = 9990 )N, IOLDSD, JTYPE, JR, 1931 $ RESULT( JR ) 1932 END IF 1933 290 CONTINUE 1934 300 CONTINUE 1935 310 CONTINUE 1936* 1937* Summary 1938* 1939 CALL SLASUM( 'SST', NOUNIT, NERRS, NTESTT ) 1940 RETURN 1941* 1942 9999 FORMAT( ' SCHKST: ', A, ' returned INFO=', I6, '.', / 9X, 'N=', 1943 $ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' ) 1944* 1945 9998 FORMAT( / 1X, A3, ' -- Real Symmetric eigenvalue problem' ) 1946 9997 FORMAT( ' Matrix types (see SCHKST for details): ' ) 1947* 1948 9996 FORMAT( / ' Special Matrices:', 1949 $ / ' 1=Zero matrix. ', 1950 $ ' 5=Diagonal: clustered entries.', 1951 $ / ' 2=Identity matrix. ', 1952 $ ' 6=Diagonal: large, evenly spaced.', 1953 $ / ' 3=Diagonal: evenly spaced entries. ', 1954 $ ' 7=Diagonal: small, evenly spaced.', 1955 $ / ' 4=Diagonal: geometr. spaced entries.' ) 1956 9995 FORMAT( ' Dense ', A, ' Matrices:', 1957 $ / ' 8=Evenly spaced eigenvals. ', 1958 $ ' 12=Small, evenly spaced eigenvals.', 1959 $ / ' 9=Geometrically spaced eigenvals. ', 1960 $ ' 13=Matrix with random O(1) entries.', 1961 $ / ' 10=Clustered eigenvalues. ', 1962 $ ' 14=Matrix with large random entries.', 1963 $ / ' 11=Large, evenly spaced eigenvals. ', 1964 $ ' 15=Matrix with small random entries.' ) 1965 9994 FORMAT( ' 16=Positive definite, evenly spaced eigenvalues', 1966 $ / ' 17=Positive definite, geometrically spaced eigenvlaues', 1967 $ / ' 18=Positive definite, clustered eigenvalues', 1968 $ / ' 19=Positive definite, small evenly spaced eigenvalues', 1969 $ / ' 20=Positive definite, large evenly spaced eigenvalues', 1970 $ / ' 21=Diagonally dominant tridiagonal, geometrically', 1971 $ ' spaced eigenvalues' ) 1972* 1973 9990 FORMAT( ' N=', I5, ', seed=', 4( I4, ',' ), ' type ', I2, 1974 $ ', test(', I2, ')=', G10.3 ) 1975* 1976 9988 FORMAT( / 'Test performed: see SCHKST for details.', / ) 1977* End of SCHKST 1978* 1979 END 1980