LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zdrvsg.f
Go to the documentation of this file.
1*> \brief \b ZDRVSG
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 ZDRVSG( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
12* NOUNIT, A, LDA, B, LDB, D, Z, LDZ, AB, BB, AP,
13* BP, WORK, NWORK, RWORK, LRWORK, IWORK, LIWORK,
14* RESULT, INFO )
15*
16* .. Scalar Arguments ..
17* INTEGER INFO, LDA, LDB, LDZ, LIWORK, LRWORK, NOUNIT,
18* $ NSIZES, NTYPES, NWORK
19* DOUBLE PRECISION THRESH
20* ..
21* .. Array Arguments ..
22* LOGICAL DOTYPE( * )
23* INTEGER ISEED( 4 ), IWORK( * ), NN( * )
24* DOUBLE PRECISION D( * ), RESULT( * ), RWORK( * )
25* COMPLEX*16 A( LDA, * ), AB( LDA, * ), AP( * ),
26* $ B( LDB, * ), BB( LDB, * ), BP( * ), WORK( * ),
27* $ Z( LDZ, * )
28* ..
29*
30*
31*> \par Purpose:
32* =============
33*>
34*> \verbatim
35*>
36*> ZDRVSG checks the complex Hermitian generalized eigenproblem
37*> drivers.
38*>
39*> ZHEGV computes all eigenvalues and, optionally,
40*> eigenvectors of a complex Hermitian-definite generalized
41*> eigenproblem.
42*>
43*> ZHEGVD computes all eigenvalues and, optionally,
44*> eigenvectors of a complex Hermitian-definite generalized
45*> eigenproblem using a divide and conquer algorithm.
46*>
47*> ZHEGVX computes selected eigenvalues and, optionally,
48*> eigenvectors of a complex Hermitian-definite generalized
49*> eigenproblem.
50*>
51*> ZHPGV computes all eigenvalues and, optionally,
52*> eigenvectors of a complex Hermitian-definite generalized
53*> eigenproblem in packed storage.
54*>
55*> ZHPGVD computes all eigenvalues and, optionally,
56*> eigenvectors of a complex Hermitian-definite generalized
57*> eigenproblem in packed storage using a divide and
58*> conquer algorithm.
59*>
60*> ZHPGVX computes selected eigenvalues and, optionally,
61*> eigenvectors of a complex Hermitian-definite generalized
62*> eigenproblem in packed storage.
63*>
64*> ZHBGV computes all eigenvalues and, optionally,
65*> eigenvectors of a complex Hermitian-definite banded
66*> generalized eigenproblem.
67*>
68*> ZHBGVD computes all eigenvalues and, optionally,
69*> eigenvectors of a complex Hermitian-definite banded
70*> generalized eigenproblem using a divide and conquer
71*> algorithm.
72*>
73*> ZHBGVX computes selected eigenvalues and, optionally,
74*> eigenvectors of a complex Hermitian-definite banded
75*> generalized eigenproblem.
76*>
77*> When ZDRVSG is called, a number of matrix "sizes" ("n's") and a
78*> number of matrix "types" are specified. For each size ("n")
79*> and each type of matrix, one matrix A of the given type will be
80*> generated; a random well-conditioned matrix B is also generated
81*> and the pair (A,B) is used to test the drivers.
82*>
83*> For each pair (A,B), the following tests are performed:
84*>
85*> (1) ZHEGV with ITYPE = 1 and UPLO ='U':
86*>
87*> | A Z - B Z D | / ( |A| |Z| n ulp )
88*>
89*> (2) as (1) but calling ZHPGV
90*> (3) as (1) but calling ZHBGV
91*> (4) as (1) but with UPLO = 'L'
92*> (5) as (4) but calling ZHPGV
93*> (6) as (4) but calling ZHBGV
94*>
95*> (7) ZHEGV with ITYPE = 2 and UPLO ='U':
96*>
97*> | A B Z - Z D | / ( |A| |Z| n ulp )
98*>
99*> (8) as (7) but calling ZHPGV
100*> (9) as (7) but with UPLO = 'L'
101*> (10) as (9) but calling ZHPGV
102*>
103*> (11) ZHEGV with ITYPE = 3 and UPLO ='U':
104*>
105*> | B A Z - Z D | / ( |A| |Z| n ulp )
106*>
107*> (12) as (11) but calling ZHPGV
108*> (13) as (11) but with UPLO = 'L'
109*> (14) as (13) but calling ZHPGV
110*>
111*> ZHEGVD, ZHPGVD and ZHBGVD performed the same 14 tests.
112*>
113*> ZHEGVX, ZHPGVX and ZHBGVX performed the above 14 tests with
114*> the parameter RANGE = 'A', 'N' and 'I', respectively.
115*>
116*> The "sizes" are specified by an array NN(1:NSIZES); the value of
117*> each element NN(j) specifies one size.
118*> The "types" are specified by a logical array DOTYPE( 1:NTYPES );
119*> if DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
120*> This type is used for the matrix A which has half-bandwidth KA.
121*> B is generated as a well-conditioned positive definite matrix
122*> with half-bandwidth KB (<= KA).
123*> Currently, the list of possible types for A is:
124*>
125*> (1) The zero matrix.
126*> (2) The identity matrix.
127*>
128*> (3) A diagonal matrix with evenly spaced entries
129*> 1, ..., ULP and random signs.
130*> (ULP = (first number larger than 1) - 1 )
131*> (4) A diagonal matrix with geometrically spaced entries
132*> 1, ..., ULP and random signs.
133*> (5) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
134*> and random signs.
135*>
136*> (6) Same as (4), but multiplied by SQRT( overflow threshold )
137*> (7) Same as (4), but multiplied by SQRT( underflow threshold )
138*>
139*> (8) A matrix of the form U* D U, where U is unitary and
140*> D has evenly spaced entries 1, ..., ULP with random signs
141*> on the diagonal.
142*>
143*> (9) A matrix of the form U* D U, where U is unitary and
144*> D has geometrically spaced entries 1, ..., ULP with random
145*> signs on the diagonal.
146*>
147*> (10) A matrix of the form U* D U, where U is unitary and
148*> D has "clustered" entries 1, ULP,..., ULP with random
149*> signs on the diagonal.
150*>
151*> (11) Same as (8), but multiplied by SQRT( overflow threshold )
152*> (12) Same as (8), but multiplied by SQRT( underflow threshold )
153*>
154*> (13) Hermitian matrix with random entries chosen from (-1,1).
155*> (14) Same as (13), but multiplied by SQRT( overflow threshold )
156*> (15) Same as (13), but multiplied by SQRT( underflow threshold )
157*>
158*> (16) Same as (8), but with KA = 1 and KB = 1
159*> (17) Same as (8), but with KA = 2 and KB = 1
160*> (18) Same as (8), but with KA = 2 and KB = 2
161*> (19) Same as (8), but with KA = 3 and KB = 1
162*> (20) Same as (8), but with KA = 3 and KB = 2
163*> (21) Same as (8), but with KA = 3 and KB = 3
164*> \endverbatim
165*
166* Arguments:
167* ==========
168*
169*> \verbatim
170*> NSIZES INTEGER
171*> The number of sizes of matrices to use. If it is zero,
172*> ZDRVSG does nothing. It must be at least zero.
173*> Not modified.
174*>
175*> NN INTEGER array, dimension (NSIZES)
176*> An array containing the sizes to be used for the matrices.
177*> Zero values will be skipped. The values must be at least
178*> zero.
179*> Not modified.
180*>
181*> NTYPES INTEGER
182*> The number of elements in DOTYPE. If it is zero, ZDRVSG
183*> does nothing. It must be at least zero. If it is MAXTYP+1
184*> and NSIZES is 1, then an additional type, MAXTYP+1 is
185*> defined, which is to use whatever matrix is in A. This
186*> is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
187*> DOTYPE(MAXTYP+1) is .TRUE. .
188*> Not modified.
189*>
190*> DOTYPE LOGICAL array, dimension (NTYPES)
191*> If DOTYPE(j) is .TRUE., then for each size in NN a
192*> matrix of that size and of type j will be generated.
193*> If NTYPES is smaller than the maximum number of types
194*> defined (PARAMETER MAXTYP), then types NTYPES+1 through
195*> MAXTYP will not be generated. If NTYPES is larger
196*> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
197*> will be ignored.
198*> Not modified.
199*>
200*> ISEED INTEGER array, dimension (4)
201*> On entry ISEED specifies the seed of the random number
202*> generator. The array elements should be between 0 and 4095;
203*> if not they will be reduced mod 4096. Also, ISEED(4) must
204*> be odd. The random number generator uses a linear
205*> congruential sequence limited to small integers, and so
206*> should produce machine independent random numbers. The
207*> values of ISEED are changed on exit, and can be used in the
208*> next call to ZDRVSG to continue the same random number
209*> sequence.
210*> Modified.
211*>
212*> THRESH DOUBLE PRECISION
213*> A test will count as "failed" if the "error", computed as
214*> described above, exceeds THRESH. Note that the error
215*> is scaled to be O(1), so THRESH should be a reasonably
216*> small multiple of 1, e.g., 10 or 100. In particular,
217*> it should not depend on the precision (single vs. double)
218*> or the size of the matrix. It must be at least zero.
219*> Not modified.
220*>
221*> NOUNIT INTEGER
222*> The FORTRAN unit number for printing out error messages
223*> (e.g., if a routine returns IINFO not equal to 0.)
224*> Not modified.
225*>
226*> A COMPLEX*16 array, dimension (LDA , max(NN))
227*> Used to hold the matrix whose eigenvalues are to be
228*> computed. On exit, A contains the last matrix actually
229*> used.
230*> Modified.
231*>
232*> LDA INTEGER
233*> The leading dimension of A. It must be at
234*> least 1 and at least max( NN ).
235*> Not modified.
236*>
237*> B COMPLEX*16 array, dimension (LDB , max(NN))
238*> Used to hold the Hermitian positive definite matrix for
239*> the generalized problem.
240*> On exit, B contains the last matrix actually
241*> used.
242*> Modified.
243*>
244*> LDB INTEGER
245*> The leading dimension of B. It must be at
246*> least 1 and at least max( NN ).
247*> Not modified.
248*>
249*> D DOUBLE PRECISION array, dimension (max(NN))
250*> The eigenvalues of A. On exit, the eigenvalues in D
251*> correspond with the matrix in A.
252*> Modified.
253*>
254*> Z COMPLEX*16 array, dimension (LDZ, max(NN))
255*> The matrix of eigenvectors.
256*> Modified.
257*>
258*> LDZ INTEGER
259*> The leading dimension of ZZ. It must be at least 1 and
260*> at least max( NN ).
261*> Not modified.
262*>
263*> AB COMPLEX*16 array, dimension (LDA, max(NN))
264*> Workspace.
265*> Modified.
266*>
267*> BB COMPLEX*16 array, dimension (LDB, max(NN))
268*> Workspace.
269*> Modified.
270*>
271*> AP COMPLEX*16 array, dimension (max(NN)**2)
272*> Workspace.
273*> Modified.
274*>
275*> BP COMPLEX*16 array, dimension (max(NN)**2)
276*> Workspace.
277*> Modified.
278*>
279*> WORK COMPLEX*16 array, dimension (NWORK)
280*> Workspace.
281*> Modified.
282*>
283*> NWORK INTEGER
284*> The number of entries in WORK. This must be at least
285*> 2*N + N**2 where N = max( NN(j), 2 ).
286*> Not modified.
287*>
288*> RWORK DOUBLE PRECISION array, dimension (LRWORK)
289*> Workspace.
290*> Modified.
291*>
292*> LRWORK INTEGER
293*> The number of entries in RWORK. This must be at least
294*> max( 7*N, 1 + 4*N + 2*N*lg(N) + 3*N**2 ) where
295*> N = max( NN(j) ) and lg( N ) = smallest integer k such
296*> that 2**k >= N .
297*> Not modified.
298*>
299*> IWORK INTEGER array, dimension (LIWORK))
300*> Workspace.
301*> Modified.
302*>
303*> LIWORK INTEGER
304*> The number of entries in IWORK. This must be at least
305*> 2 + 5*max( NN(j) ).
306*> Not modified.
307*>
308*> RESULT DOUBLE PRECISION array, dimension (70)
309*> The values computed by the 70 tests described above.
310*> Modified.
311*>
312*> INFO INTEGER
313*> If 0, then everything ran OK.
314*> -1: NSIZES < 0
315*> -2: Some NN(j) < 0
316*> -3: NTYPES < 0
317*> -5: THRESH < 0
318*> -9: LDA < 1 or LDA < NMAX, where NMAX is max( NN(j) ).
319*> -16: LDZ < 1 or LDZ < NMAX.
320*> -21: NWORK too small.
321*> -23: LRWORK too small.
322*> -25: LIWORK too small.
323*> If ZLATMR, CLATMS, ZHEGV, ZHPGV, ZHBGV, CHEGVD, CHPGVD,
324*> ZHPGVD, ZHEGVX, CHPGVX, ZHBGVX returns an error code,
325*> the absolute value of it is returned.
326*> Modified.
327*>
328*>-----------------------------------------------------------------------
329*>
330*> Some Local Variables and Parameters:
331*> ---- ----- --------- --- ----------
332*> ZERO, ONE Real 0 and 1.
333*> MAXTYP The number of types defined.
334*> NTEST The number of tests that have been run
335*> on this matrix.
336*> NTESTT The total number of tests for this call.
337*> NMAX Largest value in NN.
338*> NMATS The number of matrices generated so far.
339*> NERRS The number of tests which have exceeded THRESH
340*> so far (computed by DLAFTS).
341*> COND, IMODE Values to be passed to the matrix generators.
342*> ANORM Norm of A; passed to matrix generators.
343*>
344*> OVFL, UNFL Overflow and underflow thresholds.
345*> ULP, ULPINV Finest relative precision and its inverse.
346*> RTOVFL, RTUNFL Square roots of the previous 2 values.
347*> The following four arrays decode JTYPE:
348*> KTYPE(j) The general type (1-10) for type "j".
349*> KMODE(j) The MODE value to be passed to the matrix
350*> generator for type "j".
351*> KMAGN(j) The order of magnitude ( O(1),
352*> O(overflow^(1/2) ), O(underflow^(1/2) )
353*> \endverbatim
354*
355* Authors:
356* ========
357*
358*> \author Univ. of Tennessee
359*> \author Univ. of California Berkeley
360*> \author Univ. of Colorado Denver
361*> \author NAG Ltd.
362*
363*> \ingroup complex16_eig
364*
365* =====================================================================
366 SUBROUTINE zdrvsg( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
367 $ NOUNIT, A, LDA, B, LDB, D, Z, LDZ, AB, BB, AP,
368 $ BP, WORK, NWORK, RWORK, LRWORK, IWORK, LIWORK,
369 $ RESULT, INFO )
370*
371* -- LAPACK test routine --
372* -- LAPACK is a software package provided by Univ. of Tennessee, --
373* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
374*
375* .. Scalar Arguments ..
376 INTEGER INFO, LDA, LDB, LDZ, LIWORK, LRWORK, NOUNIT,
377 $ NSIZES, NTYPES, NWORK
378 DOUBLE PRECISION THRESH
379* ..
380* .. Array Arguments ..
381 LOGICAL DOTYPE( * )
382 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
383 DOUBLE PRECISION D( * ), RESULT( * ), RWORK( * )
384 COMPLEX*16 A( LDA, * ), AB( LDA, * ), AP( * ),
385 $ b( ldb, * ), bb( ldb, * ), bp( * ), work( * ),
386 $ z( ldz, * )
387* ..
388*
389* =====================================================================
390*
391* .. Parameters ..
392 DOUBLE PRECISION ZERO, ONE, TEN
393 PARAMETER ( ZERO = 0.0d+0, one = 1.0d+0, ten = 10.0d+0 )
394 COMPLEX*16 CZERO, CONE
395 parameter( czero = ( 0.0d+0, 0.0d+0 ),
396 $ cone = ( 1.0d+0, 0.0d+0 ) )
397 INTEGER MAXTYP
398 parameter( maxtyp = 21 )
399* ..
400* .. Local Scalars ..
401 LOGICAL BADNN
402 CHARACTER UPLO
403 INTEGER I, IBTYPE, IBUPLO, IINFO, IJ, IL, IMODE, ITEMP,
404 $ itype, iu, j, jcol, jsize, jtype, ka, ka9, kb,
405 $ kb9, m, mtypes, n, nerrs, nmats, nmax, ntest,
406 $ ntestt
407 DOUBLE PRECISION ABSTOL, ANINV, ANORM, COND, OVFL, RTOVFL,
408 $ RTUNFL, ULP, ULPINV, UNFL, VL, VU
409* ..
410* .. Local Arrays ..
411 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), ISEED2( 4 ),
412 $ KMAGN( MAXTYP ), KMODE( MAXTYP ),
413 $ KTYPE( MAXTYP )
414* ..
415* .. External Functions ..
416 LOGICAL LSAME
417 DOUBLE PRECISION DLAMCH, DLARND
418 EXTERNAL LSAME, DLAMCH, DLARND
419* ..
420* .. External Subroutines ..
421 EXTERNAL dlafts, dlasum, xerbla, zhbgv, zhbgvd,
424* ..
425* .. Intrinsic Functions ..
426 INTRINSIC abs, dble, max, min, sqrt
427* ..
428* .. Data statements ..
429 DATA ktype / 1, 2, 5*4, 5*5, 3*8, 6*9 /
430 DATA kmagn / 2*1, 1, 1, 1, 2, 3, 1, 1, 1, 2, 3, 1,
431 $ 2, 3, 6*1 /
432 DATA kmode / 2*0, 4, 3, 1, 4, 4, 4, 3, 1, 4, 4, 0,
433 $ 0, 0, 6*4 /
434* ..
435* .. Executable Statements ..
436*
437* 1) Check for errors
438*
439 ntestt = 0
440 info = 0
441*
442 badnn = .false.
443 nmax = 0
444 DO 10 j = 1, nsizes
445 nmax = max( nmax, nn( j ) )
446 IF( nn( j ).LT.0 )
447 $ badnn = .true.
448 10 CONTINUE
449*
450* Check for errors
451*
452 IF( nsizes.LT.0 ) THEN
453 info = -1
454 ELSE IF( badnn ) THEN
455 info = -2
456 ELSE IF( ntypes.LT.0 ) THEN
457 info = -3
458 ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
459 info = -9
460 ELSE IF( ldz.LE.1 .OR. ldz.LT.nmax ) THEN
461 info = -16
462 ELSE IF( 2*max( nmax, 2 )**2.GT.nwork ) THEN
463 info = -21
464 ELSE IF( 2*max( nmax, 2 )**2.GT.lrwork ) THEN
465 info = -23
466 ELSE IF( 2*max( nmax, 2 )**2.GT.liwork ) THEN
467 info = -25
468 END IF
469*
470 IF( info.NE.0 ) THEN
471 CALL xerbla( 'ZDRVSG', -info )
472 RETURN
473 END IF
474*
475* Quick return if possible
476*
477 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
478 $ RETURN
479*
480* More Important constants
481*
482 unfl = dlamch( 'Safe minimum' )
483 ovfl = dlamch( 'Overflow' )
484 ulp = dlamch( 'Epsilon' )*dlamch( 'Base' )
485 ulpinv = one / ulp
486 rtunfl = sqrt( unfl )
487 rtovfl = sqrt( ovfl )
488*
489 DO 20 i = 1, 4
490 iseed2( i ) = iseed( i )
491 20 CONTINUE
492*
493* Loop over sizes, types
494*
495 nerrs = 0
496 nmats = 0
497*
498 DO 650 jsize = 1, nsizes
499 n = nn( jsize )
500 aninv = one / dble( max( 1, n ) )
501*
502 IF( nsizes.NE.1 ) THEN
503 mtypes = min( maxtyp, ntypes )
504 ELSE
505 mtypes = min( maxtyp+1, ntypes )
506 END IF
507*
508 ka9 = 0
509 kb9 = 0
510 DO 640 jtype = 1, mtypes
511 IF( .NOT.dotype( jtype ) )
512 $ GO TO 640
513 nmats = nmats + 1
514 ntest = 0
515*
516 DO 30 j = 1, 4
517 ioldsd( j ) = iseed( j )
518 30 CONTINUE
519*
520* 2) Compute "A"
521*
522* Control parameters:
523*
524* KMAGN KMODE KTYPE
525* =1 O(1) clustered 1 zero
526* =2 large clustered 2 identity
527* =3 small exponential (none)
528* =4 arithmetic diagonal, w/ eigenvalues
529* =5 random log hermitian, w/ eigenvalues
530* =6 random (none)
531* =7 random diagonal
532* =8 random hermitian
533* =9 banded, w/ eigenvalues
534*
535 IF( mtypes.GT.maxtyp )
536 $ GO TO 90
537*
538 itype = ktype( jtype )
539 imode = kmode( jtype )
540*
541* Compute norm
542*
543 GO TO ( 40, 50, 60 )kmagn( jtype )
544*
545 40 CONTINUE
546 anorm = one
547 GO TO 70
548*
549 50 CONTINUE
550 anorm = ( rtovfl*ulp )*aninv
551 GO TO 70
552*
553 60 CONTINUE
554 anorm = rtunfl*n*ulpinv
555 GO TO 70
556*
557 70 CONTINUE
558*
559 iinfo = 0
560 cond = ulpinv
561*
562* Special Matrices -- Identity & Jordan block
563*
564 IF( itype.EQ.1 ) THEN
565*
566* Zero
567*
568 ka = 0
569 kb = 0
570 CALL zlaset( 'Full', lda, n, czero, czero, a, lda )
571*
572 ELSE IF( itype.EQ.2 ) THEN
573*
574* Identity
575*
576 ka = 0
577 kb = 0
578 CALL zlaset( 'Full', lda, n, czero, czero, a, lda )
579 DO 80 jcol = 1, n
580 a( jcol, jcol ) = anorm
581 80 CONTINUE
582*
583 ELSE IF( itype.EQ.4 ) THEN
584*
585* Diagonal Matrix, [Eigen]values Specified
586*
587 ka = 0
588 kb = 0
589 CALL zlatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
590 $ anorm, 0, 0, 'N', a, lda, work, iinfo )
591*
592 ELSE IF( itype.EQ.5 ) THEN
593*
594* Hermitian, eigenvalues specified
595*
596 ka = max( 0, n-1 )
597 kb = ka
598 CALL zlatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
599 $ anorm, n, n, 'N', a, lda, work, iinfo )
600*
601 ELSE IF( itype.EQ.7 ) THEN
602*
603* Diagonal, random eigenvalues
604*
605 ka = 0
606 kb = 0
607 CALL zlatmr( n, n, 'S', iseed, 'H', work, 6, one, cone,
608 $ 'T', 'N', work( n+1 ), 1, one,
609 $ work( 2*n+1 ), 1, one, 'N', idumma, 0, 0,
610 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
611*
612 ELSE IF( itype.EQ.8 ) THEN
613*
614* Hermitian, random eigenvalues
615*
616 ka = max( 0, n-1 )
617 kb = ka
618 CALL zlatmr( n, n, 'S', iseed, 'H', work, 6, one, cone,
619 $ 'T', 'N', work( n+1 ), 1, one,
620 $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
621 $ zero, anorm, 'NO', a, lda, iwork, iinfo )
622*
623 ELSE IF( itype.EQ.9 ) THEN
624*
625* Hermitian banded, eigenvalues specified
626*
627* The following values are used for the half-bandwidths:
628*
629* ka = 1 kb = 1
630* ka = 2 kb = 1
631* ka = 2 kb = 2
632* ka = 3 kb = 1
633* ka = 3 kb = 2
634* ka = 3 kb = 3
635*
636 kb9 = kb9 + 1
637 IF( kb9.GT.ka9 ) THEN
638 ka9 = ka9 + 1
639 kb9 = 1
640 END IF
641 ka = max( 0, min( n-1, ka9 ) )
642 kb = max( 0, min( n-1, kb9 ) )
643 CALL zlatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
644 $ anorm, ka, ka, 'N', a, lda, work, iinfo )
645*
646 ELSE
647*
648 iinfo = 1
649 END IF
650*
651 IF( iinfo.NE.0 ) THEN
652 WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
653 $ ioldsd
654 info = abs( iinfo )
655 RETURN
656 END IF
657*
658 90 CONTINUE
659*
660 abstol = unfl + unfl
661 IF( n.LE.1 ) THEN
662 il = 1
663 iu = n
664 ELSE
665 il = 1 + int( ( n-1 )*dlarnd( 1, iseed2 ) )
666 iu = 1 + int( ( n-1 )*dlarnd( 1, iseed2 ) )
667 IF( il.GT.iu ) THEN
668 itemp = il
669 il = iu
670 iu = itemp
671 END IF
672 END IF
673*
674* 3) Call ZHEGV, ZHPGV, ZHBGV, CHEGVD, CHPGVD, CHBGVD,
675* ZHEGVX, ZHPGVX and ZHBGVX, do tests.
676*
677* loop over the three generalized problems
678* IBTYPE = 1: A*x = (lambda)*B*x
679* IBTYPE = 2: A*B*x = (lambda)*x
680* IBTYPE = 3: B*A*x = (lambda)*x
681*
682 DO 630 ibtype = 1, 3
683*
684* loop over the setting UPLO
685*
686 DO 620 ibuplo = 1, 2
687 IF( ibuplo.EQ.1 )
688 $ uplo = 'U'
689 IF( ibuplo.EQ.2 )
690 $ uplo = 'L'
691*
692* Generate random well-conditioned positive definite
693* matrix B, of bandwidth not greater than that of A.
694*
695 CALL zlatms( n, n, 'U', iseed, 'P', rwork, 5, ten,
696 $ one, kb, kb, uplo, b, ldb, work( n+1 ),
697 $ iinfo )
698*
699* Test ZHEGV
700*
701 ntest = ntest + 1
702*
703 CALL zlacpy( ' ', n, n, a, lda, z, ldz )
704 CALL zlacpy( uplo, n, n, b, ldb, bb, ldb )
705*
706 CALL zhegv( ibtype, 'V', uplo, n, z, ldz, bb, ldb, d,
707 $ work, nwork, rwork, iinfo )
708 IF( iinfo.NE.0 ) THEN
709 WRITE( nounit, fmt = 9999 )'ZHEGV(V,' // uplo //
710 $ ')', iinfo, n, jtype, ioldsd
711 info = abs( iinfo )
712 IF( iinfo.LT.0 ) THEN
713 RETURN
714 ELSE
715 result( ntest ) = ulpinv
716 GO TO 100
717 END IF
718 END IF
719*
720* Do Test
721*
722 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
723 $ ldz, d, work, rwork, result( ntest ) )
724*
725* Test ZHEGVD
726*
727 ntest = ntest + 1
728*
729 CALL zlacpy( ' ', n, n, a, lda, z, ldz )
730 CALL zlacpy( uplo, n, n, b, ldb, bb, ldb )
731*
732 CALL zhegvd( ibtype, 'V', uplo, n, z, ldz, bb, ldb, d,
733 $ work, nwork, rwork, lrwork, iwork,
734 $ liwork, iinfo )
735 IF( iinfo.NE.0 ) THEN
736 WRITE( nounit, fmt = 9999 )'ZHEGVD(V,' // uplo //
737 $ ')', iinfo, n, jtype, ioldsd
738 info = abs( iinfo )
739 IF( iinfo.LT.0 ) THEN
740 RETURN
741 ELSE
742 result( ntest ) = ulpinv
743 GO TO 100
744 END IF
745 END IF
746*
747* Do Test
748*
749 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
750 $ ldz, d, work, rwork, result( ntest ) )
751*
752* Test ZHEGVX
753*
754 ntest = ntest + 1
755*
756 CALL zlacpy( ' ', n, n, a, lda, ab, lda )
757 CALL zlacpy( uplo, n, n, b, ldb, bb, ldb )
758*
759 CALL zhegvx( ibtype, 'V', 'A', uplo, n, ab, lda, bb,
760 $ ldb, vl, vu, il, iu, abstol, m, d, z,
761 $ ldz, work, nwork, rwork, iwork( n+1 ),
762 $ iwork, iinfo )
763 IF( iinfo.NE.0 ) THEN
764 WRITE( nounit, fmt = 9999 )'ZHEGVX(V,A' // uplo //
765 $ ')', iinfo, n, jtype, ioldsd
766 info = abs( iinfo )
767 IF( iinfo.LT.0 ) THEN
768 RETURN
769 ELSE
770 result( ntest ) = ulpinv
771 GO TO 100
772 END IF
773 END IF
774*
775* Do Test
776*
777 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
778 $ ldz, d, work, rwork, result( ntest ) )
779*
780 ntest = ntest + 1
781*
782 CALL zlacpy( ' ', n, n, a, lda, ab, lda )
783 CALL zlacpy( uplo, n, n, b, ldb, bb, ldb )
784*
785* since we do not know the exact eigenvalues of this
786* eigenpair, we just set VL and VU as constants.
787* It is quite possible that there are no eigenvalues
788* in this interval.
789*
790 vl = zero
791 vu = anorm
792 CALL zhegvx( ibtype, 'V', 'V', uplo, n, ab, lda, bb,
793 $ ldb, vl, vu, il, iu, abstol, m, d, z,
794 $ ldz, work, nwork, rwork, iwork( n+1 ),
795 $ iwork, iinfo )
796 IF( iinfo.NE.0 ) THEN
797 WRITE( nounit, fmt = 9999 )'ZHEGVX(V,V,' //
798 $ uplo // ')', iinfo, n, jtype, ioldsd
799 info = abs( iinfo )
800 IF( iinfo.LT.0 ) THEN
801 RETURN
802 ELSE
803 result( ntest ) = ulpinv
804 GO TO 100
805 END IF
806 END IF
807*
808* Do Test
809*
810 CALL zsgt01( ibtype, uplo, n, m, a, lda, b, ldb, z,
811 $ ldz, d, work, rwork, result( ntest ) )
812*
813 ntest = ntest + 1
814*
815 CALL zlacpy( ' ', n, n, a, lda, ab, lda )
816 CALL zlacpy( uplo, n, n, b, ldb, bb, ldb )
817*
818 CALL zhegvx( ibtype, 'V', 'I', uplo, n, ab, lda, bb,
819 $ ldb, vl, vu, il, iu, abstol, m, d, z,
820 $ ldz, work, nwork, rwork, iwork( n+1 ),
821 $ iwork, iinfo )
822 IF( iinfo.NE.0 ) THEN
823 WRITE( nounit, fmt = 9999 )'ZHEGVX(V,I,' //
824 $ uplo // ')', iinfo, n, jtype, ioldsd
825 info = abs( iinfo )
826 IF( iinfo.LT.0 ) THEN
827 RETURN
828 ELSE
829 result( ntest ) = ulpinv
830 GO TO 100
831 END IF
832 END IF
833*
834* Do Test
835*
836 CALL zsgt01( ibtype, uplo, n, m, a, lda, b, ldb, z,
837 $ ldz, d, work, rwork, result( ntest ) )
838*
839 100 CONTINUE
840*
841* Test ZHPGV
842*
843 ntest = ntest + 1
844*
845* Copy the matrices into packed storage.
846*
847 IF( lsame( uplo, 'U' ) ) THEN
848 ij = 1
849 DO 120 j = 1, n
850 DO 110 i = 1, j
851 ap( ij ) = a( i, j )
852 bp( ij ) = b( i, j )
853 ij = ij + 1
854 110 CONTINUE
855 120 CONTINUE
856 ELSE
857 ij = 1
858 DO 140 j = 1, n
859 DO 130 i = j, n
860 ap( ij ) = a( i, j )
861 bp( ij ) = b( i, j )
862 ij = ij + 1
863 130 CONTINUE
864 140 CONTINUE
865 END IF
866*
867 CALL zhpgv( ibtype, 'V', uplo, n, ap, bp, d, z, ldz,
868 $ work, rwork, iinfo )
869 IF( iinfo.NE.0 ) THEN
870 WRITE( nounit, fmt = 9999 )'ZHPGV(V,' // uplo //
871 $ ')', iinfo, n, jtype, ioldsd
872 info = abs( iinfo )
873 IF( iinfo.LT.0 ) THEN
874 RETURN
875 ELSE
876 result( ntest ) = ulpinv
877 GO TO 310
878 END IF
879 END IF
880*
881* Do Test
882*
883 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
884 $ ldz, d, work, rwork, result( ntest ) )
885*
886* Test ZHPGVD
887*
888 ntest = ntest + 1
889*
890* Copy the matrices into packed storage.
891*
892 IF( lsame( uplo, 'U' ) ) THEN
893 ij = 1
894 DO 160 j = 1, n
895 DO 150 i = 1, j
896 ap( ij ) = a( i, j )
897 bp( ij ) = b( i, j )
898 ij = ij + 1
899 150 CONTINUE
900 160 CONTINUE
901 ELSE
902 ij = 1
903 DO 180 j = 1, n
904 DO 170 i = j, n
905 ap( ij ) = a( i, j )
906 bp( ij ) = b( i, j )
907 ij = ij + 1
908 170 CONTINUE
909 180 CONTINUE
910 END IF
911*
912 CALL zhpgvd( ibtype, 'V', uplo, n, ap, bp, d, z, ldz,
913 $ work, nwork, rwork, lrwork, iwork,
914 $ liwork, iinfo )
915 IF( iinfo.NE.0 ) THEN
916 WRITE( nounit, fmt = 9999 )'ZHPGVD(V,' // uplo //
917 $ ')', iinfo, n, jtype, ioldsd
918 info = abs( iinfo )
919 IF( iinfo.LT.0 ) THEN
920 RETURN
921 ELSE
922 result( ntest ) = ulpinv
923 GO TO 310
924 END IF
925 END IF
926*
927* Do Test
928*
929 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
930 $ ldz, d, work, rwork, result( ntest ) )
931*
932* Test ZHPGVX
933*
934 ntest = ntest + 1
935*
936* Copy the matrices into packed storage.
937*
938 IF( lsame( uplo, 'U' ) ) THEN
939 ij = 1
940 DO 200 j = 1, n
941 DO 190 i = 1, j
942 ap( ij ) = a( i, j )
943 bp( ij ) = b( i, j )
944 ij = ij + 1
945 190 CONTINUE
946 200 CONTINUE
947 ELSE
948 ij = 1
949 DO 220 j = 1, n
950 DO 210 i = j, n
951 ap( ij ) = a( i, j )
952 bp( ij ) = b( i, j )
953 ij = ij + 1
954 210 CONTINUE
955 220 CONTINUE
956 END IF
957*
958 CALL zhpgvx( ibtype, 'V', 'A', uplo, n, ap, bp, vl,
959 $ vu, il, iu, abstol, m, d, z, ldz, work,
960 $ rwork, iwork( n+1 ), iwork, info )
961 IF( iinfo.NE.0 ) THEN
962 WRITE( nounit, fmt = 9999 )'ZHPGVX(V,A' // uplo //
963 $ ')', iinfo, n, jtype, ioldsd
964 info = abs( iinfo )
965 IF( iinfo.LT.0 ) THEN
966 RETURN
967 ELSE
968 result( ntest ) = ulpinv
969 GO TO 310
970 END IF
971 END IF
972*
973* Do Test
974*
975 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
976 $ ldz, d, work, rwork, result( ntest ) )
977*
978 ntest = ntest + 1
979*
980* Copy the matrices into packed storage.
981*
982 IF( lsame( uplo, 'U' ) ) THEN
983 ij = 1
984 DO 240 j = 1, n
985 DO 230 i = 1, j
986 ap( ij ) = a( i, j )
987 bp( ij ) = b( i, j )
988 ij = ij + 1
989 230 CONTINUE
990 240 CONTINUE
991 ELSE
992 ij = 1
993 DO 260 j = 1, n
994 DO 250 i = j, n
995 ap( ij ) = a( i, j )
996 bp( ij ) = b( i, j )
997 ij = ij + 1
998 250 CONTINUE
999 260 CONTINUE
1000 END IF
1001*
1002 vl = zero
1003 vu = anorm
1004 CALL zhpgvx( ibtype, 'V', 'V', uplo, n, ap, bp, vl,
1005 $ vu, il, iu, abstol, m, d, z, ldz, work,
1006 $ rwork, iwork( n+1 ), iwork, info )
1007 IF( iinfo.NE.0 ) THEN
1008 WRITE( nounit, fmt = 9999 )'ZHPGVX(V,V' // uplo //
1009 $ ')', iinfo, n, jtype, ioldsd
1010 info = abs( iinfo )
1011 IF( iinfo.LT.0 ) THEN
1012 RETURN
1013 ELSE
1014 result( ntest ) = ulpinv
1015 GO TO 310
1016 END IF
1017 END IF
1018*
1019* Do Test
1020*
1021 CALL zsgt01( ibtype, uplo, n, m, a, lda, b, ldb, z,
1022 $ ldz, d, work, rwork, result( ntest ) )
1023*
1024 ntest = ntest + 1
1025*
1026* Copy the matrices into packed storage.
1027*
1028 IF( lsame( uplo, 'U' ) ) THEN
1029 ij = 1
1030 DO 280 j = 1, n
1031 DO 270 i = 1, j
1032 ap( ij ) = a( i, j )
1033 bp( ij ) = b( i, j )
1034 ij = ij + 1
1035 270 CONTINUE
1036 280 CONTINUE
1037 ELSE
1038 ij = 1
1039 DO 300 j = 1, n
1040 DO 290 i = j, n
1041 ap( ij ) = a( i, j )
1042 bp( ij ) = b( i, j )
1043 ij = ij + 1
1044 290 CONTINUE
1045 300 CONTINUE
1046 END IF
1047*
1048 CALL zhpgvx( ibtype, 'V', 'I', uplo, n, ap, bp, vl,
1049 $ vu, il, iu, abstol, m, d, z, ldz, work,
1050 $ rwork, iwork( n+1 ), iwork, info )
1051 IF( iinfo.NE.0 ) THEN
1052 WRITE( nounit, fmt = 9999 )'ZHPGVX(V,I' // uplo //
1053 $ ')', iinfo, n, jtype, ioldsd
1054 info = abs( iinfo )
1055 IF( iinfo.LT.0 ) THEN
1056 RETURN
1057 ELSE
1058 result( ntest ) = ulpinv
1059 GO TO 310
1060 END IF
1061 END IF
1062*
1063* Do Test
1064*
1065 CALL zsgt01( ibtype, uplo, n, m, a, lda, b, ldb, z,
1066 $ ldz, d, work, rwork, result( ntest ) )
1067*
1068 310 CONTINUE
1069*
1070 IF( ibtype.EQ.1 ) THEN
1071*
1072* TEST ZHBGV
1073*
1074 ntest = ntest + 1
1075*
1076* Copy the matrices into band storage.
1077*
1078 IF( lsame( uplo, 'U' ) ) THEN
1079 DO 340 j = 1, n
1080 DO 320 i = max( 1, j-ka ), j
1081 ab( ka+1+i-j, j ) = a( i, j )
1082 320 CONTINUE
1083 DO 330 i = max( 1, j-kb ), j
1084 bb( kb+1+i-j, j ) = b( i, j )
1085 330 CONTINUE
1086 340 CONTINUE
1087 ELSE
1088 DO 370 j = 1, n
1089 DO 350 i = j, min( n, j+ka )
1090 ab( 1+i-j, j ) = a( i, j )
1091 350 CONTINUE
1092 DO 360 i = j, min( n, j+kb )
1093 bb( 1+i-j, j ) = b( i, j )
1094 360 CONTINUE
1095 370 CONTINUE
1096 END IF
1097*
1098 CALL zhbgv( 'V', uplo, n, ka, kb, ab, lda, bb, ldb,
1099 $ d, z, ldz, work, rwork, iinfo )
1100 IF( iinfo.NE.0 ) THEN
1101 WRITE( nounit, fmt = 9999 )'ZHBGV(V,' //
1102 $ uplo // ')', iinfo, n, jtype, ioldsd
1103 info = abs( iinfo )
1104 IF( iinfo.LT.0 ) THEN
1105 RETURN
1106 ELSE
1107 result( ntest ) = ulpinv
1108 GO TO 620
1109 END IF
1110 END IF
1111*
1112* Do Test
1113*
1114 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
1115 $ ldz, d, work, rwork, result( ntest ) )
1116*
1117* TEST ZHBGVD
1118*
1119 ntest = ntest + 1
1120*
1121* Copy the matrices into band storage.
1122*
1123 IF( lsame( uplo, 'U' ) ) THEN
1124 DO 400 j = 1, n
1125 DO 380 i = max( 1, j-ka ), j
1126 ab( ka+1+i-j, j ) = a( i, j )
1127 380 CONTINUE
1128 DO 390 i = max( 1, j-kb ), j
1129 bb( kb+1+i-j, j ) = b( i, j )
1130 390 CONTINUE
1131 400 CONTINUE
1132 ELSE
1133 DO 430 j = 1, n
1134 DO 410 i = j, min( n, j+ka )
1135 ab( 1+i-j, j ) = a( i, j )
1136 410 CONTINUE
1137 DO 420 i = j, min( n, j+kb )
1138 bb( 1+i-j, j ) = b( i, j )
1139 420 CONTINUE
1140 430 CONTINUE
1141 END IF
1142*
1143 CALL zhbgvd( 'V', uplo, n, ka, kb, ab, lda, bb,
1144 $ ldb, d, z, ldz, work, nwork, rwork,
1145 $ lrwork, iwork, liwork, iinfo )
1146 IF( iinfo.NE.0 ) THEN
1147 WRITE( nounit, fmt = 9999 )'ZHBGVD(V,' //
1148 $ uplo // ')', iinfo, n, jtype, ioldsd
1149 info = abs( iinfo )
1150 IF( iinfo.LT.0 ) THEN
1151 RETURN
1152 ELSE
1153 result( ntest ) = ulpinv
1154 GO TO 620
1155 END IF
1156 END IF
1157*
1158* Do Test
1159*
1160 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
1161 $ ldz, d, work, rwork, result( ntest ) )
1162*
1163* Test ZHBGVX
1164*
1165 ntest = ntest + 1
1166*
1167* Copy the matrices into band storage.
1168*
1169 IF( lsame( uplo, 'U' ) ) THEN
1170 DO 460 j = 1, n
1171 DO 440 i = max( 1, j-ka ), j
1172 ab( ka+1+i-j, j ) = a( i, j )
1173 440 CONTINUE
1174 DO 450 i = max( 1, j-kb ), j
1175 bb( kb+1+i-j, j ) = b( i, j )
1176 450 CONTINUE
1177 460 CONTINUE
1178 ELSE
1179 DO 490 j = 1, n
1180 DO 470 i = j, min( n, j+ka )
1181 ab( 1+i-j, j ) = a( i, j )
1182 470 CONTINUE
1183 DO 480 i = j, min( n, j+kb )
1184 bb( 1+i-j, j ) = b( i, j )
1185 480 CONTINUE
1186 490 CONTINUE
1187 END IF
1188*
1189 CALL zhbgvx( 'V', 'A', uplo, n, ka, kb, ab, lda,
1190 $ bb, ldb, bp, max( 1, n ), vl, vu, il,
1191 $ iu, abstol, m, d, z, ldz, work, rwork,
1192 $ iwork( n+1 ), iwork, iinfo )
1193 IF( iinfo.NE.0 ) THEN
1194 WRITE( nounit, fmt = 9999 )'ZHBGVX(V,A' //
1195 $ uplo // ')', iinfo, n, jtype, ioldsd
1196 info = abs( iinfo )
1197 IF( iinfo.LT.0 ) THEN
1198 RETURN
1199 ELSE
1200 result( ntest ) = ulpinv
1201 GO TO 620
1202 END IF
1203 END IF
1204*
1205* Do Test
1206*
1207 CALL zsgt01( ibtype, uplo, n, n, a, lda, b, ldb, z,
1208 $ ldz, d, work, rwork, result( ntest ) )
1209*
1210 ntest = ntest + 1
1211*
1212* Copy the matrices into band storage.
1213*
1214 IF( lsame( uplo, 'U' ) ) THEN
1215 DO 520 j = 1, n
1216 DO 500 i = max( 1, j-ka ), j
1217 ab( ka+1+i-j, j ) = a( i, j )
1218 500 CONTINUE
1219 DO 510 i = max( 1, j-kb ), j
1220 bb( kb+1+i-j, j ) = b( i, j )
1221 510 CONTINUE
1222 520 CONTINUE
1223 ELSE
1224 DO 550 j = 1, n
1225 DO 530 i = j, min( n, j+ka )
1226 ab( 1+i-j, j ) = a( i, j )
1227 530 CONTINUE
1228 DO 540 i = j, min( n, j+kb )
1229 bb( 1+i-j, j ) = b( i, j )
1230 540 CONTINUE
1231 550 CONTINUE
1232 END IF
1233*
1234 vl = zero
1235 vu = anorm
1236 CALL zhbgvx( 'V', 'V', uplo, n, ka, kb, ab, lda,
1237 $ bb, ldb, bp, max( 1, n ), vl, vu, il,
1238 $ iu, abstol, m, d, z, ldz, work, rwork,
1239 $ iwork( n+1 ), iwork, iinfo )
1240 IF( iinfo.NE.0 ) THEN
1241 WRITE( nounit, fmt = 9999 )'ZHBGVX(V,V' //
1242 $ uplo // ')', iinfo, n, jtype, ioldsd
1243 info = abs( iinfo )
1244 IF( iinfo.LT.0 ) THEN
1245 RETURN
1246 ELSE
1247 result( ntest ) = ulpinv
1248 GO TO 620
1249 END IF
1250 END IF
1251*
1252* Do Test
1253*
1254 CALL zsgt01( ibtype, uplo, n, m, a, lda, b, ldb, z,
1255 $ ldz, d, work, rwork, result( ntest ) )
1256*
1257 ntest = ntest + 1
1258*
1259* Copy the matrices into band storage.
1260*
1261 IF( lsame( uplo, 'U' ) ) THEN
1262 DO 580 j = 1, n
1263 DO 560 i = max( 1, j-ka ), j
1264 ab( ka+1+i-j, j ) = a( i, j )
1265 560 CONTINUE
1266 DO 570 i = max( 1, j-kb ), j
1267 bb( kb+1+i-j, j ) = b( i, j )
1268 570 CONTINUE
1269 580 CONTINUE
1270 ELSE
1271 DO 610 j = 1, n
1272 DO 590 i = j, min( n, j+ka )
1273 ab( 1+i-j, j ) = a( i, j )
1274 590 CONTINUE
1275 DO 600 i = j, min( n, j+kb )
1276 bb( 1+i-j, j ) = b( i, j )
1277 600 CONTINUE
1278 610 CONTINUE
1279 END IF
1280*
1281 CALL zhbgvx( 'V', 'I', uplo, n, ka, kb, ab, lda,
1282 $ bb, ldb, bp, max( 1, n ), vl, vu, il,
1283 $ iu, abstol, m, d, z, ldz, work, rwork,
1284 $ iwork( n+1 ), iwork, iinfo )
1285 IF( iinfo.NE.0 ) THEN
1286 WRITE( nounit, fmt = 9999 )'ZHBGVX(V,I' //
1287 $ uplo // ')', iinfo, n, jtype, ioldsd
1288 info = abs( iinfo )
1289 IF( iinfo.LT.0 ) THEN
1290 RETURN
1291 ELSE
1292 result( ntest ) = ulpinv
1293 GO TO 620
1294 END IF
1295 END IF
1296*
1297* Do Test
1298*
1299 CALL zsgt01( ibtype, uplo, n, m, a, lda, b, ldb, z,
1300 $ ldz, d, work, rwork, result( ntest ) )
1301*
1302 END IF
1303*
1304 620 CONTINUE
1305 630 CONTINUE
1306*
1307* End of Loop -- Check for RESULT(j) > THRESH
1308*
1309 ntestt = ntestt + ntest
1310 CALL dlafts( 'ZSG', n, n, jtype, ntest, result, ioldsd,
1311 $ thresh, nounit, nerrs )
1312 640 CONTINUE
1313 650 CONTINUE
1314*
1315* Summary
1316*
1317 CALL dlasum( 'ZSG', nounit, nerrs, ntestt )
1318*
1319 RETURN
1320*
1321 9999 FORMAT( ' ZDRVSG: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1322 $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1323*
1324* End of ZDRVSG
1325*
1326 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dlafts(type, m, n, imat, ntests, result, iseed, thresh, iounit, ie)
DLAFTS
Definition dlafts.f:99
subroutine dlasum(type, iounit, ie, nrun)
DLASUM
Definition dlasum.f:43
subroutine zhbgv(jobz, uplo, n, ka, kb, ab, ldab, bb, ldbb, w, z, ldz, work, rwork, info)
ZHBGV
Definition zhbgv.f:183
subroutine zhbgvd(jobz, uplo, n, ka, kb, ab, ldab, bb, ldbb, w, z, ldz, work, lwork, rwork, lrwork, iwork, liwork, info)
ZHBGVD
Definition zhbgvd.f:246
subroutine zhbgvx(jobz, range, uplo, n, ka, kb, ab, ldab, bb, ldbb, q, ldq, vl, vu, il, iu, abstol, m, w, z, ldz, work, rwork, iwork, ifail, info)
ZHBGVX
Definition zhbgvx.f:300
subroutine zhegv(itype, jobz, uplo, n, a, lda, b, ldb, w, work, lwork, rwork, info)
ZHEGV
Definition zhegv.f:181
subroutine zhegvd(itype, jobz, uplo, n, a, lda, b, ldb, w, work, lwork, rwork, lrwork, iwork, liwork, info)
ZHEGVD
Definition zhegvd.f:243
subroutine zhegvx(itype, jobz, range, uplo, n, a, lda, b, ldb, vl, vu, il, iu, abstol, m, w, z, ldz, work, lwork, rwork, iwork, ifail, info)
ZHEGVX
Definition zhegvx.f:307
subroutine zhpgv(itype, jobz, uplo, n, ap, bp, w, z, ldz, work, rwork, info)
ZHPGV
Definition zhpgv.f:165
subroutine zhpgvd(itype, jobz, uplo, n, ap, bp, w, z, ldz, work, lwork, rwork, lrwork, iwork, liwork, info)
ZHPGVD
Definition zhpgvd.f:225
subroutine zhpgvx(itype, jobz, range, uplo, n, ap, bp, vl, vu, il, iu, abstol, m, w, z, ldz, work, rwork, iwork, ifail, info)
ZHPGVX
Definition zhpgvx.f:277
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 zdrvsg(nsizes, nn, ntypes, dotype, iseed, thresh, nounit, a, lda, b, ldb, d, z, ldz, ab, bb, ap, bp, work, nwork, rwork, lrwork, iwork, liwork, result, info)
ZDRVSG
Definition zdrvsg.f:370
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 zlatms(m, n, dist, iseed, sym, d, mode, cond, dmax, kl, ku, pack, a, lda, work, info)
ZLATMS
Definition zlatms.f:332
subroutine zsgt01(itype, uplo, n, m, a, lda, b, ldb, z, ldz, d, work, rwork, result)
ZSGT01
Definition zsgt01.f:152