OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
zchkst.f
Go to the documentation of this file.
1*> \brief \b ZCHKST
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 ZCHKST( 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* DOUBLE PRECISION THRESH
21* ..
22* .. Array Arguments ..
23* LOGICAL DOTYPE( * )
24* INTEGER ISEED( 4 ), IWORK( * ), NN( * )
25* DOUBLE PRECISION D1( * ), D2( * ), D3( * ), D4( * ), D5( * ),
26* $ RESULT( * ), RWORK( * ), SD( * ), SE( * ),
27* $ WA1( * ), WA2( * ), WA3( * ), WR( * )
28* COMPLEX*16 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*> ZCHKST checks the Hermitian eigenvalue problem routines.
39*>
40*> ZHETRD factors A as U S U* , where * means conjugate transpose,
41*> S is real symmetric tridiagonal, and U is unitary.
42*> ZHETRD can use either just the lower or just the upper triangle
43*> of A; ZCHKST 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*> ZHPTRD does the same as ZHETRD, except that A and V are stored
49*> in "packed" format.
50*>
51*> ZUNGTR constructs the matrix U from the contents of V and TAU.
52*>
53*> ZUPGTR constructs the matrix U from the contents of VP and TAU.
54*>
55*> ZSTEQR 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*> DSTERF computes D3, the matrix of eigenvalues, by the
61*> PWK method, which does not yield eigenvectors.
62*>
63*> ZPTEQR 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*> DSTEBZ 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*> ZSTEIN computes Y, the eigenvectors of S, given the
75*> eigenvalues.
76*>
77*> ZSTEDC 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 ZHETRD/ZUNGTR or ZHPTRD/ZUPGTR ('V' option). It may
82*> also just compute eigenvalues ('N' option).
83*>
84*> ZSTEMR 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). ZSTEMR
87*> uses the Relatively Robust Representation whenever possible.
88*>
89*> When ZCHKST 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 ) ZHETRD( UPLO='U', ... )
96*>
97*> (2) | I - UV* | / ( n ulp ) ZUNGTR( UPLO='U', ... )
98*>
99*> (3) | A - V S V* | / ( |A| n ulp ) ZHETRD( UPLO='L', ... )
100*>
101*> (4) | I - UV* | / ( n ulp ) ZUNGTR( UPLO='L', ... )
102*>
103*> (5-8) Same as 1-4, but for ZHPTRD and ZUPGTR.
104*>
105*> (9) | S - Z D Z* | / ( |S| n ulp ) ZSTEQR('V',...)
106*>
107*> (10) | I - ZZ* | / ( n ulp ) ZSTEQR('V',...)
108*>
109*> (11) | D1 - D2 | / ( |D1| ulp ) ZSTEQR('N',...)
110*>
111*> (12) | D1 - D3 | / ( |D1| ulp ) DSTERF
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*> DSTECH)
117*>
118*> For S positive definite,
119*>
120*> (14) | S - Z4 D4 Z4* | / ( |S| n ulp ) ZPTEQR('V',...)
121*>
122*> (15) | I - Z4 Z4* | / ( n ulp ) ZPTEQR('V',...)
123*>
124*> (16) | D4 - D5 | / ( 100 |D4| ulp ) ZPTEQR('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*> DSTEBZ( 'A', 'E', ...)
132*>
133*> (18) | WA1 - D3 | / ( |D3| ulp ) DSTEBZ( '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*> DSTEBZ( 'I', 'E', ...)
140*>
141*> (20) | S - Y WA1 Y* | / ( |S| n ulp ) DSTEBZ, ZSTEIN
142*>
143*> (21) | I - Y Y* | / ( n ulp ) DSTEBZ, ZSTEIN
144*>
145*> (22) | S - Z D Z* | / ( |S| n ulp ) ZSTEDC('I')
146*>
147*> (23) | I - ZZ* | / ( n ulp ) ZSTEDC('I')
148*>
149*> (24) | S - Z D Z* | / ( |S| n ulp ) ZSTEDC('V')
150*>
151*> (25) | I - ZZ* | / ( n ulp ) ZSTEDC('V')
152*>
153*> (26) | D1 - D2 | / ( |D1| ulp ) ZSTEDC('V') and
154*> ZSTEDC('N')
155*>
156*> Test 27 is disabled at the moment because ZSTEMR 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*> ZSTEMR('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*> ZSTEMR('V', 'I')
168*>
169*> Tests 29 through 34 are disable at present because ZSTEMR
170*> does not handle partial spectrum requests.
171*>
172*> (29) | S - Z D Z* | / ( |S| n ulp ) ZSTEMR('V', 'I')
173*>
174*> (30) | I - ZZ* | / ( n ulp ) ZSTEMR('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*> ZSTEMR('N', 'I') vs. CSTEMR('V', 'I')
181*>
182*> (32) | S - Z D Z* | / ( |S| n ulp ) ZSTEMR('V', 'V')
183*>
184*> (33) | I - ZZ* | / ( n ulp ) ZSTEMR('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*> ZSTEMR('N', 'V') vs. CSTEMR('V', 'V')
191*>
192*> (35) | S - Z D Z* | / ( |S| n ulp ) ZSTEMR('V', 'A')
193*>
194*> (36) | I - ZZ* | / ( n ulp ) ZSTEMR('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*> ZSTEMR('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*> ZCHKST 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, ZCHKST
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 ZCHKST to continue the same random number
301*> sequence.
302*> \endverbatim
303*>
304*> \param[in] THRESH
305*> \verbatim
306*> THRESH is DOUBLE PRECISION
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*16 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*16 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 DOUBLE PRECISION array of
348*> dimension( max(NN) )
349*> The diagonal of the tridiagonal matrix computed by ZHETRD.
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 DOUBLE PRECISION array of
357*> dimension( max(NN) )
358*> The off-diagonal of the tridiagonal matrix computed by
359*> ZHETRD. 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 DOUBLE PRECISION array of
366*> dimension( max(NN) )
367*> The eigenvalues of A, as computed by ZSTEQR 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 DOUBLE PRECISION array of
375*> dimension( max(NN) )
376*> The eigenvalues of A, as computed by ZSTEQR 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 DOUBLE PRECISION array of
384*> dimension( max(NN) )
385*> The eigenvalues of A, as computed by DSTERF. 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 DOUBLE PRECISION array of
392*> dimension( max(NN) )
393*> The eigenvalues of A, as computed by ZPTEQR(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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DSTEBZ.
414*> \endverbatim
415*>
416*> \param[out] WA2
417*> \verbatim
418*> WA2 is DOUBLE PRECISION 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 DSTEBZ.
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 DOUBLE PRECISION 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 DSTEBZ.
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*16 array of
450*> dimension( LDU, max(NN) ).
451*> The unitary matrix computed by ZHETRD + ZUNGTR.
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*16 array of
464*> dimension( LDU, max(NN) ).
465*> The Housholder vectors computed by ZHETRD 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 ZHETRD, 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 ZUNGTR, 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*16 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*16 array of
485*> dimension( max(NN) )
486*> The Householder factors computed by ZHETRD in reducing A
487*> to tridiagonal form.
488*> \endverbatim
489*>
490*> \param[out] Z
491*> \verbatim
492*> Z is COMPLEX*16 array of
493*> dimension( LDU, max(NN) ).
494*> The unitary matrix of eigenvectors computed by ZSTEQR,
495*> ZPTEQR, and ZSTEIN.
496*> \endverbatim
497*>
498*> \param[out] WORK
499*> \verbatim
500*> WORK is COMPLEX*16 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 ZLATMR, CLATMS, ZHETRD, ZUNGTR, ZSTEQR, DSTERF,
557*> or ZUNMC2 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 complex16_eig
597*
598* =====================================================================
599 SUBROUTINE zchkst( 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 DOUBLE PRECISION THRESH
613* ..
614* .. Array Arguments ..
615 LOGICAL DOTYPE( * )
616 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
617 DOUBLE PRECISION D1( * ), D2( * ), D3( * ), D4( * ), D5( * ),
618 $ RESULT( * ), RWORK( * ), SD( * ), SE( * ),
619 $ wa1( * ), wa2( * ), wa3( * ), wr( * )
620 COMPLEX*16 A( LDA, * ), AP( * ), TAU( * ), U( LDU, * ),
621 $ v( ldu, * ), vp( * ), work( * ), z( ldu, * )
622* ..
623*
624* =====================================================================
625*
626* .. Parameters ..
627 DOUBLE PRECISION ZERO, ONE, TWO, EIGHT, TEN, HUN
628 PARAMETER ( ZERO = 0.0d0, one = 1.0d0, two = 2.0d0,
629 $ eight = 8.0d0, ten = 10.0d0, hun = 100.0d0 )
630 COMPLEX*16 CZERO, CONE
631 parameter( czero = ( 0.0d+0, 0.0d+0 ),
632 $ cone = ( 1.0d+0, 0.0d+0 ) )
633 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION DUMMA( 1 )
658* ..
659* .. External Functions ..
660 INTEGER ILAENV
661 DOUBLE PRECISION DLAMCH, DLARND, DSXT1
662 EXTERNAL ILAENV, DLAMCH, DLARND, DSXT1
663* ..
664* .. External Subroutines ..
665 EXTERNAL dcopy, dlabad, dlasum, dstebz, dstech, dsterf,
666 $ xerbla, zcopy, zhet21, zhetrd, zhpt21, zhptrd,
667 $ zlacpy, zlaset, zlatmr, zlatms, zpteqr, zstedc,
668 $ zstemr, zstein, zsteqr, zstt21, zstt22, zungtr,
669 $ zupgtr
670* ..
671* .. Intrinsic Functions ..
672 INTRINSIC abs, dble, dconjg, int, log, max, min, 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, 'ZHETRD', '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( 'ZCHKST', -info )
724 RETURN
725 END IF
726*
727* Quick return if possible
728*
729 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
730 $ RETURN
731*
732* More Important constants
733*
734 unfl = dlamch( 'Safe minimum' )
735 ovfl = one / unfl
736 CALL dlabad( unfl, ovfl )
737 ulp = dlamch( 'Epsilon' )*dlamch( '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 IF( n.GT.0 ) THEN
754 lgn = int( log( dble( n ) ) / log( two ) )
755 IF( 2**lgn.LT.n )
756 $ lgn = lgn + 1
757 IF( 2**lgn.LT.n )
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 / dble( max( 1, n ) )
769*
770 IF( nsizes.NE.1 ) 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 IF( .NOT.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 IF( mtypes.GT.maxtyp )
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 zlaset( 'Full', lda, n, czero, czero, a, lda )
827 iinfo = 0
828 IF( jtype.LE.15 ) 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 IF( itype.EQ.1 ) THEN
839 iinfo = 0
840*
841 ELSE IF( itype.EQ.2 ) THEN
842*
843* Identity
844*
845 DO 80 jc = 1, n
846 a( jc, jc ) = anorm
847 80 CONTINUE
848*
849 ELSE IF( itype.EQ.4 ) THEN
850*
851* Diagonal Matrix, [Eigen]values Specified
852*
853 CALL zlatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
854 $ anorm, 0, 0, 'N', a, lda, work, iinfo )
855*
856*
857 ELSE IF( itype.EQ.5 ) THEN
858*
859* Hermitian, eigenvalues specified
860*
861 CALL zlatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
862 $ anorm, n, n, 'N', a, lda, work, iinfo )
863*
864 ELSE IF( itype.EQ.7 ) THEN
865*
866* Diagonal, random eigenvalues
867*
868 CALL zlatmr( 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 ELSE IF( itype.EQ.8 ) THEN
874*
875* Hermitian, random eigenvalues
876*
877 CALL zlatmr( 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 ELSE IF( itype.EQ.9 ) THEN
883*
884* Positive definite, eigenvalues specified.
885*
886 CALL zlatms( n, n, 'S', iseed, 'P', rwork, imode, cond,
887 $ anorm, n, n, 'N', a, lda, work, iinfo )
888*
889 ELSE IF( itype.EQ.10 ) THEN
890*
891* Positive definite tridiagonal, eigenvalues specified.
892*
893 CALL zlatms( 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 IF( temp1.GT.half*temp2 ) THEN
899 a( i-1, i ) = a( i-1, i )*
900 $ ( half*temp2 / ( unfl+temp1 ) )
901 a( i, i-1 ) = dconjg( a( i-1, i ) )
902 END IF
903 90 CONTINUE
904*
905 ELSE
906*
907 iinfo = 1
908 END IF
909*
910 IF( iinfo.NE.0 ) 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 ZHETRD and ZUNGTR to compute S and U from
920* upper triangle.
921*
922 CALL zlacpy( 'U', n, n, a, lda, v, ldu )
923*
924 ntest = 1
925 CALL zhetrd( 'U', n, v, ldu, sd, se, tau, work, lwork,
926 $ iinfo )
927*
928 IF( iinfo.NE.0 ) THEN
929 WRITE( nounit, fmt = 9999 )'ZHETRD(U)', iinfo, n, jtype,
930 $ ioldsd
931 info = abs( iinfo )
932 IF( iinfo.LT.0 ) THEN
933 RETURN
934 ELSE
935 result( 1 ) = ulpinv
936 GO TO 280
937 END IF
938 END IF
939*
940 CALL zlacpy( 'U', n, n, v, ldu, u, ldu )
941*
942 ntest = 2
943 CALL zungtr( 'U', n, u, ldu, tau, work, lwork, iinfo )
944 IF( iinfo.NE.0 ) THEN
945 WRITE( nounit, fmt = 9999 )'ZUNGTR(U)', iinfo, n, jtype,
946 $ ioldsd
947 info = abs( iinfo )
948 IF( iinfo.LT.0 ) 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 zhet21( 2, 'Upper', n, 1, a, lda, sd, se, u, ldu, v,
959 $ ldu, tau, work, rwork, result( 1 ) )
960 CALL zhet21( 3, 'Upper', n, 1, a, lda, sd, se, u, ldu, v,
961 $ ldu, tau, work, rwork, result( 2 ) )
962*
963* Call ZHETRD and ZUNGTR to compute S and U from
964* lower triangle, do tests.
965*
966 CALL zlacpy( 'L', n, n, a, lda, v, ldu )
967*
968 ntest = 3
969 CALL zhetrd( 'L', n, v, ldu, sd, se, tau, work, lwork,
970 $ iinfo )
971*
972 IF( iinfo.NE.0 ) THEN
973 WRITE( nounit, fmt = 9999 )'ZHETRD(L)', iinfo, n, jtype,
974 $ ioldsd
975 info = abs( iinfo )
976 IF( iinfo.LT.0 ) THEN
977 RETURN
978 ELSE
979 result( 3 ) = ulpinv
980 GO TO 280
981 END IF
982 END IF
983*
984 CALL zlacpy( 'L', n, n, v, ldu, u, ldu )
985*
986 ntest = 4
987 CALL zungtr( 'L', n, u, ldu, tau, work, lwork, iinfo )
988 IF( iinfo.NE.0 ) THEN
989 WRITE( nounit, fmt = 9999 )'zungtr(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 ZHET21( 2, 'lower', N, 1, A, LDA, SD, SE, U, LDU, V,
1001 $ LDU, TAU, WORK, RWORK, RESULT( 3 ) )
1002 CALL ZHET21( 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 ZHPTRD and ZUPGTR to compute S and U from AP
1016*
1017 CALL ZCOPY( NAP, AP, 1, VP, 1 )
1018*
1019 NTEST = 5
1020 CALL ZHPTRD( 'u', N, VP, SD, SE, TAU, IINFO )
1021*
1022.NE. IF( IINFO0 ) THEN
1023 WRITE( NOUNIT, FMT = 9999 )'zhptrd(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 ZUPGTR( 'u', N, VP, TAU, U, LDU, WORK, IINFO )
1036.NE. IF( IINFO0 ) THEN
1037 WRITE( NOUNIT, FMT = 9999 )'zupgtr(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 ZHPT21( 2, 'upper', N, 1, AP, SD, SE, U, LDU, VP, TAU,
1051 $ WORK, RWORK, RESULT( 5 ) )
1052 CALL ZHPT21( 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 ZHPTRD and ZUPGTR to compute S and U from AP
1066*
1067 CALL ZCOPY( NAP, AP, 1, VP, 1 )
1068*
1069 NTEST = 7
1070 CALL ZHPTRD( 'l', n, vp, sd, se, tau, iinfo )
1071*
1072 IF( iinfo.NE.0 ) THEN
1073 WRITE( nounit, fmt = 9999 )'ZHPTRD(L)', iinfo, n, jtype,
1074 $ ioldsd
1075 info = abs( iinfo )
1076 IF( iinfo.LT.0 ) 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 zupgtr( 'L', n, vp, tau, u, ldu, work, iinfo )
1086 IF( iinfo.NE.0 ) THEN
1087 WRITE( nounit, fmt = 9999 )'ZUPGTR(L)', iinfo, n, jtype,
1088 $ ioldsd
1089 info = abs( iinfo )
1090 IF( iinfo.LT.0 ) THEN
1091 RETURN
1092 ELSE
1093 result( 8 ) = ulpinv
1094 GO TO 280
1095 END IF
1096 END IF
1097*
1098 CALL zhpt21( 2, 'Lower', n, 1, ap, sd, se, u, ldu, vp, tau,
1099 $ work, rwork, result( 7 ) )
1100 CALL zhpt21( 3, 'Lower', n, 1, ap, sd, se, u, ldu, vp, tau,
1101 $ work, rwork, result( 8 ) )
1102*
1103* Call ZSTEQR to compute D1, D2, and Z, do tests.
1104*
1105* Compute D1 and Z
1106*
1107 CALL dcopy( n, sd, 1, d1, 1 )
1108 IF( n.GT.0 )
1109 $ CALL dcopy( n-1, se, 1, rwork, 1 )
1110 CALL zlaset( 'Full', n, n, czero, cone, z, ldu )
1111*
1112 ntest = 9
1113 CALL zsteqr( 'V', n, d1, rwork, z, ldu, rwork( n+1 ),
1114 $ iinfo )
1115 IF( iinfo.NE.0 ) THEN
1116 WRITE( nounit, fmt = 9999 )'ZSTEQR(V)', iinfo, n, jtype,
1117 $ ioldsd
1118 info = abs( iinfo )
1119 IF( iinfo.LT.0 ) 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 dcopy( n, sd, 1, d2, 1 )
1130 IF( n.GT.0 )
1131 $ CALL dcopy( n-1, se, 1, rwork, 1 )
1132*
1133 ntest = 11
1134 CALL zsteqr( 'N', n, d2, rwork, work, ldu, rwork( n+1 ),
1135 $ iinfo )
1136 IF( iinfo.NE.0 ) THEN
1137 WRITE( nounit, fmt = 9999 )'ZSTEQR(N)', iinfo, n, jtype,
1138 $ ioldsd
1139 info = abs( iinfo )
1140 IF( iinfo.LT.0 ) 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 dcopy( n, sd, 1, d3, 1 )
1151 IF( n.GT.0 )
1152 $ CALL dcopy( n-1, se, 1, rwork, 1 )
1153*
1154 ntest = 12
1155 CALL dsterf( n, d3, rwork, iinfo )
1156 IF( iinfo.NE.0 ) THEN
1157 WRITE( nounit, fmt = 9999 )'DSTERF', iinfo, n, jtype,
1158 $ ioldsd
1159 info = abs( iinfo )
1160 IF( iinfo.LT.0 ) 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 zstt21( 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 dstech( n, sd, se, d1, temp1, rwork, iinfo )
1198 IF( iinfo.EQ.0 )
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 ZPTEQR
1207* and do tests 14, 15, and 16 .
1208*
1209 IF( jtype.GT.15 ) THEN
1210*
1211* Compute D4 and Z4
1212*
1213 CALL dcopy( n, sd, 1, d4, 1 )
1214 IF( n.GT.0 )
1215 $ CALL dcopy( n-1, se, 1, rwork, 1 )
1216 CALL zlaset( 'Full', n, n, czero, cone, z, ldu )
1217*
1218 ntest = 14
1219 CALL zpteqr( 'V', n, d4, rwork, z, ldu, rwork( n+1 ),
1220 $ iinfo )
1221 IF( iinfo.NE.0 ) THEN
1222 WRITE( nounit, fmt = 9999 )'ZPTEQR(V)', iinfo, n,
1223 $ jtype, ioldsd
1224 info = abs( iinfo )
1225 IF( iinfo.LT.0 ) 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 zstt21( n, 0, sd, se, d4, dumma, z, ldu, work,
1236 $ rwork, result( 14 ) )
1237*
1238* Compute D5
1239*
1240 CALL dcopy( n, sd, 1, d5, 1 )
1241 IF( n.GT.0 )
1242 $ CALL dcopy( n-1, se, 1, rwork, 1 )
1243*
1244 ntest = 16
1245 CALL zpteqr( 'N', n, d5, rwork, z, ldu, rwork( n+1 ),
1246 $ iinfo )
1247 IF( iinfo.NE.0 ) THEN
1248 WRITE( nounit, fmt = 9999 )'ZPTEQR(N)', iinfo, n,
1249 $ jtype, ioldsd
1250 info = abs( iinfo )
1251 IF( iinfo.LT.0 ) 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 DSTEBZ 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 IF( jtype.EQ.21 ) THEN
1286 ntest = 17
1287 abstol = unfl + unfl
1288 CALL dstebz( '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 IF( iinfo.NE.0 ) THEN
1292 WRITE( nounit, fmt = 9999 )'DSTEBZ(A,rel)', iinfo, n,
1293 $ jtype, ioldsd
1294 info = abs( iinfo )
1295 IF( iinfo.LT.0 ) 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 dstebz( '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 IF( iinfo.NE.0 ) THEN
1327 WRITE( nounit, fmt = 9999 )'DSTEBZ(A)', iinfo, n, jtype,
1328 $ ioldsd
1329 info = abs( iinfo )
1330 IF( iinfo.LT.0 ) 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 IF( n.LE.1 ) THEN
1354 il = 1
1355 iu = n
1356 ELSE
1357 il = 1 + ( n-1 )*int( dlarnd( 1, iseed2 ) )
1358 iu = 1 + ( n-1 )*int( dlarnd( 1, iseed2 ) )
1359 IF( iu.LT.il ) THEN
1360 itemp = iu
1361 iu = il
1362 il = itemp
1363 END IF
1364 END IF
1365*
1366 CALL dstebz( '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 IF( iinfo.NE.0 ) THEN
1370 WRITE( nounit, fmt = 9999 )'DSTEBZ(I)', iinfo, n, jtype,
1371 $ ioldsd
1372 info = abs( iinfo )
1373 IF( iinfo.LT.0 ) 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 IF( n.GT.0 ) THEN
1385 IF( il.NE.1 ) 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 IF( iu.NE.n ) 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 dstebz( '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 IF( iinfo.NE.0 ) THEN
1408 WRITE( nounit, fmt = 9999 )'dstebz(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 = DSXT1( 1, WA2, M2, WA3, M3, ABSTOL, ULP, UNFL )
1427 TEMP2 = DSXT1( 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 ZSTEIN to compute eigenvectors corresponding to
1437* eigenvalues in WA1. (First call DSTEBZ again, to make sure
1438* it returns these eigenvalues in the correct order.)
1439*
1440 NTEST = 21
1441 CALL DSTEBZ( '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 )'dstebz(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 ZSTEIN( 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 )'zstein', 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 ZSTT21( N, 0, SD, SE, WA1, DUMMA, Z, LDU, WORK, RWORK,
1476 $ RESULT( 20 ) )
1477*
1478* Call ZSTEDC(I) to compute D1 and Z, do tests.
1479*
1480* Compute D1 and Z
1481*
1482 INDE = 1
1483 INDRWK = INDE + N
1484 CALL DCOPY( N, SD, 1, D1, 1 )
1485.GT. IF( N0 )
1486 $ CALL DCOPY( N-1, SE, 1, RWORK( INDE ), 1 )
1487 CALL ZLASET( 'full', N, N, CZERO, CONE, Z, LDU )
1488*
1489 NTEST = 22
1490 CALL ZSTEDC( '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 )'zstedc(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 ZSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, RWORK,
1507 $ RESULT( 22 ) )
1508*
1509* Call ZSTEDC(V) to compute D1 and Z, do tests.
1510*
1511* Compute D1 and Z
1512*
1513 CALL DCOPY( N, SD, 1, D1, 1 )
1514.GT. IF( N0 )
1515 $ CALL DCOPY( N-1, SE, 1, RWORK( INDE ), 1 )
1516 CALL ZLASET( 'full', N, N, CZERO, CONE, Z, LDU )
1517*
1518 NTEST = 24
1519 CALL ZSTEDC( '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 )'zstedc(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 ZSTT21( N, 0, SD, SE, D1, DUMMA, Z, LDU, WORK, RWORK,
1536 $ RESULT( 24 ) )
1537*
1538* Call ZSTEDC(N) to compute D2, do tests.
1539*
1540* Compute D2
1541*
1542 CALL DCOPY( N, SD, 1, D2, 1 )
1543.GT. IF( N0 )
1544 $ CALL DCOPY( N-1, SE, 1, RWORK( INDE ), 1 )
1545 CALL ZLASET( 'full', N, N, CZERO, CONE, Z, LDU )
1546*
1547 NTEST = 26
1548 CALL ZSTEDC( '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 )'zstedc(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 ZSTEMR if IEEE compliant
1575*
1576 IF( ILAENV( 10, 'zstemr', 'va.EQ..AND.', 1, 0, 0, 0 )1
1577 $ ILAENV( 11, 'zstemr', 'va.EQ.', 1, 0, 0, 0 )1 ) THEN
1578*
1579* Call ZSTEMR, 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 ZSTEMR( '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 )'zstemr(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( DLARND( 1, ISEED2 ) )
1621 IU = 1 + ( N-1 )*INT( DLARND( 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 ZSTEMR( '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 )'zstemr(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 ZSTEMR(V,I) to compute D1 and Z, do tests.
1670*
1671* Compute D1 and Z
1672*
1673 CALL DCOPY( N, SD, 1, D5, 1 )
1674.GT. IF( N0 )
1675 $ CALL DCOPY( N-1, SE, 1, RWORK, 1 )
1676 CALL ZLASET( 'full', N, N, CZERO, CONE, Z, LDU )
1677*
1678 IF( CRANGE ) THEN
1679 NTEST = 29
1680 IL = 1 + ( N-1 )*INT( DLARND( 1, ISEED2 ) )
1681 IU = 1 + ( N-1 )*INT( DLARND( 1, ISEED2 ) )
1682.LT. IF( IUIL ) THEN
1683 ITEMP = IU
1684 IU = IL
1685 IL = ITEMP
1686 END IF
1687 CALL ZSTEMR( '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 )'zstemr(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 ZSTEMR to compute D2, do tests.
1707*
1708* Compute D2
1709*
1710 CALL DCOPY( N, SD, 1, D5, 1 )
1711.GT. IF( N0 )
1712 $ CALL DCOPY( N-1, SE, 1, RWORK, 1 )
1713*
1714 NTEST = 31
1715 CALL ZSTEMR( '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 )'zstemr(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 ZSTEMR(V,V) to compute D1 and Z, do tests.
1747*
1748* Compute D1 and Z
1749*
1750 CALL DCOPY( N, SD, 1, D5, 1 )
1751.GT. IF( N0 )
1752 $ CALL DCOPY( N-1, SE, 1, RWORK, 1 )
1753 CALL ZLASET( '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 ZSTEMR( 'v', 'v', N, D5, RWORK, VL, VU, IL, IU,
1780 $ M, D1, Z, LDU, M, 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 )'zstemr(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 ZSTT22( N, M, 0, SD, SE, D1, DUMMA, Z, LDU, WORK,
1798 $ M, RWORK, RESULT( 32 ) )
1799*
1800* Call ZSTEMR to compute D2, do tests.
1801*
1802* Compute D2
1803*
1804 CALL DCOPY( N, SD, 1, D5, 1 )
1805.GT. IF( N0 )
1806 $ CALL DCOPY( N-1, SE, 1, RWORK, 1 )
1807*
1808 NTEST = 34
1809 CALL ZSTEMR( '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 IF( iinfo.NE.0 ) THEN
1814 WRITE( nounit, fmt = 9999 )'ZSTEMR(N,V)', iinfo,
1815 $ n, jtype, ioldsd
1816 info = abs( iinfo )
1817 IF( iinfo.LT.0 ) 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 ZSTEMR(V,A) to compute D1 and Z, do tests.
1849*
1850* Compute D1 and Z
1851*
1852 CALL dcopy( n, sd, 1, d5, 1 )
1853 IF( n.GT.0 )
1854 $ CALL dcopy( n-1, se, 1, rwork, 1 )
1855*
1856 ntest = 35
1857*
1858 CALL zstemr( '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 IF( iinfo.NE.0 ) THEN
1863 WRITE( nounit, fmt = 9999 )'ZSTEMR(V,A)', iinfo, n,
1864 $ jtype, ioldsd
1865 info = abs( iinfo )
1866 IF( iinfo.LT.0 ) 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 zstt22( n, m, 0, sd, se, d1, dumma, z, ldu, work, m,
1877 $ rwork, result( 35 ) )
1878*
1879* Call ZSTEMR to compute D2, do tests.
1880*
1881* Compute D2
1882*
1883 CALL dcopy( n, sd, 1, d5, 1 )
1884 IF( n.GT.0 )
1885 $ CALL dcopy( n-1, se, 1, rwork, 1 )
1886*
1887 ntest = 37
1888 CALL zstemr( '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 IF( iinfo.NE.0 ) THEN
1893 WRITE( nounit, fmt = 9999 )'ZSTEMR(N,A)', iinfo, n,
1894 $ jtype, ioldsd
1895 info = abs( iinfo )
1896 IF( iinfo.LT.0 ) 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 IF( result( jr ).GE.thresh ) THEN
1928*
1929* If this is the first test to fail,
1930* print a header to the data file.
1931*
1932 IF( nerrs.EQ.0 ) THEN
1933 WRITE( nounit, fmt = 9998 )'ZST'
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 IF( result( jr ).LT.10000.0d0 ) 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 dlasum( 'ZST', nounit, nerrs, ntestt )
1959 RETURN
1960*
1961 9999 FORMAT( ' ZCHKST: ', 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 ZCHKST 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, d10.3 )
1996*
1997 9987 FORMAT( / 'Test performed: see ZCHKST for details.', / )
1998* End of ZCHKST
1999*
2000 END
dlabad
subroutine dlabad(small, large)
DLABAD
Definition dlabad.f:74
dstebz
subroutine dstebz(range, order, n, vl, vu, il, iu, abstol, d, e, m, nsplit, w, iblock, isplit, work, iwork, info)
DSTEBZ
Definition dstebz.f:273
dsterf
subroutine dsterf(n, d, e, info)
DSTERF
Definition dsterf.f:86
xerbla
subroutine xerbla(srname, info)
XERBLA
Definition xerbla.f:60
zhetrd
subroutine zhetrd(uplo, n, a, lda, d, e, tau, work, lwork, info)
ZHETRD
Definition zhetrd.f:192
zlacpy
subroutine zlacpy(uplo, m, n, a, lda, b, ldb)
ZLACPY copies all or part of one two-dimensional array to another.
Definition zlacpy.f:103
zlaset
subroutine zlaset(uplo, m, n, alpha, beta, a, lda)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition zlaset.f:106
zupgtr
subroutine zupgtr(uplo, n, ap, tau, q, ldq, work, info)
ZUPGTR
Definition zupgtr.f:114
zstein
subroutine zstein(n, d, e, m, w, iblock, isplit, z, ldz, work, iwork, ifail, info)
ZSTEIN
Definition zstein.f:182
zsteqr
subroutine zsteqr(compz, n, d, e, z, ldz, work, info)
ZSTEQR
Definition zsteqr.f:132
zstemr
subroutine zstemr(jobz, range, n, d, e, vl, vu, il, iu, m, w, z, ldz, nzc, isuppz, tryrac, work, lwork, iwork, liwork, info)
ZSTEMR
Definition zstemr.f:338
zungtr
subroutine zungtr(uplo, n, a, lda, tau, work, lwork, info)
ZUNGTR
Definition zungtr.f:123
zstedc
subroutine zstedc(compz, n, d, e, z, ldz, work, lwork, rwork, lrwork, iwork, liwork, info)
ZSTEDC
Definition zstedc.f:212
zhptrd
subroutine zhptrd(uplo, n, ap, d, e, tau, info)
ZHPTRD
Definition zhptrd.f:151
zpteqr
subroutine zpteqr(compz, n, d, e, z, ldz, work, info)
ZPTEQR
Definition zpteqr.f:145
zcopy
subroutine zcopy(n, zx, incx, zy, incy)
ZCOPY
Definition zcopy.f:81
zhet21
subroutine zhet21(itype, uplo, n, kband, a, lda, d, e, u, ldu, v, ldv, tau, work, rwork, result)
ZHET21
Definition zhet21.f:214
zstt22
subroutine zstt22(n, m, kband, ad, ae, sd, se, u, ldu, work, ldwork, rwork, result)
ZSTT22
Definition zstt22.f:145
zhpt21
subroutine zhpt21(itype, uplo, n, kband, ap, d, e, u, ldu, vp, tau, work, rwork, result)
ZHPT21
Definition zhpt21.f:228
zchkst
subroutine zchkst(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)
ZCHKST
Definition zchkst.f:604
zstt21
subroutine zstt21(n, kband, ad, ae, sd, se, u, ldu, work, rwork, result)
ZSTT21
Definition zstt21.f:133
zlatmr
subroutine zlatmr(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)
ZLATMR
Definition zlatmr.f:490
zlatms
subroutine zlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
ZLATMS
Definition zlatms.f:332
dcopy
subroutine dcopy(n, dx, incx, dy, incy)
DCOPY
Definition dcopy.f:82
dlasum
subroutine dlasum(type, iounit, ie, nrun)
DLASUM
Definition dlasum.f:43
dstech
subroutine dstech(n, a, b, eig, tol, work, info)
DSTECH
Definition dstech.f:101
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21