LAPACK  3.9.1
LAPACK: Linear Algebra PACKage
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,
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  IF( iinfo.LT.0 ) 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  IF( iinfo.NE.0 ) THEN
1023  WRITE( nounit, fmt = 9999 )'ZHPTRD(U)', iinfo, n, jtype,
1024  $ ioldsd
1025  info = abs( iinfo )
1026  IF( iinfo.LT.0 ) 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  IF( iinfo.NE.0 ) THEN
1037  WRITE( nounit, fmt = 9999 )'ZUPGTR(U)', iinfo, n, jtype,
1038  $ ioldsd
1039  info = abs( iinfo )
1040  IF( iinfo.LT.0 ) 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  IF( iinfo.LT.0 ) THEN
1412  RETURN
1413  ELSE
1414  result( 19 ) = ulpinv
1415  GO TO 280
1416  END IF
1417  END IF
1418 *
1419  IF( m3.EQ.0 .AND. n.NE.0 ) 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  IF( n.GT.0 ) 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  IF( iinfo.NE.0 ) THEN
1445  WRITE( nounit, fmt = 9999 )'DSTEBZ(A,B)', iinfo, n,
1446  $ jtype, ioldsd
1447  info = abs( iinfo )
1448  IF( iinfo.LT.0 ) 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  IF( iinfo.NE.0 ) THEN
1461  WRITE( nounit, fmt = 9999 )'ZSTEIN', iinfo, n, jtype,
1462  $ ioldsd
1463  info = abs( iinfo )
1464  IF( iinfo.LT.0 ) 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  IF( n.GT.0 )
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  IF( iinfo.NE.0 ) THEN
1493  WRITE( nounit, fmt = 9999 )'ZSTEDC(I)', iinfo, n, jtype,
1494  $ ioldsd
1495  info = abs( iinfo )
1496  IF( iinfo.LT.0 ) 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  IF( n.GT.0 )
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  IF( iinfo.NE.0 ) THEN
1522  WRITE( nounit, fmt = 9999 )'ZSTEDC(V)', iinfo, n, jtype,
1523  $ ioldsd
1524  info = abs( iinfo )
1525  IF( iinfo.LT.0 ) 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  IF( n.GT.0 )
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  IF( iinfo.NE.0 ) THEN
1551  WRITE( nounit, fmt = 9999 )'ZSTEDC(N)', iinfo, n, jtype,
1552  $ ioldsd
1553  info = abs( iinfo )
1554  IF( iinfo.LT.0 ) 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', 1, 0, 0, 0 ).EQ.1 .AND.
1577  $ ilaenv( 11, 'ZSTEMR', 'VA', 1, 0, 0, 0 ).EQ.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  IF( jtype.EQ.21 .AND. 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  IF( iinfo.NE.0 ) THEN
1596  WRITE( nounit, fmt = 9999 )'ZSTEMR(V,A,rel)',
1597  $ iinfo, n, jtype, ioldsd
1598  info = abs( iinfo )
1599  IF( iinfo.LT.0 ) 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  IF( iu.LT.il ) 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  IF( iinfo.NE.0 ) THEN
1637  WRITE( nounit, fmt = 9999 )'ZSTEMR(V,I,rel)',
1638  $ iinfo, n, jtype, ioldsd
1639  info = abs( iinfo )
1640  IF( iinfo.LT.0 ) 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  IF( n.GT.0 )
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  IF( iu.LT.il ) 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  IF( iinfo.NE.0 ) THEN
1692  WRITE( nounit, fmt = 9999 )'ZSTEMR(V,I)', iinfo,
1693  $ n, jtype, ioldsd
1694  info = abs( iinfo )
1695  IF( iinfo.LT.0 ) 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  IF( n.GT.0 )
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  IF( iinfo.NE.0 ) THEN
1720  WRITE( nounit, fmt = 9999 )'ZSTEMR(N,I)', iinfo,
1721  $ n, jtype, ioldsd
1722  info = abs( iinfo )
1723  IF( iinfo.LT.0 ) 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  IF( n.GT.0 )
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  IF( n.GT.0 ) THEN
1758  IF( il.NE.1 ) 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  IF( iu.NE.n ) 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  IF( iinfo.NE.0 ) THEN
1784  WRITE( nounit, fmt = 9999 )'ZSTEMR(V,V)', iinfo,
1785  $ n, jtype, ioldsd
1786  info = abs( iinfo )
1787  IF( iinfo.LT.0 ) 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  IF( n.GT.0 )
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
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
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
subroutine dsterf(N, D, E, INFO)
DSTERF
Definition: dsterf.f:86
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:81
subroutine zstt22(N, M, KBAND, AD, AE, SD, SE, U, LDU, WORK, LDWORK, RWORK, RESULT)
ZSTT22
Definition: zstt22.f:145
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
subroutine zstt21(N, KBAND, AD, AE, SD, SE, U, LDU, WORK, RWORK, RESULT)
ZSTT21
Definition: zstt21.f:133
subroutine zhet21(ITYPE, UPLO, N, KBAND, A, LDA, D, E, U, LDU, V, LDV, TAU, WORK, RWORK, RESULT)
ZHET21
Definition: zhet21.f:214
subroutine zhpt21(ITYPE, UPLO, N, KBAND, AP, D, E, U, LDU, VP, TAU, WORK, RWORK, RESULT)
ZHPT21
Definition: zhpt21.f:228
subroutine zlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
ZLATMS
Definition: zlatms.f:332
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
subroutine zhetrd(UPLO, N, A, LDA, D, E, TAU, WORK, LWORK, INFO)
ZHETRD
Definition: zhetrd.f:192
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
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
subroutine zsteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
ZSTEQR
Definition: zsteqr.f:132
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
subroutine zungtr(UPLO, N, A, LDA, TAU, WORK, LWORK, INFO)
ZUNGTR
Definition: zungtr.f:123
subroutine zhptrd(UPLO, N, AP, D, E, TAU, INFO)
ZHPTRD
Definition: zhptrd.f:151
subroutine zstein(N, D, E, M, W, IBLOCK, ISPLIT, Z, LDZ, WORK, IWORK, IFAIL, INFO)
ZSTEIN
Definition: zstein.f:182
subroutine zstedc(COMPZ, N, D, E, Z, LDZ, WORK, LWORK, RWORK, LRWORK, IWORK, LIWORK, INFO)
ZSTEDC
Definition: zstedc.f:212
subroutine zupgtr(UPLO, N, AP, TAU, Q, LDQ, WORK, INFO)
ZUPGTR
Definition: zupgtr.f:114
subroutine zpteqr(COMPZ, N, D, E, Z, LDZ, WORK, INFO)
ZPTEQR
Definition: zpteqr.f:145
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:82
subroutine dlasum(TYPE, IOUNIT, IE, NRUN)
DLASUM
Definition: dlasum.f:43
subroutine dstech(N, A, B, EIG, TOL, WORK, INFO)
DSTECH
Definition: dstech.f:101