 |
LAPACK 3.3.1
Linear Algebra PACKage
|
|
LAPACK 3.3.1
Linear Algebra PACKage
|
00001 SUBROUTINE SCHKST( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, 00002 $ NOUNIT, A, LDA, AP, SD, SE, D1, D2, D3, D4, D5, 00003 $ WA1, WA2, WA3, WR, U, LDU, V, VP, TAU, Z, WORK, 00004 $ LWORK, IWORK, LIWORK, RESULT, INFO ) 00005 IMPLICIT NONE 00006 * 00007 * -- LAPACK test routine (version 3.1) -- 00008 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd.. 00009 * November 2006 00010 * 00011 * .. Scalar Arguments .. 00012 INTEGER INFO, LDA, LDU, LIWORK, LWORK, NOUNIT, NSIZES, 00013 $ NTYPES 00014 REAL THRESH 00015 * .. 00016 * .. Array Arguments .. 00017 LOGICAL DOTYPE( * ) 00018 INTEGER ISEED( 4 ), IWORK( * ), NN( * ) 00019 REAL A( LDA, * ), AP( * ), D1( * ), D2( * ), 00020 $ D3( * ), D4( * ), D5( * ), RESULT( * ), 00021 $ SD( * ), SE( * ), TAU( * ), U( LDU, * ), 00022 $ V( LDU, * ), VP( * ), WA1( * ), WA2( * ), 00023 $ WA3( * ), WORK( * ), WR( * ), Z( LDU, * ) 00024 * .. 00025 * 00026 * Purpose 00027 * ======= 00028 * 00029 * SCHKST checks the symmetric eigenvalue problem routines. 00030 * 00031 * SSYTRD factors A as U S U' , where ' means transpose, 00032 * S is symmetric tridiagonal, and U is orthogonal. 00033 * SSYTRD can use either just the lower or just the upper triangle 00034 * of A; SCHKST checks both cases. 00035 * U is represented as a product of Householder 00036 * transformations, whose vectors are stored in the first 00037 * n-1 columns of V, and whose scale factors are in TAU. 00038 * 00039 * SSPTRD does the same as SSYTRD, except that A and V are stored 00040 * in "packed" format. 00041 * 00042 * SORGTR constructs the matrix U from the contents of V and TAU. 00043 * 00044 * SOPGTR constructs the matrix U from the contents of VP and TAU. 00045 * 00046 * SSTEQR factors S as Z D1 Z' , where Z is the orthogonal 00047 * matrix of eigenvectors and D1 is a diagonal matrix with 00048 * the eigenvalues on the diagonal. D2 is the matrix of 00049 * eigenvalues computed when Z is not computed. 00050 * 00051 * SSTERF computes D3, the matrix of eigenvalues, by the 00052 * PWK method, which does not yield eigenvectors. 00053 * 00054 * SPTEQR factors S as Z4 D4 Z4' , for a 00055 * symmetric positive definite tridiagonal matrix. 00056 * D5 is the matrix of eigenvalues computed when Z is not 00057 * computed. 00058 * 00059 * SSTEBZ computes selected eigenvalues. WA1, WA2, and 00060 * WA3 will denote eigenvalues computed to high 00061 * absolute accuracy, with different range options. 00062 * WR will denote eigenvalues computed to high relative 00063 * accuracy. 00064 * 00065 * SSTEIN computes Y, the eigenvectors of S, given the 00066 * eigenvalues. 00067 * 00068 * SSTEDC factors S as Z D1 Z' , where Z is the orthogonal 00069 * matrix of eigenvectors and D1 is a diagonal matrix with 00070 * the eigenvalues on the diagonal ('I' option). It may also 00071 * update an input orthogonal matrix, usually the output 00072 * from SSYTRD/SORGTR or SSPTRD/SOPGTR ('V' option). It may 00073 * also just compute eigenvalues ('N' option). 00074 * 00075 * SSTEMR factors S as Z D1 Z' , where Z is the orthogonal 00076 * matrix of eigenvectors and D1 is a diagonal matrix with 00077 * the eigenvalues on the diagonal ('I' option). SSTEMR 00078 * uses the Relatively Robust Representation whenever possible. 00079 * 00080 * When SCHKST is called, a number of matrix "sizes" ("n's") and a 00081 * number of matrix "types" are specified. For each size ("n") 00082 * and each type of matrix, one matrix will be generated and used 00083 * to test the symmetric eigenroutines. For each matrix, a number 00084 * of tests will be performed: 00085 * 00086 * (1) | A - V S V' | / ( |A| n ulp ) SSYTRD( UPLO='U', ... ) 00087 * 00088 * (2) | I - UV' | / ( n ulp ) SORGTR( UPLO='U', ... ) 00089 * 00090 * (3) | A - V S V' | / ( |A| n ulp ) SSYTRD( UPLO='L', ... ) 00091 * 00092 * (4) | I - UV' | / ( n ulp ) SORGTR( UPLO='L', ... ) 00093 * 00094 * (5-8) Same as 1-4, but for SSPTRD and SOPGTR. 00095 * 00096 * (9) | S - Z D Z' | / ( |S| n ulp ) SSTEQR('V',...) 00097 * 00098 * (10) | I - ZZ' | / ( n ulp ) SSTEQR('V',...) 00099 * 00100 * (11) | D1 - D2 | / ( |D1| ulp ) SSTEQR('N',...) 00101 * 00102 * (12) | D1 - D3 | / ( |D1| ulp ) SSTERF 00103 * 00104 * (13) 0 if the true eigenvalues (computed by sturm count) 00105 * of S are within THRESH of 00106 * those in D1. 2*THRESH if they are not. (Tested using 00107 * SSTECH) 00108 * 00109 * For S positive definite, 00110 * 00111 * (14) | S - Z4 D4 Z4' | / ( |S| n ulp ) SPTEQR('V',...) 00112 * 00113 * (15) | I - Z4 Z4' | / ( n ulp ) SPTEQR('V',...) 00114 * 00115 * (16) | D4 - D5 | / ( 100 |D4| ulp ) SPTEQR('N',...) 00116 * 00117 * When S is also diagonally dominant by the factor gamma < 1, 00118 * 00119 * (17) max | D4(i) - WR(i) | / ( |D4(i)| omega ) , 00120 * i 00121 * omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4 00122 * SSTEBZ( 'A', 'E', ...) 00123 * 00124 * (18) | WA1 - D3 | / ( |D3| ulp ) SSTEBZ( 'A', 'E', ...) 00125 * 00126 * (19) ( max { min | WA2(i)-WA3(j) | } + 00127 * i j 00128 * max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) 00129 * i j 00130 * SSTEBZ( 'I', 'E', ...) 00131 * 00132 * (20) | S - Y WA1 Y' | / ( |S| n ulp ) SSTEBZ, SSTEIN 00133 * 00134 * (21) | I - Y Y' | / ( n ulp ) SSTEBZ, SSTEIN 00135 * 00136 * (22) | S - Z D Z' | / ( |S| n ulp ) SSTEDC('I') 00137 * 00138 * (23) | I - ZZ' | / ( n ulp ) SSTEDC('I') 00139 * 00140 * (24) | S - Z D Z' | / ( |S| n ulp ) SSTEDC('V') 00141 * 00142 * (25) | I - ZZ' | / ( n ulp ) SSTEDC('V') 00143 * 00144 * (26) | D1 - D2 | / ( |D1| ulp ) SSTEDC('V') and 00145 * SSTEDC('N') 00146 * 00147 * Test 27 is disabled at the moment because SSTEMR does not 00148 * guarantee high relatvie accuracy. 00149 * 00150 * (27) max | D6(i) - WR(i) | / ( |D6(i)| omega ) , 00151 * i 00152 * omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4 00153 * SSTEMR('V', 'A') 00154 * 00155 * (28) max | D6(i) - WR(i) | / ( |D6(i)| omega ) , 00156 * i 00157 * omega = 2 (2n-1) ULP (1 + 8 gamma**2) / (1 - gamma)**4 00158 * SSTEMR('V', 'I') 00159 * 00160 * Tests 29 through 34 are disable at present because SSTEMR 00161 * does not handle partial specturm requests. 00162 * 00163 * (29) | S - Z D Z' | / ( |S| n ulp ) SSTEMR('V', 'I') 00164 * 00165 * (30) | I - ZZ' | / ( n ulp ) SSTEMR('V', 'I') 00166 * 00167 * (31) ( max { min | WA2(i)-WA3(j) | } + 00168 * i j 00169 * max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) 00170 * i j 00171 * SSTEMR('N', 'I') vs. SSTEMR('V', 'I') 00172 * 00173 * (32) | S - Z D Z' | / ( |S| n ulp ) SSTEMR('V', 'V') 00174 * 00175 * (33) | I - ZZ' | / ( n ulp ) SSTEMR('V', 'V') 00176 * 00177 * (34) ( max { min | WA2(i)-WA3(j) | } + 00178 * i j 00179 * max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) 00180 * i j 00181 * SSTEMR('N', 'V') vs. SSTEMR('V', 'V') 00182 * 00183 * (35) | S - Z D Z' | / ( |S| n ulp ) SSTEMR('V', 'A') 00184 * 00185 * (36) | I - ZZ' | / ( n ulp ) SSTEMR('V', 'A') 00186 * 00187 * (37) ( max { min | WA2(i)-WA3(j) | } + 00188 * i j 00189 * max { min | WA3(i)-WA2(j) | } ) / ( |D3| ulp ) 00190 * i j 00191 * SSTEMR('N', 'A') vs. SSTEMR('V', 'A') 00192 * 00193 * The "sizes" are specified by an array NN(1:NSIZES); the value of 00194 * each element NN(j) specifies one size. 00195 * The "types" are specified by a logical array DOTYPE( 1:NTYPES ); 00196 * if DOTYPE(j) is .TRUE., then matrix type "j" will be generated. 00197 * Currently, the list of possible types is: 00198 * 00199 * (1) The zero matrix. 00200 * (2) The identity matrix. 00201 * 00202 * (3) A diagonal matrix with evenly spaced entries 00203 * 1, ..., ULP and random signs. 00204 * (ULP = (first number larger than 1) - 1 ) 00205 * (4) A diagonal matrix with geometrically spaced entries 00206 * 1, ..., ULP and random signs. 00207 * (5) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP 00208 * and random signs. 00209 * 00210 * (6) Same as (4), but multiplied by SQRT( overflow threshold ) 00211 * (7) Same as (4), but multiplied by SQRT( underflow threshold ) 00212 * 00213 * (8) A matrix of the form U' D U, where U is orthogonal and 00214 * D has evenly spaced entries 1, ..., ULP with random signs 00215 * on the diagonal. 00216 * 00217 * (9) A matrix of the form U' D U, where U is orthogonal and 00218 * D has geometrically spaced entries 1, ..., ULP with random 00219 * signs on the diagonal. 00220 * 00221 * (10) A matrix of the form U' D U, where U is orthogonal and 00222 * D has "clustered" entries 1, ULP,..., ULP with random 00223 * signs on the diagonal. 00224 * 00225 * (11) Same as (8), but multiplied by SQRT( overflow threshold ) 00226 * (12) Same as (8), but multiplied by SQRT( underflow threshold ) 00227 * 00228 * (13) Symmetric matrix with random entries chosen from (-1,1). 00229 * (14) Same as (13), but multiplied by SQRT( overflow threshold ) 00230 * (15) Same as (13), but multiplied by SQRT( underflow threshold ) 00231 * (16) Same as (8), but diagonal elements are all positive. 00232 * (17) Same as (9), but diagonal elements are all positive. 00233 * (18) Same as (10), but diagonal elements are all positive. 00234 * (19) Same as (16), but multiplied by SQRT( overflow threshold ) 00235 * (20) Same as (16), but multiplied by SQRT( underflow threshold ) 00236 * (21) A diagonally dominant tridiagonal matrix with geometrically 00237 * spaced diagonal entries 1, ..., ULP. 00238 * 00239 * Arguments 00240 * ========= 00241 * 00242 * NSIZES (input) INTEGER 00243 * The number of sizes of matrices to use. If it is zero, 00244 * SCHKST does nothing. It must be at least zero. 00245 * 00246 * NN (input) INTEGER array, dimension (NSIZES) 00247 * An array containing the sizes to be used for the matrices. 00248 * Zero values will be skipped. The values must be at least 00249 * zero. 00250 * 00251 * NTYPES (input) INTEGER 00252 * The number of elements in DOTYPE. If it is zero, SCHKST 00253 * does nothing. It must be at least zero. If it is MAXTYP+1 00254 * and NSIZES is 1, then an additional type, MAXTYP+1 is 00255 * defined, which is to use whatever matrix is in A. This 00256 * is only useful if DOTYPE(1:MAXTYP) is .FALSE. and 00257 * DOTYPE(MAXTYP+1) is .TRUE. . 00258 * 00259 * DOTYPE (input) LOGICAL array, dimension (NTYPES) 00260 * If DOTYPE(j) is .TRUE., then for each size in NN a 00261 * matrix of that size and of type j will be generated. 00262 * If NTYPES is smaller than the maximum number of types 00263 * defined (PARAMETER MAXTYP), then types NTYPES+1 through 00264 * MAXTYP will not be generated. If NTYPES is larger 00265 * than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES) 00266 * will be ignored. 00267 * 00268 * ISEED (input/output) INTEGER array, dimension (4) 00269 * On entry ISEED specifies the seed of the random number 00270 * generator. The array elements should be between 0 and 4095; 00271 * if not they will be reduced mod 4096. Also, ISEED(4) must 00272 * be odd. The random number generator uses a linear 00273 * congruential sequence limited to small integers, and so 00274 * should produce machine independent random numbers. The 00275 * values of ISEED are changed on exit, and can be used in the 00276 * next call to SCHKST to continue the same random number 00277 * sequence. 00278 * 00279 * THRESH (input) REAL 00280 * A test will count as "failed" if the "error", computed as 00281 * described above, exceeds THRESH. Note that the error 00282 * is scaled to be O(1), so THRESH should be a reasonably 00283 * small multiple of 1, e.g., 10 or 100. In particular, 00284 * it should not depend on the precision (single vs. double) 00285 * or the size of the matrix. It must be at least zero. 00286 * 00287 * NOUNIT (input) INTEGER 00288 * The FORTRAN unit number for printing out error messages 00289 * (e.g., if a routine returns IINFO not equal to 0.) 00290 * 00291 * A (input/workspace/output) REAL array of 00292 * dimension ( LDA , max(NN) ) 00293 * Used to hold the matrix whose eigenvalues are to be 00294 * computed. On exit, A contains the last matrix actually 00295 * used. 00296 * 00297 * LDA (input) INTEGER 00298 * The leading dimension of A. It must be at 00299 * least 1 and at least max( NN ). 00300 * 00301 * AP (workspace) REAL array of 00302 * dimension( max(NN)*max(NN+1)/2 ) 00303 * The matrix A stored in packed format. 00304 * 00305 * SD (workspace/output) REAL array of 00306 * dimension( max(NN) ) 00307 * The diagonal of the tridiagonal matrix computed by SSYTRD. 00308 * On exit, SD and SE contain the tridiagonal form of the 00309 * matrix in A. 00310 * 00311 * SE (workspace/output) REAL array of 00312 * dimension( max(NN) ) 00313 * The off-diagonal of the tridiagonal matrix computed by 00314 * SSYTRD. On exit, SD and SE contain the tridiagonal form of 00315 * the matrix in A. 00316 * 00317 * D1 (workspace/output) REAL array of 00318 * dimension( max(NN) ) 00319 * The eigenvalues of A, as computed by SSTEQR simlutaneously 00320 * with Z. On exit, the eigenvalues in D1 correspond with the 00321 * matrix in A. 00322 * 00323 * D2 (workspace/output) REAL array of 00324 * dimension( max(NN) ) 00325 * The eigenvalues of A, as computed by SSTEQR if Z is not 00326 * computed. On exit, the eigenvalues in D2 correspond with 00327 * the matrix in A. 00328 * 00329 * D3 (workspace/output) REAL array of 00330 * dimension( max(NN) ) 00331 * The eigenvalues of A, as computed by SSTERF. On exit, the 00332 * eigenvalues in D3 correspond with the matrix in A. 00333 * 00334 * U (workspace/output) REAL array of 00335 * dimension( LDU, max(NN) ). 00336 * The orthogonal matrix computed by SSYTRD + SORGTR. 00337 * 00338 * LDU (input) INTEGER 00339 * The leading dimension of U, Z, and V. It must be at least 1 00340 * and at least max( NN ). 00341 * 00342 * V (workspace/output) REAL array of 00343 * dimension( LDU, max(NN) ). 00344 * The Housholder vectors computed by SSYTRD in reducing A to 00345 * tridiagonal form. The vectors computed with UPLO='U' are 00346 * in the upper triangle, and the vectors computed with UPLO='L' 00347 * are in the lower triangle. (As described in SSYTRD, the 00348 * sub- and superdiagonal are not set to 1, although the 00349 * true Householder vector has a 1 in that position. The 00350 * routines that use V, such as SORGTR, set those entries to 00351 * 1 before using them, and then restore them later.) 00352 * 00353 * VP (workspace) REAL array of 00354 * dimension( max(NN)*max(NN+1)/2 ) 00355 * The matrix V stored in packed format. 00356 * 00357 * TAU (workspace/output) REAL array of 00358 * dimension( max(NN) ) 00359 * The Householder factors computed by SSYTRD in reducing A 00360 * to tridiagonal form. 00361 * 00362 * Z (workspace/output) REAL array of 00363 * dimension( LDU, max(NN) ). 00364 * The orthogonal matrix of eigenvectors computed by SSTEQR, 00365 * SPTEQR, and SSTEIN. 00366 * 00367 * WORK (workspace/output) REAL array of 00368 * dimension( LWORK ) 00369 * 00370 * LWORK (input) INTEGER 00371 * The number of entries in WORK. This must be at least 00372 * 1 + 4 * Nmax + 2 * Nmax * lg Nmax + 3 * Nmax**2 00373 * where Nmax = max( NN(j), 2 ) and lg = log base 2. 00374 * 00375 * IWORK (workspace/output) INTEGER array, 00376 * dimension (6 + 6*Nmax + 5 * Nmax * lg Nmax ) 00377 * where Nmax = max( NN(j), 2 ) and lg = log base 2. 00378 * Workspace. 00379 * 00380 * RESULT (output) REAL array, dimension (26) 00381 * The values computed by the tests described above. 00382 * The values are currently limited to 1/ulp, to avoid 00383 * overflow. 00384 * 00385 * INFO (output) INTEGER 00386 * If 0, then everything ran OK. 00387 * -1: NSIZES < 0 00388 * -2: Some NN(j) < 0 00389 * -3: NTYPES < 0 00390 * -5: THRESH < 0 00391 * -9: LDA < 1 or LDA < NMAX, where NMAX is max( NN(j) ). 00392 * -23: LDU < 1 or LDU < NMAX. 00393 * -29: LWORK too small. 00394 * If SLATMR, SLATMS, SSYTRD, SORGTR, SSTEQR, SSTERF, 00395 * or SORMC2 returns an error code, the 00396 * absolute value of it is returned. 00397 * 00398 *----------------------------------------------------------------------- 00399 * 00400 * Some Local Variables and Parameters: 00401 * ---- ----- --------- --- ---------- 00402 * ZERO, ONE Real 0 and 1. 00403 * MAXTYP The number of types defined. 00404 * NTEST The number of tests performed, or which can 00405 * be performed so far, for the current matrix. 00406 * NTESTT The total number of tests performed so far. 00407 * NBLOCK Blocksize as returned by ENVIR. 00408 * NMAX Largest value in NN. 00409 * NMATS The number of matrices generated so far. 00410 * NERRS The number of tests which have exceeded THRESH 00411 * so far. 00412 * COND, IMODE Values to be passed to the matrix generators. 00413 * ANORM Norm of A; passed to matrix generators. 00414 * 00415 * OVFL, UNFL Overflow and underflow thresholds. 00416 * ULP, ULPINV Finest relative precision and its inverse. 00417 * RTOVFL, RTUNFL Square roots of the previous 2 values. 00418 * The following four arrays decode JTYPE: 00419 * KTYPE(j) The general type (1-10) for type "j". 00420 * KMODE(j) The MODE value to be passed to the matrix 00421 * generator for type "j". 00422 * KMAGN(j) The order of magnitude ( O(1), 00423 * O(overflow^(1/2) ), O(underflow^(1/2) ) 00424 * 00425 * ===================================================================== 00426 * 00427 * .. Parameters .. 00428 REAL ZERO, ONE, TWO, EIGHT, TEN, HUN 00429 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0, 00430 $ EIGHT = 8.0E0, TEN = 10.0E0, HUN = 100.0E0 ) 00431 REAL HALF 00432 PARAMETER ( HALF = ONE / TWO ) 00433 INTEGER MAXTYP 00434 PARAMETER ( MAXTYP = 21 ) 00435 LOGICAL SRANGE 00436 PARAMETER ( SRANGE = .FALSE. ) 00437 LOGICAL SREL 00438 PARAMETER ( SREL = .FALSE. ) 00439 * .. 00440 * .. Local Scalars .. 00441 LOGICAL BADNN, TRYRAC 00442 INTEGER I, IINFO, IL, IMODE, ITEMP, ITYPE, IU, J, JC, 00443 $ JR, JSIZE, JTYPE, LGN, LIWEDC, LOG2UI, LWEDC, 00444 $ M, M2, M3, MTYPES, N, NAP, NBLOCK, NERRS, 00445 $ NMATS, NMAX, NSPLIT, NTEST, NTESTT 00446 REAL ABSTOL, ANINV, ANORM, COND, OVFL, RTOVFL, 00447 $ RTUNFL, TEMP1, TEMP2, TEMP3, TEMP4, ULP, 00448 $ ULPINV, UNFL, VL, VU 00449 * .. 00450 * .. Local Arrays .. 00451 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ), 00452 $ KMAGN( MAXTYP ), KMODE( MAXTYP ), 00453 $ KTYPE( MAXTYP ) 00454 REAL DUMMA( 1 ) 00455 * .. 00456 * .. External Functions .. 00457 INTEGER ILAENV 00458 REAL SLAMCH, SLARND, SSXT1 00459 EXTERNAL ILAENV, SLAMCH, SLARND, SSXT1 00460 * .. 00461 * .. External Subroutines .. 00462 EXTERNAL SCOPY, SLABAD, SLACPY, SLASET, SLASUM, SLATMR, 00463 $ SLATMS, SOPGTR, SORGTR, SPTEQR, SSPT21, SSPTRD, 00464 $ SSTEBZ, SSTECH, SSTEDC, SSTEMR, SSTEIN, SSTEQR, 00465 $ SSTERF, SSTT21, SSTT22, SSYT21, SSYTRD, XERBLA 00466 * .. 00467 * .. Intrinsic Functions .. 00468 INTRINSIC ABS, INT, LOG, MAX, MIN, REAL, SQRT 00469 * .. 00470 * .. Data statements .. 00471 DATA KTYPE / 1, 2, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 8, 00472 $ 8, 8, 9, 9, 9, 9, 9, 10 / 00473 DATA KMAGN / 1, 1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1, 00474 $ 2, 3, 1, 1, 1, 2, 3, 1 / 00475 DATA KMODE / 0, 0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0, 00476 $ 0, 0, 4, 3, 1, 4, 4, 3 / 00477 * .. 00478 * .. Executable Statements .. 00479 * 00480 * Keep ftnchek happy 00481 IDUMMA( 1 ) = 1 00482 * 00483 * Check for errors 00484 * 00485 NTESTT = 0 00486 INFO = 0 00487 * 00488 * Important constants 00489 * 00490 BADNN = .FALSE. 00491 TRYRAC = .TRUE. 00492 NMAX = 1 00493 DO 10 J = 1, NSIZES 00494 NMAX = MAX( NMAX, NN( J ) ) 00495 IF( NN( J ).LT.0 ) 00496 $ BADNN = .TRUE. 00497 10 CONTINUE 00498 * 00499 NBLOCK = ILAENV( 1, 'SSYTRD', 'L', NMAX, -1, -1, -1 ) 00500 NBLOCK = MIN( NMAX, MAX( 1, NBLOCK ) ) 00501 * 00502 * Check for errors 00503 * 00504 IF( NSIZES.LT.0 ) THEN 00505 INFO = -1 00506 ELSE IF( BADNN ) THEN 00507 INFO = -2 00508 ELSE IF( NTYPES.LT.0 ) THEN 00509 INFO = -3 00510 ELSE IF( LDA.LT.NMAX ) THEN 00511 INFO = -9 00512 ELSE IF( LDU.LT.NMAX ) THEN 00513 INFO = -23 00514 ELSE IF( 2*MAX( 2, NMAX )**2.GT.LWORK ) THEN 00515 INFO = -29 00516 END IF 00517 * 00518 IF( INFO.NE.0 ) THEN 00519 CALL XERBLA( 'SCHKST', -INFO ) 00520 RETURN 00521 END IF 00522 * 00523 * Quick return if possible 00524 * 00525 IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 ) 00526 $ RETURN 00527 * 00528 * More Important constants 00529 * 00530 UNFL = SLAMCH( 'Safe minimum' ) 00531 OVFL = ONE / UNFL 00532 CALL SLABAD( UNFL, OVFL ) 00533 ULP = SLAMCH( 'Epsilon' )*SLAMCH( 'Base' ) 00534 ULPINV = ONE / ULP 00535 LOG2UI = INT( LOG( ULPINV ) / LOG( TWO ) ) 00536 RTUNFL = SQRT( UNFL ) 00537 RTOVFL = SQRT( OVFL ) 00538 * 00539 * Loop over sizes, types 00540 * 00541 DO 20 I = 1, 4 00542 ISEED2( I ) = ISEED( I ) 00543 20 CONTINUE 00544 NERRS = 0 00545 NMATS = 0 00546 * 00547 DO 310 JSIZE = 1, NSIZES 00548 N = NN( JSIZE ) 00549 IF( N.GT.0 ) THEN 00550 LGN = INT( LOG( REAL( N ) ) / LOG( TWO ) ) 00551 IF( 2**LGN.LT.N ) 00552 $ LGN = LGN + 1 00553 IF( 2**LGN.LT.N ) 00554 $ LGN = LGN + 1 00555 LWEDC = 1 + 4*N + 2*N*LGN + 3*N**2 00556 LIWEDC = 6 + 6*N + 5*N*LGN 00557 ELSE 00558 LWEDC = 8 00559 LIWEDC = 12 00560 END IF 00561 NAP = ( N*( N+1 ) ) / 2 00562 ANINV = ONE / REAL( MAX( 1, N ) ) 00563 * 00564 IF( NSIZES.NE.1 ) THEN 00565 MTYPES = MIN( MAXTYP, NTYPES ) 00566 ELSE 00567 MTYPES = MIN( MAXTYP+1, NTYPES ) 00568 END IF 00569 * 00570 DO 300 JTYPE = 1, MTYPES 00571 IF( .NOT.DOTYPE( JTYPE ) ) 00572 $ GO TO 300 00573 NMATS = NMATS + 1 00574 NTEST = 0 00575 * 00576 DO 30 J = 1, 4 00577 IOLDSD( J ) = ISEED( J ) 00578 30 CONTINUE 00579 * 00580 * Compute "A" 00581 * 00582 * Control parameters: 00583 * 00584 * KMAGN KMODE KTYPE 00585 * =1 O(1) clustered 1 zero 00586 * =2 large clustered 2 identity 00587 * =3 small exponential (none) 00588 * =4 arithmetic diagonal, (w/ eigenvalues) 00589 * =5 random log symmetric, w/ eigenvalues 00590 * =6 random (none) 00591 * =7 random diagonal 00592 * =8 random symmetric 00593 * =9 positive definite 00594 * =10 diagonally dominant tridiagonal 00595 * 00596 IF( MTYPES.GT.MAXTYP ) 00597 $ GO TO 100 00598 * 00599 ITYPE = KTYPE( JTYPE ) 00600 IMODE = KMODE( JTYPE ) 00601 * 00602 * Compute norm 00603 * 00604 GO TO ( 40, 50, 60 )KMAGN( JTYPE ) 00605 * 00606 40 CONTINUE 00607 ANORM = ONE 00608 GO TO 70 00609 * 00610 50 CONTINUE 00611 ANORM = ( RTOVFL*ULP )*ANINV 00612 GO TO 70 00613 * 00614 60 CONTINUE 00615 ANORM = RTUNFL*N*ULPINV 00616 GO TO 70 00617 * 00618 70 CONTINUE 00619 * 00620 CALL SLASET( 'Full', LDA, N, ZERO, ZERO, A, LDA ) 00621 IINFO = 0 00622 IF( JTYPE.LE.15 ) THEN 00623 COND = ULPINV 00624 ELSE 00625 COND = ULPINV*ANINV / TEN 00626 END IF 00627 * 00628 * Special Matrices -- Identity & Jordan block 00629 * 00630 * Zero 00631 * 00632 IF( ITYPE.EQ.1 ) THEN 00633 IINFO = 0 00634 * 00635 ELSE IF( ITYPE.EQ.2 ) THEN 00636 * 00637 * Identity 00638 * 00639 DO 80 JC = 1, N 00640 A( JC, JC ) = ANORM 00641 80 CONTINUE 00642 * 00643 ELSE IF( ITYPE.EQ.4 ) THEN 00644 * 00645 * Diagonal Matrix, [Eigen]values Specified 00646 * 00647 CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND, 00648 $ ANORM, 0, 0, 'N', A, LDA, WORK( N+1 ), 00649 $ IINFO ) 00650 * 00651 * 00652 ELSE IF( ITYPE.EQ.5 ) THEN 00653 * 00654 * Symmetric, eigenvalues specified 00655 * 00656 CALL SLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND, 00657 $ ANORM, N, N, 'N', A, LDA, WORK( N+1 ), 00658 $ IINFO ) 00659 * 00660 ELSE IF( ITYPE.EQ.7 ) THEN 00661 * 00662 * Diagonal, random eigenvalues 00663 * 00664 CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE, 00665 $ 'T', 'N', WORK( N+1 ), 1, ONE, 00666 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0, 00667 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO ) 00668 * 00669 ELSE IF( ITYPE.EQ.8 ) THEN 00670 * 00671 * Symmetric, random eigenvalues 00672 * 00673 CALL SLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE, 00674 $ 'T', 'N', WORK( N+1 ), 1, ONE, 00675 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N, 00676 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO ) 00677 * 00678 ELSE IF( ITYPE.EQ.9 ) THEN 00679 * 00680 * Positive definite, eigenvalues specified. 00681 * 00682 CALL SLATMS( N, N, 'S', ISEED, 'P', WORK, IMODE, COND, 00683 $ ANORM, N, N, 'N', A, LDA, WORK( N+1 ), 00684 $ IINFO ) 00685 * 00686 ELSE IF( ITYPE.EQ.10 ) THEN 00687 * 00688 * Positive definite tridiagonal, eigenvalues specified. 00689 * 00690 CALL SLATMS( N, N, 'S', ISEED, 'P', WORK, IMODE, COND, 00691 $ ANORM, 1, 1, 'N', A, LDA, WORK( N+1 ), 00692 $ IINFO ) 00693 DO 90 I = 2, N 00694 TEMP1 = ABS( A( I-1, I ) ) / 00695 $ SQRT( ABS( A( I-1, I-1 )*A( I, I ) ) ) 00696 IF( TEMP1.GT.HALF ) THEN 00697 A( I-1, I ) = HALF*SQRT( ABS( A( I-1, I-1 )*A( I, 00698 $ I ) ) ) 00699 A( I, I-1 ) = A( I-1, I ) 00700 END IF 00701 90 CONTINUE 00702 * 00703 ELSE 00704 * 00705 IINFO = 1 00706 END IF 00707 * 00708 IF( IINFO.NE.0 ) THEN 00709 WRITE( NOUNIT, FMT = 9999 )'Generator', IINFO, N, JTYPE, 00710 $ IOLDSD 00711 INFO = ABS( IINFO ) 00712 RETURN 00713 END IF 00714 * 00715 100 CONTINUE 00716 * 00717 * Call SSYTRD and SORGTR to compute S and U from 00718 * upper triangle. 00719 * 00720 CALL SLACPY( 'U', N, N, A, LDA, V, LDU ) 00721 * 00722 NTEST = 1 00723 CALL SSYTRD( 'U', N, V, LDU, SD, SE, TAU, WORK, LWORK, 00724 $ IINFO ) 00725 * 00726 IF( IINFO.NE.0 ) THEN 00727 WRITE( NOUNIT, FMT = 9999 )'SSYTRD(U)', IINFO, N, JTYPE, 00728 $ IOLDSD 00729 INFO = ABS( IINFO ) 00730 IF( IINFO.LT.0 ) THEN 00731 RETURN 00732 ELSE 00733 RESULT( 1 ) = ULPINV 00734 GO TO 280 00735 END IF 00736 END IF 00737 * 00738 CALL SLACPY( 'U', N, N, V, LDU, U, LDU ) 00739 * 00740 NTEST = 2 00741 CALL SORGTR( 'U', N, U, LDU, TAU, WORK, LWORK, IINFO ) 00742 IF( IINFO.NE.0 ) THEN 00743 WRITE( NOUNIT, FMT = 9999 )'SORGTR(U)', IINFO, N, JTYPE, 00744 $ IOLDSD 00745 INFO = ABS( IINFO ) 00746 IF( IINFO.LT.0 ) THEN 00747 RETURN 00748 ELSE 00749 RESULT( 2 ) = ULPINV 00750 GO TO 280 00751 END IF 00752 END IF 00753 * 00754 * Do tests 1 and 2 00755 * 00756 CALL SSYT21( 2, 'Upper', N, 1, A, LDA, SD, SE, U, LDU, V, 00757 $ LDU, TAU, WORK, RESULT( 1 ) ) 00758 CALL SSYT21( 3, 'Upper', N, 1, A, LDA, SD, SE, U, LDU, V, 00759 $ LDU, TAU, WORK, RESULT( 2 ) ) 00760 * 00761 * Call SSYTRD and SORGTR to compute S and U from 00762 * lower triangle, do tests. 00763 * 00764 CALL SLACPY( 'L', N, N, A, LDA, V, LDU ) 00765 * 00766 NTEST = 3 00767 CALL SSYTRD( 'L', N, V, LDU, SD, SE, TAU, WORK, LWORK, 00768 $ IINFO ) 00769 * 00770 IF( IINFO.NE.0 ) THEN 00771 WRITE( NOUNIT, FMT = 9999 )'SSYTRD(L)', IINFO, N, JTYPE, 00772 $ IOLDSD 00773 INFO = ABS( IINFO ) 00774 IF( IINFO.LT.0 ) THEN 00775 RETURN 00776 ELSE 00777 RESULT( 3 ) = ULPINV 00778 GO TO 280 00779 END IF 00780 END IF 00781 * 00782 CALL SLACPY( 'L', N, N, V, LDU, U, LDU ) 00783 * 00784 NTEST = 4 00785 CALL SORGTR( 'L', N, U, LDU, TAU, WORK, LWORK, IINFO ) 00786 IF( IINFO.NE.0 ) THEN 00787 WRITE( NOUNIT, FMT = 9999 )'SORGTR(L)', IINFO, N, JTYPE, 00788 $ IOLDSD 00789 INFO = ABS( IINFO ) 00790 IF( IINFO.LT.0 ) THEN 00791 RETURN 00792 ELSE 00793 RESULT( 4 ) = ULPINV 00794 GO TO 280 00795 END IF 00796 END IF 00797 * 00798 CALL SSYT21( 2, 'Lower', N, 1, A, LDA, SD, SE, U, LDU, V, 00799 $ LDU, TAU, WORK, RESULT( 3 ) ) 00800 CALL SSYT21( 3, 'Lower', N, 1, A, LDA, SD, SE, U, LDU, V, 00801 $ LDU, TAU, WORK, RESULT( 4 ) ) 00802 * 00803 * Store the upper triangle of A in AP 00804 * 00805 I = 0 00806 DO 120 JC = 1, N 00807 DO 110 JR = 1, JC 00808 I = I + 1 00809 AP( I ) = A( JR, JC ) 00810 110 CONTINUE 00811 120 CONTINUE 00812 * 00813 * Call SSPTRD and SOPGTR to compute S and U from AP 00814 * 00815 CALL SCOPY( NAP, AP, 1, VP, 1 ) 00816 * 00817 NTEST = 5 00818 CALL SSPTRD( 'U', N, VP, SD, SE, TAU, IINFO ) 00819 * 00820 IF( IINFO.NE.0 ) THEN 00821 WRITE( NOUNIT, FMT = 9999 )'SSPTRD(U)', IINFO, N, JTYPE, 00822 $ IOLDSD 00823 INFO = ABS( IINFO ) 00824 IF( IINFO.LT.0 ) THEN 00825 RETURN 00826 ELSE 00827 RESULT( 5 ) = ULPINV 00828 GO TO 280 00829 END IF 00830 END IF 00831 * 00832 NTEST = 6 00833 CALL SOPGTR( 'U', N, VP, TAU, U, LDU, WORK, IINFO ) 00834 IF( IINFO.NE.0 ) THEN 00835 WRITE( NOUNIT, FMT = 9999 )'SOPGTR(U)', IINFO, N, JTYPE, 00836 $ IOLDSD 00837 INFO = ABS( IINFO ) 00838 IF( IINFO.LT.0 ) THEN 00839 RETURN 00840 ELSE 00841 RESULT( 6 ) = ULPINV 00842 GO TO 280 00843 END IF 00844 END IF 00845 * 00846 * Do tests 5 and 6 00847 * 00848 CALL SSPT21( 2, 'Upper', N, 1, AP, SD, SE, U, LDU, VP, TAU, 00849 $ WORK, RESULT( 5 ) ) 00850 CALL SSPT21( 3, 'Upper', N, 1, AP, SD, SE, U, LDU, VP, TAU, 00851 $ WORK, RESULT( 6 ) ) 00852 * 00853 * Store the lower triangle of A in AP 00854 * 00855 I = 0 00856 DO 140 JC = 1, N 00857 DO 130 JR = JC, N 00858 I = I + 1 00859 AP( I ) = A( JR, JC ) 00860 130 CONTINUE 00861 140 CONTINUE 00862 * 00863 * Call SSPTRD and SOPGTR to compute S and U from AP 00864 * 00865 CALL SCOPY( NAP, AP, 1, VP, 1 ) 00866 * 00867 NTEST = 7 00868 CALL SSPTRD( 'L', N, VP, SD, SE, TAU, IINFO ) 00869 * 00870 IF( IINFO.NE.0 ) THEN 00871 WRITE( NOUNIT, FMT = 9999 )'SSPTRD(L)', IINFO, N, JTYPE, 00872 $ IOLDSD 00873 INFO = ABS( IINFO ) 00874 IF( IINFO.LT.0 ) THEN 00875 RETURN 00876 ELSE 00877 RESULT( 7 ) = ULPINV 00878 GO TO 280 00879 END IF 00880 END IF 00881 * 00882 NTEST = 8 00883 CALL SOPGTR( 'L', N, VP, TAU, U, LDU, WORK, IINFO ) 00884 IF( IINFO.NE.0 ) THEN 00885 WRITE( NOUNIT, FMT = 9999 )'SOPGTR(L)', IINFO, N, JTYPE, 00886 $ IOLDSD 00887 INFO = ABS( IINFO ) 00888 IF( IINFO.LT.0 ) THEN 00889 RETURN 00890 ELSE 00891 RESULT( 8 ) = ULPINV 00892 GO TO 280 00893 END IF 00894 END IF 00895 * 00896 CALL SSPT21( 2, 'Lower', N, 1, AP, SD, SE, U, LDU, VP, TAU, 00897 $ WORK, RESULT( 7 ) ) 00898 CALL SSPT21( 3, 'Lower', N, 1, AP, SD, SE, U, LDU, VP, TAU, 00899 $ WORK, RESULT( 8 ) ) 00900 * 00901 * Call SSTEQR to compute D1, D2, and Z, do tests. 00902 * 00903 * Compute D1 and Z 00904 * 00905 CALL SCOPY( N, SD, 1, D1, 1 ) 00906 IF( N.GT.0 ) 00907 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 00908 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 00909 * 00910 NTEST = 9 00911 CALL SSTEQR( 'V', N, D1, WORK, Z, LDU, WORK( N+1 ), IINFO ) 00912 IF( IINFO.NE.0 ) THEN 00913 WRITE( NOUNIT, FMT = 9999 )'SSTEQR(V)', IINFO, N, JTYPE, 00914 $ IOLDSD 00915 INFO = ABS( IINFO ) 00916 IF( IINFO.LT.0 ) THEN 00917 RETURN 00918 ELSE 00919 RESULT( 9 ) = ULPINV 00920 GO TO 280 00921 END IF 00922 END IF 00923 * 00924 * Compute D2 00925 * 00926 CALL SCOPY( N, SD, 1, D2, 1 ) 00927 IF( N.GT.0 ) 00928 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 00929 * 00930 NTEST = 11 00931 CALL SSTEQR( 'N', N, D2, WORK, WORK( N+1 ), LDU, 00932 $ WORK( N+1 ), IINFO ) 00933 IF( IINFO.NE.0 ) THEN 00934 WRITE( NOUNIT, FMT = 9999 )'SSTEQR(N)', IINFO, N, JTYPE, 00935 $ IOLDSD 00936 INFO = ABS( IINFO ) 00937 IF( IINFO.LT.0 ) THEN 00938 RETURN 00939 ELSE 00940 RESULT( 11 ) = ULPINV 00941 GO TO 280 00942 END IF 00943 END IF 00944 * 00945 * Compute D3 (using PWK method) 00946 * 00947 CALL SCOPY( N, SD, 1, D3, 1 ) 00948 IF( N.GT.0 ) 00949 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 00950 * 00951 NTEST = 12 00952 CALL SSTERF( N, D3, WORK, IINFO ) 00953 IF( IINFO.NE.0 ) THEN 00954 WRITE( NOUNIT, FMT = 9999 )'SSTERF', IINFO, N, JTYPE, 00955 $ IOLDSD 00956 INFO = ABS( IINFO ) 00957 IF( IINFO.LT.0 ) THEN 00958 RETURN 00959 ELSE 00960 RESULT( 12 ) = ULPINV 00961 GO TO 280 00962 END IF 00963 END IF 00964 * 00965 * Do Tests 9 and 10 00966 * 00967 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 00968 $ RESULT( 9 ) ) 00969 * 00970 * Do Tests 11 and 12 00971 * 00972 TEMP1 = ZERO 00973 TEMP2 = ZERO 00974 TEMP3 = ZERO 00975 TEMP4 = ZERO 00976 * 00977 DO 150 J = 1, N 00978 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 00979 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 00980 TEMP3 = MAX( TEMP3, ABS( D1( J ) ), ABS( D3( J ) ) ) 00981 TEMP4 = MAX( TEMP4, ABS( D1( J )-D3( J ) ) ) 00982 150 CONTINUE 00983 * 00984 RESULT( 11 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 00985 RESULT( 12 ) = TEMP4 / MAX( UNFL, ULP*MAX( TEMP3, TEMP4 ) ) 00986 * 00987 * Do Test 13 -- Sturm Sequence Test of Eigenvalues 00988 * Go up by factors of two until it succeeds 00989 * 00990 NTEST = 13 00991 TEMP1 = THRESH*( HALF-ULP ) 00992 * 00993 DO 160 J = 0, LOG2UI 00994 CALL SSTECH( N, SD, SE, D1, TEMP1, WORK, IINFO ) 00995 IF( IINFO.EQ.0 ) 00996 $ GO TO 170 00997 TEMP1 = TEMP1*TWO 00998 160 CONTINUE 00999 * 01000 170 CONTINUE 01001 RESULT( 13 ) = TEMP1 01002 * 01003 * For positive definite matrices ( JTYPE.GT.15 ) call SPTEQR 01004 * and do tests 14, 15, and 16 . 01005 * 01006 IF( JTYPE.GT.15 ) THEN 01007 * 01008 * Compute D4 and Z4 01009 * 01010 CALL SCOPY( N, SD, 1, D4, 1 ) 01011 IF( N.GT.0 ) 01012 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01013 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 01014 * 01015 NTEST = 14 01016 CALL SPTEQR( 'V', N, D4, WORK, Z, LDU, WORK( N+1 ), 01017 $ IINFO ) 01018 IF( IINFO.NE.0 ) THEN 01019 WRITE( NOUNIT, FMT = 9999 )'SPTEQR(V)', IINFO, N, 01020 $ JTYPE, IOLDSD 01021 INFO = ABS( IINFO ) 01022 IF( IINFO.LT.0 ) THEN 01023 RETURN 01024 ELSE 01025 RESULT( 14 ) = ULPINV 01026 GO TO 280 01027 END IF 01028 END IF 01029 * 01030 * Do Tests 14 and 15 01031 * 01032 CALL SSTT21( N, 0, SD, SE, D4, DUMMA, Z, LDU, WORK, 01033 $ RESULT( 14 ) ) 01034 * 01035 * Compute D5 01036 * 01037 CALL SCOPY( N, SD, 1, D5, 1 ) 01038 IF( N.GT.0 ) 01039 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01040 * 01041 NTEST = 16 01042 CALL SPTEQR( 'N', N, D5, WORK, Z, LDU, WORK( N+1 ), 01043 $ IINFO ) 01044 IF( IINFO.NE.0 ) THEN 01045 WRITE( NOUNIT, FMT = 9999 )'SPTEQR(N)', IINFO, N, 01046 $ JTYPE, IOLDSD 01047 INFO = ABS( IINFO ) 01048 IF( IINFO.LT.0 ) THEN 01049 RETURN 01050 ELSE 01051 RESULT( 16 ) = ULPINV 01052 GO TO 280 01053 END IF 01054 END IF 01055 * 01056 * Do Test 16 01057 * 01058 TEMP1 = ZERO 01059 TEMP2 = ZERO 01060 DO 180 J = 1, N 01061 TEMP1 = MAX( TEMP1, ABS( D4( J ) ), ABS( D5( J ) ) ) 01062 TEMP2 = MAX( TEMP2, ABS( D4( J )-D5( J ) ) ) 01063 180 CONTINUE 01064 * 01065 RESULT( 16 ) = TEMP2 / MAX( UNFL, 01066 $ HUN*ULP*MAX( TEMP1, TEMP2 ) ) 01067 ELSE 01068 RESULT( 14 ) = ZERO 01069 RESULT( 15 ) = ZERO 01070 RESULT( 16 ) = ZERO 01071 END IF 01072 * 01073 * Call SSTEBZ with different options and do tests 17-18. 01074 * 01075 * If S is positive definite and diagonally dominant, 01076 * ask for all eigenvalues with high relative accuracy. 01077 * 01078 VL = ZERO 01079 VU = ZERO 01080 IL = 0 01081 IU = 0 01082 IF( JTYPE.EQ.21 ) THEN 01083 NTEST = 17 01084 ABSTOL = UNFL + UNFL 01085 CALL SSTEBZ( 'A', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 01086 $ M, NSPLIT, WR, IWORK( 1 ), IWORK( N+1 ), 01087 $ WORK, IWORK( 2*N+1 ), IINFO ) 01088 IF( IINFO.NE.0 ) THEN 01089 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A,rel)', IINFO, N, 01090 $ JTYPE, IOLDSD 01091 INFO = ABS( IINFO ) 01092 IF( IINFO.LT.0 ) THEN 01093 RETURN 01094 ELSE 01095 RESULT( 17 ) = ULPINV 01096 GO TO 280 01097 END IF 01098 END IF 01099 * 01100 * Do test 17 01101 * 01102 TEMP2 = TWO*( TWO*N-ONE )*ULP*( ONE+EIGHT*HALF**2 ) / 01103 $ ( ONE-HALF )**4 01104 * 01105 TEMP1 = ZERO 01106 DO 190 J = 1, N 01107 TEMP1 = MAX( TEMP1, ABS( D4( J )-WR( N-J+1 ) ) / 01108 $ ( ABSTOL+ABS( D4( J ) ) ) ) 01109 190 CONTINUE 01110 * 01111 RESULT( 17 ) = TEMP1 / TEMP2 01112 ELSE 01113 RESULT( 17 ) = ZERO 01114 END IF 01115 * 01116 * Now ask for all eigenvalues with high absolute accuracy. 01117 * 01118 NTEST = 18 01119 ABSTOL = UNFL + UNFL 01120 CALL SSTEBZ( 'A', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, M, 01121 $ NSPLIT, WA1, IWORK( 1 ), IWORK( N+1 ), WORK, 01122 $ IWORK( 2*N+1 ), IINFO ) 01123 IF( IINFO.NE.0 ) THEN 01124 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A)', IINFO, N, JTYPE, 01125 $ IOLDSD 01126 INFO = ABS( IINFO ) 01127 IF( IINFO.LT.0 ) THEN 01128 RETURN 01129 ELSE 01130 RESULT( 18 ) = ULPINV 01131 GO TO 280 01132 END IF 01133 END IF 01134 * 01135 * Do test 18 01136 * 01137 TEMP1 = ZERO 01138 TEMP2 = ZERO 01139 DO 200 J = 1, N 01140 TEMP1 = MAX( TEMP1, ABS( D3( J ) ), ABS( WA1( J ) ) ) 01141 TEMP2 = MAX( TEMP2, ABS( D3( J )-WA1( J ) ) ) 01142 200 CONTINUE 01143 * 01144 RESULT( 18 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 01145 * 01146 * Choose random values for IL and IU, and ask for the 01147 * IL-th through IU-th eigenvalues. 01148 * 01149 NTEST = 19 01150 IF( N.LE.1 ) THEN 01151 IL = 1 01152 IU = N 01153 ELSE 01154 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 01155 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 01156 IF( IU.LT.IL ) THEN 01157 ITEMP = IU 01158 IU = IL 01159 IL = ITEMP 01160 END IF 01161 END IF 01162 * 01163 CALL SSTEBZ( 'I', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 01164 $ M2, NSPLIT, WA2, IWORK( 1 ), IWORK( N+1 ), 01165 $ WORK, IWORK( 2*N+1 ), IINFO ) 01166 IF( IINFO.NE.0 ) THEN 01167 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(I)', IINFO, N, JTYPE, 01168 $ IOLDSD 01169 INFO = ABS( IINFO ) 01170 IF( IINFO.LT.0 ) THEN 01171 RETURN 01172 ELSE 01173 RESULT( 19 ) = ULPINV 01174 GO TO 280 01175 END IF 01176 END IF 01177 * 01178 * Determine the values VL and VU of the IL-th and IU-th 01179 * eigenvalues and ask for all eigenvalues in this range. 01180 * 01181 IF( N.GT.0 ) THEN 01182 IF( IL.NE.1 ) THEN 01183 VL = WA1( IL ) - MAX( HALF*( WA1( IL )-WA1( IL-1 ) ), 01184 $ ULP*ANORM, TWO*RTUNFL ) 01185 ELSE 01186 VL = WA1( 1 ) - MAX( HALF*( WA1( N )-WA1( 1 ) ), 01187 $ ULP*ANORM, TWO*RTUNFL ) 01188 END IF 01189 IF( IU.NE.N ) THEN 01190 VU = WA1( IU ) + MAX( HALF*( WA1( IU+1 )-WA1( IU ) ), 01191 $ ULP*ANORM, TWO*RTUNFL ) 01192 ELSE 01193 VU = WA1( N ) + MAX( HALF*( WA1( N )-WA1( 1 ) ), 01194 $ ULP*ANORM, TWO*RTUNFL ) 01195 END IF 01196 ELSE 01197 VL = ZERO 01198 VU = ONE 01199 END IF 01200 * 01201 CALL SSTEBZ( 'V', 'E', N, VL, VU, IL, IU, ABSTOL, SD, SE, 01202 $ M3, NSPLIT, WA3, IWORK( 1 ), IWORK( N+1 ), 01203 $ WORK, IWORK( 2*N+1 ), IINFO ) 01204 IF( IINFO.NE.0 ) THEN 01205 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(V)', IINFO, N, JTYPE, 01206 $ IOLDSD 01207 INFO = ABS( IINFO ) 01208 IF( IINFO.LT.0 ) THEN 01209 RETURN 01210 ELSE 01211 RESULT( 19 ) = ULPINV 01212 GO TO 280 01213 END IF 01214 END IF 01215 * 01216 IF( M3.EQ.0 .AND. N.NE.0 ) THEN 01217 RESULT( 19 ) = ULPINV 01218 GO TO 280 01219 END IF 01220 * 01221 * Do test 19 01222 * 01223 TEMP1 = SSXT1( 1, WA2, M2, WA3, M3, ABSTOL, ULP, UNFL ) 01224 TEMP2 = SSXT1( 1, WA3, M3, WA2, M2, ABSTOL, ULP, UNFL ) 01225 IF( N.GT.0 ) THEN 01226 TEMP3 = MAX( ABS( WA1( N ) ), ABS( WA1( 1 ) ) ) 01227 ELSE 01228 TEMP3 = ZERO 01229 END IF 01230 * 01231 RESULT( 19 ) = ( TEMP1+TEMP2 ) / MAX( UNFL, TEMP3*ULP ) 01232 * 01233 * Call SSTEIN to compute eigenvectors corresponding to 01234 * eigenvalues in WA1. (First call SSTEBZ again, to make sure 01235 * it returns these eigenvalues in the correct order.) 01236 * 01237 NTEST = 21 01238 CALL SSTEBZ( 'A', 'B', N, VL, VU, IL, IU, ABSTOL, SD, SE, M, 01239 $ NSPLIT, WA1, IWORK( 1 ), IWORK( N+1 ), WORK, 01240 $ IWORK( 2*N+1 ), IINFO ) 01241 IF( IINFO.NE.0 ) THEN 01242 WRITE( NOUNIT, FMT = 9999 )'SSTEBZ(A,B)', IINFO, N, 01243 $ JTYPE, IOLDSD 01244 INFO = ABS( IINFO ) 01245 IF( IINFO.LT.0 ) THEN 01246 RETURN 01247 ELSE 01248 RESULT( 20 ) = ULPINV 01249 RESULT( 21 ) = ULPINV 01250 GO TO 280 01251 END IF 01252 END IF 01253 * 01254 CALL SSTEIN( N, SD, SE, M, WA1, IWORK( 1 ), IWORK( N+1 ), Z, 01255 $ LDU, WORK, IWORK( 2*N+1 ), IWORK( 3*N+1 ), 01256 $ IINFO ) 01257 IF( IINFO.NE.0 ) THEN 01258 WRITE( NOUNIT, FMT = 9999 )'SSTEIN', IINFO, N, JTYPE, 01259 $ IOLDSD 01260 INFO = ABS( IINFO ) 01261 IF( IINFO.LT.0 ) THEN 01262 RETURN 01263 ELSE 01264 RESULT( 20 ) = ULPINV 01265 RESULT( 21 ) = ULPINV 01266 GO TO 280 01267 END IF 01268 END IF 01269 * 01270 * Do tests 20 and 21 01271 * 01272 CALL SSTT21( N, 0, SD, SE, WA1, DUMMA, Z, LDU, WORK, 01273 $ RESULT( 20 ) ) 01274 * 01275 * Call SSTEDC(I) to compute D1 and Z, do tests. 01276 * 01277 * Compute D1 and Z 01278 * 01279 CALL SCOPY( N, SD, 1, D1, 1 ) 01280 IF( N.GT.0 ) 01281 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01282 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 01283 * 01284 NTEST = 22 01285 CALL SSTEDC( 'I', N, D1, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 01286 $ IWORK, LIWEDC, IINFO ) 01287 IF( IINFO.NE.0 ) THEN 01288 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(I)', IINFO, N, JTYPE, 01289 $ IOLDSD 01290 INFO = ABS( IINFO ) 01291 IF( IINFO.LT.0 ) THEN 01292 RETURN 01293 ELSE 01294 RESULT( 22 ) = ULPINV 01295 GO TO 280 01296 END IF 01297 END IF 01298 * 01299 * Do Tests 22 and 23 01300 * 01301 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 01302 $ RESULT( 22 ) ) 01303 * 01304 * Call SSTEDC(V) to compute D1 and Z, do tests. 01305 * 01306 * Compute D1 and Z 01307 * 01308 CALL SCOPY( N, SD, 1, D1, 1 ) 01309 IF( N.GT.0 ) 01310 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01311 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 01312 * 01313 NTEST = 24 01314 CALL SSTEDC( 'V', N, D1, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 01315 $ IWORK, LIWEDC, IINFO ) 01316 IF( IINFO.NE.0 ) THEN 01317 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(V)', IINFO, N, JTYPE, 01318 $ IOLDSD 01319 INFO = ABS( IINFO ) 01320 IF( IINFO.LT.0 ) THEN 01321 RETURN 01322 ELSE 01323 RESULT( 24 ) = ULPINV 01324 GO TO 280 01325 END IF 01326 END IF 01327 * 01328 * Do Tests 24 and 25 01329 * 01330 CALL SSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 01331 $ RESULT( 24 ) ) 01332 * 01333 * Call SSTEDC(N) to compute D2, do tests. 01334 * 01335 * Compute D2 01336 * 01337 CALL SCOPY( N, SD, 1, D2, 1 ) 01338 IF( N.GT.0 ) 01339 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01340 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 01341 * 01342 NTEST = 26 01343 CALL SSTEDC( 'N', N, D2, WORK, Z, LDU, WORK( N+1 ), LWEDC-N, 01344 $ IWORK, LIWEDC, IINFO ) 01345 IF( IINFO.NE.0 ) THEN 01346 WRITE( NOUNIT, FMT = 9999 )'SSTEDC(N)', IINFO, N, JTYPE, 01347 $ IOLDSD 01348 INFO = ABS( IINFO ) 01349 IF( IINFO.LT.0 ) THEN 01350 RETURN 01351 ELSE 01352 RESULT( 26 ) = ULPINV 01353 GO TO 280 01354 END IF 01355 END IF 01356 * 01357 * Do Test 26 01358 * 01359 TEMP1 = ZERO 01360 TEMP2 = ZERO 01361 * 01362 DO 210 J = 1, N 01363 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 01364 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 01365 210 CONTINUE 01366 * 01367 RESULT( 26 ) = TEMP2 / MAX( UNFL, ULP*MAX( TEMP1, TEMP2 ) ) 01368 * 01369 * Only test SSTEMR if IEEE compliant 01370 * 01371 IF( ILAENV( 10, 'SSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 .AND. 01372 $ ILAENV( 11, 'SSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 ) THEN 01373 * 01374 * Call SSTEMR, do test 27 (relative eigenvalue accuracy) 01375 * 01376 * If S is positive definite and diagonally dominant, 01377 * ask for all eigenvalues with high relative accuracy. 01378 * 01379 VL = ZERO 01380 VU = ZERO 01381 IL = 0 01382 IU = 0 01383 IF( JTYPE.EQ.21 .AND. SREL ) THEN 01384 NTEST = 27 01385 ABSTOL = UNFL + UNFL 01386 CALL SSTEMR( 'V', 'A', N, SD, SE, VL, VU, IL, IU, 01387 $ M, WR, Z, LDU, N, IWORK( 1 ), TRYRAC, 01388 $ WORK, LWORK, IWORK( 2*N+1 ), LWORK-2*N, 01389 $ IINFO ) 01390 IF( IINFO.NE.0 ) THEN 01391 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,A,rel)', 01392 $ IINFO, N, JTYPE, IOLDSD 01393 INFO = ABS( IINFO ) 01394 IF( IINFO.LT.0 ) THEN 01395 RETURN 01396 ELSE 01397 RESULT( 27 ) = ULPINV 01398 GO TO 270 01399 END IF 01400 END IF 01401 * 01402 * Do test 27 01403 * 01404 TEMP2 = TWO*( TWO*N-ONE )*ULP*( ONE+EIGHT*HALF**2 ) / 01405 $ ( ONE-HALF )**4 01406 * 01407 TEMP1 = ZERO 01408 DO 220 J = 1, N 01409 TEMP1 = MAX( TEMP1, ABS( D4( J )-WR( N-J+1 ) ) / 01410 $ ( ABSTOL+ABS( D4( J ) ) ) ) 01411 220 CONTINUE 01412 * 01413 RESULT( 27 ) = TEMP1 / TEMP2 01414 * 01415 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 01416 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 01417 IF( IU.LT.IL ) THEN 01418 ITEMP = IU 01419 IU = IL 01420 IL = ITEMP 01421 END IF 01422 * 01423 IF( SRANGE ) THEN 01424 NTEST = 28 01425 ABSTOL = UNFL + UNFL 01426 CALL SSTEMR( 'V', 'I', N, SD, SE, VL, VU, IL, IU, 01427 $ M, WR, Z, LDU, N, IWORK( 1 ), TRYRAC, 01428 $ WORK, LWORK, IWORK( 2*N+1 ), 01429 $ LWORK-2*N, IINFO ) 01430 * 01431 IF( IINFO.NE.0 ) THEN 01432 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,I,rel)', 01433 $ IINFO, N, JTYPE, IOLDSD 01434 INFO = ABS( IINFO ) 01435 IF( IINFO.LT.0 ) THEN 01436 RETURN 01437 ELSE 01438 RESULT( 28 ) = ULPINV 01439 GO TO 270 01440 END IF 01441 END IF 01442 * 01443 * 01444 * Do test 28 01445 * 01446 TEMP2 = TWO*( TWO*N-ONE )*ULP* 01447 $ ( ONE+EIGHT*HALF**2 ) / ( ONE-HALF )**4 01448 * 01449 TEMP1 = ZERO 01450 DO 230 J = IL, IU 01451 TEMP1 = MAX( TEMP1, ABS( WR( J-IL+1 )-D4( N-J+ 01452 $ 1 ) ) / ( ABSTOL+ABS( WR( J-IL+1 ) ) ) ) 01453 230 CONTINUE 01454 * 01455 RESULT( 28 ) = TEMP1 / TEMP2 01456 ELSE 01457 RESULT( 28 ) = ZERO 01458 END IF 01459 ELSE 01460 RESULT( 27 ) = ZERO 01461 RESULT( 28 ) = ZERO 01462 END IF 01463 * 01464 * Call SSTEMR(V,I) to compute D1 and Z, do tests. 01465 * 01466 * Compute D1 and Z 01467 * 01468 CALL SCOPY( N, SD, 1, D5, 1 ) 01469 IF( N.GT.0 ) 01470 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01471 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 01472 * 01473 IF( SRANGE ) THEN 01474 NTEST = 29 01475 IL = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 01476 IU = 1 + ( N-1 )*INT( SLARND( 1, ISEED2 ) ) 01477 IF( IU.LT.IL ) THEN 01478 ITEMP = IU 01479 IU = IL 01480 IL = ITEMP 01481 END IF 01482 CALL SSTEMR( 'V', 'I', N, D5, WORK, VL, VU, IL, IU, 01483 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 01484 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 01485 $ LIWORK-2*N, IINFO ) 01486 IF( IINFO.NE.0 ) THEN 01487 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,I)', IINFO, 01488 $ N, JTYPE, IOLDSD 01489 INFO = ABS( IINFO ) 01490 IF( IINFO.LT.0 ) THEN 01491 RETURN 01492 ELSE 01493 RESULT( 29 ) = ULPINV 01494 GO TO 280 01495 END IF 01496 END IF 01497 * 01498 * Do Tests 29 and 30 01499 * 01500 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 01501 $ M, RESULT( 29 ) ) 01502 * 01503 * Call SSTEMR to compute D2, do tests. 01504 * 01505 * Compute D2 01506 * 01507 CALL SCOPY( N, SD, 1, D5, 1 ) 01508 IF( N.GT.0 ) 01509 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01510 * 01511 NTEST = 31 01512 CALL SSTEMR( 'N', 'I', N, D5, WORK, VL, VU, IL, IU, 01513 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 01514 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 01515 $ LIWORK-2*N, IINFO ) 01516 IF( IINFO.NE.0 ) THEN 01517 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,I)', IINFO, 01518 $ N, JTYPE, IOLDSD 01519 INFO = ABS( IINFO ) 01520 IF( IINFO.LT.0 ) THEN 01521 RETURN 01522 ELSE 01523 RESULT( 31 ) = ULPINV 01524 GO TO 280 01525 END IF 01526 END IF 01527 * 01528 * Do Test 31 01529 * 01530 TEMP1 = ZERO 01531 TEMP2 = ZERO 01532 * 01533 DO 240 J = 1, IU - IL + 1 01534 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), 01535 $ ABS( D2( J ) ) ) 01536 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 01537 240 CONTINUE 01538 * 01539 RESULT( 31 ) = TEMP2 / MAX( UNFL, 01540 $ ULP*MAX( TEMP1, TEMP2 ) ) 01541 * 01542 * 01543 * Call SSTEMR(V,V) to compute D1 and Z, do tests. 01544 * 01545 * Compute D1 and Z 01546 * 01547 CALL SCOPY( N, SD, 1, D5, 1 ) 01548 IF( N.GT.0 ) 01549 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01550 CALL SLASET( 'Full', N, N, ZERO, ONE, Z, LDU ) 01551 * 01552 NTEST = 32 01553 * 01554 IF( N.GT.0 ) THEN 01555 IF( IL.NE.1 ) THEN 01556 VL = D2( IL ) - MAX( HALF* 01557 $ ( D2( IL )-D2( IL-1 ) ), ULP*ANORM, 01558 $ TWO*RTUNFL ) 01559 ELSE 01560 VL = D2( 1 ) - MAX( HALF*( D2( N )-D2( 1 ) ), 01561 $ ULP*ANORM, TWO*RTUNFL ) 01562 END IF 01563 IF( IU.NE.N ) THEN 01564 VU = D2( IU ) + MAX( HALF* 01565 $ ( D2( IU+1 )-D2( IU ) ), ULP*ANORM, 01566 $ TWO*RTUNFL ) 01567 ELSE 01568 VU = D2( N ) + MAX( HALF*( D2( N )-D2( 1 ) ), 01569 $ ULP*ANORM, TWO*RTUNFL ) 01570 END IF 01571 ELSE 01572 VL = ZERO 01573 VU = ONE 01574 END IF 01575 * 01576 CALL SSTEMR( 'V', 'V', N, D5, WORK, VL, VU, IL, IU, 01577 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 01578 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 01579 $ LIWORK-2*N, IINFO ) 01580 IF( IINFO.NE.0 ) THEN 01581 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,V)', IINFO, 01582 $ N, JTYPE, IOLDSD 01583 INFO = ABS( IINFO ) 01584 IF( IINFO.LT.0 ) THEN 01585 RETURN 01586 ELSE 01587 RESULT( 32 ) = ULPINV 01588 GO TO 280 01589 END IF 01590 END IF 01591 * 01592 * Do Tests 32 and 33 01593 * 01594 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, 01595 $ M, RESULT( 32 ) ) 01596 * 01597 * Call SSTEMR to compute D2, do tests. 01598 * 01599 * Compute D2 01600 * 01601 CALL SCOPY( N, SD, 1, D5, 1 ) 01602 IF( N.GT.0 ) 01603 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01604 * 01605 NTEST = 34 01606 CALL SSTEMR( 'N', 'V', N, D5, WORK, VL, VU, IL, IU, 01607 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 01608 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 01609 $ LIWORK-2*N, IINFO ) 01610 IF( IINFO.NE.0 ) THEN 01611 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,V)', IINFO, 01612 $ N, JTYPE, IOLDSD 01613 INFO = ABS( IINFO ) 01614 IF( IINFO.LT.0 ) THEN 01615 RETURN 01616 ELSE 01617 RESULT( 34 ) = ULPINV 01618 GO TO 280 01619 END IF 01620 END IF 01621 * 01622 * Do Test 34 01623 * 01624 TEMP1 = ZERO 01625 TEMP2 = ZERO 01626 * 01627 DO 250 J = 1, IU - IL + 1 01628 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), 01629 $ ABS( D2( J ) ) ) 01630 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 01631 250 CONTINUE 01632 * 01633 RESULT( 34 ) = TEMP2 / MAX( UNFL, 01634 $ ULP*MAX( TEMP1, TEMP2 ) ) 01635 ELSE 01636 RESULT( 29 ) = ZERO 01637 RESULT( 30 ) = ZERO 01638 RESULT( 31 ) = ZERO 01639 RESULT( 32 ) = ZERO 01640 RESULT( 33 ) = ZERO 01641 RESULT( 34 ) = ZERO 01642 END IF 01643 * 01644 * 01645 * Call SSTEMR(V,A) to compute D1 and Z, do tests. 01646 * 01647 * Compute D1 and Z 01648 * 01649 CALL SCOPY( N, SD, 1, D5, 1 ) 01650 IF( N.GT.0 ) 01651 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01652 * 01653 NTEST = 35 01654 * 01655 CALL SSTEMR( 'V', 'A', N, D5, WORK, VL, VU, IL, IU, 01656 $ M, D1, Z, LDU, N, IWORK( 1 ), TRYRAC, 01657 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 01658 $ LIWORK-2*N, IINFO ) 01659 IF( IINFO.NE.0 ) THEN 01660 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(V,A)', IINFO, N, 01661 $ JTYPE, IOLDSD 01662 INFO = ABS( IINFO ) 01663 IF( IINFO.LT.0 ) THEN 01664 RETURN 01665 ELSE 01666 RESULT( 35 ) = ULPINV 01667 GO TO 280 01668 END IF 01669 END IF 01670 * 01671 * Do Tests 35 and 36 01672 * 01673 CALL SSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, M, 01674 $ RESULT( 35 ) ) 01675 * 01676 * Call SSTEMR to compute D2, do tests. 01677 * 01678 * Compute D2 01679 * 01680 CALL SCOPY( N, SD, 1, D5, 1 ) 01681 IF( N.GT.0 ) 01682 $ CALL SCOPY( N-1, SE, 1, WORK, 1 ) 01683 * 01684 NTEST = 37 01685 CALL SSTEMR( 'N', 'A', N, D5, WORK, VL, VU, IL, IU, 01686 $ M, D2, Z, LDU, N, IWORK( 1 ), TRYRAC, 01687 $ WORK( N+1 ), LWORK-N, IWORK( 2*N+1 ), 01688 $ LIWORK-2*N, IINFO ) 01689 IF( IINFO.NE.0 ) THEN 01690 WRITE( NOUNIT, FMT = 9999 )'SSTEMR(N,A)', IINFO, N, 01691 $ JTYPE, IOLDSD 01692 INFO = ABS( IINFO ) 01693 IF( IINFO.LT.0 ) THEN 01694 RETURN 01695 ELSE 01696 RESULT( 37 ) = ULPINV 01697 GO TO 280 01698 END IF 01699 END IF 01700 * 01701 * Do Test 34 01702 * 01703 TEMP1 = ZERO 01704 TEMP2 = ZERO 01705 * 01706 DO 260 J = 1, N 01707 TEMP1 = MAX( TEMP1, ABS( D1( J ) ), ABS( D2( J ) ) ) 01708 TEMP2 = MAX( TEMP2, ABS( D1( J )-D2( J ) ) ) 01709 260 CONTINUE 01710 * 01711 RESULT( 37 ) = TEMP2 / MAX( UNFL, 01712 $ ULP*MAX( TEMP1, TEMP2 ) ) 01713 END IF 01714 270 CONTINUE 01715 280 CONTINUE 01716 NTESTT = NTESTT + NTEST 01717 * 01718 * End of Loop -- Check for RESULT(j) > THRESH 01719 * 01720 * 01721 * Print out tests which fail. 01722 * 01723 DO 290 JR = 1, NTEST 01724 IF( RESULT( JR ).GE.THRESH ) THEN 01725 * 01726 * If this is the first test to fail, 01727 * print a header to the data file. 01728 * 01729 IF( NERRS.EQ.0 ) THEN 01730 WRITE( NOUNIT, FMT = 9998 )'SST' 01731 WRITE( NOUNIT, FMT = 9997 ) 01732 WRITE( NOUNIT, FMT = 9996 ) 01733 WRITE( NOUNIT, FMT = 9995 )'Symmetric' 01734 WRITE( NOUNIT, FMT = 9994 ) 01735 * 01736 * Tests performed 01737 * 01738 WRITE( NOUNIT, FMT = 9988 ) 01739 END IF 01740 NERRS = NERRS + 1 01741 WRITE( NOUNIT, FMT = 9990 )N, IOLDSD, JTYPE, JR, 01742 $ RESULT( JR ) 01743 END IF 01744 290 CONTINUE 01745 300 CONTINUE 01746 310 CONTINUE 01747 * 01748 * Summary 01749 * 01750 CALL SLASUM( 'SST', NOUNIT, NERRS, NTESTT ) 01751 RETURN 01752 * 01753 9999 FORMAT( ' SCHKST: ', A, ' returned INFO=', I6, '.', / 9X, 'N=', 01754 $ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' ) 01755 * 01756 9998 FORMAT( / 1X, A3, ' -- Real Symmetric eigenvalue problem' ) 01757 9997 FORMAT( ' Matrix types (see SCHKST for details): ' ) 01758 * 01759 9996 FORMAT( / ' Special Matrices:', 01760 $ / ' 1=Zero matrix. ', 01761 $ ' 5=Diagonal: clustered entries.', 01762 $ / ' 2=Identity matrix. ', 01763 $ ' 6=Diagonal: large, evenly spaced.', 01764 $ / ' 3=Diagonal: evenly spaced entries. ', 01765 $ ' 7=Diagonal: small, evenly spaced.', 01766 $ / ' 4=Diagonal: geometr. spaced entries.' ) 01767 9995 FORMAT( ' Dense ', A, ' Matrices:', 01768 $ / ' 8=Evenly spaced eigenvals. ', 01769 $ ' 12=Small, evenly spaced eigenvals.', 01770 $ / ' 9=Geometrically spaced eigenvals. ', 01771 $ ' 13=Matrix with random O(1) entries.', 01772 $ / ' 10=Clustered eigenvalues. ', 01773 $ ' 14=Matrix with large random entries.', 01774 $ / ' 11=Large, evenly spaced eigenvals. ', 01775 $ ' 15=Matrix with small random entries.' ) 01776 9994 FORMAT( ' 16=Positive definite, evenly spaced eigenvalues', 01777 $ / ' 17=Positive definite, geometrically spaced eigenvlaues', 01778 $ / ' 18=Positive definite, clustered eigenvalues', 01779 $ / ' 19=Positive definite, small evenly spaced eigenvalues', 01780 $ / ' 20=Positive definite, large evenly spaced eigenvalues', 01781 $ / ' 21=Diagonally dominant tridiagonal, geometrically', 01782 $ ' spaced eigenvalues' ) 01783 * 01784 9993 FORMAT( / ' Tests performed: ', 01785 $ '(S is Tridiag, D is diagonal, U and Z are ', A, ',', / 20X, 01786 $ A, ', W is a diagonal matrix of eigenvalues,', / 20X, 01787 $ ' V is U represented by Householder vectors, and', / 20X, 01788 $ ' Y is a matrix of eigenvectors of S.)', 01789 $ / ' SSYTRD, UPLO=''U'':', / ' 1= | A - V S V', A1, 01790 $ ' | / ( |A| n ulp ) ', ' 2= | I - U V', A1, 01791 $ ' | / ( n ulp )', / ' SSYTRD, UPLO=''L'':', 01792 $ / ' 3= | A - V S V', A1, ' | / ( |A| n ulp ) ', 01793 $ ' 4= | I - U V', A1, ' | / ( n ulp )' ) 01794 9992 FORMAT( ' SSPTRD, UPLO=''U'':', / ' 5= | A - V S V', A1, 01795 $ ' | / ( |A| n ulp ) ', ' 6= | I - U V', A1, 01796 $ ' | / ( n ulp )', / ' SSPTRD, UPLO=''L'':', 01797 $ / ' 7= | A - V S V', A1, ' | / ( |A| n ulp ) ', 01798 $ ' 8= | I - U V', A1, ' | / ( n ulp )', 01799 $ / ' 9= | S - Z D Z', A1, ' | / ( |S| n ulp ) ', 01800 $ ' 10= | I - Z Z', A1, ' | / ( n ulp )', 01801 $ / ' 11= |D(with Z) - D(w/o Z)| / (|D| ulp) ', 01802 $ ' 12= | D(PWK) - D(QR) | / (|D| ulp)', 01803 $ / ' 13= Sturm sequence test on W ' ) 01804 9991 FORMAT( ' 14= | S - Z4 D4 Z4', A1, ' | / (|S| n ulp)', 01805 $ / ' 15= | I - Z4 Z4', A1, ' | / (n ulp ) ', 01806 $ ' 16= | D4 - D5 | / ( 100 |D4| ulp ) ', 01807 $ / ' 17= max | D4(i) - WR(i) | / ( |D4(i)| (2n-1) ulp )', 01808 $ / ' 18= | WA1 - D3 | / ( |D3| ulp )', 01809 $ / ' 19= max | WA2(i) - WA3(ii) | / ( |D3| ulp )', 01810 $ / ' 20= | S - Y WA1 Y', A1, ' | / ( |S| n ulp )', 01811 $ / ' 21= | I - Y Y', A1, ' | / ( n ulp )' ) 01812 9990 FORMAT( ' N=', I5, ', seed=', 4( I4, ',' ), ' type ', I2, 01813 $ ', test(', I2, ')=', G10.3 ) 01814 9989 FORMAT( ' 22= | S - Z D Z', A1, '| / ( |S| n ulp ) for SSTEDC(I)', 01815 $ / ' 23= | I - Z Z', A1, '| / ( n ulp ) for SSTEDC(I)', 01816 $ / ' 24= | S - Z D Z', A1, '| / ( |S| n ulp ) for SSTEDC(V)', 01817 $ / ' 25= | I - Z Z', A1, '| / ( n ulp ) for SSTEDC(V)', 01818 $ / ' 26= | D1(SSTEDC(V)) - D2(SSTEDC(N)) | / ( |D1| ulp )' ) 01819 * 01820 9988 FORMAT( / 'Test performed: see SCHKST for details.', / ) 01821 * End of SCHKST 01822 * 01823 END
1.7.3 
