LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
schkst.f
Go to the documentation of this file.
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 simultaneously
367*> with Z. On exit, the eigenvalues in D1 correspond with the
368*> matrix in A.
369*> \endverbatim
370*>
371*> \param[out] D2
372*> \verbatim
373*> D2 is REAL array of
374*> dimension( max(NN) )
375*> The eigenvalues of A, as computed by SSTEQR if Z is not
376*> computed. On exit, the eigenvalues in D2 correspond with
377*> the matrix in A.
378*> \endverbatim
379*>
380*> \param[out] D3
381*> \verbatim
382*> D3 is REAL array of
383*> dimension( max(NN) )
384*> The eigenvalues of A, as computed by SSTERF. On exit, the
385*> eigenvalues in D3 correspond with the matrix in A.
386*> \endverbatim
387*>
388*> \param[out] D4
389*> \verbatim
390*> D4 is REAL array of
391*> dimension( max(NN) )
392*> The eigenvalues of A, as computed by SPTEQR(V).
393*> ZPTEQR factors S as Z4 D4 Z4*
394*> On exit, the eigenvalues in D4 correspond with the matrix in A.
395*> \endverbatim
396*>
397*> \param[out] D5
398*> \verbatim
399*> D5 is REAL array of
400*> dimension( max(NN) )
401*> The eigenvalues of A, as computed by SPTEQR(N)
402*> when Z is not computed. On exit, the
403*> eigenvalues in D4 correspond with the matrix in A.
404*> \endverbatim
405*>
406*> \param[out] WA1
407*> \verbatim
408*> WA1 is REAL array of
409*> dimension( max(NN) )
410*> All eigenvalues of A, computed to high
411*> absolute accuracy, with different range options.
412*> as computed by SSTEBZ.
413*> \endverbatim
414*>
415*> \param[out] WA2
416*> \verbatim
417*> WA2 is REAL array of
418*> dimension( max(NN) )
419*> Selected eigenvalues of A, computed to high
420*> absolute accuracy, with different range options.
421*> as computed by SSTEBZ.
422*> Choose random values for IL and IU, and ask for the
423*> IL-th through IU-th eigenvalues.
424*> \endverbatim
425*>
426*> \param[out] WA3
427*> \verbatim
428*> WA3 is REAL array of
429*> dimension( max(NN) )
430*> Selected eigenvalues of A, computed to high
431*> absolute accuracy, with different range options.
432*> as computed by SSTEBZ.
433*> Determine the values VL and VU of the IL-th and IU-th
434*> eigenvalues and ask for all eigenvalues in this range.
435*> \endverbatim
436*>
437*> \param[out] WR
438*> \verbatim
439*> WR is REAL array of
440*> dimension( max(NN) )
441*> All eigenvalues of A, computed to high
442*> absolute accuracy, with different options.
443*> as computed by SSTEBZ.
444*> \endverbatim
445*>
446*> \param[out] U
447*> \verbatim
448*> U is REAL array of
449*> dimension( LDU, max(NN) ).
450*> The orthogonal matrix computed by SSYTRD + SORGTR.
451*> \endverbatim
452*>
453*> \param[in] LDU
454*> \verbatim
455*> LDU is INTEGER
456*> The leading dimension of U, Z, and V. It must be at least 1
457*> and at least max( NN ).
458*> \endverbatim
459*>
460*> \param[out] V
461*> \verbatim
462*> V is REAL array of
463*> dimension( LDU, max(NN) ).
464*> The Housholder vectors computed by SSYTRD in reducing A to
465*> tridiagonal form. The vectors computed with UPLO='U' are
466*> in the upper triangle, and the vectors computed with UPLO='L'
467*> are in the lower triangle. (As described in SSYTRD, the
468*> sub- and superdiagonal are not set to 1, although the
469*> true Householder vector has a 1 in that position. The
470*> routines that use V, such as SORGTR, set those entries to
471*> 1 before using them, and then restore them later.)
472*> \endverbatim
473*>
474*> \param[out] VP
475*> \verbatim
476*> VP is REAL array of
477*> dimension( max(NN)*max(NN+1)/2 )
478*> The matrix V stored in packed format.
479*> \endverbatim
480*>
481*> \param[out] TAU
482*> \verbatim
483*> TAU is REAL array of
484*> dimension( max(NN) )
485*> The Householder factors computed by SSYTRD in reducing A
486*> to tridiagonal form.
487*> \endverbatim
488*>
489*> \param[out] Z
490*> \verbatim
491*> Z is REAL array of
492*> dimension( LDU, max(NN) ).
493*> The orthogonal matrix of eigenvectors computed by SSTEQR,
494*> SPTEQR, and SSTEIN.
495*> \endverbatim
496*>
497*> \param[out] WORK
498*> \verbatim
499*> WORK is REAL array of
500*> dimension( LWORK )
501*> \endverbatim
502*>
503*> \param[in] LWORK
504*> \verbatim
505*> LWORK is INTEGER
506*> The number of entries in WORK. This must be at least
507*> 1 + 4 * Nmax + 2 * Nmax * lg Nmax + 3 * Nmax**2
508*> where Nmax = max( NN(j), 2 ) and lg = log base 2.
509*> \endverbatim
510*>
511*> \param[out] IWORK
512*> \verbatim
513*> IWORK is INTEGER array,
514*> Workspace.
515*> \endverbatim
516*>
517*> \param[out] LIWORK
518*> \verbatim
519*> LIWORK is INTEGER
520*> The number of entries in IWORK. This must be at least
521*> 6 + 6*Nmax + 5 * Nmax * lg Nmax
522*> where Nmax = max( NN(j), 2 ) and lg = log base 2.
523*> \endverbatim
524*>
525*> \param[out] RESULT
526*> \verbatim
527*> RESULT is REAL array, dimension (26)
528*> The values computed by the tests described above.
529*> The values are currently limited to 1/ulp, to avoid
530*> overflow.
531*> \endverbatim
532*>
533*> \param[out] INFO
534*> \verbatim
535*> INFO is INTEGER
536*> If 0, then everything ran OK.
537*> -1: NSIZES < 0
538*> -2: Some NN(j) < 0
539*> -3: NTYPES < 0
540*> -5: THRESH < 0
541*> -9: LDA < 1 or LDA < NMAX, where NMAX is max( NN(j) ).
542*> -23: LDU < 1 or LDU < NMAX.
543*> -29: LWORK too small.
544*> If SLATMR, SLATMS, SSYTRD, SORGTR, SSTEQR, SSTERF,
545*> or SORMC2 returns an error code, the
546*> absolute value of it is returned.
547*>
548*>-----------------------------------------------------------------------
549*>
550*> Some Local Variables and Parameters:
551*> ---- ----- --------- --- ----------
552*> ZERO, ONE Real 0 and 1.
553*> MAXTYP The number of types defined.
554*> NTEST The number of tests performed, or which can
555*> be performed so far, for the current matrix.
556*> NTESTT The total number of tests performed so far.
557*> NBLOCK Blocksize as returned by ENVIR.
558*> NMAX Largest value in NN.
559*> NMATS The number of matrices generated so far.
560*> NERRS The number of tests which have exceeded THRESH
561*> so far.
562*> COND, IMODE Values to be passed to the matrix generators.
563*> ANORM Norm of A; passed to matrix generators.
564*>
565*> OVFL, UNFL Overflow and underflow thresholds.
566*> ULP, ULPINV Finest relative precision and its inverse.
567*> RTOVFL, RTUNFL Square roots of the previous 2 values.
568*> The following four arrays decode JTYPE:
569*> KTYPE(j) The general type (1-10) for type "j".
570*> KMODE(j) The MODE value to be passed to the matrix
571*> generator for type "j".
572*> KMAGN(j) The order of magnitude ( O(1),
573*> O(overflow^(1/2) ), O(underflow^(1/2) )
574*> \endverbatim
575*
576* Authors:
577* ========
578*
579*> \author Univ. of Tennessee
580*> \author Univ. of California Berkeley
581*> \author Univ. of Colorado Denver
582*> \author NAG Ltd.
583*
584*> \ingroup single_eig
585*
586* =====================================================================
587 SUBROUTINE schkst( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
588 $ NOUNIT, A, LDA, AP, SD, SE, D1, D2, D3, D4, D5,
589 $ WA1, WA2, WA3, WR, U, LDU, V, VP, TAU, Z, WORK,
590 $ LWORK, IWORK, LIWORK, RESULT, INFO )
591*
592* -- LAPACK test routine --
593* -- LAPACK is a software package provided by Univ. of Tennessee, --
594* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
595*
596* .. Scalar Arguments ..
597 INTEGER INFO, LDA, LDU, LIWORK, LWORK, NOUNIT, NSIZES,
598 $ NTYPES
599 REAL THRESH
600* ..
601* .. Array Arguments ..
602 LOGICAL DOTYPE( * )
603 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
604 REAL A( LDA, * ), AP( * ), D1( * ), D2( * ),
605 $ d3( * ), d4( * ), d5( * ), result( * ),
606 $ sd( * ), se( * ), tau( * ), u( ldu, * ),
607 $ v( ldu, * ), vp( * ), wa1( * ), wa2( * ),
608 $ wa3( * ), work( * ), wr( * ), z( ldu, * )
609* ..
610*
611* =====================================================================
612*
613* .. Parameters ..
614 REAL ZERO, ONE, TWO, EIGHT, TEN, HUN
615 PARAMETER ( ZERO = 0.0e0, one = 1.0e0, two = 2.0e0,
616 $ eight = 8.0e0, ten = 10.0e0, hun = 100.0e0 )
617 REAL HALF
618 parameter( half = one / two )
619 INTEGER MAXTYP
620 parameter( maxtyp = 21 )
621 LOGICAL SRANGE
622 parameter( srange = .false. )
623 LOGICAL SREL
624 parameter( srel = .false. )
625* ..
626* .. Local Scalars ..
627 LOGICAL BADNN, TRYRAC
628 INTEGER I, IINFO, IL, IMODE, ITEMP, ITYPE, IU, J, JC,
629 $ JR, JSIZE, JTYPE, LGN, LIWEDC, LOG2UI, LWEDC,
630 $ m, m2, m3, mtypes, n, nap, nblock, nerrs,
631 $ nmats, nmax, nsplit, ntest, ntestt
632 REAL ABSTOL, ANINV, ANORM, COND, OVFL, RTOVFL,
633 $ RTUNFL, TEMP1, TEMP2, TEMP3, TEMP4, ULP,
634 $ ULPINV, UNFL, VL, VU
635* ..
636* .. Local Arrays ..
637 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ),
638 $ KMAGN( MAXTYP ), KMODE( MAXTYP ),
639 $ KTYPE( MAXTYP )
640 REAL DUMMA( 1 )
641* ..
642* .. External Functions ..
643 INTEGER ILAENV
644 REAL SLAMCH, SLARND, SSXT1
645 EXTERNAL ILAENV, SLAMCH, SLARND, SSXT1
646* ..
647* .. External Subroutines ..
648 EXTERNAL scopy, slacpy, slaset, slasum, slatmr, slatms,
649 $ sopgtr, sorgtr, spteqr, sspt21, ssptrd, sstebz,
650 $ sstech, sstedc, sstemr, sstein, ssteqr, ssterf,
651 $ sstt21, sstt22, ssyt21, ssytrd, xerbla
652* ..
653* .. Intrinsic Functions ..
654 INTRINSIC abs, int, log, max, min, real, sqrt
655* ..
656* .. Data statements ..
657 DATA ktype / 1, 2, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 8,
658 $ 8, 8, 9, 9, 9, 9, 9, 10 /
659 DATA kmagn / 1, 1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1,
660 $ 2, 3, 1, 1, 1, 2, 3, 1 /
661 DATA kmode / 0, 0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
662 $ 0, 0, 4, 3, 1, 4, 4, 3 /
663* ..
664* .. Executable Statements ..
665*
666* Keep ftnchek happy
667 idumma( 1 ) = 1
668*
669* Check for errors
670*
671 ntestt = 0
672 info = 0
673*
674* Important constants
675*
676 badnn = .false.
677 tryrac = .true.
678 nmax = 1
679 DO 10 j = 1, nsizes
680 nmax = max( nmax, nn( j ) )
681 IF( nn( j ).LT.0 )
682 $ badnn = .true.
683 10 CONTINUE
684*
685 nblock = ilaenv( 1, 'SSYTRD', 'L', nmax, -1, -1, -1 )
686 nblock = min( nmax, max( 1, nblock ) )
687*
688* Check for errors
689*
690 IF( nsizes.LT.0 ) THEN
691 info = -1
692 ELSE IF( badnn ) THEN
693 info = -2
694 ELSE IF( ntypes.LT.0 ) THEN
695 info = -3
696 ELSE IF( lda.LT.nmax ) THEN
697 info = -9
698 ELSE IF( ldu.LT.nmax ) THEN
699 info = -23
700 ELSE IF( 2*max( 2, nmax )**2.GT.lwork ) THEN
701 info = -29
702 END IF
703*
704 IF( info.NE.0 ) THEN
705 CALL xerbla( 'SCHKST', -info )
706 RETURN
707 END IF
708*
709* Quick return if possible
710*
711 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
712 $ RETURN
713*
714* More Important constants
715*
716 unfl = slamch( 'Safe minimum' )
717 ovfl = one / unfl
718 ulp = slamch( 'Epsilon' )*slamch( 'Base' )
719 ulpinv = one / ulp
720 log2ui = int( log( ulpinv ) / log( two ) )
721 rtunfl = sqrt( unfl )
722 rtovfl = sqrt( ovfl )
723*
724* Loop over sizes, types
725*
726 DO 20 i = 1, 4
727 iseed2( i ) = iseed( i )
728 20 CONTINUE
729 nerrs = 0
730 nmats = 0
731*
732 DO 310 jsize = 1, nsizes
733 n = nn( jsize )
734 IF( n.GT.0 ) THEN
735 lgn = int( log( real( n ) ) / log( two ) )
736 IF( 2**lgn.LT.n )
737 $ lgn = lgn + 1
738 IF( 2**lgn.LT.n )
739 $ lgn = lgn + 1
740 lwedc = 1 + 4*n + 2*n*lgn + 4*n**2
741 liwedc = 6 + 6*n + 5*n*lgn
742 ELSE
743 lwedc = 8
744 liwedc = 12
745 END IF
746 nap = ( n*( n+1 ) ) / 2
747 aninv = one / real( max( 1, n ) )
748*
749 IF( nsizes.NE.1 ) THEN
750 mtypes = min( maxtyp, ntypes )
751 ELSE
752 mtypes = min( maxtyp+1, ntypes )
753 END IF
754*
755 DO 300 jtype = 1, mtypes
756 IF( .NOT.dotype( jtype ) )
757 $ GO TO 300
758 nmats = nmats + 1
759 ntest = 0
760*
761 DO 30 j = 1, 4
762 ioldsd( j ) = iseed( j )
763 30 CONTINUE
764*
765* Compute "A"
766*
767* Control parameters:
768*
769* KMAGN KMODE KTYPE
770* =1 O(1) clustered 1 zero
771* =2 large clustered 2 identity
772* =3 small exponential (none)
773* =4 arithmetic diagonal, (w/ eigenvalues)
774* =5 random log symmetric, w/ eigenvalues
775* =6 random (none)
776* =7 random diagonal
777* =8 random symmetric
778* =9 positive definite
779* =10 diagonally dominant tridiagonal
780*
781 IF( mtypes.GT.maxtyp )
782 $ GO TO 100
783*
784 itype = ktype( jtype )
785 imode = kmode( jtype )
786*
787* Compute norm
788*
789 GO TO ( 40, 50, 60 )kmagn( jtype )
790*
791 40 CONTINUE
792 anorm = one
793 GO TO 70
794*
795 50 CONTINUE
796 anorm = ( rtovfl*ulp )*aninv
797 GO TO 70
798*
799 60 CONTINUE
800 anorm = rtunfl*n*ulpinv
801 GO TO 70
802*
803 70 CONTINUE
804*
805 CALL slaset( 'Full', lda, n, zero, zero, a, lda )
806 iinfo = 0
807 IF( jtype.LE.15 ) THEN
808 cond = ulpinv
809 ELSE
810 cond = ulpinv*aninv / ten
811 END IF
812*
813* Special Matrices -- Identity & Jordan block
814*
815* Zero
816*
817 IF( itype.EQ.1 ) THEN
818 iinfo = 0
819*
820 ELSE IF( itype.EQ.2 ) THEN
821*
822* Identity
823*
824 DO 80 jc = 1, n
825 a( jc, jc ) = anorm
826 80 CONTINUE
827*
828 ELSE IF( itype.EQ.4 ) THEN
829*
830* Diagonal Matrix, [Eigen]values Specified
831*
832 CALL slatms( n, n, 'S', iseed, 'S', work, imode, cond,
833 $ anorm, 0, 0, 'N', a, lda, work( n+1 ),
834 $ iinfo )
835*
836*
837 ELSE IF( itype.EQ.5 ) THEN
838*
839* Symmetric, eigenvalues specified
840*
841 CALL slatms( n, n, 'S', iseed, 'S', work, imode, cond,
842 $ anorm, n, n, 'N', a, lda, work( n+1 ),
843 $ iinfo )
844*
845 ELSE IF( itype.EQ.7 ) THEN
846*
847* Diagonal, random eigenvalues
848*
849 CALL slatmr( n, n, 'S', iseed, 'S', work, 6, one, one,
850 $ 'T', 'N', work( n+1 ), 1, one,
851 $ work( 2*n+1 ), 1, one, 'N', idumma, 0, 0,
852 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
853*
854 ELSE IF( itype.EQ.8 ) THEN
855*
856* Symmetric, random eigenvalues
857*
858 CALL slatmr( n, n, 'S', iseed, 'S', work, 6, one, one,
859 $ 'T', 'N', work( n+1 ), 1, one,
860 $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
861 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
862*
863 ELSE IF( itype.EQ.9 ) THEN
864*
865* Positive definite, eigenvalues specified.
866*
867 CALL slatms( n, n, 'S', iseed, 'P', work, imode, cond,
868 $ anorm, n, n, 'N', a, lda, work( n+1 ),
869 $ iinfo )
870*
871 ELSE IF( itype.EQ.10 ) THEN
872*
873* Positive definite tridiagonal, eigenvalues specified.
874*
875 CALL slatms( n, n, 'S', iseed, 'P', work, imode, cond,
876 $ anorm, 1, 1, 'N', a, lda, work( n+1 ),
877 $ iinfo )
878 DO 90 i = 2, n
879 temp1 = abs( a( i-1, i ) ) /
880 $ sqrt( abs( a( i-1, i-1 )*a( i, i ) ) )
881 IF( temp1.GT.half ) THEN
882 a( i-1, i ) = half*sqrt( abs( a( i-1, i-1 )*a( i,
883 $ i ) ) )
884 a( i, i-1 ) = a( i-1, i )
885 END IF
886 90 CONTINUE
887*
888 ELSE
889*
890 iinfo = 1
891 END IF
892*
893 IF( iinfo.NE.0 ) THEN
894 WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
895 $ ioldsd
896 info = abs( iinfo )
897 RETURN
898 END IF
899*
900 100 CONTINUE
901*
902* Call SSYTRD and SORGTR to compute S and U from
903* upper triangle.
904*
905 CALL slacpy( 'U', n, n, a, lda, v, ldu )
906*
907 ntest = 1
908 CALL ssytrd( 'U', n, v, ldu, sd, se, tau, work, lwork,
909 $ iinfo )
910*
911 IF( iinfo.NE.0 ) THEN
912 WRITE( nounit, fmt = 9999 )'SSYTRD(U)', iinfo, n, jtype,
913 $ ioldsd
914 info = abs( iinfo )
915 IF( iinfo.LT.0 ) THEN
916 RETURN
917 ELSE
918 result( 1 ) = ulpinv
919 GO TO 280
920 END IF
921 END IF
922*
923 CALL slacpy( 'U', n, n, v, ldu, u, ldu )
924*
925 ntest = 2
926 CALL sorgtr( 'U', n, u, ldu, tau, work, lwork, iinfo )
927 IF( iinfo.NE.0 ) THEN
928 WRITE( nounit, fmt = 9999 )'SORGTR(U)', iinfo, n, jtype,
929 $ ioldsd
930 info = abs( iinfo )
931 IF( iinfo.LT.0 ) THEN
932 RETURN
933 ELSE
934 result( 2 ) = ulpinv
935 GO TO 280
936 END IF
937 END IF
938*
939* Do tests 1 and 2
940*
941 CALL ssyt21( 2, 'Upper', n, 1, a, lda, sd, se, u, ldu, v,
942 $ ldu, tau, work, result( 1 ) )
943 CALL ssyt21( 3, 'Upper', n, 1, a, lda, sd, se, u, ldu, v,
944 $ ldu, tau, work, result( 2 ) )
945*
946* Call SSYTRD and SORGTR to compute S and U from
947* lower triangle, do tests.
948*
949 CALL slacpy( 'L', n, n, a, lda, v, ldu )
950*
951 ntest = 3
952 CALL ssytrd( 'L', n, v, ldu, sd, se, tau, work, lwork,
953 $ iinfo )
954*
955 IF( iinfo.NE.0 ) THEN
956 WRITE( nounit, fmt = 9999 )'SSYTRD(L)', iinfo, n, jtype,
957 $ ioldsd
958 info = abs( iinfo )
959 IF( iinfo.LT.0 ) THEN
960 RETURN
961 ELSE
962 result( 3 ) = ulpinv
963 GO TO 280
964 END IF
965 END IF
966*
967 CALL slacpy( 'L', n, n, v, ldu, u, ldu )
968*
969 ntest = 4
970 CALL sorgtr( 'L', n, u, ldu, tau, work, lwork, iinfo )
971 IF( iinfo.NE.0 ) THEN
972 WRITE( nounit, fmt = 9999 )'SORGTR(L)', iinfo, n, jtype,
973 $ ioldsd
974 info = abs( iinfo )
975 IF( iinfo.LT.0 ) THEN
976 RETURN
977 ELSE
978 result( 4 ) = ulpinv
979 GO TO 280
980 END IF
981 END IF
982*
983 CALL ssyt21( 2, 'Lower', n, 1, a, lda, sd, se, u, ldu, v,
984 $ ldu, tau, work, result( 3 ) )
985 CALL ssyt21( 3, 'Lower', n, 1, a, lda, sd, se, u, ldu, v,
986 $ ldu, tau, work, result( 4 ) )
987*
988* Store the upper triangle of A in AP
989*
990 i = 0
991 DO 120 jc = 1, n
992 DO 110 jr = 1, jc
993 i = i + 1
994 ap( i ) = a( jr, jc )
995 110 CONTINUE
996 120 CONTINUE
997*
998* Call SSPTRD and SOPGTR to compute S and U from AP
999*
1000 CALL scopy( nap, ap, 1, vp, 1 )
1001*
1002 ntest = 5
1003 CALL ssptrd( 'U', n, vp, sd, se, tau, iinfo )
1004*
1005 IF( iinfo.NE.0 ) THEN
1006 WRITE( nounit, fmt = 9999 )'SSPTRD(U)', iinfo, n, jtype,
1007 $ ioldsd
1008 info = abs( iinfo )
1009 IF( iinfo.LT.0 ) THEN
1010 RETURN
1011 ELSE
1012 result( 5 ) = ulpinv
1013 GO TO 280
1014 END IF
1015 END IF
1016*
1017 ntest = 6
1018 CALL sopgtr( 'U', n, vp, tau, u, ldu, work, iinfo )
1019 IF( iinfo.NE.0 ) THEN
1020 WRITE( nounit, fmt = 9999 )'SOPGTR(U)', iinfo, n, jtype,
1021 $ ioldsd
1022 info = abs( iinfo )
1023 IF( iinfo.LT.0 ) THEN
1024 RETURN
1025 ELSE
1026 result( 6 ) = ulpinv
1027 GO TO 280
1028 END IF
1029 END IF
1030*
1031* Do tests 5 and 6
1032*
1033 CALL sspt21( 2, 'Upper', n, 1, ap, sd, se, u, ldu, vp, tau,
1034 $ work, result( 5 ) )
1035 CALL sspt21( 3, 'Upper', n, 1, ap, sd, se, u, ldu, vp, tau,
1036 $ work, result( 6 ) )
1037*
1038* Store the lower triangle of A in AP
1039*
1040 i = 0
1041 DO 140 jc = 1, n
1042 DO 130 jr = jc, n
1043 i = i + 1
1044 ap( i ) = a( jr, jc )
1045 130 CONTINUE
1046 140 CONTINUE
1047*
1048* Call SSPTRD and SOPGTR to compute S and U from AP
1049*
1050 CALL scopy( nap, ap, 1, vp, 1 )
1051*
1052 ntest = 7
1053 CALL ssptrd( 'L', n, vp, sd, se, tau, iinfo )
1054*
1055 IF( iinfo.NE.0 ) THEN
1056 WRITE( nounit, fmt = 9999 )'SSPTRD(L)', iinfo, n, jtype,
1057 $ ioldsd
1058 info = abs( iinfo )
1059 IF( iinfo.LT.0 ) THEN
1060 RETURN
1061 ELSE
1062 result( 7 ) = ulpinv
1063 GO TO 280
1064 END IF
1065 END IF
1066*
1067 ntest = 8
1068 CALL sopgtr( 'L', n, vp, tau, u, ldu, work, iinfo )
1069 IF( iinfo.NE.0 ) THEN
1070 WRITE( nounit, fmt = 9999 )'SOPGTR(L)', iinfo, n, jtype,
1071 $ ioldsd
1072 info = abs( iinfo )
1073 IF( iinfo.LT.0 ) THEN
1074 RETURN
1075 ELSE
1076 result( 8 ) = ulpinv
1077 GO TO 280
1078 END IF
1079 END IF
1080*
1081 CALL sspt21( 2, 'Lower', n, 1, ap, sd, se, u, ldu, vp, tau,
1082 $ work, result( 7 ) )
1083 CALL sspt21( 3, 'Lower', n, 1, ap, sd, se, u, ldu, vp, tau,
1084 $ work, result( 8 ) )
1085*
1086* Call SSTEQR to compute D1, D2, and Z, do tests.
1087*
1088* Compute D1 and Z
1089*
1090 CALL scopy( n, sd, 1, d1, 1 )
1091 IF( n.GT.0 )
1092 $ CALL scopy( n-1, se, 1, work, 1 )
1093 CALL slaset( 'Full', n, n, zero, one, z, ldu )
1094*
1095 ntest = 9
1096 CALL ssteqr( 'V', n, d1, work, z, ldu, work( n+1 ), iinfo )
1097 IF( iinfo.NE.0 ) THEN
1098 WRITE( nounit, fmt = 9999 )'SSTEQR(V)', iinfo, n, jtype,
1099 $ ioldsd
1100 info = abs( iinfo )
1101 IF( iinfo.LT.0 ) THEN
1102 RETURN
1103 ELSE
1104 result( 9 ) = ulpinv
1105 GO TO 280
1106 END IF
1107 END IF
1108*
1109* Compute D2
1110*
1111 CALL scopy( n, sd, 1, d2, 1 )
1112 IF( n.GT.0 )
1113 $ CALL scopy( n-1, se, 1, work, 1 )
1114*
1115 ntest = 11
1116 CALL ssteqr( 'N', n, d2, work, work( n+1 ), ldu,
1117 $ work( n+1 ), iinfo )
1118 IF( iinfo.NE.0 ) THEN
1119 WRITE( nounit, fmt = 9999 )'SSTEQR(N)', iinfo, n, jtype,
1120 $ ioldsd
1121 info = abs( iinfo )
1122 IF( iinfo.LT.0 ) THEN
1123 RETURN
1124 ELSE
1125 result( 11 ) = ulpinv
1126 GO TO 280
1127 END IF
1128 END IF
1129*
1130* Compute D3 (using PWK method)
1131*
1132 CALL scopy( n, sd, 1, d3, 1 )
1133 IF( n.GT.0 )
1134 $ CALL scopy( n-1, se, 1, work, 1 )
1135*
1136 ntest = 12
1137 CALL ssterf( n, d3, work, iinfo )
1138 IF( iinfo.NE.0 ) THEN
1139 WRITE( nounit, fmt = 9999 )'SSTERF', iinfo, n, jtype,
1140 $ ioldsd
1141 info = abs( iinfo )
1142 IF( iinfo.LT.0 ) THEN
1143 RETURN
1144 ELSE
1145 result( 12 ) = ulpinv
1146 GO TO 280
1147 END IF
1148 END IF
1149*
1150* Do Tests 9 and 10
1151*
1152 CALL sstt21( n, 0, sd, se, d1, dumma, z, ldu, work,
1153 $ result( 9 ) )
1154*
1155* Do Tests 11 and 12
1156*
1157 temp1 = zero
1158 temp2 = zero
1159 temp3 = zero
1160 temp4 = zero
1161*
1162 DO 150 j = 1, n
1163 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1164 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1165 temp3 = max( temp3, abs( d1( j ) ), abs( d3( j ) ) )
1166 temp4 = max( temp4, abs( d1( j )-d3( j ) ) )
1167 150 CONTINUE
1168*
1169 result( 11 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1170 result( 12 ) = temp4 / max( unfl, ulp*max( temp3, temp4 ) )
1171*
1172* Do Test 13 -- Sturm Sequence Test of Eigenvalues
1173* Go up by factors of two until it succeeds
1174*
1175 ntest = 13
1176 temp1 = thresh*( half-ulp )
1177*
1178 DO 160 j = 0, log2ui
1179 CALL sstech( n, sd, se, d1, temp1, work, iinfo )
1180 IF( iinfo.EQ.0 )
1181 $ GO TO 170
1182 temp1 = temp1*two
1183 160 CONTINUE
1184*
1185 170 CONTINUE
1186 result( 13 ) = temp1
1187*
1188* For positive definite matrices ( JTYPE.GT.15 ) call SPTEQR
1189* and do tests 14, 15, and 16 .
1190*
1191 IF( jtype.GT.15 ) THEN
1192*
1193* Compute D4 and Z4
1194*
1195 CALL scopy( n, sd, 1, d4, 1 )
1196 IF( n.GT.0 )
1197 $ CALL scopy( n-1, se, 1, work, 1 )
1198 CALL slaset( 'Full', n, n, zero, one, z, ldu )
1199*
1200 ntest = 14
1201 CALL spteqr( 'V', n, d4, work, z, ldu, work( n+1 ),
1202 $ iinfo )
1203 IF( iinfo.NE.0 ) THEN
1204 WRITE( nounit, fmt = 9999 )'SPTEQR(V)', iinfo, n,
1205 $ jtype, ioldsd
1206 info = abs( iinfo )
1207 IF( iinfo.LT.0 ) THEN
1208 RETURN
1209 ELSE
1210 result( 14 ) = ulpinv
1211 GO TO 280
1212 END IF
1213 END IF
1214*
1215* Do Tests 14 and 15
1216*
1217 CALL sstt21( n, 0, sd, se, d4, dumma, z, ldu, work,
1218 $ result( 14 ) )
1219*
1220* Compute D5
1221*
1222 CALL scopy( n, sd, 1, d5, 1 )
1223 IF( n.GT.0 )
1224 $ CALL scopy( n-1, se, 1, work, 1 )
1225*
1226 ntest = 16
1227 CALL spteqr( 'N', n, d5, work, z, ldu, work( n+1 ),
1228 $ iinfo )
1229 IF( iinfo.NE.0 ) THEN
1230 WRITE( nounit, fmt = 9999 )'SPTEQR(N)', iinfo, n,
1231 $ jtype, ioldsd
1232 info = abs( iinfo )
1233 IF( iinfo.LT.0 ) THEN
1234 RETURN
1235 ELSE
1236 result( 16 ) = ulpinv
1237 GO TO 280
1238 END IF
1239 END IF
1240*
1241* Do Test 16
1242*
1243 temp1 = zero
1244 temp2 = zero
1245 DO 180 j = 1, n
1246 temp1 = max( temp1, abs( d4( j ) ), abs( d5( j ) ) )
1247 temp2 = max( temp2, abs( d4( j )-d5( j ) ) )
1248 180 CONTINUE
1249*
1250 result( 16 ) = temp2 / max( unfl,
1251 $ hun*ulp*max( temp1, temp2 ) )
1252 ELSE
1253 result( 14 ) = zero
1254 result( 15 ) = zero
1255 result( 16 ) = zero
1256 END IF
1257*
1258* Call SSTEBZ with different options and do tests 17-18.
1259*
1260* If S is positive definite and diagonally dominant,
1261* ask for all eigenvalues with high relative accuracy.
1262*
1263 vl = zero
1264 vu = zero
1265 il = 0
1266 iu = 0
1267 IF( jtype.EQ.21 ) THEN
1268 ntest = 17
1269 abstol = unfl + unfl
1270 CALL sstebz( 'A', 'E', n, vl, vu, il, iu, abstol, sd, se,
1271 $ m, nsplit, wr, iwork( 1 ), iwork( n+1 ),
1272 $ work, iwork( 2*n+1 ), iinfo )
1273 IF( iinfo.NE.0 ) THEN
1274 WRITE( nounit, fmt = 9999 )'SSTEBZ(A,rel)', iinfo, n,
1275 $ jtype, ioldsd
1276 info = abs( iinfo )
1277 IF( iinfo.LT.0 ) THEN
1278 RETURN
1279 ELSE
1280 result( 17 ) = ulpinv
1281 GO TO 280
1282 END IF
1283 END IF
1284*
1285* Do test 17
1286*
1287 temp2 = two*( two*n-one )*ulp*( one+eight*half**2 ) /
1288 $ ( one-half )**4
1289*
1290 temp1 = zero
1291 DO 190 j = 1, n
1292 temp1 = max( temp1, abs( d4( j )-wr( n-j+1 ) ) /
1293 $ ( abstol+abs( d4( j ) ) ) )
1294 190 CONTINUE
1295*
1296 result( 17 ) = temp1 / temp2
1297 ELSE
1298 result( 17 ) = zero
1299 END IF
1300*
1301* Now ask for all eigenvalues with high absolute accuracy.
1302*
1303 ntest = 18
1304 abstol = unfl + unfl
1305 CALL sstebz( 'A', 'E', n, vl, vu, il, iu, abstol, sd, se, m,
1306 $ nsplit, wa1, iwork( 1 ), iwork( n+1 ), work,
1307 $ iwork( 2*n+1 ), iinfo )
1308 IF( iinfo.NE.0 ) THEN
1309 WRITE( nounit, fmt = 9999 )'SSTEBZ(A)', iinfo, n, jtype,
1310 $ ioldsd
1311 info = abs( iinfo )
1312 IF( iinfo.LT.0 ) THEN
1313 RETURN
1314 ELSE
1315 result( 18 ) = ulpinv
1316 GO TO 280
1317 END IF
1318 END IF
1319*
1320* Do test 18
1321*
1322 temp1 = zero
1323 temp2 = zero
1324 DO 200 j = 1, n
1325 temp1 = max( temp1, abs( d3( j ) ), abs( wa1( j ) ) )
1326 temp2 = max( temp2, abs( d3( j )-wa1( j ) ) )
1327 200 CONTINUE
1328*
1329 result( 18 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1330*
1331* Choose random values for IL and IU, and ask for the
1332* IL-th through IU-th eigenvalues.
1333*
1334 ntest = 19
1335 IF( n.LE.1 ) THEN
1336 il = 1
1337 iu = n
1338 ELSE
1339 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1340 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1341 IF( iu.LT.il ) THEN
1342 itemp = iu
1343 iu = il
1344 il = itemp
1345 END IF
1346 END IF
1347*
1348 CALL sstebz( 'I', 'E', n, vl, vu, il, iu, abstol, sd, se,
1349 $ m2, nsplit, wa2, iwork( 1 ), iwork( n+1 ),
1350 $ work, iwork( 2*n+1 ), iinfo )
1351 IF( iinfo.NE.0 ) THEN
1352 WRITE( nounit, fmt = 9999 )'SSTEBZ(I)', iinfo, n, jtype,
1353 $ ioldsd
1354 info = abs( iinfo )
1355 IF( iinfo.LT.0 ) THEN
1356 RETURN
1357 ELSE
1358 result( 19 ) = ulpinv
1359 GO TO 280
1360 END IF
1361 END IF
1362*
1363* Determine the values VL and VU of the IL-th and IU-th
1364* eigenvalues and ask for all eigenvalues in this range.
1365*
1366 IF( n.GT.0 ) THEN
1367 IF( il.NE.1 ) THEN
1368 vl = wa1( il ) - max( half*( wa1( il )-wa1( il-1 ) ),
1369 $ ulp*anorm, two*rtunfl )
1370 ELSE
1371 vl = wa1( 1 ) - max( half*( wa1( n )-wa1( 1 ) ),
1372 $ ulp*anorm, two*rtunfl )
1373 END IF
1374 IF( iu.NE.n ) THEN
1375 vu = wa1( iu ) + max( half*( wa1( iu+1 )-wa1( iu ) ),
1376 $ ulp*anorm, two*rtunfl )
1377 ELSE
1378 vu = wa1( n ) + max( half*( wa1( n )-wa1( 1 ) ),
1379 $ ulp*anorm, two*rtunfl )
1380 END IF
1381 ELSE
1382 vl = zero
1383 vu = one
1384 END IF
1385*
1386 CALL sstebz( 'V', 'E', n, vl, vu, il, iu, abstol, sd, se,
1387 $ m3, nsplit, wa3, iwork( 1 ), iwork( n+1 ),
1388 $ work, iwork( 2*n+1 ), iinfo )
1389 IF( iinfo.NE.0 ) THEN
1390 WRITE( nounit, fmt = 9999 )'SSTEBZ(V)', iinfo, n, jtype,
1391 $ ioldsd
1392 info = abs( iinfo )
1393 IF( iinfo.LT.0 ) THEN
1394 RETURN
1395 ELSE
1396 result( 19 ) = ulpinv
1397 GO TO 280
1398 END IF
1399 END IF
1400*
1401 IF( m3.EQ.0 .AND. n.NE.0 ) THEN
1402 result( 19 ) = ulpinv
1403 GO TO 280
1404 END IF
1405*
1406* Do test 19
1407*
1408 temp1 = ssxt1( 1, wa2, m2, wa3, m3, abstol, ulp, unfl )
1409 temp2 = ssxt1( 1, wa3, m3, wa2, m2, abstol, ulp, unfl )
1410 IF( n.GT.0 ) THEN
1411 temp3 = max( abs( wa1( n ) ), abs( wa1( 1 ) ) )
1412 ELSE
1413 temp3 = zero
1414 END IF
1415*
1416 result( 19 ) = ( temp1+temp2 ) / max( unfl, temp3*ulp )
1417*
1418* Call SSTEIN to compute eigenvectors corresponding to
1419* eigenvalues in WA1. (First call SSTEBZ again, to make sure
1420* it returns these eigenvalues in the correct order.)
1421*
1422 ntest = 21
1423 CALL sstebz( 'A', 'B', n, vl, vu, il, iu, abstol, sd, se, m,
1424 $ nsplit, wa1, iwork( 1 ), iwork( n+1 ), work,
1425 $ iwork( 2*n+1 ), iinfo )
1426 IF( iinfo.NE.0 ) THEN
1427 WRITE( nounit, fmt = 9999 )'SSTEBZ(A,B)', iinfo, n,
1428 $ jtype, ioldsd
1429 info = abs( iinfo )
1430 IF( iinfo.LT.0 ) THEN
1431 RETURN
1432 ELSE
1433 result( 20 ) = ulpinv
1434 result( 21 ) = ulpinv
1435 GO TO 280
1436 END IF
1437 END IF
1438*
1439 CALL sstein( n, sd, se, m, wa1, iwork( 1 ), iwork( n+1 ), z,
1440 $ ldu, work, iwork( 2*n+1 ), iwork( 3*n+1 ),
1441 $ iinfo )
1442 IF( iinfo.NE.0 ) THEN
1443 WRITE( nounit, fmt = 9999 )'SSTEIN', iinfo, n, jtype,
1444 $ ioldsd
1445 info = abs( iinfo )
1446 IF( iinfo.LT.0 ) THEN
1447 RETURN
1448 ELSE
1449 result( 20 ) = ulpinv
1450 result( 21 ) = ulpinv
1451 GO TO 280
1452 END IF
1453 END IF
1454*
1455* Do tests 20 and 21
1456*
1457 CALL sstt21( n, 0, sd, se, wa1, dumma, z, ldu, work,
1458 $ result( 20 ) )
1459*
1460* Call SSTEDC(I) to compute D1 and Z, do tests.
1461*
1462* Compute D1 and Z
1463*
1464 CALL scopy( n, sd, 1, d1, 1 )
1465 IF( n.GT.0 )
1466 $ CALL scopy( n-1, se, 1, work, 1 )
1467 CALL slaset( 'Full', n, n, zero, one, z, ldu )
1468*
1469 ntest = 22
1470 CALL sstedc( 'I', n, d1, work, z, ldu, work( n+1 ), lwedc-n,
1471 $ iwork, liwedc, iinfo )
1472 IF( iinfo.NE.0 ) THEN
1473 WRITE( nounit, fmt = 9999 )'SSTEDC(I)', iinfo, n, jtype,
1474 $ ioldsd
1475 info = abs( iinfo )
1476 IF( iinfo.LT.0 ) THEN
1477 RETURN
1478 ELSE
1479 result( 22 ) = ulpinv
1480 GO TO 280
1481 END IF
1482 END IF
1483*
1484* Do Tests 22 and 23
1485*
1486 CALL sstt21( n, 0, sd, se, d1, dumma, z, ldu, work,
1487 $ result( 22 ) )
1488*
1489* Call SSTEDC(V) to compute D1 and Z, do tests.
1490*
1491* Compute D1 and Z
1492*
1493 CALL scopy( n, sd, 1, d1, 1 )
1494 IF( n.GT.0 )
1495 $ CALL scopy( n-1, se, 1, work, 1 )
1496 CALL slaset( 'Full', n, n, zero, one, z, ldu )
1497*
1498 ntest = 24
1499 CALL sstedc( 'V', n, d1, work, z, ldu, work( n+1 ), lwedc-n,
1500 $ iwork, liwedc, iinfo )
1501 IF( iinfo.NE.0 ) THEN
1502 WRITE( nounit, fmt = 9999 )'SSTEDC(V)', iinfo, n, jtype,
1503 $ ioldsd
1504 info = abs( iinfo )
1505 IF( iinfo.LT.0 ) THEN
1506 RETURN
1507 ELSE
1508 result( 24 ) = ulpinv
1509 GO TO 280
1510 END IF
1511 END IF
1512*
1513* Do Tests 24 and 25
1514*
1515 CALL sstt21( n, 0, sd, se, d1, dumma, z, ldu, work,
1516 $ result( 24 ) )
1517*
1518* Call SSTEDC(N) to compute D2, do tests.
1519*
1520* Compute D2
1521*
1522 CALL scopy( n, sd, 1, d2, 1 )
1523 IF( n.GT.0 )
1524 $ CALL scopy( n-1, se, 1, work, 1 )
1525 CALL slaset( 'Full', n, n, zero, one, z, ldu )
1526*
1527 ntest = 26
1528 CALL sstedc( 'N', n, d2, work, z, ldu, work( n+1 ), lwedc-n,
1529 $ iwork, liwedc, iinfo )
1530 IF( iinfo.NE.0 ) THEN
1531 WRITE( nounit, fmt = 9999 )'SSTEDC(N)', iinfo, n, jtype,
1532 $ ioldsd
1533 info = abs( iinfo )
1534 IF( iinfo.LT.0 ) THEN
1535 RETURN
1536 ELSE
1537 result( 26 ) = ulpinv
1538 GO TO 280
1539 END IF
1540 END IF
1541*
1542* Do Test 26
1543*
1544 temp1 = zero
1545 temp2 = zero
1546*
1547 DO 210 j = 1, n
1548 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1549 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1550 210 CONTINUE
1551*
1552 result( 26 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1553*
1554* Only test SSTEMR if IEEE compliant
1555*
1556 IF( ilaenv( 10, 'SSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 .AND.
1557 $ ilaenv( 11, 'SSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 ) THEN
1558*
1559* Call SSTEMR, do test 27 (relative eigenvalue accuracy)
1560*
1561* If S is positive definite and diagonally dominant,
1562* ask for all eigenvalues with high relative accuracy.
1563*
1564 vl = zero
1565 vu = zero
1566 il = 0
1567 iu = 0
1568 IF( jtype.EQ.21 .AND. srel ) THEN
1569 ntest = 27
1570 abstol = unfl + unfl
1571 CALL sstemr( 'V', 'A', n, sd, se, vl, vu, il, iu,
1572 $ m, wr, z, ldu, n, iwork( 1 ), tryrac,
1573 $ work, lwork, iwork( 2*n+1 ), lwork-2*n,
1574 $ iinfo )
1575 IF( iinfo.NE.0 ) THEN
1576 WRITE( nounit, fmt = 9999 )'SSTEMR(V,A,rel)',
1577 $ iinfo, n, jtype, ioldsd
1578 info = abs( iinfo )
1579 IF( iinfo.LT.0 ) THEN
1580 RETURN
1581 ELSE
1582 result( 27 ) = ulpinv
1583 GO TO 270
1584 END IF
1585 END IF
1586*
1587* Do test 27
1588*
1589 temp2 = two*( two*n-one )*ulp*( one+eight*half**2 ) /
1590 $ ( one-half )**4
1591*
1592 temp1 = zero
1593 DO 220 j = 1, n
1594 temp1 = max( temp1, abs( d4( j )-wr( n-j+1 ) ) /
1595 $ ( abstol+abs( d4( j ) ) ) )
1596 220 CONTINUE
1597*
1598 result( 27 ) = temp1 / temp2
1599*
1600 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1601 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1602 IF( iu.LT.il ) THEN
1603 itemp = iu
1604 iu = il
1605 il = itemp
1606 END IF
1607*
1608 IF( srange ) THEN
1609 ntest = 28
1610 abstol = unfl + unfl
1611 CALL sstemr( 'V', 'I', n, sd, se, vl, vu, il, iu,
1612 $ m, wr, z, ldu, n, iwork( 1 ), tryrac,
1613 $ work, lwork, iwork( 2*n+1 ),
1614 $ lwork-2*n, iinfo )
1615*
1616 IF( iinfo.NE.0 ) THEN
1617 WRITE( nounit, fmt = 9999 )'SSTEMR(V,I,rel)',
1618 $ iinfo, n, jtype, ioldsd
1619 info = abs( iinfo )
1620 IF( iinfo.LT.0 ) THEN
1621 RETURN
1622 ELSE
1623 result( 28 ) = ulpinv
1624 GO TO 270
1625 END IF
1626 END IF
1627*
1628*
1629* Do test 28
1630*
1631 temp2 = two*( two*n-one )*ulp*
1632 $ ( one+eight*half**2 ) / ( one-half )**4
1633*
1634 temp1 = zero
1635 DO 230 j = il, iu
1636 temp1 = max( temp1, abs( wr( j-il+1 )-d4( n-j+
1637 $ 1 ) ) / ( abstol+abs( wr( j-il+1 ) ) ) )
1638 230 CONTINUE
1639*
1640 result( 28 ) = temp1 / temp2
1641 ELSE
1642 result( 28 ) = zero
1643 END IF
1644 ELSE
1645 result( 27 ) = zero
1646 result( 28 ) = zero
1647 END IF
1648*
1649* Call SSTEMR(V,I) to compute D1 and Z, do tests.
1650*
1651* Compute D1 and Z
1652*
1653 CALL scopy( n, sd, 1, d5, 1 )
1654 IF( n.GT.0 )
1655 $ CALL scopy( n-1, se, 1, work, 1 )
1656 CALL slaset( 'Full', n, n, zero, one, z, ldu )
1657*
1658 IF( srange ) THEN
1659 ntest = 29
1660 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1661 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1662 IF( iu.LT.il ) THEN
1663 itemp = iu
1664 iu = il
1665 il = itemp
1666 END IF
1667 CALL sstemr( 'V', 'I', n, d5, work, vl, vu, il, iu,
1668 $ m, d1, z, ldu, n, iwork( 1 ), tryrac,
1669 $ work( n+1 ), lwork-n, iwork( 2*n+1 ),
1670 $ liwork-2*n, iinfo )
1671 IF( iinfo.NE.0 ) THEN
1672 WRITE( nounit, fmt = 9999 )'SSTEMR(V,I)', iinfo,
1673 $ n, jtype, ioldsd
1674 info = abs( iinfo )
1675 IF( iinfo.LT.0 ) THEN
1676 RETURN
1677 ELSE
1678 result( 29 ) = ulpinv
1679 GO TO 280
1680 END IF
1681 END IF
1682*
1683* Do Tests 29 and 30
1684*
1685 CALL sstt22( n, m, 0, sd, se, d1, dumma, z, ldu, work,
1686 $ m, result( 29 ) )
1687*
1688* Call SSTEMR to compute D2, do tests.
1689*
1690* Compute D2
1691*
1692 CALL scopy( n, sd, 1, d5, 1 )
1693 IF( n.GT.0 )
1694 $ CALL scopy( n-1, se, 1, work, 1 )
1695*
1696 ntest = 31
1697 CALL sstemr( 'N', 'I', n, d5, work, vl, vu, il, iu,
1698 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1699 $ work( n+1 ), lwork-n, iwork( 2*n+1 ),
1700 $ liwork-2*n, iinfo )
1701 IF( iinfo.NE.0 ) THEN
1702 WRITE( nounit, fmt = 9999 )'SSTEMR(N,I)', iinfo,
1703 $ n, jtype, ioldsd
1704 info = abs( iinfo )
1705 IF( iinfo.LT.0 ) THEN
1706 RETURN
1707 ELSE
1708 result( 31 ) = ulpinv
1709 GO TO 280
1710 END IF
1711 END IF
1712*
1713* Do Test 31
1714*
1715 temp1 = zero
1716 temp2 = zero
1717*
1718 DO 240 j = 1, iu - il + 1
1719 temp1 = max( temp1, abs( d1( j ) ),
1720 $ abs( d2( j ) ) )
1721 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1722 240 CONTINUE
1723*
1724 result( 31 ) = temp2 / max( unfl,
1725 $ ulp*max( temp1, temp2 ) )
1726*
1727*
1728* Call SSTEMR(V,V) to compute D1 and Z, do tests.
1729*
1730* Compute D1 and Z
1731*
1732 CALL scopy( n, sd, 1, d5, 1 )
1733 IF( n.GT.0 )
1734 $ CALL scopy( n-1, se, 1, work, 1 )
1735 CALL slaset( 'Full', n, n, zero, one, z, ldu )
1736*
1737 ntest = 32
1738*
1739 IF( n.GT.0 ) THEN
1740 IF( il.NE.1 ) THEN
1741 vl = d2( il ) - max( half*
1742 $ ( d2( il )-d2( il-1 ) ), ulp*anorm,
1743 $ two*rtunfl )
1744 ELSE
1745 vl = d2( 1 ) - max( half*( d2( n )-d2( 1 ) ),
1746 $ ulp*anorm, two*rtunfl )
1747 END IF
1748 IF( iu.NE.n ) THEN
1749 vu = d2( iu ) + max( half*
1750 $ ( d2( iu+1 )-d2( iu ) ), ulp*anorm,
1751 $ two*rtunfl )
1752 ELSE
1753 vu = d2( n ) + max( half*( d2( n )-d2( 1 ) ),
1754 $ ulp*anorm, two*rtunfl )
1755 END IF
1756 ELSE
1757 vl = zero
1758 vu = one
1759 END IF
1760*
1761 CALL sstemr( 'V', 'V', n, d5, work, vl, vu, il, iu,
1762 $ m, d1, z, ldu, n, iwork( 1 ), tryrac,
1763 $ work( n+1 ), lwork-n, iwork( 2*n+1 ),
1764 $ liwork-2*n, iinfo )
1765 IF( iinfo.NE.0 ) THEN
1766 WRITE( nounit, fmt = 9999 )'SSTEMR(V,V)', iinfo,
1767 $ n, jtype, ioldsd
1768 info = abs( iinfo )
1769 IF( iinfo.LT.0 ) THEN
1770 RETURN
1771 ELSE
1772 result( 32 ) = ulpinv
1773 GO TO 280
1774 END IF
1775 END IF
1776*
1777* Do Tests 32 and 33
1778*
1779 CALL sstt22( n, m, 0, sd, se, d1, dumma, z, ldu, work,
1780 $ m, result( 32 ) )
1781*
1782* Call SSTEMR to compute D2, do tests.
1783*
1784* Compute D2
1785*
1786 CALL scopy( n, sd, 1, d5, 1 )
1787 IF( n.GT.0 )
1788 $ CALL scopy( n-1, se, 1, work, 1 )
1789*
1790 ntest = 34
1791 CALL sstemr( 'N', 'V', n, d5, work, vl, vu, il, iu,
1792 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1793 $ work( n+1 ), lwork-n, iwork( 2*n+1 ),
1794 $ liwork-2*n, iinfo )
1795 IF( iinfo.NE.0 ) THEN
1796 WRITE( nounit, fmt = 9999 )'SSTEMR(N,V)', iinfo,
1797 $ n, jtype, ioldsd
1798 info = abs( iinfo )
1799 IF( iinfo.LT.0 ) THEN
1800 RETURN
1801 ELSE
1802 result( 34 ) = ulpinv
1803 GO TO 280
1804 END IF
1805 END IF
1806*
1807* Do Test 34
1808*
1809 temp1 = zero
1810 temp2 = zero
1811*
1812 DO 250 j = 1, iu - il + 1
1813 temp1 = max( temp1, abs( d1( j ) ),
1814 $ abs( d2( j ) ) )
1815 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1816 250 CONTINUE
1817*
1818 result( 34 ) = temp2 / max( unfl,
1819 $ ulp*max( temp1, temp2 ) )
1820 ELSE
1821 result( 29 ) = zero
1822 result( 30 ) = zero
1823 result( 31 ) = zero
1824 result( 32 ) = zero
1825 result( 33 ) = zero
1826 result( 34 ) = zero
1827 END IF
1828*
1829*
1830* Call SSTEMR(V,A) to compute D1 and Z, do tests.
1831*
1832* Compute D1 and Z
1833*
1834 CALL scopy( n, sd, 1, d5, 1 )
1835 IF( n.GT.0 )
1836 $ CALL scopy( n-1, se, 1, work, 1 )
1837*
1838 ntest = 35
1839*
1840 CALL sstemr( 'V', 'A', n, d5, work, vl, vu, il, iu,
1841 $ m, d1, z, ldu, n, iwork( 1 ), tryrac,
1842 $ work( n+1 ), lwork-n, iwork( 2*n+1 ),
1843 $ liwork-2*n, iinfo )
1844 IF( iinfo.NE.0 ) THEN
1845 WRITE( nounit, fmt = 9999 )'SSTEMR(V,A)', iinfo, n,
1846 $ jtype, ioldsd
1847 info = abs( iinfo )
1848 IF( iinfo.LT.0 ) THEN
1849 RETURN
1850 ELSE
1851 result( 35 ) = ulpinv
1852 GO TO 280
1853 END IF
1854 END IF
1855*
1856* Do Tests 35 and 36
1857*
1858 CALL sstt22( n, m, 0, sd, se, d1, dumma, z, ldu, work, m,
1859 $ result( 35 ) )
1860*
1861* Call SSTEMR to compute D2, do tests.
1862*
1863* Compute D2
1864*
1865 CALL scopy( n, sd, 1, d5, 1 )
1866 IF( n.GT.0 )
1867 $ CALL scopy( n-1, se, 1, work, 1 )
1868*
1869 ntest = 37
1870 CALL sstemr( 'N', 'A', n, d5, work, vl, vu, il, iu,
1871 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1872 $ work( n+1 ), lwork-n, iwork( 2*n+1 ),
1873 $ liwork-2*n, iinfo )
1874 IF( iinfo.NE.0 ) THEN
1875 WRITE( nounit, fmt = 9999 )'SSTEMR(N,A)', iinfo, n,
1876 $ jtype, ioldsd
1877 info = abs( iinfo )
1878 IF( iinfo.LT.0 ) THEN
1879 RETURN
1880 ELSE
1881 result( 37 ) = ulpinv
1882 GO TO 280
1883 END IF
1884 END IF
1885*
1886* Do Test 34
1887*
1888 temp1 = zero
1889 temp2 = zero
1890*
1891 DO 260 j = 1, n
1892 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1893 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1894 260 CONTINUE
1895*
1896 result( 37 ) = temp2 / max( unfl,
1897 $ ulp*max( temp1, temp2 ) )
1898 END IF
1899 270 CONTINUE
1900 280 CONTINUE
1901 ntestt = ntestt + ntest
1902*
1903* End of Loop -- Check for RESULT(j) > THRESH
1904*
1905*
1906* Print out tests which fail.
1907*
1908 DO 290 jr = 1, ntest
1909 IF( result( jr ).GE.thresh ) THEN
1910*
1911* If this is the first test to fail,
1912* print a header to the data file.
1913*
1914 IF( nerrs.EQ.0 ) THEN
1915 WRITE( nounit, fmt = 9998 )'SST'
1916 WRITE( nounit, fmt = 9997 )
1917 WRITE( nounit, fmt = 9996 )
1918 WRITE( nounit, fmt = 9995 )'Symmetric'
1919 WRITE( nounit, fmt = 9994 )
1920*
1921* Tests performed
1922*
1923 WRITE( nounit, fmt = 9988 )
1924 END IF
1925 nerrs = nerrs + 1
1926 WRITE( nounit, fmt = 9990 )n, ioldsd, jtype, jr,
1927 $ result( jr )
1928 END IF
1929 290 CONTINUE
1930 300 CONTINUE
1931 310 CONTINUE
1932*
1933* Summary
1934*
1935 CALL slasum( 'SST', nounit, nerrs, ntestt )
1936 RETURN
1937*
1938 9999 FORMAT( ' SCHKST: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1939 $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1940*
1941 9998 FORMAT( / 1x, a3, ' -- Real Symmetric eigenvalue problem' )
1942 9997 FORMAT( ' Matrix types (see SCHKST for details): ' )
1943*
1944 9996 FORMAT( / ' Special Matrices:',
1945 $ / ' 1=Zero matrix. ',
1946 $ ' 5=Diagonal: clustered entries.',
1947 $ / ' 2=Identity matrix. ',
1948 $ ' 6=Diagonal: large, evenly spaced.',
1949 $ / ' 3=Diagonal: evenly spaced entries. ',
1950 $ ' 7=Diagonal: small, evenly spaced.',
1951 $ / ' 4=Diagonal: geometr. spaced entries.' )
1952 9995 FORMAT( ' Dense ', a, ' Matrices:',
1953 $ / ' 8=Evenly spaced eigenvals. ',
1954 $ ' 12=Small, evenly spaced eigenvals.',
1955 $ / ' 9=Geometrically spaced eigenvals. ',
1956 $ ' 13=Matrix with random O(1) entries.',
1957 $ / ' 10=Clustered eigenvalues. ',
1958 $ ' 14=Matrix with large random entries.',
1959 $ / ' 11=Large, evenly spaced eigenvals. ',
1960 $ ' 15=Matrix with small random entries.' )
1961 9994 FORMAT( ' 16=Positive definite, evenly spaced eigenvalues',
1962 $ / ' 17=Positive definite, geometrically spaced eigenvlaues',
1963 $ / ' 18=Positive definite, clustered eigenvalues',
1964 $ / ' 19=Positive definite, small evenly spaced eigenvalues',
1965 $ / ' 20=Positive definite, large evenly spaced eigenvalues',
1966 $ / ' 21=Diagonally dominant tridiagonal, geometrically',
1967 $ ' spaced eigenvalues' )
1968*
1969 9990 FORMAT( ' N=', i5, ', seed=', 4( i4, ',' ), ' type ', i2,
1970 $ ', test(', i2, ')=', g10.3 )
1971*
1972 9988 FORMAT( / 'Test performed: see SCHKST for details.', / )
1973* End of SCHKST
1974*
1975 END
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
scopy
subroutine scopy(n, sx, incx, sy, incy)
SCOPY
Definition scopy.f:82
ssytrd
subroutine ssytrd(uplo, n, a, lda, d, e, tau, work, lwork, info)
SSYTRD
Definition ssytrd.f:191
ssptrd
subroutine ssptrd(uplo, n, ap, d, e, tau, info)
SSPTRD
Definition ssptrd.f:148
slacpy
subroutine slacpy(uplo, m, n, a, lda, b, ldb)
SLACPY copies all or part of one two-dimensional array to another.
Definition slacpy.f:101
slaset
subroutine slaset(uplo, m, n, alpha, beta, a, lda)
SLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition slaset.f:108
spteqr
subroutine spteqr(compz, n, d, e, z, ldz, work, info)
SPTEQR
Definition spteqr.f:143
sstebz
subroutine sstebz(range, order, n, vl, vu, il, iu, abstol, d, e, m, nsplit, w, iblock, isplit, work, iwork, info)
SSTEBZ
Definition sstebz.f:272
sstedc
subroutine sstedc(compz, n, d, e, z, ldz, work, lwork, iwork, liwork, info)
SSTEDC
Definition sstedc.f:180
sstein
subroutine sstein(n, d, e, m, w, iblock, isplit, z, ldz, work, iwork, ifail, info)
SSTEIN
Definition sstein.f:172
sstemr
subroutine sstemr(jobz, range, n, d, e, vl, vu, il, iu, m, w, z, ldz, nzc, isuppz, tryrac, work, lwork, iwork, liwork, info)
SSTEMR
Definition sstemr.f:320
ssteqr
subroutine ssteqr(compz, n, d, e, z, ldz, work, info)
SSTEQR
Definition ssteqr.f:129
ssterf
subroutine ssterf(n, d, e, info)
SSTERF
Definition ssterf.f:84
sorgtr
subroutine sorgtr(uplo, n, a, lda, tau, work, lwork, info)
SORGTR
Definition sorgtr.f:121
sopgtr
subroutine sopgtr(uplo, n, ap, tau, q, ldq, work, info)
SOPGTR
Definition sopgtr.f:112
schkst
subroutine schkst(nsizes, nn, ntypes, dotype, iseed, thresh, nounit, a, lda, ap, sd, se, d1, d2, d3, d4, d5, wa1, wa2, wa3, wr, u, ldu, v, vp, tau, z, work, lwork, iwork, liwork, result, info)
SCHKST
Definition schkst.f:591
slasum
subroutine slasum(type, iounit, ie, nrun)
SLASUM
Definition slasum.f:41
slatmr
subroutine slatmr(m, n, dist, iseed, sym, d, mode, cond, dmax, rsign, grade, dl, model, condl, dr, moder, condr, pivtng, ipivot, kl, ku, sparse, anorm, pack, a, lda, iwork, info)
SLATMR
Definition slatmr.f:471
slatms
subroutine slatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
SLATMS
Definition slatms.f:321
sspt21
subroutine sspt21(itype, uplo, n, kband, ap, d, e, u, ldu, vp, tau, work, result)
SSPT21
Definition sspt21.f:221
sstech
subroutine sstech(n, a, b, eig, tol, work, info)
SSTECH
Definition sstech.f:101
sstt21
subroutine sstt21(n, kband, ad, ae, sd, se, u, ldu, work, result)
SSTT21
Definition sstt21.f:127
sstt22
subroutine sstt22(n, m, kband, ad, ae, sd, se, u, ldu, work, ldwork, result)
SSTT22
Definition sstt22.f:139
ssyt21
subroutine ssyt21(itype, uplo, n, kband, a, lda, d, e, u, ldu, v, ldv, tau, work, result)
SSYT21
Definition ssyt21.f:207