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