LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
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 simultaneously
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,
669* ..
670* .. Intrinsic Functions ..
671 INTRINSIC abs, conjg, int, log, max, min, real, sqrt
672* ..
673* .. Data statements ..
674 DATA ktype / 1, 2, 4, 4, 4, 4, 4, 5, 5, 5, 5, 5, 8,
675 $ 8, 8, 9, 9, 9, 9, 9, 10 /
676 DATA kmagn / 1, 1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1,
677 $ 2, 3, 1, 1, 1, 2, 3, 1 /
678 DATA kmode / 0, 0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
679 $ 0, 0, 4, 3, 1, 4, 4, 3 /
680* ..
681* .. Executable Statements ..
682*
683* Keep ftnchek happy
684 idumma( 1 ) = 1
685*
686* Check for errors
687*
688 ntestt = 0
689 info = 0
690*
691* Important constants
692*
693 badnn = .false.
694 tryrac = .true.
695 nmax = 1
696 DO 10 j = 1, nsizes
697 nmax = max( nmax, nn( j ) )
698 IF( nn( j ).LT.0 )
699 $ badnn = .true.
700 10 CONTINUE
701*
702 nblock = ilaenv( 1, 'CHETRD', 'L', nmax, -1, -1, -1 )
703 nblock = min( nmax, max( 1, nblock ) )
704*
705* Check for errors
706*
707 IF( nsizes.LT.0 ) THEN
708 info = -1
709 ELSE IF( badnn ) THEN
710 info = -2
711 ELSE IF( ntypes.LT.0 ) THEN
712 info = -3
713 ELSE IF( lda.LT.nmax ) THEN
714 info = -9
715 ELSE IF( ldu.LT.nmax ) THEN
716 info = -23
717 ELSE IF( 2*max( 2, nmax )**2.GT.lwork ) THEN
718 info = -29
719 END IF
720*
721 IF( info.NE.0 ) THEN
722 CALL xerbla( 'CCHKST', -info )
723 RETURN
724 END IF
725*
726* Quick return if possible
727*
728 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
729 $ RETURN
730*
731* More Important constants
732*
733 unfl = slamch( 'Safe minimum' )
734 ovfl = one / unfl
735 ulp = slamch( 'Epsilon' )*slamch( 'Base' )
736 ulpinv = one / ulp
737 log2ui = int( log( ulpinv ) / log( two ) )
738 rtunfl = sqrt( unfl )
739 rtovfl = sqrt( ovfl )
740*
741* Loop over sizes, types
742*
743 DO 20 i = 1, 4
744 iseed2( i ) = iseed( i )
745 20 CONTINUE
746 nerrs = 0
747 nmats = 0
748*
749 DO 310 jsize = 1, nsizes
750 n = nn( jsize )
751 IF( n.GT.0 ) THEN
752 lgn = int( log( real( n ) ) / log( two ) )
753 IF( 2**lgn.LT.n )
754 $ lgn = lgn + 1
755 IF( 2**lgn.LT.n )
756 $ lgn = lgn + 1
757 lwedc = 1 + 4*n + 2*n*lgn + 4*n**2
758 lrwedc = 1 + 3*n + 2*n*lgn + 4*n**2
759 liwedc = 6 + 6*n + 5*n*lgn
760 ELSE
761 lwedc = 8
762 lrwedc = 7
763 liwedc = 12
764 END IF
765 nap = ( n*( n+1 ) ) / 2
766 aninv = one / real( max( 1, n ) )
767*
768 IF( nsizes.NE.1 ) THEN
769 mtypes = min( maxtyp, ntypes )
770 ELSE
771 mtypes = min( maxtyp+1, ntypes )
772 END IF
773*
774 DO 300 jtype = 1, mtypes
775 IF( .NOT.dotype( jtype ) )
776 $ GO TO 300
777 nmats = nmats + 1
778 ntest = 0
779*
780 DO 30 j = 1, 4
781 ioldsd( j ) = iseed( j )
782 30 CONTINUE
783*
784* Compute "A"
785*
786* Control parameters:
787*
788* KMAGN KMODE KTYPE
789* =1 O(1) clustered 1 zero
790* =2 large clustered 2 identity
791* =3 small exponential (none)
792* =4 arithmetic diagonal, (w/ eigenvalues)
793* =5 random log Hermitian, w/ eigenvalues
794* =6 random (none)
795* =7 random diagonal
796* =8 random Hermitian
797* =9 positive definite
798* =10 diagonally dominant tridiagonal
799*
800 IF( mtypes.GT.maxtyp )
801 $ GO TO 100
802*
803 itype = ktype( jtype )
804 imode = kmode( jtype )
805*
806* Compute norm
807*
808 GO TO ( 40, 50, 60 )kmagn( jtype )
809*
810 40 CONTINUE
811 anorm = one
812 GO TO 70
813*
814 50 CONTINUE
815 anorm = ( rtovfl*ulp )*aninv
816 GO TO 70
817*
818 60 CONTINUE
819 anorm = rtunfl*n*ulpinv
820 GO TO 70
821*
822 70 CONTINUE
823*
824 CALL claset( 'Full', lda, n, czero, czero, a, lda )
825 iinfo = 0
826 IF( jtype.LE.15 ) THEN
827 cond = ulpinv
828 ELSE
829 cond = ulpinv*aninv / ten
830 END IF
831*
832* Special Matrices -- Identity & Jordan block
833*
834* Zero
835*
836 IF( itype.EQ.1 ) THEN
837 iinfo = 0
838*
839 ELSE IF( itype.EQ.2 ) THEN
840*
841* Identity
842*
843 DO 80 jc = 1, n
844 a( jc, jc ) = anorm
845 80 CONTINUE
846*
847 ELSE IF( itype.EQ.4 ) THEN
848*
849* Diagonal Matrix, [Eigen]values Specified
850*
851 CALL clatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
852 $ anorm, 0, 0, 'N', a, lda, work, iinfo )
853*
854*
855 ELSE IF( itype.EQ.5 ) THEN
856*
857* Hermitian, eigenvalues specified
858*
859 CALL clatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
860 $ anorm, n, n, 'N', a, lda, work, iinfo )
861*
862 ELSE IF( itype.EQ.7 ) THEN
863*
864* Diagonal, random eigenvalues
865*
866 CALL clatmr( n, n, 'S', iseed, 'H', work, 6, one, cone,
867 $ 'T', 'N', work( n+1 ), 1, one,
868 $ work( 2*n+1 ), 1, one, 'N', idumma, 0, 0,
869 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
870*
871 ELSE IF( itype.EQ.8 ) THEN
872*
873* Hermitian, random eigenvalues
874*
875 CALL clatmr( n, n, 'S', iseed, 'H', work, 6, one, cone,
876 $ 'T', 'N', work( n+1 ), 1, one,
877 $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
878 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
879*
880 ELSE IF( itype.EQ.9 ) THEN
881*
882* Positive definite, eigenvalues specified.
883*
884 CALL clatms( n, n, 'S', iseed, 'P', rwork, imode, cond,
885 $ anorm, n, n, 'N', a, lda, work, iinfo )
886*
887 ELSE IF( itype.EQ.10 ) THEN
888*
889* Positive definite tridiagonal, eigenvalues specified.
890*
891 CALL clatms( n, n, 'S', iseed, 'P', rwork, imode, cond,
892 $ anorm, 1, 1, 'N', a, lda, work, iinfo )
893 DO 90 i = 2, n
894 temp1 = abs( a( i-1, i ) )
895 temp2 = sqrt( abs( a( i-1, i-1 )*a( i, i ) ) )
896 IF( temp1.GT.half*temp2 ) THEN
897 a( i-1, i ) = a( i-1, i )*
898 $ ( half*temp2 / ( unfl+temp1 ) )
899 a( i, i-1 ) = conjg( a( i-1, i ) )
900 END IF
901 90 CONTINUE
902*
903 ELSE
904*
905 iinfo = 1
906 END IF
907*
908 IF( iinfo.NE.0 ) THEN
909 WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
910 $ ioldsd
911 info = abs( iinfo )
912 RETURN
913 END IF
914*
915 100 CONTINUE
916*
917* Call CHETRD and CUNGTR to compute S and U from
918* upper triangle.
919*
920 CALL clacpy( 'U', n, n, a, lda, v, ldu )
921*
922 ntest = 1
923 CALL chetrd( 'U', n, v, ldu, sd, se, tau, work, lwork,
924 $ iinfo )
925*
926 IF( iinfo.NE.0 ) THEN
927 WRITE( nounit, fmt = 9999 )'CHETRD(U)', iinfo, n, jtype,
928 $ ioldsd
929 info = abs( iinfo )
930 IF( iinfo.LT.0 ) THEN
931 RETURN
932 ELSE
933 result( 1 ) = ulpinv
934 GO TO 280
935 END IF
936 END IF
937*
938 CALL clacpy( 'U', n, n, v, ldu, u, ldu )
939*
940 ntest = 2
941 CALL cungtr( 'U', n, u, ldu, tau, work, lwork, iinfo )
942 IF( iinfo.NE.0 ) THEN
943 WRITE( nounit, fmt = 9999 )'CUNGTR(U)', iinfo, n, jtype,
944 $ ioldsd
945 info = abs( iinfo )
946 IF( iinfo.LT.0 ) THEN
947 RETURN
948 ELSE
949 result( 2 ) = ulpinv
950 GO TO 280
951 END IF
952 END IF
953*
954* Do tests 1 and 2
955*
956 CALL chet21( 2, 'Upper', n, 1, a, lda, sd, se, u, ldu, v,
957 $ ldu, tau, work, rwork, result( 1 ) )
958 CALL chet21( 3, 'Upper', n, 1, a, lda, sd, se, u, ldu, v,
959 $ ldu, tau, work, rwork, result( 2 ) )
960*
961* Call CHETRD and CUNGTR to compute S and U from
962* lower triangle, do tests.
963*
964 CALL clacpy( 'L', n, n, a, lda, v, ldu )
965*
966 ntest = 3
967 CALL chetrd( 'L', n, v, ldu, sd, se, tau, work, lwork,
968 $ iinfo )
969*
970 IF( iinfo.NE.0 ) THEN
971 WRITE( nounit, fmt = 9999 )'CHETRD(L)', iinfo, n, jtype,
972 $ ioldsd
973 info = abs( iinfo )
974 IF( iinfo.LT.0 ) THEN
975 RETURN
976 ELSE
977 result( 3 ) = ulpinv
978 GO TO 280
979 END IF
980 END IF
981*
982 CALL clacpy( 'L', n, n, v, ldu, u, ldu )
983*
984 ntest = 4
985 CALL cungtr( 'L', n, u, ldu, tau, work, lwork, iinfo )
986 IF( iinfo.NE.0 ) THEN
987 WRITE( nounit, fmt = 9999 )'CUNGTR(L)', iinfo, n, jtype,
988 $ ioldsd
989 info = abs( iinfo )
990 IF( iinfo.LT.0 ) THEN
991 RETURN
992 ELSE
993 result( 4 ) = ulpinv
994 GO TO 280
995 END IF
996 END IF
997*
998 CALL chet21( 2, 'Lower', n, 1, a, lda, sd, se, u, ldu, v,
999 $ ldu, tau, work, rwork, result( 3 ) )
1000 CALL chet21( 3, 'Lower', n, 1, a, lda, sd, se, u, ldu, v,
1001 $ ldu, tau, work, rwork, result( 4 ) )
1002*
1003* Store the upper triangle of A in AP
1004*
1005 i = 0
1006 DO 120 jc = 1, n
1007 DO 110 jr = 1, jc
1008 i = i + 1
1009 ap( i ) = a( jr, jc )
1010 110 CONTINUE
1011 120 CONTINUE
1012*
1013* Call CHPTRD and CUPGTR to compute S and U from AP
1014*
1015 CALL ccopy( nap, ap, 1, vp, 1 )
1016*
1017 ntest = 5
1018 CALL chptrd( 'U', n, vp, sd, se, tau, iinfo )
1019*
1020 IF( iinfo.NE.0 ) THEN
1021 WRITE( nounit, fmt = 9999 )'CHPTRD(U)', iinfo, n, jtype,
1022 $ ioldsd
1023 info = abs( iinfo )
1024 IF( iinfo.LT.0 ) THEN
1025 RETURN
1026 ELSE
1027 result( 5 ) = ulpinv
1028 GO TO 280
1029 END IF
1030 END IF
1031*
1032 ntest = 6
1033 CALL cupgtr( 'U', n, vp, tau, u, ldu, work, iinfo )
1034 IF( iinfo.NE.0 ) THEN
1035 WRITE( nounit, fmt = 9999 )'CUPGTR(U)', iinfo, n, jtype,
1036 $ ioldsd
1037 info = abs( iinfo )
1038 IF( iinfo.LT.0 ) THEN
1039 RETURN
1040 ELSE
1041 result( 6 ) = ulpinv
1042 GO TO 280
1043 END IF
1044 END IF
1045*
1046* Do tests 5 and 6
1047*
1048 CALL chpt21( 2, 'Upper', n, 1, ap, sd, se, u, ldu, vp, tau,
1049 $ work, rwork, result( 5 ) )
1050 CALL chpt21( 3, 'Upper', n, 1, ap, sd, se, u, ldu, vp, tau,
1051 $ work, rwork, result( 6 ) )
1052*
1053* Store the lower triangle of A in AP
1054*
1055 i = 0
1056 DO 140 jc = 1, n
1057 DO 130 jr = jc, n
1058 i = i + 1
1059 ap( i ) = a( jr, jc )
1060 130 CONTINUE
1061 140 CONTINUE
1062*
1063* Call CHPTRD and CUPGTR to compute S and U from AP
1064*
1065 CALL ccopy( nap, ap, 1, vp, 1 )
1066*
1067 ntest = 7
1068 CALL chptrd( 'L', n, vp, sd, se, tau, iinfo )
1069*
1070 IF( iinfo.NE.0 ) THEN
1071 WRITE( nounit, fmt = 9999 )'CHPTRD(L)', iinfo, n, jtype,
1072 $ ioldsd
1073 info = abs( iinfo )
1074 IF( iinfo.LT.0 ) THEN
1075 RETURN
1076 ELSE
1077 result( 7 ) = ulpinv
1078 GO TO 280
1079 END IF
1080 END IF
1081*
1082 ntest = 8
1083 CALL cupgtr( 'L', n, vp, tau, u, ldu, work, iinfo )
1084 IF( iinfo.NE.0 ) THEN
1085 WRITE( nounit, fmt = 9999 )'CUPGTR(L)', iinfo, n, jtype,
1086 $ ioldsd
1087 info = abs( iinfo )
1088 IF( iinfo.LT.0 ) THEN
1089 RETURN
1090 ELSE
1091 result( 8 ) = ulpinv
1092 GO TO 280
1093 END IF
1094 END IF
1095*
1096 CALL chpt21( 2, 'Lower', n, 1, ap, sd, se, u, ldu, vp, tau,
1097 $ work, rwork, result( 7 ) )
1098 CALL chpt21( 3, 'Lower', n, 1, ap, sd, se, u, ldu, vp, tau,
1099 $ work, rwork, result( 8 ) )
1100*
1101* Call CSTEQR to compute D1, D2, and Z, do tests.
1102*
1103* Compute D1 and Z
1104*
1105 CALL scopy( n, sd, 1, d1, 1 )
1106 IF( n.GT.0 )
1107 $ CALL scopy( n-1, se, 1, rwork, 1 )
1108 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1109*
1110 ntest = 9
1111 CALL csteqr( 'V', n, d1, rwork, z, ldu, rwork( n+1 ),
1112 $ iinfo )
1113 IF( iinfo.NE.0 ) THEN
1114 WRITE( nounit, fmt = 9999 )'CSTEQR(V)', iinfo, n, jtype,
1115 $ ioldsd
1116 info = abs( iinfo )
1117 IF( iinfo.LT.0 ) THEN
1118 RETURN
1119 ELSE
1120 result( 9 ) = ulpinv
1121 GO TO 280
1122 END IF
1123 END IF
1124*
1125* Compute D2
1126*
1127 CALL scopy( n, sd, 1, d2, 1 )
1128 IF( n.GT.0 )
1129 $ CALL scopy( n-1, se, 1, rwork, 1 )
1130*
1131 ntest = 11
1132 CALL csteqr( 'N', n, d2, rwork, work, ldu, rwork( n+1 ),
1133 $ iinfo )
1134 IF( iinfo.NE.0 ) THEN
1135 WRITE( nounit, fmt = 9999 )'CSTEQR(N)', iinfo, n, jtype,
1136 $ ioldsd
1137 info = abs( iinfo )
1138 IF( iinfo.LT.0 ) THEN
1139 RETURN
1140 ELSE
1141 result( 11 ) = ulpinv
1142 GO TO 280
1143 END IF
1144 END IF
1145*
1146* Compute D3 (using PWK method)
1147*
1148 CALL scopy( n, sd, 1, d3, 1 )
1149 IF( n.GT.0 )
1150 $ CALL scopy( n-1, se, 1, rwork, 1 )
1151*
1152 ntest = 12
1153 CALL ssterf( n, d3, rwork, iinfo )
1154 IF( iinfo.NE.0 ) THEN
1155 WRITE( nounit, fmt = 9999 )'SSTERF', iinfo, n, jtype,
1156 $ ioldsd
1157 info = abs( iinfo )
1158 IF( iinfo.LT.0 ) THEN
1159 RETURN
1160 ELSE
1161 result( 12 ) = ulpinv
1162 GO TO 280
1163 END IF
1164 END IF
1165*
1166* Do Tests 9 and 10
1167*
1168 CALL cstt21( n, 0, sd, se, d1, dumma, z, ldu, work, rwork,
1169 $ result( 9 ) )
1170*
1171* Do Tests 11 and 12
1172*
1173 temp1 = zero
1174 temp2 = zero
1175 temp3 = zero
1176 temp4 = zero
1177*
1178 DO 150 j = 1, n
1179 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1180 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1181 temp3 = max( temp3, abs( d1( j ) ), abs( d3( j ) ) )
1182 temp4 = max( temp4, abs( d1( j )-d3( j ) ) )
1183 150 CONTINUE
1184*
1185 result( 11 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1186 result( 12 ) = temp4 / max( unfl, ulp*max( temp3, temp4 ) )
1187*
1188* Do Test 13 -- Sturm Sequence Test of Eigenvalues
1189* Go up by factors of two until it succeeds
1190*
1191 ntest = 13
1192 temp1 = thresh*( half-ulp )
1193*
1194 DO 160 j = 0, log2ui
1195 CALL sstech( n, sd, se, d1, temp1, rwork, iinfo )
1196 IF( iinfo.EQ.0 )
1197 $ GO TO 170
1198 temp1 = temp1*two
1199 160 CONTINUE
1200*
1201 170 CONTINUE
1202 result( 13 ) = temp1
1203*
1204* For positive definite matrices ( JTYPE.GT.15 ) call CPTEQR
1205* and do tests 14, 15, and 16 .
1206*
1207 IF( jtype.GT.15 ) THEN
1208*
1209* Compute D4 and Z4
1210*
1211 CALL scopy( n, sd, 1, d4, 1 )
1212 IF( n.GT.0 )
1213 $ CALL scopy( n-1, se, 1, rwork, 1 )
1214 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1215*
1216 ntest = 14
1217 CALL cpteqr( 'V', n, d4, rwork, z, ldu, rwork( n+1 ),
1218 $ iinfo )
1219 IF( iinfo.NE.0 ) THEN
1220 WRITE( nounit, fmt = 9999 )'CPTEQR(V)', iinfo, n,
1221 $ jtype, ioldsd
1222 info = abs( iinfo )
1223 IF( iinfo.LT.0 ) THEN
1224 RETURN
1225 ELSE
1226 result( 14 ) = ulpinv
1227 GO TO 280
1228 END IF
1229 END IF
1230*
1231* Do Tests 14 and 15
1232*
1233 CALL cstt21( n, 0, sd, se, d4, dumma, z, ldu, work,
1234 $ rwork, result( 14 ) )
1235*
1236* Compute D5
1237*
1238 CALL scopy( n, sd, 1, d5, 1 )
1239 IF( n.GT.0 )
1240 $ CALL scopy( n-1, se, 1, rwork, 1 )
1241*
1242 ntest = 16
1243 CALL cpteqr( 'N', n, d5, rwork, z, ldu, rwork( n+1 ),
1244 $ iinfo )
1245 IF( iinfo.NE.0 ) THEN
1246 WRITE( nounit, fmt = 9999 )'CPTEQR(N)', iinfo, n,
1247 $ jtype, ioldsd
1248 info = abs( iinfo )
1249 IF( iinfo.LT.0 ) THEN
1250 RETURN
1251 ELSE
1252 result( 16 ) = ulpinv
1253 GO TO 280
1254 END IF
1255 END IF
1256*
1257* Do Test 16
1258*
1259 temp1 = zero
1260 temp2 = zero
1261 DO 180 j = 1, n
1262 temp1 = max( temp1, abs( d4( j ) ), abs( d5( j ) ) )
1263 temp2 = max( temp2, abs( d4( j )-d5( j ) ) )
1264 180 CONTINUE
1265*
1266 result( 16 ) = temp2 / max( unfl,
1267 $ hun*ulp*max( temp1, temp2 ) )
1268 ELSE
1269 result( 14 ) = zero
1270 result( 15 ) = zero
1271 result( 16 ) = zero
1272 END IF
1273*
1274* Call SSTEBZ with different options and do tests 17-18.
1275*
1276* If S is positive definite and diagonally dominant,
1277* ask for all eigenvalues with high relative accuracy.
1278*
1279 vl = zero
1280 vu = zero
1281 il = 0
1282 iu = 0
1283 IF( jtype.EQ.21 ) THEN
1284 ntest = 17
1285 abstol = unfl + unfl
1286 CALL sstebz( 'A', 'E', n, vl, vu, il, iu, abstol, sd, se,
1287 $ m, nsplit, wr, iwork( 1 ), iwork( n+1 ),
1288 $ rwork, iwork( 2*n+1 ), iinfo )
1289 IF( iinfo.NE.0 ) THEN
1290 WRITE( nounit, fmt = 9999 )'SSTEBZ(A,rel)', iinfo, n,
1291 $ jtype, ioldsd
1292 info = abs( iinfo )
1293 IF( iinfo.LT.0 ) THEN
1294 RETURN
1295 ELSE
1296 result( 17 ) = ulpinv
1297 GO TO 280
1298 END IF
1299 END IF
1300*
1301* Do test 17
1302*
1303 temp2 = two*( two*n-one )*ulp*( one+eight*half**2 ) /
1304 $ ( one-half )**4
1305*
1306 temp1 = zero
1307 DO 190 j = 1, n
1308 temp1 = max( temp1, abs( d4( j )-wr( n-j+1 ) ) /
1309 $ ( abstol+abs( d4( j ) ) ) )
1310 190 CONTINUE
1311*
1312 result( 17 ) = temp1 / temp2
1313 ELSE
1314 result( 17 ) = zero
1315 END IF
1316*
1317* Now ask for all eigenvalues with high absolute accuracy.
1318*
1319 ntest = 18
1320 abstol = unfl + unfl
1321 CALL sstebz( 'A', 'E', n, vl, vu, il, iu, abstol, sd, se, m,
1322 $ nsplit, wa1, iwork( 1 ), iwork( n+1 ), rwork,
1323 $ iwork( 2*n+1 ), iinfo )
1324 IF( iinfo.NE.0 ) THEN
1325 WRITE( nounit, fmt = 9999 )'SSTEBZ(A)', iinfo, n, jtype,
1326 $ ioldsd
1327 info = abs( iinfo )
1328 IF( iinfo.LT.0 ) THEN
1329 RETURN
1330 ELSE
1331 result( 18 ) = ulpinv
1332 GO TO 280
1333 END IF
1334 END IF
1335*
1336* Do test 18
1337*
1338 temp1 = zero
1339 temp2 = zero
1340 DO 200 j = 1, n
1341 temp1 = max( temp1, abs( d3( j ) ), abs( wa1( j ) ) )
1342 temp2 = max( temp2, abs( d3( j )-wa1( j ) ) )
1343 200 CONTINUE
1344*
1345 result( 18 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1346*
1347* Choose random values for IL and IU, and ask for the
1348* IL-th through IU-th eigenvalues.
1349*
1350 ntest = 19
1351 IF( n.LE.1 ) THEN
1352 il = 1
1353 iu = n
1354 ELSE
1355 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1356 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1357 IF( iu.LT.il ) THEN
1358 itemp = iu
1359 iu = il
1360 il = itemp
1361 END IF
1362 END IF
1363*
1364 CALL sstebz( 'I', 'E', n, vl, vu, il, iu, abstol, sd, se,
1365 $ m2, nsplit, wa2, iwork( 1 ), iwork( n+1 ),
1366 $ rwork, iwork( 2*n+1 ), iinfo )
1367 IF( iinfo.NE.0 ) THEN
1368 WRITE( nounit, fmt = 9999 )'SSTEBZ(I)', iinfo, n, jtype,
1369 $ ioldsd
1370 info = abs( iinfo )
1371 IF( iinfo.LT.0 ) THEN
1372 RETURN
1373 ELSE
1374 result( 19 ) = ulpinv
1375 GO TO 280
1376 END IF
1377 END IF
1378*
1379* Determine the values VL and VU of the IL-th and IU-th
1380* eigenvalues and ask for all eigenvalues in this range.
1381*
1382 IF( n.GT.0 ) THEN
1383 IF( il.NE.1 ) THEN
1384 vl = wa1( il ) - max( half*( wa1( il )-wa1( il-1 ) ),
1385 $ ulp*anorm, two*rtunfl )
1386 ELSE
1387 vl = wa1( 1 ) - max( half*( wa1( n )-wa1( 1 ) ),
1388 $ ulp*anorm, two*rtunfl )
1389 END IF
1390 IF( iu.NE.n ) THEN
1391 vu = wa1( iu ) + max( half*( wa1( iu+1 )-wa1( iu ) ),
1392 $ ulp*anorm, two*rtunfl )
1393 ELSE
1394 vu = wa1( n ) + max( half*( wa1( n )-wa1( 1 ) ),
1395 $ ulp*anorm, two*rtunfl )
1396 END IF
1397 ELSE
1398 vl = zero
1399 vu = one
1400 END IF
1401*
1402 CALL sstebz( 'V', 'E', n, vl, vu, il, iu, abstol, sd, se,
1403 $ m3, nsplit, wa3, iwork( 1 ), iwork( n+1 ),
1404 $ rwork, iwork( 2*n+1 ), iinfo )
1405 IF( iinfo.NE.0 ) THEN
1406 WRITE( nounit, fmt = 9999 )'SSTEBZ(V)', iinfo, n, jtype,
1407 $ ioldsd
1408 info = abs( iinfo )
1409 IF( iinfo.LT.0 ) THEN
1410 RETURN
1411 ELSE
1412 result( 19 ) = ulpinv
1413 GO TO 280
1414 END IF
1415 END IF
1416*
1417 IF( m3.EQ.0 .AND. n.NE.0 ) THEN
1418 result( 19 ) = ulpinv
1419 GO TO 280
1420 END IF
1421*
1422* Do test 19
1423*
1424 temp1 = ssxt1( 1, wa2, m2, wa3, m3, abstol, ulp, unfl )
1425 temp2 = ssxt1( 1, wa3, m3, wa2, m2, abstol, ulp, unfl )
1426 IF( n.GT.0 ) THEN
1427 temp3 = max( abs( wa1( n ) ), abs( wa1( 1 ) ) )
1428 ELSE
1429 temp3 = zero
1430 END IF
1431*
1432 result( 19 ) = ( temp1+temp2 ) / max( unfl, temp3*ulp )
1433*
1434* Call CSTEIN to compute eigenvectors corresponding to
1435* eigenvalues in WA1. (First call SSTEBZ again, to make sure
1436* it returns these eigenvalues in the correct order.)
1437*
1438 ntest = 21
1439 CALL sstebz( 'A', 'B', n, vl, vu, il, iu, abstol, sd, se, m,
1440 $ nsplit, wa1, iwork( 1 ), iwork( n+1 ), rwork,
1441 $ iwork( 2*n+1 ), iinfo )
1442 IF( iinfo.NE.0 ) THEN
1443 WRITE( nounit, fmt = 9999 )'SSTEBZ(A,B)', iinfo, n,
1444 $ jtype, 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 CALL cstein( n, sd, se, m, wa1, iwork( 1 ), iwork( n+1 ), z,
1456 $ ldu, rwork, iwork( 2*n+1 ), iwork( 3*n+1 ),
1457 $ iinfo )
1458 IF( iinfo.NE.0 ) THEN
1459 WRITE( nounit, fmt = 9999 )'CSTEIN', iinfo, n, jtype,
1460 $ ioldsd
1461 info = abs( iinfo )
1462 IF( iinfo.LT.0 ) THEN
1463 RETURN
1464 ELSE
1465 result( 20 ) = ulpinv
1466 result( 21 ) = ulpinv
1467 GO TO 280
1468 END IF
1469 END IF
1470*
1471* Do tests 20 and 21
1472*
1473 CALL cstt21( n, 0, sd, se, wa1, dumma, z, ldu, work, rwork,
1474 $ result( 20 ) )
1475*
1476* Call CSTEDC(I) to compute D1 and Z, do tests.
1477*
1478* Compute D1 and Z
1479*
1480 inde = 1
1481 indrwk = inde + n
1482 CALL scopy( n, sd, 1, d1, 1 )
1483 IF( n.GT.0 )
1484 $ CALL scopy( n-1, se, 1, rwork( inde ), 1 )
1485 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1486*
1487 ntest = 22
1488 CALL cstedc( 'I', n, d1, rwork( inde ), z, ldu, work, lwedc,
1489 $ rwork( indrwk ), lrwedc, iwork, liwedc, iinfo )
1490 IF( iinfo.NE.0 ) THEN
1491 WRITE( nounit, fmt = 9999 )'CSTEDC(I)', iinfo, n, jtype,
1492 $ ioldsd
1493 info = abs( iinfo )
1494 IF( iinfo.LT.0 ) THEN
1495 RETURN
1496 ELSE
1497 result( 22 ) = ulpinv
1498 GO TO 280
1499 END IF
1500 END IF
1501*
1502* Do Tests 22 and 23
1503*
1504 CALL cstt21( n, 0, sd, se, d1, dumma, z, ldu, work, rwork,
1505 $ result( 22 ) )
1506*
1507* Call CSTEDC(V) to compute D1 and Z, do tests.
1508*
1509* Compute D1 and Z
1510*
1511 CALL scopy( n, sd, 1, d1, 1 )
1512 IF( n.GT.0 )
1513 $ CALL scopy( n-1, se, 1, rwork( inde ), 1 )
1514 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1515*
1516 ntest = 24
1517 CALL cstedc( 'V', n, d1, rwork( inde ), z, ldu, work, lwedc,
1518 $ rwork( indrwk ), lrwedc, iwork, liwedc, iinfo )
1519 IF( iinfo.NE.0 ) THEN
1520 WRITE( nounit, fmt = 9999 )'CSTEDC(V)', iinfo, n, jtype,
1521 $ ioldsd
1522 info = abs( iinfo )
1523 IF( iinfo.LT.0 ) THEN
1524 RETURN
1525 ELSE
1526 result( 24 ) = ulpinv
1527 GO TO 280
1528 END IF
1529 END IF
1530*
1531* Do Tests 24 and 25
1532*
1533 CALL cstt21( n, 0, sd, se, d1, dumma, z, ldu, work, rwork,
1534 $ result( 24 ) )
1535*
1536* Call CSTEDC(N) to compute D2, do tests.
1537*
1538* Compute D2
1539*
1540 CALL scopy( n, sd, 1, d2, 1 )
1541 IF( n.GT.0 )
1542 $ CALL scopy( n-1, se, 1, rwork( inde ), 1 )
1543 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1544*
1545 ntest = 26
1546 CALL cstedc( 'N', n, d2, rwork( inde ), z, ldu, work, lwedc,
1547 $ rwork( indrwk ), lrwedc, iwork, liwedc, iinfo )
1548 IF( iinfo.NE.0 ) THEN
1549 WRITE( nounit, fmt = 9999 )'CSTEDC(N)', iinfo, n, jtype,
1550 $ ioldsd
1551 info = abs( iinfo )
1552 IF( iinfo.LT.0 ) THEN
1553 RETURN
1554 ELSE
1555 result( 26 ) = ulpinv
1556 GO TO 280
1557 END IF
1558 END IF
1559*
1560* Do Test 26
1561*
1562 temp1 = zero
1563 temp2 = zero
1564*
1565 DO 210 j = 1, n
1566 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1567 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1568 210 CONTINUE
1569*
1570 result( 26 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
1571*
1572* Only test CSTEMR if IEEE compliant
1573*
1574 IF( ilaenv( 10, 'CSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 .AND.
1575 $ ilaenv( 11, 'CSTEMR', 'VA', 1, 0, 0, 0 ).EQ.1 ) THEN
1576*
1577* Call CSTEMR, do test 27 (relative eigenvalue accuracy)
1578*
1579* If S is positive definite and diagonally dominant,
1580* ask for all eigenvalues with high relative accuracy.
1581*
1582 vl = zero
1583 vu = zero
1584 il = 0
1585 iu = 0
1586 IF( jtype.EQ.21 .AND. crel ) THEN
1587 ntest = 27
1588 abstol = unfl + unfl
1589 CALL cstemr( 'V', 'A', n, sd, se, vl, vu, il, iu,
1590 $ m, wr, z, ldu, n, iwork( 1 ), tryrac,
1591 $ rwork, lrwork, iwork( 2*n+1 ), lwork-2*n,
1592 $ iinfo )
1593 IF( iinfo.NE.0 ) THEN
1594 WRITE( nounit, fmt = 9999 )'CSTEMR(V,A,rel)',
1595 $ iinfo, n, jtype, ioldsd
1596 info = abs( iinfo )
1597 IF( iinfo.LT.0 ) THEN
1598 RETURN
1599 ELSE
1600 result( 27 ) = ulpinv
1601 GO TO 270
1602 END IF
1603 END IF
1604*
1605* Do test 27
1606*
1607 temp2 = two*( two*n-one )*ulp*( one+eight*half**2 ) /
1608 $ ( one-half )**4
1609*
1610 temp1 = zero
1611 DO 220 j = 1, n
1612 temp1 = max( temp1, abs( d4( j )-wr( n-j+1 ) ) /
1613 $ ( abstol+abs( d4( j ) ) ) )
1614 220 CONTINUE
1615*
1616 result( 27 ) = temp1 / temp2
1617*
1618 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1619 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1620 IF( iu.LT.il ) THEN
1621 itemp = iu
1622 iu = il
1623 il = itemp
1624 END IF
1625*
1626 IF( crange ) THEN
1627 ntest = 28
1628 abstol = unfl + unfl
1629 CALL cstemr( 'V', 'I', n, sd, se, vl, vu, il, iu,
1630 $ m, wr, z, ldu, n, iwork( 1 ), tryrac,
1631 $ rwork, lrwork, iwork( 2*n+1 ),
1632 $ lwork-2*n, iinfo )
1633*
1634 IF( iinfo.NE.0 ) THEN
1635 WRITE( nounit, fmt = 9999 )'CSTEMR(V,I,rel)',
1636 $ iinfo, n, jtype, ioldsd
1637 info = abs( iinfo )
1638 IF( iinfo.LT.0 ) THEN
1639 RETURN
1640 ELSE
1641 result( 28 ) = ulpinv
1642 GO TO 270
1643 END IF
1644 END IF
1645*
1646*
1647* Do test 28
1648*
1649 temp2 = two*( two*n-one )*ulp*
1650 $ ( one+eight*half**2 ) / ( one-half )**4
1651*
1652 temp1 = zero
1653 DO 230 j = il, iu
1654 temp1 = max( temp1, abs( wr( j-il+1 )-d4( n-j+
1655 $ 1 ) ) / ( abstol+abs( wr( j-il+1 ) ) ) )
1656 230 CONTINUE
1657*
1658 result( 28 ) = temp1 / temp2
1659 ELSE
1660 result( 28 ) = zero
1661 END IF
1662 ELSE
1663 result( 27 ) = zero
1664 result( 28 ) = zero
1665 END IF
1666*
1667* Call CSTEMR(V,I) to compute D1 and Z, do tests.
1668*
1669* Compute D1 and Z
1670*
1671 CALL scopy( n, sd, 1, d5, 1 )
1672 IF( n.GT.0 )
1673 $ CALL scopy( n-1, se, 1, rwork, 1 )
1674 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1675*
1676 IF( crange ) THEN
1677 ntest = 29
1678 il = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1679 iu = 1 + ( n-1 )*int( slarnd( 1, iseed2 ) )
1680 IF( iu.LT.il ) THEN
1681 itemp = iu
1682 iu = il
1683 il = itemp
1684 END IF
1685 CALL cstemr( 'V', 'I', n, d5, rwork, vl, vu, il, iu,
1686 $ m, d1, z, ldu, n, iwork( 1 ), tryrac,
1687 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1688 $ liwork-2*n, iinfo )
1689 IF( iinfo.NE.0 ) THEN
1690 WRITE( nounit, fmt = 9999 )'CSTEMR(V,I)', iinfo,
1691 $ n, jtype, ioldsd
1692 info = abs( iinfo )
1693 IF( iinfo.LT.0 ) THEN
1694 RETURN
1695 ELSE
1696 result( 29 ) = ulpinv
1697 GO TO 280
1698 END IF
1699 END IF
1700*
1701* Do Tests 29 and 30
1702*
1703*
1704* Call CSTEMR to compute D2, do tests.
1705*
1706* Compute D2
1707*
1708 CALL scopy( n, sd, 1, d5, 1 )
1709 IF( n.GT.0 )
1710 $ CALL scopy( n-1, se, 1, rwork, 1 )
1711*
1712 ntest = 31
1713 CALL cstemr( 'N', 'I', n, d5, rwork, vl, vu, il, iu,
1714 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1715 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1716 $ liwork-2*n, iinfo )
1717 IF( iinfo.NE.0 ) THEN
1718 WRITE( nounit, fmt = 9999 )'CSTEMR(N,I)', iinfo,
1719 $ n, jtype, ioldsd
1720 info = abs( iinfo )
1721 IF( iinfo.LT.0 ) THEN
1722 RETURN
1723 ELSE
1724 result( 31 ) = ulpinv
1725 GO TO 280
1726 END IF
1727 END IF
1728*
1729* Do Test 31
1730*
1731 temp1 = zero
1732 temp2 = zero
1733*
1734 DO 240 j = 1, iu - il + 1
1735 temp1 = max( temp1, abs( d1( j ) ),
1736 $ abs( d2( j ) ) )
1737 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1738 240 CONTINUE
1739*
1740 result( 31 ) = temp2 / max( unfl,
1741 $ ulp*max( temp1, temp2 ) )
1742*
1743*
1744* Call CSTEMR(V,V) to compute D1 and Z, do tests.
1745*
1746* Compute D1 and Z
1747*
1748 CALL scopy( n, sd, 1, d5, 1 )
1749 IF( n.GT.0 )
1750 $ CALL scopy( n-1, se, 1, rwork, 1 )
1751 CALL claset( 'Full', n, n, czero, cone, z, ldu )
1752*
1753 ntest = 32
1754*
1755 IF( n.GT.0 ) THEN
1756 IF( il.NE.1 ) THEN
1757 vl = d2( il ) - max( half*
1758 $ ( d2( il )-d2( il-1 ) ), ulp*anorm,
1759 $ two*rtunfl )
1760 ELSE
1761 vl = d2( 1 ) - max( half*( d2( n )-d2( 1 ) ),
1762 $ ulp*anorm, two*rtunfl )
1763 END IF
1764 IF( iu.NE.n ) THEN
1765 vu = d2( iu ) + max( half*
1766 $ ( d2( iu+1 )-d2( iu ) ), ulp*anorm,
1767 $ two*rtunfl )
1768 ELSE
1769 vu = d2( n ) + max( half*( d2( n )-d2( 1 ) ),
1770 $ ulp*anorm, two*rtunfl )
1771 END IF
1772 ELSE
1773 vl = zero
1774 vu = one
1775 END IF
1776*
1777 CALL cstemr( 'V', 'V', n, d5, rwork, vl, vu, il, iu,
1778 $ m, d1, z, ldu, n, iwork( 1 ), tryrac,
1779 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1780 $ liwork-2*n, iinfo )
1781 IF( iinfo.NE.0 ) THEN
1782 WRITE( nounit, fmt = 9999 )'CSTEMR(V,V)', iinfo,
1783 $ n, jtype, ioldsd
1784 info = abs( iinfo )
1785 IF( iinfo.LT.0 ) THEN
1786 RETURN
1787 ELSE
1788 result( 32 ) = ulpinv
1789 GO TO 280
1790 END IF
1791 END IF
1792*
1793* Do Tests 32 and 33
1794*
1795 CALL cstt22( n, m, 0, sd, se, d1, dumma, z, ldu, work,
1796 $ m, rwork, result( 32 ) )
1797*
1798* Call CSTEMR to compute D2, do tests.
1799*
1800* Compute D2
1801*
1802 CALL scopy( n, sd, 1, d5, 1 )
1803 IF( n.GT.0 )
1804 $ CALL scopy( n-1, se, 1, rwork, 1 )
1805*
1806 ntest = 34
1807 CALL cstemr( 'N', 'V', n, d5, rwork, vl, vu, il, iu,
1808 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1809 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1810 $ liwork-2*n, iinfo )
1811 IF( iinfo.NE.0 ) THEN
1812 WRITE( nounit, fmt = 9999 )'CSTEMR(N,V)', iinfo,
1813 $ n, jtype, ioldsd
1814 info = abs( iinfo )
1815 IF( iinfo.LT.0 ) THEN
1816 RETURN
1817 ELSE
1818 result( 34 ) = ulpinv
1819 GO TO 280
1820 END IF
1821 END IF
1822*
1823* Do Test 34
1824*
1825 temp1 = zero
1826 temp2 = zero
1827*
1828 DO 250 j = 1, iu - il + 1
1829 temp1 = max( temp1, abs( d1( j ) ),
1830 $ abs( d2( j ) ) )
1831 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1832 250 CONTINUE
1833*
1834 result( 34 ) = temp2 / max( unfl,
1835 $ ulp*max( temp1, temp2 ) )
1836 ELSE
1837 result( 29 ) = zero
1838 result( 30 ) = zero
1839 result( 31 ) = zero
1840 result( 32 ) = zero
1841 result( 33 ) = zero
1842 result( 34 ) = zero
1843 END IF
1844*
1845*
1846* Call CSTEMR(V,A) to compute D1 and Z, do tests.
1847*
1848* Compute D1 and Z
1849*
1850 CALL scopy( n, sd, 1, d5, 1 )
1851 IF( n.GT.0 )
1852 $ CALL scopy( n-1, se, 1, rwork, 1 )
1853*
1854 ntest = 35
1855*
1856 CALL cstemr( 'V', 'A', n, d5, rwork, vl, vu, il, iu,
1857 $ m, d1, z, ldu, n, iwork( 1 ), tryrac,
1858 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1859 $ liwork-2*n, iinfo )
1860 IF( iinfo.NE.0 ) THEN
1861 WRITE( nounit, fmt = 9999 )'CSTEMR(V,A)', iinfo, n,
1862 $ jtype, ioldsd
1863 info = abs( iinfo )
1864 IF( iinfo.LT.0 ) THEN
1865 RETURN
1866 ELSE
1867 result( 35 ) = ulpinv
1868 GO TO 280
1869 END IF
1870 END IF
1871*
1872* Do Tests 35 and 36
1873*
1874 CALL cstt22( n, m, 0, sd, se, d1, dumma, z, ldu, work, m,
1875 $ rwork, result( 35 ) )
1876*
1877* Call CSTEMR to compute D2, do tests.
1878*
1879* Compute D2
1880*
1881 CALL scopy( n, sd, 1, d5, 1 )
1882 IF( n.GT.0 )
1883 $ CALL scopy( n-1, se, 1, rwork, 1 )
1884*
1885 ntest = 37
1886 CALL cstemr( 'N', 'A', n, d5, rwork, vl, vu, il, iu,
1887 $ m, d2, z, ldu, n, iwork( 1 ), tryrac,
1888 $ rwork( n+1 ), lrwork-n, iwork( 2*n+1 ),
1889 $ liwork-2*n, iinfo )
1890 IF( iinfo.NE.0 ) THEN
1891 WRITE( nounit, fmt = 9999 )'CSTEMR(N,A)', iinfo, n,
1892 $ jtype, ioldsd
1893 info = abs( iinfo )
1894 IF( iinfo.LT.0 ) THEN
1895 RETURN
1896 ELSE
1897 result( 37 ) = ulpinv
1898 GO TO 280
1899 END IF
1900 END IF
1901*
1902* Do Test 34
1903*
1904 temp1 = zero
1905 temp2 = zero
1906*
1907 DO 260 j = 1, n
1908 temp1 = max( temp1, abs( d1( j ) ), abs( d2( j ) ) )
1909 temp2 = max( temp2, abs( d1( j )-d2( j ) ) )
1910 260 CONTINUE
1911*
1912 result( 37 ) = temp2 / max( unfl,
1913 $ ulp*max( temp1, temp2 ) )
1914 END IF
1915 270 CONTINUE
1916 280 CONTINUE
1917 ntestt = ntestt + ntest
1918*
1919* End of Loop -- Check for RESULT(j) > THRESH
1920*
1921*
1922* Print out tests which fail.
1923*
1924 DO 290 jr = 1, ntest
1925 IF( result( jr ).GE.thresh ) THEN
1926*
1927* If this is the first test to fail,
1928* print a header to the data file.
1929*
1930 IF( nerrs.EQ.0 ) THEN
1931 WRITE( nounit, fmt = 9998 )'CST'
1932 WRITE( nounit, fmt = 9997 )
1933 WRITE( nounit, fmt = 9996 )
1934 WRITE( nounit, fmt = 9995 )'Hermitian'
1935 WRITE( nounit, fmt = 9994 )
1936*
1937* Tests performed
1938*
1939 WRITE( nounit, fmt = 9987 )
1940 END IF
1941 nerrs = nerrs + 1
1942 IF( result( jr ).LT.10000.0e0 ) THEN
1943 WRITE( nounit, fmt = 9989 )n, jtype, ioldsd, jr,
1944 $ result( jr )
1945 ELSE
1946 WRITE( nounit, fmt = 9988 )n, jtype, ioldsd, jr,
1947 $ result( jr )
1948 END IF
1949 END IF
1950 290 CONTINUE
1951 300 CONTINUE
1952 310 CONTINUE
1953*
1954* Summary
1955*
1956 CALL slasum( 'CST', nounit, nerrs, ntestt )
1957 RETURN
1958*
1959 9999 FORMAT( ' CCHKST: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1960 $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1961*
1962 9998 FORMAT( / 1x, a3, ' -- Complex Hermitian eigenvalue problem' )
1963 9997 FORMAT( ' Matrix types (see CCHKST for details): ' )
1964*
1965 9996 FORMAT( / ' Special Matrices:',
1966 $ / ' 1=Zero matrix. ',
1967 $ ' 5=Diagonal: clustered entries.',
1968 $ / ' 2=Identity matrix. ',
1969 $ ' 6=Diagonal: large, evenly spaced.',
1970 $ / ' 3=Diagonal: evenly spaced entries. ',
1971 $ ' 7=Diagonal: small, evenly spaced.',
1972 $ / ' 4=Diagonal: geometr. spaced entries.' )
1973 9995 FORMAT( ' Dense ', a, ' Matrices:',
1974 $ / ' 8=Evenly spaced eigenvals. ',
1975 $ ' 12=Small, evenly spaced eigenvals.',
1976 $ / ' 9=Geometrically spaced eigenvals. ',
1977 $ ' 13=Matrix with random O(1) entries.',
1978 $ / ' 10=Clustered eigenvalues. ',
1979 $ ' 14=Matrix with large random entries.',
1980 $ / ' 11=Large, evenly spaced eigenvals. ',
1981 $ ' 15=Matrix with small random entries.' )
1982 9994 FORMAT( ' 16=Positive definite, evenly spaced eigenvalues',
1983 $ / ' 17=Positive definite, geometrically spaced eigenvlaues',
1984 $ / ' 18=Positive definite, clustered eigenvalues',
1985 $ / ' 19=Positive definite, small evenly spaced eigenvalues',
1986 $ / ' 20=Positive definite, large evenly spaced eigenvalues',
1987 $ / ' 21=Diagonally dominant tridiagonal, geometrically',
1988 $ ' spaced eigenvalues' )
1989*
1990 9989 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
1991 $ 4( i4, ',' ), ' result ', i3, ' is', 0p, f8.2 )
1992 9988 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
1993 $ 4( i4, ',' ), ' result ', i3, ' is', 1p, e10.3 )
1994*
1995 9987 FORMAT( / 'Test performed: see CCHKST for details.', / )
1996* End of CCHKST
1997*
1998 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
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
subroutine chet21(itype, uplo, n, kband, a, lda, d, e, u, ldu, v, ldv, tau, work, rwork, result)
CHET21
Definition chet21.f:214
subroutine chpt21(itype, uplo, n, kband, ap, d, e, u, ldu, vp, tau, work, rwork, result)
CHPT21
Definition chpt21.f:228
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
subroutine clatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
CLATMS
Definition clatms.f:332
subroutine cstt21(n, kband, ad, ae, sd, se, u, ldu, work, rwork, result)
CSTT21
Definition cstt21.f:133
subroutine cstt22(n, m, kband, ad, ae, sd, se, u, ldu, work, ldwork, rwork, result)
CSTT22
Definition cstt22.f:145
subroutine scopy(n, sx, incx, sy, incy)
SCOPY
Definition scopy.f:82
subroutine ccopy(n, cx, incx, cy, incy)
CCOPY
Definition ccopy.f:81
subroutine chetrd(uplo, n, a, lda, d, e, tau, work, lwork, info)
CHETRD
Definition chetrd.f:192
subroutine chptrd(uplo, n, ap, d, e, tau, info)
CHPTRD
Definition chptrd.f:151
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
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
subroutine cpteqr(compz, n, d, e, z, ldz, work, info)
CPTEQR
Definition cpteqr.f:145
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
subroutine cstedc(compz, n, d, e, z, ldz, work, lwork, rwork, lrwork, iwork, liwork, info)
CSTEDC
Definition cstedc.f:206
subroutine cstein(n, d, e, m, w, iblock, isplit, z, ldz, work, iwork, ifail, info)
CSTEIN
Definition cstein.f:182
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:339
subroutine csteqr(compz, n, d, e, z, ldz, work, info)
CSTEQR
Definition csteqr.f:132
subroutine ssterf(n, d, e, info)
SSTERF
Definition ssterf.f:86
subroutine cungtr(uplo, n, a, lda, tau, work, lwork, info)
CUNGTR
Definition cungtr.f:123
subroutine cupgtr(uplo, n, ap, tau, q, ldq, work, info)
CUPGTR
Definition cupgtr.f:114
subroutine slasum(type, iounit, ie, nrun)
SLASUM
Definition slasum.f:41
subroutine sstech(n, a, b, eig, tol, work, info)
SSTECH
Definition sstech.f:101