LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zchkgg.f
Go to the documentation of this file.
1*> \brief \b ZCHKGG
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 ZCHKGG( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
12* TSTDIF, THRSHN, NOUNIT, A, LDA, B, H, T, S1,
13* S2, P1, P2, U, LDU, V, Q, Z, ALPHA1, BETA1,
14* ALPHA3, BETA3, EVECTL, EVECTR, WORK, LWORK,
15* RWORK, LLWORK, RESULT, INFO )
16*
17* .. Scalar Arguments ..
18* LOGICAL TSTDIF
19* INTEGER INFO, LDA, LDU, LWORK, NOUNIT, NSIZES, NTYPES
20* DOUBLE PRECISION THRESH, THRSHN
21* ..
22* .. Array Arguments ..
23* LOGICAL DOTYPE( * ), LLWORK( * )
24* INTEGER ISEED( 4 ), NN( * )
25* DOUBLE PRECISION RESULT( 15 ), RWORK( * )
26* COMPLEX*16 A( LDA, * ), ALPHA1( * ), ALPHA3( * ),
27* $ B( LDA, * ), BETA1( * ), BETA3( * ),
28* $ EVECTL( LDU, * ), EVECTR( LDU, * ),
29* $ H( LDA, * ), P1( LDA, * ), P2( LDA, * ),
30* $ Q( LDU, * ), S1( LDA, * ), S2( LDA, * ),
31* $ T( LDA, * ), U( LDU, * ), V( LDU, * ),
32* $ WORK( * ), Z( LDU, * )
33* ..
34*
35*
36*> \par Purpose:
37* =============
38*>
39*> \verbatim
40*>
41*> ZCHKGG checks the nonsymmetric generalized eigenvalue problem
42*> routines.
43*> H H H
44*> ZGGHRD factors A and B as U H V and U T V , where means conjugate
45*> transpose, H is hessenberg, T is triangular and U and V are unitary.
46*>
47*> H H
48*> ZHGEQZ factors H and T as Q S Z and Q P Z , where P and S are upper
49*> triangular and Q and Z are unitary. It also computes the generalized
50*> eigenvalues (alpha(1),beta(1)),...,(alpha(n),beta(n)), where
51*> alpha(j)=S(j,j) and beta(j)=P(j,j) -- thus, w(j) = alpha(j)/beta(j)
52*> is a root of the generalized eigenvalue problem
53*>
54*> det( A - w(j) B ) = 0
55*>
56*> and m(j) = beta(j)/alpha(j) is a root of the essentially equivalent
57*> problem
58*>
59*> det( m(j) A - B ) = 0
60*>
61*> ZTGEVC computes the matrix L of left eigenvectors and the matrix R
62*> of right eigenvectors for the matrix pair ( S, P ). In the
63*> description below, l and r are left and right eigenvectors
64*> corresponding to the generalized eigenvalues (alpha,beta).
65*>
66*> When ZCHKGG is called, a number of matrix "sizes" ("n's") and a
67*> number of matrix "types" are specified. For each size ("n")
68*> and each type of matrix, one matrix will be generated and used
69*> to test the nonsymmetric eigenroutines. For each matrix, 13
70*> tests will be performed. The first twelve "test ratios" should be
71*> small -- O(1). They will be compared with the threshold THRESH:
72*>
73*> H
74*> (1) | A - U H V | / ( |A| n ulp )
75*>
76*> H
77*> (2) | B - U T V | / ( |B| n ulp )
78*>
79*> H
80*> (3) | I - UU | / ( n ulp )
81*>
82*> H
83*> (4) | I - VV | / ( n ulp )
84*>
85*> H
86*> (5) | H - Q S Z | / ( |H| n ulp )
87*>
88*> H
89*> (6) | T - Q P Z | / ( |T| n ulp )
90*>
91*> H
92*> (7) | I - QQ | / ( n ulp )
93*>
94*> H
95*> (8) | I - ZZ | / ( n ulp )
96*>
97*> (9) max over all left eigenvalue/-vector pairs (beta/alpha,l) of
98*> H
99*> | (beta A - alpha B) l | / ( ulp max( |beta A|, |alpha B| ) )
100*>
101*> (10) max over all left eigenvalue/-vector pairs (beta/alpha,l') of
102*> H
103*> | (beta H - alpha T) l' | / ( ulp max( |beta H|, |alpha T| ) )
104*>
105*> where the eigenvectors l' are the result of passing Q to
106*> DTGEVC and back transforming (JOB='B').
107*>
108*> (11) max over all right eigenvalue/-vector pairs (beta/alpha,r) of
109*>
110*> | (beta A - alpha B) r | / ( ulp max( |beta A|, |alpha B| ) )
111*>
112*> (12) max over all right eigenvalue/-vector pairs (beta/alpha,r') of
113*>
114*> | (beta H - alpha T) r' | / ( ulp max( |beta H|, |alpha T| ) )
115*>
116*> where the eigenvectors r' are the result of passing Z to
117*> DTGEVC and back transforming (JOB='B').
118*>
119*> The last three test ratios will usually be small, but there is no
120*> mathematical requirement that they be so. They are therefore
121*> compared with THRESH only if TSTDIF is .TRUE.
122*>
123*> (13) | S(Q,Z computed) - S(Q,Z not computed) | / ( |S| ulp )
124*>
125*> (14) | P(Q,Z computed) - P(Q,Z not computed) | / ( |P| ulp )
126*>
127*> (15) max( |alpha(Q,Z computed) - alpha(Q,Z not computed)|/|S| ,
128*> |beta(Q,Z computed) - beta(Q,Z not computed)|/|P| ) / ulp
129*>
130*> In addition, the normalization of L and R are checked, and compared
131*> with the threshold THRSHN.
132*>
133*> Test Matrices
134*> ---- --------
135*>
136*> The sizes of the test matrices are specified by an array
137*> NN(1:NSIZES); the value of each element NN(j) specifies one size.
138*> The "types" are specified by a logical array DOTYPE( 1:NTYPES ); if
139*> DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
140*> Currently, the list of possible types is:
141*>
142*> (1) ( 0, 0 ) (a pair of zero matrices)
143*>
144*> (2) ( I, 0 ) (an identity and a zero matrix)
145*>
146*> (3) ( 0, I ) (an identity and a zero matrix)
147*>
148*> (4) ( I, I ) (a pair of identity matrices)
149*>
150*> t t
151*> (5) ( J , J ) (a pair of transposed Jordan blocks)
152*>
153*> t ( I 0 )
154*> (6) ( X, Y ) where X = ( J 0 ) and Y = ( t )
155*> ( 0 I ) ( 0 J )
156*> and I is a k x k identity and J a (k+1)x(k+1)
157*> Jordan block; k=(N-1)/2
158*>
159*> (7) ( D, I ) where D is P*D1, P is a random unitary diagonal
160*> matrix (i.e., with random magnitude 1 entries
161*> on the diagonal), and D1=diag( 0, 1,..., N-1 )
162*> (i.e., a diagonal matrix with D1(1,1)=0,
163*> D1(2,2)=1, ..., D1(N,N)=N-1.)
164*> (8) ( I, D )
165*>
166*> (9) ( big*D, small*I ) where "big" is near overflow and small=1/big
167*>
168*> (10) ( small*D, big*I )
169*>
170*> (11) ( big*I, small*D )
171*>
172*> (12) ( small*I, big*D )
173*>
174*> (13) ( big*D, big*I )
175*>
176*> (14) ( small*D, small*I )
177*>
178*> (15) ( D1, D2 ) where D1=P*diag( 0, 0, 1, ..., N-3, 0 ) and
179*> D2=Q*diag( 0, N-3, N-4,..., 1, 0, 0 ), and
180*> P and Q are random unitary diagonal matrices.
181*> t t
182*> (16) U ( J , J ) V where U and V are random unitary matrices.
183*>
184*> (17) U ( T1, T2 ) V where T1 and T2 are upper triangular matrices
185*> with random O(1) entries above the diagonal
186*> and diagonal entries diag(T1) =
187*> P*( 0, 0, 1, ..., N-3, 0 ) and diag(T2) =
188*> Q*( 0, N-3, N-4,..., 1, 0, 0 )
189*>
190*> (18) U ( T1, T2 ) V diag(T1) = ( 0, 0, 1, 1, s, ..., s, 0 )
191*> diag(T2) = ( 0, 1, 0, 1,..., 1, 0 )
192*> s = machine precision.
193*>
194*> (19) U ( T1, T2 ) V diag(T1)=( 0,0,1,1, 1-d, ..., 1-(N-5)*d=s, 0 )
195*> diag(T2) = ( 0, 1, 0, 1, ..., 1, 0 )
196*>
197*> N-5
198*> (20) U ( T1, T2 ) V diag(T1)=( 0, 0, 1, 1, a, ..., a =s, 0 )
199*> diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )
200*>
201*> (21) U ( T1, T2 ) V diag(T1)=( 0, 0, 1, r1, r2, ..., r(N-4), 0 )
202*> diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )
203*> where r1,..., r(N-4) are random.
204*>
205*> (22) U ( big*T1, small*T2 ) V diag(T1) = P*( 0, 0, 1, ..., N-3, 0 )
206*> diag(T2) = ( 0, 1, ..., 1, 0, 0 )
207*>
208*> (23) U ( small*T1, big*T2 ) V diag(T1) = P*( 0, 0, 1, ..., N-3, 0 )
209*> diag(T2) = ( 0, 1, ..., 1, 0, 0 )
210*>
211*> (24) U ( small*T1, small*T2 ) V diag(T1) = P*( 0, 0, 1, ..., N-3, 0 )
212*> diag(T2) = ( 0, 1, ..., 1, 0, 0 )
213*>
214*> (25) U ( big*T1, big*T2 ) V diag(T1) = P*( 0, 0, 1, ..., N-3, 0 )
215*> diag(T2) = ( 0, 1, ..., 1, 0, 0 )
216*>
217*> (26) U ( T1, T2 ) V where T1 and T2 are random upper-triangular
218*> matrices.
219*> \endverbatim
220*
221* Arguments:
222* ==========
223*
224*> \param[in] NSIZES
225*> \verbatim
226*> NSIZES is INTEGER
227*> The number of sizes of matrices to use. If it is zero,
228*> ZCHKGG does nothing. It must be at least zero.
229*> \endverbatim
230*>
231*> \param[in] NN
232*> \verbatim
233*> NN is INTEGER array, dimension (NSIZES)
234*> An array containing the sizes to be used for the matrices.
235*> Zero values will be skipped. The values must be at least
236*> zero.
237*> \endverbatim
238*>
239*> \param[in] NTYPES
240*> \verbatim
241*> NTYPES is INTEGER
242*> The number of elements in DOTYPE. If it is zero, ZCHKGG
243*> does nothing. It must be at least zero. If it is MAXTYP+1
244*> and NSIZES is 1, then an additional type, MAXTYP+1 is
245*> defined, which is to use whatever matrix is in A. This
246*> is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
247*> DOTYPE(MAXTYP+1) is .TRUE. .
248*> \endverbatim
249*>
250*> \param[in] DOTYPE
251*> \verbatim
252*> DOTYPE is LOGICAL array, dimension (NTYPES)
253*> If DOTYPE(j) is .TRUE., then for each size in NN a
254*> matrix of that size and of type j will be generated.
255*> If NTYPES is smaller than the maximum number of types
256*> defined (PARAMETER MAXTYP), then types NTYPES+1 through
257*> MAXTYP will not be generated. If NTYPES is larger
258*> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
259*> will be ignored.
260*> \endverbatim
261*>
262*> \param[in,out] ISEED
263*> \verbatim
264*> ISEED is INTEGER array, dimension (4)
265*> On entry ISEED specifies the seed of the random number
266*> generator. The array elements should be between 0 and 4095;
267*> if not they will be reduced mod 4096. Also, ISEED(4) must
268*> be odd. The random number generator uses a linear
269*> congruential sequence limited to small integers, and so
270*> should produce machine independent random numbers. The
271*> values of ISEED are changed on exit, and can be used in the
272*> next call to ZCHKGG to continue the same random number
273*> sequence.
274*> \endverbatim
275*>
276*> \param[in] THRESH
277*> \verbatim
278*> THRESH is DOUBLE PRECISION
279*> A test will count as "failed" if the "error", computed as
280*> described above, exceeds THRESH. Note that the error
281*> is scaled to be O(1), so THRESH should be a reasonably
282*> small multiple of 1, e.g., 10 or 100. In particular,
283*> it should not depend on the precision (single vs. double)
284*> or the size of the matrix. It must be at least zero.
285*> \endverbatim
286*>
287*> \param[in] TSTDIF
288*> \verbatim
289*> TSTDIF is LOGICAL
290*> Specifies whether test ratios 13-15 will be computed and
291*> compared with THRESH.
292*> = .FALSE.: Only test ratios 1-12 will be computed and tested.
293*> Ratios 13-15 will be set to zero.
294*> = .TRUE.: All the test ratios 1-15 will be computed and
295*> tested.
296*> \endverbatim
297*>
298*> \param[in] THRSHN
299*> \verbatim
300*> THRSHN is DOUBLE PRECISION
301*> Threshold for reporting eigenvector normalization error.
302*> If the normalization of any eigenvector differs from 1 by
303*> more than THRSHN*ulp, then a special error message will be
304*> printed. (This is handled separately from the other tests,
305*> since only a compiler or programming error should cause an
306*> error message, at least if THRSHN is at least 5--10.)
307*> \endverbatim
308*>
309*> \param[in] NOUNIT
310*> \verbatim
311*> NOUNIT is INTEGER
312*> The FORTRAN unit number for printing out error messages
313*> (e.g., if a routine returns IINFO not equal to 0.)
314*> \endverbatim
315*>
316*> \param[in,out] A
317*> \verbatim
318*> A is COMPLEX*16 array, dimension (LDA, max(NN))
319*> Used to hold the original A matrix. Used as input only
320*> if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
321*> DOTYPE(MAXTYP+1)=.TRUE.
322*> \endverbatim
323*>
324*> \param[in] LDA
325*> \verbatim
326*> LDA is INTEGER
327*> The leading dimension of A, B, H, T, S1, P1, S2, and P2.
328*> It must be at least 1 and at least max( NN ).
329*> \endverbatim
330*>
331*> \param[in,out] B
332*> \verbatim
333*> B is COMPLEX*16 array, dimension (LDA, max(NN))
334*> Used to hold the original B matrix. Used as input only
335*> if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
336*> DOTYPE(MAXTYP+1)=.TRUE.
337*> \endverbatim
338*>
339*> \param[out] H
340*> \verbatim
341*> H is COMPLEX*16 array, dimension (LDA, max(NN))
342*> The upper Hessenberg matrix computed from A by ZGGHRD.
343*> \endverbatim
344*>
345*> \param[out] T
346*> \verbatim
347*> T is COMPLEX*16 array, dimension (LDA, max(NN))
348*> The upper triangular matrix computed from B by ZGGHRD.
349*> \endverbatim
350*>
351*> \param[out] S1
352*> \verbatim
353*> S1 is COMPLEX*16 array, dimension (LDA, max(NN))
354*> The Schur (upper triangular) matrix computed from H by ZHGEQZ
355*> when Q and Z are also computed.
356*> \endverbatim
357*>
358*> \param[out] S2
359*> \verbatim
360*> S2 is COMPLEX*16 array, dimension (LDA, max(NN))
361*> The Schur (upper triangular) matrix computed from H by ZHGEQZ
362*> when Q and Z are not computed.
363*> \endverbatim
364*>
365*> \param[out] P1
366*> \verbatim
367*> P1 is COMPLEX*16 array, dimension (LDA, max(NN))
368*> The upper triangular matrix computed from T by ZHGEQZ
369*> when Q and Z are also computed.
370*> \endverbatim
371*>
372*> \param[out] P2
373*> \verbatim
374*> P2 is COMPLEX*16 array, dimension (LDA, max(NN))
375*> The upper triangular matrix computed from T by ZHGEQZ
376*> when Q and Z are not computed.
377*> \endverbatim
378*>
379*> \param[out] U
380*> \verbatim
381*> U is COMPLEX*16 array, dimension (LDU, max(NN))
382*> The (left) unitary matrix computed by ZGGHRD.
383*> \endverbatim
384*>
385*> \param[in] LDU
386*> \verbatim
387*> LDU is INTEGER
388*> The leading dimension of U, V, Q, Z, EVECTL, and EVEZTR. It
389*> must be at least 1 and at least max( NN ).
390*> \endverbatim
391*>
392*> \param[out] V
393*> \verbatim
394*> V is COMPLEX*16 array, dimension (LDU, max(NN))
395*> The (right) unitary matrix computed by ZGGHRD.
396*> \endverbatim
397*>
398*> \param[out] Q
399*> \verbatim
400*> Q is COMPLEX*16 array, dimension (LDU, max(NN))
401*> The (left) unitary matrix computed by ZHGEQZ.
402*> \endverbatim
403*>
404*> \param[out] Z
405*> \verbatim
406*> Z is COMPLEX*16 array, dimension (LDU, max(NN))
407*> The (left) unitary matrix computed by ZHGEQZ.
408*> \endverbatim
409*>
410*> \param[out] ALPHA1
411*> \verbatim
412*> ALPHA1 is COMPLEX*16 array, dimension (max(NN))
413*> \endverbatim
414*>
415*> \param[out] BETA1
416*> \verbatim
417*> BETA1 is COMPLEX*16 array, dimension (max(NN))
418*> The generalized eigenvalues of (A,B) computed by ZHGEQZ
419*> when Q, Z, and the full Schur matrices are computed.
420*> \endverbatim
421*>
422*> \param[out] ALPHA3
423*> \verbatim
424*> ALPHA3 is COMPLEX*16 array, dimension (max(NN))
425*> \endverbatim
426*>
427*> \param[out] BETA3
428*> \verbatim
429*> BETA3 is COMPLEX*16 array, dimension (max(NN))
430*> The generalized eigenvalues of (A,B) computed by ZHGEQZ
431*> when neither Q, Z, nor the Schur matrices are computed.
432*> \endverbatim
433*>
434*> \param[out] EVECTL
435*> \verbatim
436*> EVECTL is COMPLEX*16 array, dimension (LDU, max(NN))
437*> The (lower triangular) left eigenvector matrix for the
438*> matrices in S1 and P1.
439*> \endverbatim
440*>
441*> \param[out] EVECTR
442*> \verbatim
443*> EVECTR is COMPLEX*16 array, dimension (LDU, max(NN))
444*> The (upper triangular) right eigenvector matrix for the
445*> matrices in S1 and P1.
446*> \endverbatim
447*>
448*> \param[out] WORK
449*> \verbatim
450*> WORK is COMPLEX*16 array, dimension (LWORK)
451*> \endverbatim
452*>
453*> \param[in] LWORK
454*> \verbatim
455*> LWORK is INTEGER
456*> The number of entries in WORK. This must be at least
457*> max( 4*N, 2 * N**2, 1 ), for all N=NN(j).
458*> \endverbatim
459*>
460*> \param[out] RWORK
461*> \verbatim
462*> RWORK is DOUBLE PRECISION array, dimension (2*max(NN))
463*> \endverbatim
464*>
465*> \param[out] LLWORK
466*> \verbatim
467*> LLWORK is LOGICAL array, dimension (max(NN))
468*> \endverbatim
469*>
470*> \param[out] RESULT
471*> \verbatim
472*> RESULT is DOUBLE PRECISION array, dimension (15)
473*> The values computed by the tests described above.
474*> The values are currently limited to 1/ulp, to avoid
475*> overflow.
476*> \endverbatim
477*>
478*> \param[out] INFO
479*> \verbatim
480*> INFO is INTEGER
481*> = 0: successful exit.
482*> < 0: if INFO = -i, the i-th argument had an illegal value.
483*> > 0: A routine returned an error code. INFO is the
484*> absolute value of the INFO value returned.
485*> \endverbatim
486*
487* Authors:
488* ========
489*
490*> \author Univ. of Tennessee
491*> \author Univ. of California Berkeley
492*> \author Univ. of Colorado Denver
493*> \author NAG Ltd.
494*
495*> \ingroup complex16_eig
496*
497* =====================================================================
498 SUBROUTINE zchkgg( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
499 $ TSTDIF, THRSHN, NOUNIT, A, LDA, B, H, T, S1,
500 $ S2, P1, P2, U, LDU, V, Q, Z, ALPHA1, BETA1,
501 $ ALPHA3, BETA3, EVECTL, EVECTR, WORK, LWORK,
502 $ RWORK, LLWORK, RESULT, INFO )
503*
504* -- LAPACK test routine --
505* -- LAPACK is a software package provided by Univ. of Tennessee, --
506* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
507*
508* .. Scalar Arguments ..
509 LOGICAL TSTDIF
510 INTEGER INFO, LDA, LDU, LWORK, NOUNIT, NSIZES, NTYPES
511 DOUBLE PRECISION THRESH, THRSHN
512* ..
513* .. Array Arguments ..
514 LOGICAL DOTYPE( * ), LLWORK( * )
515 INTEGER ISEED( 4 ), NN( * )
516 DOUBLE PRECISION RESULT( 15 ), RWORK( * )
517 COMPLEX*16 A( LDA, * ), ALPHA1( * ), ALPHA3( * ),
518 $ b( lda, * ), beta1( * ), beta3( * ),
519 $ evectl( ldu, * ), evectr( ldu, * ),
520 $ h( lda, * ), p1( lda, * ), p2( lda, * ),
521 $ q( ldu, * ), s1( lda, * ), s2( lda, * ),
522 $ t( lda, * ), u( ldu, * ), v( ldu, * ),
523 $ work( * ), z( ldu, * )
524* ..
525*
526* =====================================================================
527*
528* .. Parameters ..
529 DOUBLE PRECISION ZERO, ONE
530 PARAMETER ( ZERO = 0.0d+0, one = 1.0d+0 )
531 COMPLEX*16 CZERO, CONE
532 PARAMETER ( CZERO = ( 0.0d+0, 0.0d+0 ),
533 $ cone = ( 1.0d+0, 0.0d+0 ) )
534 INTEGER MAXTYP
535 parameter( maxtyp = 26 )
536* ..
537* .. Local Scalars ..
538 LOGICAL BADNN
539 INTEGER I1, IADD, IINFO, IN, J, JC, JR, JSIZE, JTYPE,
540 $ LWKOPT, MTYPES, N, N1, NERRS, NMATS, NMAX,
541 $ NTEST, NTESTT
542 DOUBLE PRECISION ANORM, BNORM, SAFMAX, SAFMIN, TEMP1, TEMP2,
543 $ ulp, ulpinv
544 COMPLEX*16 CTEMP
545* ..
546* .. Local Arrays ..
547 LOGICAL LASIGN( MAXTYP ), LBSIGN( MAXTYP )
548 INTEGER IOLDSD( 4 ), KADD( 6 ), KAMAGN( MAXTYP ),
549 $ KATYPE( MAXTYP ), KAZERO( MAXTYP ),
550 $ KBMAGN( MAXTYP ), KBTYPE( MAXTYP ),
551 $ kbzero( maxtyp ), kclass( maxtyp ),
552 $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
553 DOUBLE PRECISION DUMMA( 4 ), RMAGN( 0: 3 )
554 COMPLEX*16 CDUMMA( 4 )
555* ..
556* .. External Functions ..
557 DOUBLE PRECISION DLAMCH, ZLANGE
558 COMPLEX*16 ZLARND
559 EXTERNAL DLAMCH, ZLANGE, ZLARND
560* ..
561* .. External Subroutines ..
562 EXTERNAL dlabad, dlasum, xerbla, zgeqr2, zget51, zget52,
564 $ ztgevc, zunm2r
565* ..
566* .. Intrinsic Functions ..
567 INTRINSIC abs, dble, dconjg, max, min, sign
568* ..
569* .. Data statements ..
570 DATA kclass / 15*1, 10*2, 1*3 /
571 DATA kz1 / 0, 1, 2, 1, 3, 3 /
572 DATA kz2 / 0, 0, 1, 2, 1, 1 /
573 DATA kadd / 0, 0, 0, 0, 3, 2 /
574 DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
575 $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
576 DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
577 $ 1, 1, -4, 2, -4, 8*8, 0 /
578 DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
579 $ 4*5, 4*3, 1 /
580 DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
581 $ 4*6, 4*4, 1 /
582 DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
583 $ 2, 1 /
584 DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
585 $ 2, 1 /
586 DATA ktrian / 16*0, 10*1 /
587 DATA lasign / 6*.false., .true., .false., 2*.true.,
588 $ 2*.false., 3*.true., .false., .true.,
589 $ 3*.false., 5*.true., .false. /
590 DATA lbsign / 7*.false., .true., 2*.false.,
591 $ 2*.true., 2*.false., .true., .false., .true.,
592 $ 9*.false. /
593* ..
594* .. Executable Statements ..
595*
596* Check for errors
597*
598 info = 0
599*
600 badnn = .false.
601 nmax = 1
602 DO 10 j = 1, nsizes
603 nmax = max( nmax, nn( j ) )
604 IF( nn( j ).LT.0 )
605 $ badnn = .true.
606 10 CONTINUE
607*
608 lwkopt = max( 2*nmax*nmax, 4*nmax, 1 )
609*
610* Check for errors
611*
612 IF( nsizes.LT.0 ) THEN
613 info = -1
614 ELSE IF( badnn ) THEN
615 info = -2
616 ELSE IF( ntypes.LT.0 ) THEN
617 info = -3
618 ELSE IF( thresh.LT.zero ) THEN
619 info = -6
620 ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
621 info = -10
622 ELSE IF( ldu.LE.1 .OR. ldu.LT.nmax ) THEN
623 info = -19
624 ELSE IF( lwkopt.GT.lwork ) THEN
625 info = -30
626 END IF
627*
628 IF( info.NE.0 ) THEN
629 CALL xerbla( 'ZCHKGG', -info )
630 RETURN
631 END IF
632*
633* Quick return if possible
634*
635 IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
636 $ RETURN
637*
638 safmin = dlamch( 'Safe minimum' )
639 ulp = dlamch( 'Epsilon' )*dlamch( 'Base' )
640 safmin = safmin / ulp
641 safmax = one / safmin
642 CALL dlabad( safmin, safmax )
643 ulpinv = one / ulp
644*
645* The values RMAGN(2:3) depend on N, see below.
646*
647 rmagn( 0 ) = zero
648 rmagn( 1 ) = one
649*
650* Loop over sizes, types
651*
652 ntestt = 0
653 nerrs = 0
654 nmats = 0
655*
656 DO 240 jsize = 1, nsizes
657 n = nn( jsize )
658 n1 = max( 1, n )
659 rmagn( 2 ) = safmax*ulp / dble( n1 )
660 rmagn( 3 ) = safmin*ulpinv*n1
661*
662 IF( nsizes.NE.1 ) THEN
663 mtypes = min( maxtyp, ntypes )
664 ELSE
665 mtypes = min( maxtyp+1, ntypes )
666 END IF
667*
668 DO 230 jtype = 1, mtypes
669 IF( .NOT.dotype( jtype ) )
670 $ GO TO 230
671 nmats = nmats + 1
672 ntest = 0
673*
674* Save ISEED in case of an error.
675*
676 DO 20 j = 1, 4
677 ioldsd( j ) = iseed( j )
678 20 CONTINUE
679*
680* Initialize RESULT
681*
682 DO 30 j = 1, 15
683 result( j ) = zero
684 30 CONTINUE
685*
686* Compute A and B
687*
688* Description of control parameters:
689*
690* KZLASS: =1 means w/o rotation, =2 means w/ rotation,
691* =3 means random.
692* KATYPE: the "type" to be passed to ZLATM4 for computing A.
693* KAZERO: the pattern of zeros on the diagonal for A:
694* =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ),
695* =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ),
696* =6: ( 0, 1, 0, xxx, 0 ). (xxx means a string of
697* non-zero entries.)
698* KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1),
699* =2: large, =3: small.
700* LASIGN: .TRUE. if the diagonal elements of A are to be
701* multiplied by a random magnitude 1 number.
702* KBTYPE, KBZERO, KBMAGN, LBSIGN: the same, but for B.
703* KTRIAN: =0: don't fill in the upper triangle, =1: do.
704* KZ1, KZ2, KADD: used to implement KAZERO and KBZERO.
705* RMAGN: used to implement KAMAGN and KBMAGN.
706*
707 IF( mtypes.GT.maxtyp )
708 $ GO TO 110
709 iinfo = 0
710 IF( kclass( jtype ).LT.3 ) THEN
711*
712* Generate A (w/o rotation)
713*
714 IF( abs( katype( jtype ) ).EQ.3 ) THEN
715 in = 2*( ( n-1 ) / 2 ) + 1
716 IF( in.NE.n )
717 $ CALL zlaset( 'Full', n, n, czero, czero, a, lda )
718 ELSE
719 in = n
720 END IF
721 CALL zlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
722 $ kz2( kazero( jtype ) ), lasign( jtype ),
723 $ rmagn( kamagn( jtype ) ), ulp,
724 $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 4,
725 $ iseed, a, lda )
726 iadd = kadd( kazero( jtype ) )
727 IF( iadd.GT.0 .AND. iadd.LE.n )
728 $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
729*
730* Generate B (w/o rotation)
731*
732 IF( abs( kbtype( jtype ) ).EQ.3 ) THEN
733 in = 2*( ( n-1 ) / 2 ) + 1
734 IF( in.NE.n )
735 $ CALL zlaset( 'Full', n, n, czero, czero, b, lda )
736 ELSE
737 in = n
738 END IF
739 CALL zlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
740 $ kz2( kbzero( jtype ) ), lbsign( jtype ),
741 $ rmagn( kbmagn( jtype ) ), one,
742 $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 4,
743 $ iseed, b, lda )
744 iadd = kadd( kbzero( jtype ) )
745 IF( iadd.NE.0 )
746 $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
747*
748 IF( kclass( jtype ).EQ.2 .AND. n.GT.0 ) THEN
749*
750* Include rotations
751*
752* Generate U, V as Householder transformations times a
753* diagonal matrix. (Note that ZLARFG makes U(j,j) and
754* V(j,j) real.)
755*
756 DO 50 jc = 1, n - 1
757 DO 40 jr = jc, n
758 u( jr, jc ) = zlarnd( 3, iseed )
759 v( jr, jc ) = zlarnd( 3, iseed )
760 40 CONTINUE
761 CALL zlarfg( n+1-jc, u( jc, jc ), u( jc+1, jc ), 1,
762 $ work( jc ) )
763 work( 2*n+jc ) = sign( one, dble( u( jc, jc ) ) )
764 u( jc, jc ) = cone
765 CALL zlarfg( n+1-jc, v( jc, jc ), v( jc+1, jc ), 1,
766 $ work( n+jc ) )
767 work( 3*n+jc ) = sign( one, dble( v( jc, jc ) ) )
768 v( jc, jc ) = cone
769 50 CONTINUE
770 ctemp = zlarnd( 3, iseed )
771 u( n, n ) = cone
772 work( n ) = czero
773 work( 3*n ) = ctemp / abs( ctemp )
774 ctemp = zlarnd( 3, iseed )
775 v( n, n ) = cone
776 work( 2*n ) = czero
777 work( 4*n ) = ctemp / abs( ctemp )
778*
779* Apply the diagonal matrices
780*
781 DO 70 jc = 1, n
782 DO 60 jr = 1, n
783 a( jr, jc ) = work( 2*n+jr )*
784 $ dconjg( work( 3*n+jc ) )*
785 $ a( jr, jc )
786 b( jr, jc ) = work( 2*n+jr )*
787 $ dconjg( work( 3*n+jc ) )*
788 $ b( jr, jc )
789 60 CONTINUE
790 70 CONTINUE
791 CALL zunm2r( 'L', 'N', n, n, n-1, u, ldu, work, a,
792 $ lda, work( 2*n+1 ), iinfo )
793 IF( iinfo.NE.0 )
794 $ GO TO 100
795 CALL zunm2r( 'R', 'C', n, n, n-1, v, ldu, work( n+1 ),
796 $ a, lda, work( 2*n+1 ), iinfo )
797 IF( iinfo.NE.0 )
798 $ GO TO 100
799 CALL zunm2r( 'L', 'N', n, n, n-1, u, ldu, work, b,
800 $ lda, work( 2*n+1 ), iinfo )
801 IF( iinfo.NE.0 )
802 $ GO TO 100
803 CALL zunm2r( 'R', 'C', n, n, n-1, v, ldu, work( n+1 ),
804 $ b, lda, work( 2*n+1 ), iinfo )
805 IF( iinfo.NE.0 )
806 $ GO TO 100
807 END IF
808 ELSE
809*
810* Random matrices
811*
812 DO 90 jc = 1, n
813 DO 80 jr = 1, n
814 a( jr, jc ) = rmagn( kamagn( jtype ) )*
815 $ zlarnd( 4, iseed )
816 b( jr, jc ) = rmagn( kbmagn( jtype ) )*
817 $ zlarnd( 4, iseed )
818 80 CONTINUE
819 90 CONTINUE
820 END IF
821*
822 anorm = zlange( '1', n, n, a, lda, rwork )
823 bnorm = zlange( '1', n, n, b, lda, rwork )
824*
825 100 CONTINUE
826*
827 IF( iinfo.NE.0 ) THEN
828 WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
829 $ ioldsd
830 info = abs( iinfo )
831 RETURN
832 END IF
833*
834 110 CONTINUE
835*
836* Call ZGEQR2, ZUNM2R, and ZGGHRD to compute H, T, U, and V
837*
838 CALL zlacpy( ' ', n, n, a, lda, h, lda )
839 CALL zlacpy( ' ', n, n, b, lda, t, lda )
840 ntest = 1
841 result( 1 ) = ulpinv
842*
843 CALL zgeqr2( n, n, t, lda, work, work( n+1 ), iinfo )
844 IF( iinfo.NE.0 ) THEN
845 WRITE( nounit, fmt = 9999 )'ZGEQR2', iinfo, n, jtype,
846 $ ioldsd
847 info = abs( iinfo )
848 GO TO 210
849 END IF
850*
851 CALL zunm2r( 'L', 'C', n, n, n, t, lda, work, h, lda,
852 $ work( n+1 ), iinfo )
853 IF( iinfo.NE.0 ) THEN
854 WRITE( nounit, fmt = 9999 )'ZUNM2R', iinfo, n, jtype,
855 $ ioldsd
856 info = abs( iinfo )
857 GO TO 210
858 END IF
859*
860 CALL zlaset( 'Full', n, n, czero, cone, u, ldu )
861 CALL zunm2r( 'R', 'N', n, n, n, t, lda, work, u, ldu,
862 $ work( n+1 ), iinfo )
863 IF( iinfo.NE.0 ) THEN
864 WRITE( nounit, fmt = 9999 )'ZUNM2R', iinfo, n, jtype,
865 $ ioldsd
866 info = abs( iinfo )
867 GO TO 210
868 END IF
869*
870 CALL zgghrd( 'V', 'I', n, 1, n, h, lda, t, lda, u, ldu, v,
871 $ ldu, iinfo )
872 IF( iinfo.NE.0 ) THEN
873 WRITE( nounit, fmt = 9999 )'ZGGHRD', iinfo, n, jtype,
874 $ ioldsd
875 info = abs( iinfo )
876 GO TO 210
877 END IF
878 ntest = 4
879*
880* Do tests 1--4
881*
882 CALL zget51( 1, n, a, lda, h, lda, u, ldu, v, ldu, work,
883 $ rwork, result( 1 ) )
884 CALL zget51( 1, n, b, lda, t, lda, u, ldu, v, ldu, work,
885 $ rwork, result( 2 ) )
886 CALL zget51( 3, n, b, lda, t, lda, u, ldu, u, ldu, work,
887 $ rwork, result( 3 ) )
888 CALL zget51( 3, n, b, lda, t, lda, v, ldu, v, ldu, work,
889 $ rwork, result( 4 ) )
890*
891* Call ZHGEQZ to compute S1, P1, S2, P2, Q, and Z, do tests.
892*
893* Compute T1 and UZ
894*
895* Eigenvalues only
896*
897 CALL zlacpy( ' ', n, n, h, lda, s2, lda )
898 CALL zlacpy( ' ', n, n, t, lda, p2, lda )
899 ntest = 5
900 result( 5 ) = ulpinv
901*
902 CALL zhgeqz( 'E', 'N', 'N', n, 1, n, s2, lda, p2, lda,
903 $ alpha3, beta3, q, ldu, z, ldu, work, lwork,
904 $ rwork, iinfo )
905 IF( iinfo.NE.0 ) THEN
906 WRITE( nounit, fmt = 9999 )'ZHGEQZ(E)', iinfo, n, jtype,
907 $ ioldsd
908 info = abs( iinfo )
909 GO TO 210
910 END IF
911*
912* Eigenvalues and Full Schur Form
913*
914 CALL zlacpy( ' ', n, n, h, lda, s2, lda )
915 CALL zlacpy( ' ', n, n, t, lda, p2, lda )
916*
917 CALL zhgeqz( 'S', 'N', 'N', n, 1, n, s2, lda, p2, lda,
918 $ alpha1, beta1, q, ldu, z, ldu, work, lwork,
919 $ rwork, iinfo )
920 IF( iinfo.NE.0 ) THEN
921 WRITE( nounit, fmt = 9999 )'ZHGEQZ(S)', iinfo, n, jtype,
922 $ ioldsd
923 info = abs( iinfo )
924 GO TO 210
925 END IF
926*
927* Eigenvalues, Schur Form, and Schur Vectors
928*
929 CALL zlacpy( ' ', n, n, h, lda, s1, lda )
930 CALL zlacpy( ' ', n, n, t, lda, p1, lda )
931*
932 CALL zhgeqz( 'S', 'I', 'I', n, 1, n, s1, lda, p1, lda,
933 $ alpha1, beta1, q, ldu, z, ldu, work, lwork,
934 $ rwork, iinfo )
935 IF( iinfo.NE.0 ) THEN
936 WRITE( nounit, fmt = 9999 )'ZHGEQZ(V)', iinfo, n, jtype,
937 $ ioldsd
938 info = abs( iinfo )
939 GO TO 210
940 END IF
941*
942 ntest = 8
943*
944* Do Tests 5--8
945*
946 CALL zget51( 1, n, h, lda, s1, lda, q, ldu, z, ldu, work,
947 $ rwork, result( 5 ) )
948 CALL zget51( 1, n, t, lda, p1, lda, q, ldu, z, ldu, work,
949 $ rwork, result( 6 ) )
950 CALL zget51( 3, n, t, lda, p1, lda, q, ldu, q, ldu, work,
951 $ rwork, result( 7 ) )
952 CALL zget51( 3, n, t, lda, p1, lda, z, ldu, z, ldu, work,
953 $ rwork, result( 8 ) )
954*
955* Compute the Left and Right Eigenvectors of (S1,P1)
956*
957* 9: Compute the left eigenvector Matrix without
958* back transforming:
959*
960 ntest = 9
961 result( 9 ) = ulpinv
962*
963* To test "SELECT" option, compute half of the eigenvectors
964* in one call, and half in another
965*
966 i1 = n / 2
967 DO 120 j = 1, i1
968 llwork( j ) = .true.
969 120 CONTINUE
970 DO 130 j = i1 + 1, n
971 llwork( j ) = .false.
972 130 CONTINUE
973*
974 CALL ztgevc( 'L', 'S', llwork, n, s1, lda, p1, lda, evectl,
975 $ ldu, cdumma, ldu, n, in, work, rwork, iinfo )
976 IF( iinfo.NE.0 ) THEN
977 WRITE( nounit, fmt = 9999 )'ZTGEVC(L,S1)', iinfo, n,
978 $ jtype, ioldsd
979 info = abs( iinfo )
980 GO TO 210
981 END IF
982*
983 i1 = in
984 DO 140 j = 1, i1
985 llwork( j ) = .false.
986 140 CONTINUE
987 DO 150 j = i1 + 1, n
988 llwork( j ) = .true.
989 150 CONTINUE
990*
991 CALL ztgevc( 'L', 'S', llwork, n, s1, lda, p1, lda,
992 $ evectl( 1, i1+1 ), ldu, cdumma, ldu, n, in,
993 $ work, rwork, iinfo )
994 IF( iinfo.NE.0 ) THEN
995 WRITE( nounit, fmt = 9999 )'ZTGEVC(L,S2)', iinfo, n,
996 $ jtype, ioldsd
997 info = abs( iinfo )
998 GO TO 210
999 END IF
1000*
1001 CALL zget52( .true., n, s1, lda, p1, lda, evectl, ldu,
1002 $ alpha1, beta1, work, rwork, dumma( 1 ) )
1003 result( 9 ) = dumma( 1 )
1004 IF( dumma( 2 ).GT.thrshn ) THEN
1005 WRITE( nounit, fmt = 9998 )'Left', 'ZTGEVC(HOWMNY=S)',
1006 $ dumma( 2 ), n, jtype, ioldsd
1007 END IF
1008*
1009* 10: Compute the left eigenvector Matrix with
1010* back transforming:
1011*
1012 ntest = 10
1013 result( 10 ) = ulpinv
1014 CALL zlacpy( 'F', n, n, q, ldu, evectl, ldu )
1015 CALL ztgevc( 'L', 'B', llwork, n, s1, lda, p1, lda, evectl,
1016 $ ldu, cdumma, ldu, n, in, work, rwork, iinfo )
1017 IF( iinfo.NE.0 ) THEN
1018 WRITE( nounit, fmt = 9999 )'ZTGEVC(L,B)', iinfo, n,
1019 $ jtype, ioldsd
1020 info = abs( iinfo )
1021 GO TO 210
1022 END IF
1023*
1024 CALL zget52( .true., n, h, lda, t, lda, evectl, ldu, alpha1,
1025 $ beta1, work, rwork, dumma( 1 ) )
1026 result( 10 ) = dumma( 1 )
1027 IF( dumma( 2 ).GT.thrshn ) THEN
1028 WRITE( nounit, fmt = 9998 )'Left', 'ZTGEVC(HOWMNY=B)',
1029 $ dumma( 2 ), n, jtype, ioldsd
1030 END IF
1031*
1032* 11: Compute the right eigenvector Matrix without
1033* back transforming:
1034*
1035 ntest = 11
1036 result( 11 ) = ulpinv
1037*
1038* To test "SELECT" option, compute half of the eigenvectors
1039* in one call, and half in another
1040*
1041 i1 = n / 2
1042 DO 160 j = 1, i1
1043 llwork( j ) = .true.
1044 160 CONTINUE
1045 DO 170 j = i1 + 1, n
1046 llwork( j ) = .false.
1047 170 CONTINUE
1048*
1049 CALL ztgevc( 'R', 'S', llwork, n, s1, lda, p1, lda, cdumma,
1050 $ ldu, evectr, ldu, n, in, work, rwork, iinfo )
1051 IF( iinfo.NE.0 ) THEN
1052 WRITE( nounit, fmt = 9999 )'ZTGEVC(R,S1)', iinfo, n,
1053 $ jtype, ioldsd
1054 info = abs( iinfo )
1055 GO TO 210
1056 END IF
1057*
1058 i1 = in
1059 DO 180 j = 1, i1
1060 llwork( j ) = .false.
1061 180 CONTINUE
1062 DO 190 j = i1 + 1, n
1063 llwork( j ) = .true.
1064 190 CONTINUE
1065*
1066 CALL ztgevc( 'R', 'S', llwork, n, s1, lda, p1, lda, cdumma,
1067 $ ldu, evectr( 1, i1+1 ), ldu, n, in, work,
1068 $ rwork, iinfo )
1069 IF( iinfo.NE.0 ) THEN
1070 WRITE( nounit, fmt = 9999 )'ZTGEVC(R,S2)', iinfo, n,
1071 $ jtype, ioldsd
1072 info = abs( iinfo )
1073 GO TO 210
1074 END IF
1075*
1076 CALL zget52( .false., n, s1, lda, p1, lda, evectr, ldu,
1077 $ alpha1, beta1, work, rwork, dumma( 1 ) )
1078 result( 11 ) = dumma( 1 )
1079 IF( dumma( 2 ).GT.thresh ) THEN
1080 WRITE( nounit, fmt = 9998 )'Right', 'ZTGEVC(HOWMNY=S)',
1081 $ dumma( 2 ), n, jtype, ioldsd
1082 END IF
1083*
1084* 12: Compute the right eigenvector Matrix with
1085* back transforming:
1086*
1087 ntest = 12
1088 result( 12 ) = ulpinv
1089 CALL zlacpy( 'F', n, n, z, ldu, evectr, ldu )
1090 CALL ztgevc( 'R', 'B', llwork, n, s1, lda, p1, lda, cdumma,
1091 $ ldu, evectr, ldu, n, in, work, rwork, iinfo )
1092 IF( iinfo.NE.0 ) THEN
1093 WRITE( nounit, fmt = 9999 )'ZTGEVC(R,B)', iinfo, n,
1094 $ jtype, ioldsd
1095 info = abs( iinfo )
1096 GO TO 210
1097 END IF
1098*
1099 CALL zget52( .false., n, h, lda, t, lda, evectr, ldu,
1100 $ alpha1, beta1, work, rwork, dumma( 1 ) )
1101 result( 12 ) = dumma( 1 )
1102 IF( dumma( 2 ).GT.thresh ) THEN
1103 WRITE( nounit, fmt = 9998 )'Right', 'ZTGEVC(HOWMNY=B)',
1104 $ dumma( 2 ), n, jtype, ioldsd
1105 END IF
1106*
1107* Tests 13--15 are done only on request
1108*
1109 IF( tstdif ) THEN
1110*
1111* Do Tests 13--14
1112*
1113 CALL zget51( 2, n, s1, lda, s2, lda, q, ldu, z, ldu,
1114 $ work, rwork, result( 13 ) )
1115 CALL zget51( 2, n, p1, lda, p2, lda, q, ldu, z, ldu,
1116 $ work, rwork, result( 14 ) )
1117*
1118* Do Test 15
1119*
1120 temp1 = zero
1121 temp2 = zero
1122 DO 200 j = 1, n
1123 temp1 = max( temp1, abs( alpha1( j )-alpha3( j ) ) )
1124 temp2 = max( temp2, abs( beta1( j )-beta3( j ) ) )
1125 200 CONTINUE
1126*
1127 temp1 = temp1 / max( safmin, ulp*max( temp1, anorm ) )
1128 temp2 = temp2 / max( safmin, ulp*max( temp2, bnorm ) )
1129 result( 15 ) = max( temp1, temp2 )
1130 ntest = 15
1131 ELSE
1132 result( 13 ) = zero
1133 result( 14 ) = zero
1134 result( 15 ) = zero
1135 ntest = 12
1136 END IF
1137*
1138* End of Loop -- Check for RESULT(j) > THRESH
1139*
1140 210 CONTINUE
1141*
1142 ntestt = ntestt + ntest
1143*
1144* Print out tests which fail.
1145*
1146 DO 220 jr = 1, ntest
1147 IF( result( jr ).GE.thresh ) THEN
1148*
1149* If this is the first test to fail,
1150* print a header to the data file.
1151*
1152 IF( nerrs.EQ.0 ) THEN
1153 WRITE( nounit, fmt = 9997 )'ZGG'
1154*
1155* Matrix types
1156*
1157 WRITE( nounit, fmt = 9996 )
1158 WRITE( nounit, fmt = 9995 )
1159 WRITE( nounit, fmt = 9994 )'Unitary'
1160*
1161* Tests performed
1162*
1163 WRITE( nounit, fmt = 9993 )'unitary', '*',
1164 $ 'conjugate transpose', ( '*', j = 1, 10 )
1165*
1166 END IF
1167 nerrs = nerrs + 1
1168 IF( result( jr ).LT.10000.0d0 ) THEN
1169 WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
1170 $ result( jr )
1171 ELSE
1172 WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
1173 $ result( jr )
1174 END IF
1175 END IF
1176 220 CONTINUE
1177*
1178 230 CONTINUE
1179 240 CONTINUE
1180*
1181* Summary
1182*
1183 CALL dlasum( 'ZGG', nounit, nerrs, ntestt )
1184 RETURN
1185*
1186 9999 FORMAT( ' ZCHKGG: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1187 $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1188*
1189 9998 FORMAT( ' ZCHKGG: ', a, ' Eigenvectors from ', a, ' incorrectly ',
1190 $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
1191 $ 'N=', i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5,
1192 $ ')' )
1193*
1194 9997 FORMAT( 1x, a3, ' -- Complex Generalized eigenvalue problem' )
1195*
1196 9996 FORMAT( ' Matrix types (see ZCHKGG for details): ' )
1197*
1198 9995 FORMAT( ' Special Matrices:', 23x,
1199 $ '(J''=transposed Jordan block)',
1200 $ / ' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
1201 $ '6=(diag(J'',I), diag(I,J''))', / ' Diagonal Matrices: ( ',
1202 $ 'D=diag(0,1,2,...) )', / ' 7=(D,I) 9=(large*D, small*I',
1203 $ ') 11=(large*I, small*D) 13=(large*D, large*I)', /
1204 $ ' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
1205 $ ' 14=(small*D, small*I)', / ' 15=(D, reversed D)' )
1206 9994 FORMAT( ' Matrices Rotated by Random ', a, ' Matrices U, V:',
1207 $ / ' 16=Transposed Jordan Blocks 19=geometric ',
1208 $ 'alpha, beta=0,1', / ' 17=arithm. alpha&beta ',
1209 $ ' 20=arithmetic alpha, beta=0,1', / ' 18=clustered ',
1210 $ 'alpha, beta=0,1 21=random alpha, beta=0,1',
1211 $ / ' Large & Small Matrices:', / ' 22=(large, small) ',
1212 $ '23=(small,large) 24=(small,small) 25=(large,large)',
1213 $ / ' 26=random O(1) matrices.' )
1214*
1215 9993 FORMAT( / ' Tests performed: (H is Hessenberg, S is Schur, B, ',
1216 $ 'T, P are triangular,', / 20x, 'U, V, Q, and Z are ', a,
1217 $ ', l and r are the', / 20x,
1218 $ 'appropriate left and right eigenvectors, resp., a is',
1219 $ / 20x, 'alpha, b is beta, and ', a, ' means ', a, '.)',
1220 $ / ' 1 = | A - U H V', a,
1221 $ ' | / ( |A| n ulp ) 2 = | B - U T V', a,
1222 $ ' | / ( |B| n ulp )', / ' 3 = | I - UU', a,
1223 $ ' | / ( n ulp ) 4 = | I - VV', a,
1224 $ ' | / ( n ulp )', / ' 5 = | H - Q S Z', a,
1225 $ ' | / ( |H| n ulp )', 6x, '6 = | T - Q P Z', a,
1226 $ ' | / ( |T| n ulp )', / ' 7 = | I - QQ', a,
1227 $ ' | / ( n ulp ) 8 = | I - ZZ', a,
1228 $ ' | / ( n ulp )', / ' 9 = max | ( b S - a P )', a,
1229 $ ' l | / const. 10 = max | ( b H - a T )', a,
1230 $ ' l | / const.', /
1231 $ ' 11= max | ( b S - a P ) r | / const. 12 = max | ( b H',
1232 $ ' - a T ) r | / const.', / 1x )
1233*
1234 9992 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
1235 $ 4( i4, ',' ), ' result ', i2, ' is', 0p, f8.2 )
1236 9991 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
1237 $ 4( i4, ',' ), ' result ', i2, ' is', 1p, d10.3 )
1238*
1239* End of ZCHKGG
1240*
1241 END
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine zchkgg(NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, TSTDIF, THRSHN, NOUNIT, A, LDA, B, H, T, S1, S2, P1, P2, U, LDU, V, Q, Z, ALPHA1, BETA1, ALPHA3, BETA3, EVECTL, EVECTR, WORK, LWORK, RWORK, LLWORK, RESULT, INFO)
ZCHKGG
Definition: zchkgg.f:503
subroutine zget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHA, BETA, WORK, RWORK, RESULT)
ZGET52
Definition: zget52.f:162
subroutine zlatm4(ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND, TRIANG, IDIST, ISEED, A, LDA)
ZLATM4
Definition: zlatm4.f:171
subroutine zget51(ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK, RWORK, RESULT)
ZGET51
Definition: zget51.f:155
subroutine ztgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
ZTGEVC
Definition: ztgevc.f:219
subroutine zhgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
ZHGEQZ
Definition: zhgeqz.f:284
subroutine zgeqr2(M, N, A, LDA, TAU, WORK, INFO)
ZGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm.
Definition: zgeqr2.f:130
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 zlarfg(N, ALPHA, X, INCX, TAU)
ZLARFG generates an elementary reflector (Householder matrix).
Definition: zlarfg.f:106
subroutine zgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
ZGGHRD
Definition: zgghrd.f:204
subroutine zunm2r(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
ZUNM2R multiplies a general matrix by the unitary matrix from a QR factorization determined by cgeqrf...
Definition: zunm2r.f:159
subroutine dlasum(TYPE, IOUNIT, IE, NRUN)
DLASUM
Definition: dlasum.f:43