LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Classes Files Functions Variables Typedefs Macros
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 threshhold 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 threshhold 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 *> Threshhold 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 *> \date November 2011
496 *
497 *> \ingroup complex16_eig
498 *
499 * =====================================================================
500  SUBROUTINE zchkgg( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
501  $ tstdif, thrshn, nounit, a, lda, b, h, t, s1,
502  $ s2, p1, p2, u, ldu, v, q, z, alpha1, beta1,
503  $ alpha3, beta3, evectl, evectr, work, lwork,
504  $ rwork, llwork, result, info )
505 *
506 * -- LAPACK test routine (version 3.4.0) --
507 * -- LAPACK is a software package provided by Univ. of Tennessee, --
508 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
509 * November 2011
510 *
511 * .. Scalar Arguments ..
512  LOGICAL tstdif
513  INTEGER info, lda, ldu, lwork, nounit, nsizes, ntypes
514  DOUBLE PRECISION thresh, thrshn
515 * ..
516 * .. Array Arguments ..
517  LOGICAL dotype( * ), llwork( * )
518  INTEGER iseed( 4 ), nn( * )
519  DOUBLE PRECISION result( 15 ), rwork( * )
520  COMPLEX*16 a( lda, * ), alpha1( * ), alpha3( * ),
521  $ b( lda, * ), beta1( * ), beta3( * ),
522  $ evectl( ldu, * ), evectr( ldu, * ),
523  $ h( lda, * ), p1( lda, * ), p2( lda, * ),
524  $ q( ldu, * ), s1( lda, * ), s2( lda, * ),
525  $ t( lda, * ), u( ldu, * ), v( ldu, * ),
526  $ work( * ), z( ldu, * )
527 * ..
528 *
529 * =====================================================================
530 *
531 * .. Parameters ..
532  DOUBLE PRECISION zero, one
533  parameter( zero = 0.0d+0, one = 1.0d+0 )
534  COMPLEX*16 czero, cone
535  parameter( czero = ( 0.0d+0, 0.0d+0 ),
536  $ cone = ( 1.0d+0, 0.0d+0 ) )
537  INTEGER maxtyp
538  parameter( maxtyp = 26 )
539 * ..
540 * .. Local Scalars ..
541  LOGICAL badnn
542  INTEGER i1, iadd, iinfo, in, j, jc, jr, jsize, jtype,
543  $ lwkopt, mtypes, n, n1, nerrs, nmats, nmax,
544  $ ntest, ntestt
545  DOUBLE PRECISION anorm, bnorm, safmax, safmin, temp1, temp2,
546  $ ulp, ulpinv
547  COMPLEX*16 ctemp
548 * ..
549 * .. Local Arrays ..
550  LOGICAL lasign( maxtyp ), lbsign( maxtyp )
551  INTEGER ioldsd( 4 ), kadd( 6 ), kamagn( maxtyp ),
552  $ katype( maxtyp ), kazero( maxtyp ),
553  $ kbmagn( maxtyp ), kbtype( maxtyp ),
554  $ kbzero( maxtyp ), kclass( maxtyp ),
555  $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
556  DOUBLE PRECISION dumma( 4 ), rmagn( 0: 3 )
557  COMPLEX*16 cdumma( 4 )
558 * ..
559 * .. External Functions ..
560  DOUBLE PRECISION dlamch, zlange
561  COMPLEX*16 zlarnd
562  EXTERNAL dlamch, zlange, zlarnd
563 * ..
564 * .. External Subroutines ..
565  EXTERNAL dlabad, dlasum, xerbla, zgeqr2, zget51, zget52,
567  $ ztgevc, zunm2r
568 * ..
569 * .. Intrinsic Functions ..
570  INTRINSIC abs, dble, dconjg, max, min, sign
571 * ..
572 * .. Data statements ..
573  DATA kclass / 15*1, 10*2, 1*3 /
574  DATA kz1 / 0, 1, 2, 1, 3, 3 /
575  DATA kz2 / 0, 0, 1, 2, 1, 1 /
576  DATA kadd / 0, 0, 0, 0, 3, 2 /
577  DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
578  $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
579  DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
580  $ 1, 1, -4, 2, -4, 8*8, 0 /
581  DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
582  $ 4*5, 4*3, 1 /
583  DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
584  $ 4*6, 4*4, 1 /
585  DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
586  $ 2, 1 /
587  DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
588  $ 2, 1 /
589  DATA ktrian / 16*0, 10*1 /
590  DATA lasign / 6*.false., .true., .false., 2*.true.,
591  $ 2*.false., 3*.true., .false., .true.,
592  $ 3*.false., 5*.true., .false. /
593  DATA lbsign / 7*.false., .true., 2*.false.,
594  $ 2*.true., 2*.false., .true., .false., .true.,
595  $ 9*.false. /
596 * ..
597 * .. Executable Statements ..
598 *
599 * Check for errors
600 *
601  info = 0
602 *
603  badnn = .false.
604  nmax = 1
605  DO 10 j = 1, nsizes
606  nmax = max( nmax, nn( j ) )
607  IF( nn( j ).LT.0 )
608  $ badnn = .true.
609  10 CONTINUE
610 *
611  lwkopt = max( 2*nmax*nmax, 4*nmax, 1 )
612 *
613 * Check for errors
614 *
615  IF( nsizes.LT.0 ) THEN
616  info = -1
617  ELSE IF( badnn ) THEN
618  info = -2
619  ELSE IF( ntypes.LT.0 ) THEN
620  info = -3
621  ELSE IF( thresh.LT.zero ) THEN
622  info = -6
623  ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
624  info = -10
625  ELSE IF( ldu.LE.1 .OR. ldu.LT.nmax ) THEN
626  info = -19
627  ELSE IF( lwkopt.GT.lwork ) THEN
628  info = -30
629  END IF
630 *
631  IF( info.NE.0 ) THEN
632  CALL xerbla( 'ZCHKGG', -info )
633  RETURN
634  END IF
635 *
636 * Quick return if possible
637 *
638  IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
639  $ RETURN
640 *
641  safmin = dlamch( 'Safe minimum' )
642  ulp = dlamch( 'Epsilon' )*dlamch( 'Base' )
643  safmin = safmin / ulp
644  safmax = one / safmin
645  CALL dlabad( safmin, safmax )
646  ulpinv = one / ulp
647 *
648 * The values RMAGN(2:3) depend on N, see below.
649 *
650  rmagn( 0 ) = zero
651  rmagn( 1 ) = one
652 *
653 * Loop over sizes, types
654 *
655  ntestt = 0
656  nerrs = 0
657  nmats = 0
658 *
659  DO 240 jsize = 1, nsizes
660  n = nn( jsize )
661  n1 = max( 1, n )
662  rmagn( 2 ) = safmax*ulp / dble( n1 )
663  rmagn( 3 ) = safmin*ulpinv*n1
664 *
665  IF( nsizes.NE.1 ) THEN
666  mtypes = min( maxtyp, ntypes )
667  ELSE
668  mtypes = min( maxtyp+1, ntypes )
669  END IF
670 *
671  DO 230 jtype = 1, mtypes
672  IF( .NOT.dotype( jtype ) )
673  $ go to 230
674  nmats = nmats + 1
675  ntest = 0
676 *
677 * Save ISEED in case of an error.
678 *
679  DO 20 j = 1, 4
680  ioldsd( j ) = iseed( j )
681  20 CONTINUE
682 *
683 * Initialize RESULT
684 *
685  DO 30 j = 1, 15
686  result( j ) = zero
687  30 CONTINUE
688 *
689 * Compute A and B
690 *
691 * Description of control parameters:
692 *
693 * KZLASS: =1 means w/o rotation, =2 means w/ rotation,
694 * =3 means random.
695 * KATYPE: the "type" to be passed to ZLATM4 for computing A.
696 * KAZERO: the pattern of zeros on the diagonal for A:
697 * =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ),
698 * =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ),
699 * =6: ( 0, 1, 0, xxx, 0 ). (xxx means a string of
700 * non-zero entries.)
701 * KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1),
702 * =2: large, =3: small.
703 * LASIGN: .TRUE. if the diagonal elements of A are to be
704 * multiplied by a random magnitude 1 number.
705 * KBTYPE, KBZERO, KBMAGN, LBSIGN: the same, but for B.
706 * KTRIAN: =0: don't fill in the upper triangle, =1: do.
707 * KZ1, KZ2, KADD: used to implement KAZERO and KBZERO.
708 * RMAGN: used to implement KAMAGN and KBMAGN.
709 *
710  IF( mtypes.GT.maxtyp )
711  $ go to 110
712  iinfo = 0
713  IF( kclass( jtype ).LT.3 ) THEN
714 *
715 * Generate A (w/o rotation)
716 *
717  IF( abs( katype( jtype ) ).EQ.3 ) THEN
718  in = 2*( ( n-1 ) / 2 ) + 1
719  IF( in.NE.n )
720  $ CALL zlaset( 'Full', n, n, czero, czero, a, lda )
721  ELSE
722  in = n
723  END IF
724  CALL zlatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
725  $ kz2( kazero( jtype ) ), lasign( jtype ),
726  $ rmagn( kamagn( jtype ) ), ulp,
727  $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 4,
728  $ iseed, a, lda )
729  iadd = kadd( kazero( jtype ) )
730  IF( iadd.GT.0 .AND. iadd.LE.n )
731  $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
732 *
733 * Generate B (w/o rotation)
734 *
735  IF( abs( kbtype( jtype ) ).EQ.3 ) THEN
736  in = 2*( ( n-1 ) / 2 ) + 1
737  IF( in.NE.n )
738  $ CALL zlaset( 'Full', n, n, czero, czero, b, lda )
739  ELSE
740  in = n
741  END IF
742  CALL zlatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
743  $ kz2( kbzero( jtype ) ), lbsign( jtype ),
744  $ rmagn( kbmagn( jtype ) ), one,
745  $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 4,
746  $ iseed, b, lda )
747  iadd = kadd( kbzero( jtype ) )
748  IF( iadd.NE.0 )
749  $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
750 *
751  IF( kclass( jtype ).EQ.2 .AND. n.GT.0 ) THEN
752 *
753 * Include rotations
754 *
755 * Generate U, V as Householder transformations times a
756 * diagonal matrix. (Note that ZLARFG makes U(j,j) and
757 * V(j,j) real.)
758 *
759  DO 50 jc = 1, n - 1
760  DO 40 jr = jc, n
761  u( jr, jc ) = zlarnd( 3, iseed )
762  v( jr, jc ) = zlarnd( 3, iseed )
763  40 CONTINUE
764  CALL zlarfg( n+1-jc, u( jc, jc ), u( jc+1, jc ), 1,
765  $ work( jc ) )
766  work( 2*n+jc ) = sign( one, dble( u( jc, jc ) ) )
767  u( jc, jc ) = cone
768  CALL zlarfg( n+1-jc, v( jc, jc ), v( jc+1, jc ), 1,
769  $ work( n+jc ) )
770  work( 3*n+jc ) = sign( one, dble( v( jc, jc ) ) )
771  v( jc, jc ) = cone
772  50 CONTINUE
773  ctemp = zlarnd( 3, iseed )
774  u( n, n ) = cone
775  work( n ) = czero
776  work( 3*n ) = ctemp / abs( ctemp )
777  ctemp = zlarnd( 3, iseed )
778  v( n, n ) = cone
779  work( 2*n ) = czero
780  work( 4*n ) = ctemp / abs( ctemp )
781 *
782 * Apply the diagonal matrices
783 *
784  DO 70 jc = 1, n
785  DO 60 jr = 1, n
786  a( jr, jc ) = work( 2*n+jr )*
787  $ dconjg( work( 3*n+jc ) )*
788  $ a( jr, jc )
789  b( jr, jc ) = work( 2*n+jr )*
790  $ dconjg( work( 3*n+jc ) )*
791  $ b( jr, jc )
792  60 CONTINUE
793  70 CONTINUE
794  CALL zunm2r( 'L', 'N', n, n, n-1, u, ldu, work, a,
795  $ lda, work( 2*n+1 ), iinfo )
796  IF( iinfo.NE.0 )
797  $ go to 100
798  CALL zunm2r( 'R', 'C', n, n, n-1, v, ldu, work( n+1 ),
799  $ a, lda, work( 2*n+1 ), iinfo )
800  IF( iinfo.NE.0 )
801  $ go to 100
802  CALL zunm2r( 'L', 'N', n, n, n-1, u, ldu, work, b,
803  $ lda, work( 2*n+1 ), iinfo )
804  IF( iinfo.NE.0 )
805  $ go to 100
806  CALL zunm2r( 'R', 'C', n, n, n-1, v, ldu, work( n+1 ),
807  $ b, lda, work( 2*n+1 ), iinfo )
808  IF( iinfo.NE.0 )
809  $ go to 100
810  END IF
811  ELSE
812 *
813 * Random matrices
814 *
815  DO 90 jc = 1, n
816  DO 80 jr = 1, n
817  a( jr, jc ) = rmagn( kamagn( jtype ) )*
818  $ zlarnd( 4, iseed )
819  b( jr, jc ) = rmagn( kbmagn( jtype ) )*
820  $ zlarnd( 4, iseed )
821  80 CONTINUE
822  90 CONTINUE
823  END IF
824 *
825  anorm = zlange( '1', n, n, a, lda, rwork )
826  bnorm = zlange( '1', n, n, b, lda, rwork )
827 *
828  100 CONTINUE
829 *
830  IF( iinfo.NE.0 ) THEN
831  WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
832  $ ioldsd
833  info = abs( iinfo )
834  RETURN
835  END IF
836 *
837  110 CONTINUE
838 *
839 * Call ZGEQR2, ZUNM2R, and ZGGHRD to compute H, T, U, and V
840 *
841  CALL zlacpy( ' ', n, n, a, lda, h, lda )
842  CALL zlacpy( ' ', n, n, b, lda, t, lda )
843  ntest = 1
844  result( 1 ) = ulpinv
845 *
846  CALL zgeqr2( n, n, t, lda, work, work( n+1 ), iinfo )
847  IF( iinfo.NE.0 ) THEN
848  WRITE( nounit, fmt = 9999 )'ZGEQR2', iinfo, n, jtype,
849  $ ioldsd
850  info = abs( iinfo )
851  go to 210
852  END IF
853 *
854  CALL zunm2r( 'L', 'C', n, n, n, t, lda, work, h, lda,
855  $ work( n+1 ), iinfo )
856  IF( iinfo.NE.0 ) THEN
857  WRITE( nounit, fmt = 9999 )'ZUNM2R', iinfo, n, jtype,
858  $ ioldsd
859  info = abs( iinfo )
860  go to 210
861  END IF
862 *
863  CALL zlaset( 'Full', n, n, czero, cone, u, ldu )
864  CALL zunm2r( 'R', 'N', n, n, n, t, lda, work, u, ldu,
865  $ work( n+1 ), iinfo )
866  IF( iinfo.NE.0 ) THEN
867  WRITE( nounit, fmt = 9999 )'ZUNM2R', iinfo, n, jtype,
868  $ ioldsd
869  info = abs( iinfo )
870  go to 210
871  END IF
872 *
873  CALL zgghrd( 'V', 'I', n, 1, n, h, lda, t, lda, u, ldu, v,
874  $ ldu, iinfo )
875  IF( iinfo.NE.0 ) THEN
876  WRITE( nounit, fmt = 9999 )'ZGGHRD', iinfo, n, jtype,
877  $ ioldsd
878  info = abs( iinfo )
879  go to 210
880  END IF
881  ntest = 4
882 *
883 * Do tests 1--4
884 *
885  CALL zget51( 1, n, a, lda, h, lda, u, ldu, v, ldu, work,
886  $ rwork, result( 1 ) )
887  CALL zget51( 1, n, b, lda, t, lda, u, ldu, v, ldu, work,
888  $ rwork, result( 2 ) )
889  CALL zget51( 3, n, b, lda, t, lda, u, ldu, u, ldu, work,
890  $ rwork, result( 3 ) )
891  CALL zget51( 3, n, b, lda, t, lda, v, ldu, v, ldu, work,
892  $ rwork, result( 4 ) )
893 *
894 * Call ZHGEQZ to compute S1, P1, S2, P2, Q, and Z, do tests.
895 *
896 * Compute T1 and UZ
897 *
898 * Eigenvalues only
899 *
900  CALL zlacpy( ' ', n, n, h, lda, s2, lda )
901  CALL zlacpy( ' ', n, n, t, lda, p2, lda )
902  ntest = 5
903  result( 5 ) = ulpinv
904 *
905  CALL zhgeqz( 'E', 'N', 'N', n, 1, n, s2, lda, p2, lda,
906  $ alpha3, beta3, q, ldu, z, ldu, work, lwork,
907  $ rwork, iinfo )
908  IF( iinfo.NE.0 ) THEN
909  WRITE( nounit, fmt = 9999 )'ZHGEQZ(E)', iinfo, n, jtype,
910  $ ioldsd
911  info = abs( iinfo )
912  go to 210
913  END IF
914 *
915 * Eigenvalues and Full Schur Form
916 *
917  CALL zlacpy( ' ', n, n, h, lda, s2, lda )
918  CALL zlacpy( ' ', n, n, t, lda, p2, lda )
919 *
920  CALL zhgeqz( 'S', 'N', 'N', n, 1, n, s2, lda, p2, lda,
921  $ alpha1, beta1, q, ldu, z, ldu, work, lwork,
922  $ rwork, iinfo )
923  IF( iinfo.NE.0 ) THEN
924  WRITE( nounit, fmt = 9999 )'ZHGEQZ(S)', iinfo, n, jtype,
925  $ ioldsd
926  info = abs( iinfo )
927  go to 210
928  END IF
929 *
930 * Eigenvalues, Schur Form, and Schur Vectors
931 *
932  CALL zlacpy( ' ', n, n, h, lda, s1, lda )
933  CALL zlacpy( ' ', n, n, t, lda, p1, lda )
934 *
935  CALL zhgeqz( 'S', 'I', 'I', n, 1, n, s1, lda, p1, lda,
936  $ alpha1, beta1, q, ldu, z, ldu, work, lwork,
937  $ rwork, iinfo )
938  IF( iinfo.NE.0 ) THEN
939  WRITE( nounit, fmt = 9999 )'ZHGEQZ(V)', iinfo, n, jtype,
940  $ ioldsd
941  info = abs( iinfo )
942  go to 210
943  END IF
944 *
945  ntest = 8
946 *
947 * Do Tests 5--8
948 *
949  CALL zget51( 1, n, h, lda, s1, lda, q, ldu, z, ldu, work,
950  $ rwork, result( 5 ) )
951  CALL zget51( 1, n, t, lda, p1, lda, q, ldu, z, ldu, work,
952  $ rwork, result( 6 ) )
953  CALL zget51( 3, n, t, lda, p1, lda, q, ldu, q, ldu, work,
954  $ rwork, result( 7 ) )
955  CALL zget51( 3, n, t, lda, p1, lda, z, ldu, z, ldu, work,
956  $ rwork, result( 8 ) )
957 *
958 * Compute the Left and Right Eigenvectors of (S1,P1)
959 *
960 * 9: Compute the left eigenvector Matrix without
961 * back transforming:
962 *
963  ntest = 9
964  result( 9 ) = ulpinv
965 *
966 * To test "SELECT" option, compute half of the eigenvectors
967 * in one call, and half in another
968 *
969  i1 = n / 2
970  DO 120 j = 1, i1
971  llwork( j ) = .true.
972  120 CONTINUE
973  DO 130 j = i1 + 1, n
974  llwork( j ) = .false.
975  130 CONTINUE
976 *
977  CALL ztgevc( 'L', 'S', llwork, n, s1, lda, p1, lda, evectl,
978  $ ldu, cdumma, ldu, n, in, work, rwork, iinfo )
979  IF( iinfo.NE.0 ) THEN
980  WRITE( nounit, fmt = 9999 )'ZTGEVC(L,S1)', iinfo, n,
981  $ jtype, ioldsd
982  info = abs( iinfo )
983  go to 210
984  END IF
985 *
986  i1 = in
987  DO 140 j = 1, i1
988  llwork( j ) = .false.
989  140 CONTINUE
990  DO 150 j = i1 + 1, n
991  llwork( j ) = .true.
992  150 CONTINUE
993 *
994  CALL ztgevc( 'L', 'S', llwork, n, s1, lda, p1, lda,
995  $ evectl( 1, i1+1 ), ldu, cdumma, ldu, n, in,
996  $ work, rwork, iinfo )
997  IF( iinfo.NE.0 ) THEN
998  WRITE( nounit, fmt = 9999 )'ZTGEVC(L,S2)', iinfo, n,
999  $ jtype, ioldsd
1000  info = abs( iinfo )
1001  go to 210
1002  END IF
1003 *
1004  CALL zget52( .true., n, s1, lda, p1, lda, evectl, ldu,
1005  $ alpha1, beta1, work, rwork, dumma( 1 ) )
1006  result( 9 ) = dumma( 1 )
1007  IF( dumma( 2 ).GT.thrshn ) THEN
1008  WRITE( nounit, fmt = 9998 )'Left', 'ZTGEVC(HOWMNY=S)',
1009  $ dumma( 2 ), n, jtype, ioldsd
1010  END IF
1011 *
1012 * 10: Compute the left eigenvector Matrix with
1013 * back transforming:
1014 *
1015  ntest = 10
1016  result( 10 ) = ulpinv
1017  CALL zlacpy( 'F', n, n, q, ldu, evectl, ldu )
1018  CALL ztgevc( 'L', 'B', llwork, n, s1, lda, p1, lda, evectl,
1019  $ ldu, cdumma, ldu, n, in, work, rwork, iinfo )
1020  IF( iinfo.NE.0 ) THEN
1021  WRITE( nounit, fmt = 9999 )'ZTGEVC(L,B)', iinfo, n,
1022  $ jtype, ioldsd
1023  info = abs( iinfo )
1024  go to 210
1025  END IF
1026 *
1027  CALL zget52( .true., n, h, lda, t, lda, evectl, ldu, alpha1,
1028  $ beta1, work, rwork, dumma( 1 ) )
1029  result( 10 ) = dumma( 1 )
1030  IF( dumma( 2 ).GT.thrshn ) THEN
1031  WRITE( nounit, fmt = 9998 )'Left', 'ZTGEVC(HOWMNY=B)',
1032  $ dumma( 2 ), n, jtype, ioldsd
1033  END IF
1034 *
1035 * 11: Compute the right eigenvector Matrix without
1036 * back transforming:
1037 *
1038  ntest = 11
1039  result( 11 ) = ulpinv
1040 *
1041 * To test "SELECT" option, compute half of the eigenvectors
1042 * in one call, and half in another
1043 *
1044  i1 = n / 2
1045  DO 160 j = 1, i1
1046  llwork( j ) = .true.
1047  160 CONTINUE
1048  DO 170 j = i1 + 1, n
1049  llwork( j ) = .false.
1050  170 CONTINUE
1051 *
1052  CALL ztgevc( 'R', 'S', llwork, n, s1, lda, p1, lda, cdumma,
1053  $ ldu, evectr, ldu, n, in, work, rwork, iinfo )
1054  IF( iinfo.NE.0 ) THEN
1055  WRITE( nounit, fmt = 9999 )'ZTGEVC(R,S1)', iinfo, n,
1056  $ jtype, ioldsd
1057  info = abs( iinfo )
1058  go to 210
1059  END IF
1060 *
1061  i1 = in
1062  DO 180 j = 1, i1
1063  llwork( j ) = .false.
1064  180 CONTINUE
1065  DO 190 j = i1 + 1, n
1066  llwork( j ) = .true.
1067  190 CONTINUE
1068 *
1069  CALL ztgevc( 'R', 'S', llwork, n, s1, lda, p1, lda, cdumma,
1070  $ ldu, evectr( 1, i1+1 ), ldu, n, in, work,
1071  $ rwork, iinfo )
1072  IF( iinfo.NE.0 ) THEN
1073  WRITE( nounit, fmt = 9999 )'ZTGEVC(R,S2)', iinfo, n,
1074  $ jtype, ioldsd
1075  info = abs( iinfo )
1076  go to 210
1077  END IF
1078 *
1079  CALL zget52( .false., n, s1, lda, p1, lda, evectr, ldu,
1080  $ alpha1, beta1, work, rwork, dumma( 1 ) )
1081  result( 11 ) = dumma( 1 )
1082  IF( dumma( 2 ).GT.thresh ) THEN
1083  WRITE( nounit, fmt = 9998 )'Right', 'ZTGEVC(HOWMNY=S)',
1084  $ dumma( 2 ), n, jtype, ioldsd
1085  END IF
1086 *
1087 * 12: Compute the right eigenvector Matrix with
1088 * back transforming:
1089 *
1090  ntest = 12
1091  result( 12 ) = ulpinv
1092  CALL zlacpy( 'F', n, n, z, ldu, evectr, ldu )
1093  CALL ztgevc( 'R', 'B', llwork, n, s1, lda, p1, lda, cdumma,
1094  $ ldu, evectr, ldu, n, in, work, rwork, iinfo )
1095  IF( iinfo.NE.0 ) THEN
1096  WRITE( nounit, fmt = 9999 )'ZTGEVC(R,B)', iinfo, n,
1097  $ jtype, ioldsd
1098  info = abs( iinfo )
1099  go to 210
1100  END IF
1101 *
1102  CALL zget52( .false., n, h, lda, t, lda, evectr, ldu,
1103  $ alpha1, beta1, work, rwork, dumma( 1 ) )
1104  result( 12 ) = dumma( 1 )
1105  IF( dumma( 2 ).GT.thresh ) THEN
1106  WRITE( nounit, fmt = 9998 )'Right', 'ZTGEVC(HOWMNY=B)',
1107  $ dumma( 2 ), n, jtype, ioldsd
1108  END IF
1109 *
1110 * Tests 13--15 are done only on request
1111 *
1112  IF( tstdif ) THEN
1113 *
1114 * Do Tests 13--14
1115 *
1116  CALL zget51( 2, n, s1, lda, s2, lda, q, ldu, z, ldu,
1117  $ work, rwork, result( 13 ) )
1118  CALL zget51( 2, n, p1, lda, p2, lda, q, ldu, z, ldu,
1119  $ work, rwork, result( 14 ) )
1120 *
1121 * Do Test 15
1122 *
1123  temp1 = zero
1124  temp2 = zero
1125  DO 200 j = 1, n
1126  temp1 = max( temp1, abs( alpha1( j )-alpha3( j ) ) )
1127  temp2 = max( temp2, abs( beta1( j )-beta3( j ) ) )
1128  200 CONTINUE
1129 *
1130  temp1 = temp1 / max( safmin, ulp*max( temp1, anorm ) )
1131  temp2 = temp2 / max( safmin, ulp*max( temp2, bnorm ) )
1132  result( 15 ) = max( temp1, temp2 )
1133  ntest = 15
1134  ELSE
1135  result( 13 ) = zero
1136  result( 14 ) = zero
1137  result( 15 ) = zero
1138  ntest = 12
1139  END IF
1140 *
1141 * End of Loop -- Check for RESULT(j) > THRESH
1142 *
1143  210 CONTINUE
1144 *
1145  ntestt = ntestt + ntest
1146 *
1147 * Print out tests which fail.
1148 *
1149  DO 220 jr = 1, ntest
1150  IF( result( jr ).GE.thresh ) THEN
1151 *
1152 * If this is the first test to fail,
1153 * print a header to the data file.
1154 *
1155  IF( nerrs.EQ.0 ) THEN
1156  WRITE( nounit, fmt = 9997 )'ZGG'
1157 *
1158 * Matrix types
1159 *
1160  WRITE( nounit, fmt = 9996 )
1161  WRITE( nounit, fmt = 9995 )
1162  WRITE( nounit, fmt = 9994 )'Unitary'
1163 *
1164 * Tests performed
1165 *
1166  WRITE( nounit, fmt = 9993 )'unitary', '*',
1167  $ 'conjugate transpose', ( '*', j = 1, 10 )
1168 *
1169  END IF
1170  nerrs = nerrs + 1
1171  IF( result( jr ).LT.10000.0d0 ) THEN
1172  WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
1173  $ result( jr )
1174  ELSE
1175  WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
1176  $ result( jr )
1177  END IF
1178  END IF
1179  220 CONTINUE
1180 *
1181  230 CONTINUE
1182  240 CONTINUE
1183 *
1184 * Summary
1185 *
1186  CALL dlasum( 'ZGG', nounit, nerrs, ntestt )
1187  RETURN
1188 *
1189  9999 FORMAT( ' ZCHKGG: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1190  $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1191 *
1192  9998 FORMAT( ' ZCHKGG: ', a, ' Eigenvectors from ', a, ' incorrectly ',
1193  $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
1194  $ 'N=', i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5,
1195  $ ')' )
1196 *
1197  9997 FORMAT( 1x, a3, ' -- Complex Generalized eigenvalue problem' )
1198 *
1199  9996 FORMAT( ' Matrix types (see ZCHKGG for details): ' )
1200 *
1201  9995 FORMAT( ' Special Matrices:', 23x,
1202  $ '(J''=transposed Jordan block)',
1203  $ / ' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
1204  $ '6=(diag(J'',I), diag(I,J''))', / ' Diagonal Matrices: ( ',
1205  $ 'D=diag(0,1,2,...) )', / ' 7=(D,I) 9=(large*D, small*I',
1206  $ ') 11=(large*I, small*D) 13=(large*D, large*I)', /
1207  $ ' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
1208  $ ' 14=(small*D, small*I)', / ' 15=(D, reversed D)' )
1209  9994 FORMAT( ' Matrices Rotated by Random ', a, ' Matrices U, V:',
1210  $ / ' 16=Transposed Jordan Blocks 19=geometric ',
1211  $ 'alpha, beta=0,1', / ' 17=arithm. alpha&beta ',
1212  $ ' 20=arithmetic alpha, beta=0,1', / ' 18=clustered ',
1213  $ 'alpha, beta=0,1 21=random alpha, beta=0,1',
1214  $ / ' Large & Small Matrices:', / ' 22=(large, small) ',
1215  $ '23=(small,large) 24=(small,small) 25=(large,large)',
1216  $ / ' 26=random O(1) matrices.' )
1217 *
1218  9993 FORMAT( / ' Tests performed: (H is Hessenberg, S is Schur, B, ',
1219  $ 'T, P are triangular,', / 20x, 'U, V, Q, and Z are ', a,
1220  $ ', l and r are the', / 20x,
1221  $ 'appropriate left and right eigenvectors, resp., a is',
1222  $ / 20x, 'alpha, b is beta, and ', a, ' means ', a, '.)',
1223  $ / ' 1 = | A - U H V', a,
1224  $ ' | / ( |A| n ulp ) 2 = | B - U T V', a,
1225  $ ' | / ( |B| n ulp )', / ' 3 = | I - UU', a,
1226  $ ' | / ( n ulp ) 4 = | I - VV', a,
1227  $ ' | / ( n ulp )', / ' 5 = | H - Q S Z', a,
1228  $ ' | / ( |H| n ulp )', 6x, '6 = | T - Q P Z', a,
1229  $ ' | / ( |T| n ulp )', / ' 7 = | I - QQ', a,
1230  $ ' | / ( n ulp ) 8 = | I - ZZ', a,
1231  $ ' | / ( n ulp )', / ' 9 = max | ( b S - a P )', a,
1232  $ ' l | / const. 10 = max | ( b H - a T )', a,
1233  $ ' l | / const.', /
1234  $ ' 11= max | ( b S - a P ) r | / const. 12 = max | ( b H',
1235  $ ' - a T ) r | / const.', / 1x )
1236 *
1237  9992 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
1238  $ 4( i4, ',' ), ' result ', i2, ' is', 0p, f8.2 )
1239  9991 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
1240  $ 4( i4, ',' ), ' result ', i2, ' is', 1p, d10.3 )
1241 *
1242 * End of ZCHKGG
1243 *
1244  END