LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
cchkgg.f
Go to the documentation of this file.
1 *> \brief \b CCHKGG
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 CCHKGG( 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 * REAL THRESH, THRSHN
21 * ..
22 * .. Array Arguments ..
23 * LOGICAL DOTYPE( * ), LLWORK( * )
24 * INTEGER ISEED( 4 ), NN( * )
25 * REAL RESULT( 15 ), RWORK( * )
26 * COMPLEX 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 *> CCHKGG checks the nonsymmetric generalized eigenvalue problem
42 *> routines.
43 *> H H H
44 *> CGGHRD 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 *> CHGEQZ 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 *> CTGEVC 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 CCHKGG 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 *> STGEVC 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 *> STGEVC 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 *> CCHKGG 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, CCHKGG
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 CCHKGG to continue the same random number
273 *> sequence.
274 *> \endverbatim
275 *>
276 *> \param[in] THRESH
277 *> \verbatim
278 *> THRESH is REAL
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 REAL
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 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 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 array, dimension (LDA, max(NN))
342 *> The upper Hessenberg matrix computed from A by CGGHRD.
343 *> \endverbatim
344 *>
345 *> \param[out] T
346 *> \verbatim
347 *> T is COMPLEX array, dimension (LDA, max(NN))
348 *> The upper triangular matrix computed from B by CGGHRD.
349 *> \endverbatim
350 *>
351 *> \param[out] S1
352 *> \verbatim
353 *> S1 is COMPLEX array, dimension (LDA, max(NN))
354 *> The Schur (upper triangular) matrix computed from H by CHGEQZ
355 *> when Q and Z are also computed.
356 *> \endverbatim
357 *>
358 *> \param[out] S2
359 *> \verbatim
360 *> S2 is COMPLEX array, dimension (LDA, max(NN))
361 *> The Schur (upper triangular) matrix computed from H by CHGEQZ
362 *> when Q and Z are not computed.
363 *> \endverbatim
364 *>
365 *> \param[out] P1
366 *> \verbatim
367 *> P1 is COMPLEX array, dimension (LDA, max(NN))
368 *> The upper triangular matrix computed from T by CHGEQZ
369 *> when Q and Z are also computed.
370 *> \endverbatim
371 *>
372 *> \param[out] P2
373 *> \verbatim
374 *> P2 is COMPLEX array, dimension (LDA, max(NN))
375 *> The upper triangular matrix computed from T by CHGEQZ
376 *> when Q and Z are not computed.
377 *> \endverbatim
378 *>
379 *> \param[out] U
380 *> \verbatim
381 *> U is COMPLEX array, dimension (LDU, max(NN))
382 *> The (left) unitary matrix computed by CGGHRD.
383 *> \endverbatim
384 *>
385 *> \param[in] LDU
386 *> \verbatim
387 *> LDU is INTEGER
388 *> The leading dimension of U, V, Q, Z, EVECTL, and EVECTR. 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 array, dimension (LDU, max(NN))
395 *> The (right) unitary matrix computed by CGGHRD.
396 *> \endverbatim
397 *>
398 *> \param[out] Q
399 *> \verbatim
400 *> Q is COMPLEX array, dimension (LDU, max(NN))
401 *> The (left) unitary matrix computed by CHGEQZ.
402 *> \endverbatim
403 *>
404 *> \param[out] Z
405 *> \verbatim
406 *> Z is COMPLEX array, dimension (LDU, max(NN))
407 *> The (left) unitary matrix computed by CHGEQZ.
408 *> \endverbatim
409 *>
410 *> \param[out] ALPHA1
411 *> \verbatim
412 *> ALPHA1 is COMPLEX array, dimension (max(NN))
413 *> \endverbatim
414 *>
415 *> \param[out] BETA1
416 *> \verbatim
417 *> BETA1 is COMPLEX array, dimension (max(NN))
418 *> The generalized eigenvalues of (A,B) computed by CHGEQZ
419 *> when Q, Z, and the full Schur matrices are computed.
420 *> \endverbatim
421 *>
422 *> \param[out] ALPHA3
423 *> \verbatim
424 *> ALPHA3 is COMPLEX array, dimension (max(NN))
425 *> \endverbatim
426 *>
427 *> \param[out] BETA3
428 *> \verbatim
429 *> BETA3 is COMPLEX array, dimension (max(NN))
430 *> The generalized eigenvalues of (A,B) computed by CHGEQZ
431 *> when neither Q, Z, nor the Schur matrices are computed.
432 *> \endverbatim
433 *>
434 *> \param[out] EVECTL
435 *> \verbatim
436 *> EVECTL is COMPLEX 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 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 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 REAL 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 REAL 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 June 2016
496 *
497 *> \ingroup complex_eig
498 *
499 * =====================================================================
500  SUBROUTINE cchkgg( 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.6.1) --
507 * -- LAPACK is a software package provided by Univ. of Tennessee, --
508 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
509 * June 2016
510 *
511 * .. Scalar Arguments ..
512  LOGICAL TSTDIF
513  INTEGER INFO, LDA, LDU, LWORK, NOUNIT, NSIZES, NTYPES
514  REAL THRESH, THRSHN
515 * ..
516 * .. Array Arguments ..
517  LOGICAL DOTYPE( * ), LLWORK( * )
518  INTEGER ISEED( 4 ), NN( * )
519  REAL RESULT( 15 ), RWORK( * )
520  COMPLEX 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  REAL ZERO, ONE
533  parameter ( zero = 0.0e+0, one = 1.0e+0 )
534  COMPLEX CZERO, CONE
535  parameter ( czero = ( 0.0e+0, 0.0e+0 ),
536  $ cone = ( 1.0e+0, 0.0e+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  REAL ANORM, BNORM, SAFMAX, SAFMIN, TEMP1, TEMP2,
546  $ ulp, ulpinv
547  COMPLEX 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  REAL DUMMA( 4 ), RMAGN( 0: 3 )
557  COMPLEX CDUMMA( 4 )
558 * ..
559 * .. External Functions ..
560  REAL CLANGE, SLAMCH
561  COMPLEX CLARND
562  EXTERNAL clange, slamch, clarnd
563 * ..
564 * .. External Subroutines ..
565  EXTERNAL cgeqr2, cget51, cget52, cgghrd, chgeqz, clacpy,
567  $ slasum, xerbla
568 * ..
569 * .. Intrinsic Functions ..
570  INTRINSIC abs, conjg, max, min, REAL, 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( 'CCHKGG', -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 = slamch( 'Safe minimum' )
642  ulp = slamch( 'Epsilon' )*slamch( 'Base' )
643  safmin = safmin / ulp
644  safmax = one / safmin
645  CALL slabad( 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 / REAL( 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 * KCLASS: =1 means w/o rotation, =2 means w/ rotation,
694 * =3 means random.
695 * KATYPE: the "type" to be passed to CLATM4 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 claset( 'Full', n, n, czero, czero, a, lda )
721  ELSE
722  in = n
723  END IF
724  CALL clatm4( 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 claset( 'Full', n, n, czero, czero, b, lda )
739  ELSE
740  in = n
741  END IF
742  CALL clatm4( 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 CLARFG 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 ) = clarnd( 3, iseed )
762  v( jr, jc ) = clarnd( 3, iseed )
763  40 CONTINUE
764  CALL clarfg( n+1-jc, u( jc, jc ), u( jc+1, jc ), 1,
765  $ work( jc ) )
766  work( 2*n+jc ) = sign( one, REAL( U( JC, JC ) ) )
767  u( jc, jc ) = cone
768  CALL clarfg( n+1-jc, v( jc, jc ), v( jc+1, jc ), 1,
769  $ work( n+jc ) )
770  work( 3*n+jc ) = sign( one, REAL( V( JC, JC ) ) )
771  v( jc, jc ) = cone
772  50 CONTINUE
773  ctemp = clarnd( 3, iseed )
774  u( n, n ) = cone
775  work( n ) = czero
776  work( 3*n ) = ctemp / abs( ctemp )
777  ctemp = clarnd( 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  $ conjg( work( 3*n+jc ) )*
788  $ a( jr, jc )
789  b( jr, jc ) = work( 2*n+jr )*
790  $ conjg( work( 3*n+jc ) )*
791  $ b( jr, jc )
792  60 CONTINUE
793  70 CONTINUE
794  CALL cunm2r( '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 cunm2r( '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 cunm2r( '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 cunm2r( '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  $ clarnd( 4, iseed )
819  b( jr, jc ) = rmagn( kbmagn( jtype ) )*
820  $ clarnd( 4, iseed )
821  80 CONTINUE
822  90 CONTINUE
823  END IF
824 *
825  anorm = clange( '1', n, n, a, lda, rwork )
826  bnorm = clange( '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 CGEQR2, CUNM2R, and CGGHRD to compute H, T, U, and V
840 *
841  CALL clacpy( ' ', n, n, a, lda, h, lda )
842  CALL clacpy( ' ', n, n, b, lda, t, lda )
843  ntest = 1
844  result( 1 ) = ulpinv
845 *
846  CALL cgeqr2( n, n, t, lda, work, work( n+1 ), iinfo )
847  IF( iinfo.NE.0 ) THEN
848  WRITE( nounit, fmt = 9999 )'CGEQR2', iinfo, n, jtype,
849  $ ioldsd
850  info = abs( iinfo )
851  GO TO 210
852  END IF
853 *
854  CALL cunm2r( '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 )'CUNM2R', iinfo, n, jtype,
858  $ ioldsd
859  info = abs( iinfo )
860  GO TO 210
861  END IF
862 *
863  CALL claset( 'Full', n, n, czero, cone, u, ldu )
864  CALL cunm2r( '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 )'CUNM2R', iinfo, n, jtype,
868  $ ioldsd
869  info = abs( iinfo )
870  GO TO 210
871  END IF
872 *
873  CALL cgghrd( '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 )'CGGHRD', 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 cget51( 1, n, a, lda, h, lda, u, ldu, v, ldu, work,
886  $ rwork, result( 1 ) )
887  CALL cget51( 1, n, b, lda, t, lda, u, ldu, v, ldu, work,
888  $ rwork, result( 2 ) )
889  CALL cget51( 3, n, b, lda, t, lda, u, ldu, u, ldu, work,
890  $ rwork, result( 3 ) )
891  CALL cget51( 3, n, b, lda, t, lda, v, ldu, v, ldu, work,
892  $ rwork, result( 4 ) )
893 *
894 * Call CHGEQZ to compute S1, P1, S2, P2, Q, and Z, do tests.
895 *
896 * Compute T1 and UZ
897 *
898 * Eigenvalues only
899 *
900  CALL clacpy( ' ', n, n, h, lda, s2, lda )
901  CALL clacpy( ' ', n, n, t, lda, p2, lda )
902  ntest = 5
903  result( 5 ) = ulpinv
904 *
905  CALL chgeqz( '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 )'CHGEQZ(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 clacpy( ' ', n, n, h, lda, s2, lda )
918  CALL clacpy( ' ', n, n, t, lda, p2, lda )
919 *
920  CALL chgeqz( '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 )'CHGEQZ(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 clacpy( ' ', n, n, h, lda, s1, lda )
933  CALL clacpy( ' ', n, n, t, lda, p1, lda )
934 *
935  CALL chgeqz( '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 )'CHGEQZ(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 cget51( 1, n, h, lda, s1, lda, q, ldu, z, ldu, work,
950  $ rwork, result( 5 ) )
951  CALL cget51( 1, n, t, lda, p1, lda, q, ldu, z, ldu, work,
952  $ rwork, result( 6 ) )
953  CALL cget51( 3, n, t, lda, p1, lda, q, ldu, q, ldu, work,
954  $ rwork, result( 7 ) )
955  CALL cget51( 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 ctgevc( '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 )'CTGEVC(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 ctgevc( '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 )'CTGEVC(L,S2)', iinfo, n,
999  $ jtype, ioldsd
1000  info = abs( iinfo )
1001  GO TO 210
1002  END IF
1003 *
1004  CALL cget52( .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', 'CTGEVC(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 clacpy( 'F', n, n, q, ldu, evectl, ldu )
1018  CALL ctgevc( '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 )'CTGEVC(L,B)', iinfo, n,
1022  $ jtype, ioldsd
1023  info = abs( iinfo )
1024  GO TO 210
1025  END IF
1026 *
1027  CALL cget52( .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', 'CTGEVC(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 ctgevc( '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 )'CTGEVC(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 ctgevc( '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 )'CTGEVC(R,S2)', iinfo, n,
1074  $ jtype, ioldsd
1075  info = abs( iinfo )
1076  GO TO 210
1077  END IF
1078 *
1079  CALL cget52( .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', 'CTGEVC(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 clacpy( 'F', n, n, z, ldu, evectr, ldu )
1093  CALL ctgevc( '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 )'CTGEVC(R,B)', iinfo, n,
1097  $ jtype, ioldsd
1098  info = abs( iinfo )
1099  GO TO 210
1100  END IF
1101 *
1102  CALL cget52( .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', 'CTGEVC(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 cget51( 2, n, s1, lda, s2, lda, q, ldu, z, ldu,
1117  $ work, rwork, result( 13 ) )
1118  CALL cget51( 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 )'CGG'
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.0 ) 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 slasum( 'CGG', nounit, nerrs, ntestt )
1187  RETURN
1188 *
1189  9999 FORMAT( ' CCHKGG: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1190  $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1191 *
1192  9998 FORMAT( ' CCHKGG: ', 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 CCHKGG 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, e10.3 )
1241 *
1242 * End of CCHKGG
1243 *
1244  END
subroutine cgeqr2(M, N, A, LDA, TAU, WORK, INFO)
CGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm...
Definition: cgeqr2.f:123
subroutine cchkgg(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)
CCHKGG
Definition: cchkgg.f:505
subroutine cunm2r(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, INFO)
CUNM2R multiplies a general matrix by the unitary matrix from a QR factorization determined by cgeqrf...
Definition: cunm2r.f:161
subroutine slabad(SMALL, LARGE)
SLABAD
Definition: slabad.f:76
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
CGGHRD
Definition: cgghrd.f:206
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
subroutine cget51(ITYPE, N, A, LDA, B, LDB, U, LDU, V, LDV, WORK, RWORK, RESULT)
CGET51
Definition: cget51.f:156
subroutine cget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHA, BETA, WORK, RWORK, RESULT)
CGET52
Definition: cget52.f:163
subroutine clatm4(ITYPE, N, NZ1, NZ2, RSIGN, AMAGN, RCOND, TRIANG, IDIST, ISEED, A, LDA)
CLATM4
Definition: clatm4.f:173
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine ctgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTGEVC
Definition: ctgevc.f:221
subroutine slasum(TYPE, IOUNIT, IE, NRUN)
SLASUM
Definition: slasum.f:42
subroutine clarfg(N, ALPHA, X, INCX, TAU)
CLARFG generates an elementary reflector (Householder matrix).
Definition: clarfg.f:108