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