LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Classes Files Functions Variables Typedefs Macros
cdrvgg.f
Go to the documentation of this file.
1 *> \brief \b CDRVGG
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 CDRVGG( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
12 * THRSHN, NOUNIT, A, LDA, B, S, T, S2, T2, Q,
13 * LDQ, Z, ALPHA1, BETA1, ALPHA2, BETA2, VL, VR,
14 * WORK, LWORK, RWORK, RESULT, INFO )
15 *
16 * .. Scalar Arguments ..
17 * INTEGER INFO, LDA, LDQ, LWORK, NOUNIT, NSIZES, NTYPES
18 * REAL THRESH, THRSHN
19 * ..
20 * .. Array Arguments ..
21 *
22 *
23 *> \par Purpose:
24 * =============
25 *>
26 *> \verbatim
27 *>
28 *> CDRVGG checks the nonsymmetric generalized eigenvalue driver
29 *> routines.
30 *> T T T
31 *> CGEGS factors A and B as Q S Z and Q T Z , where means
32 *> transpose, T is upper triangular, S is in generalized Schur form
33 *> (upper triangular), and Q and Z are unitary. It also
34 *> computes the generalized eigenvalues (alpha(1),beta(1)), ...,
35 *> (alpha(n),beta(n)), where alpha(j)=S(j,j) and beta(j)=T(j,j) --
36 *> thus, w(j) = alpha(j)/beta(j) is a root of the generalized
37 *> eigenvalue problem
38 *>
39 *> det( A - w(j) B ) = 0
40 *>
41 *> and m(j) = beta(j)/alpha(j) is a root of the essentially equivalent
42 *> problem
43 *>
44 *> det( m(j) A - B ) = 0
45 *>
46 *> CGEGV computes the generalized eigenvalues (alpha(1),beta(1)), ...,
47 *> (alpha(n),beta(n)), the matrix L whose columns contain the
48 *> generalized left eigenvectors l, and the matrix R whose columns
49 *> contain the generalized right eigenvectors r for the pair (A,B).
50 *>
51 *> When CDRVGG is called, a number of matrix "sizes" ("n's") and a
52 *> number of matrix "types" are specified. For each size ("n")
53 *> and each type of matrix, one matrix will be generated and used
54 *> to test the nonsymmetric eigenroutines. For each matrix, 7
55 *> tests will be performed and compared with the threshhold THRESH:
56 *>
57 *> Results from CGEGS:
58 *>
59 *> H
60 *> (1) | A - Q S Z | / ( |A| n ulp )
61 *>
62 *> H
63 *> (2) | B - Q T Z | / ( |B| n ulp )
64 *>
65 *> H
66 *> (3) | I - QQ | / ( n ulp )
67 *>
68 *> H
69 *> (4) | I - ZZ | / ( n ulp )
70 *>
71 *> (5) maximum over j of D(j) where:
72 *>
73 *> |alpha(j) - S(j,j)| |beta(j) - T(j,j)|
74 *> D(j) = ------------------------ + -----------------------
75 *> max(|alpha(j)|,|S(j,j)|) max(|beta(j)|,|T(j,j)|)
76 *>
77 *> Results from CGEGV:
78 *>
79 *> (6) max over all left eigenvalue/-vector pairs (beta/alpha,l) of
80 *>
81 *> | l**H * (beta A - alpha B) | / ( ulp max( |beta A|, |alpha B| ) )
82 *>
83 *> where l**H is the conjugate tranpose of l.
84 *>
85 *> (7) max over all right eigenvalue/-vector pairs (beta/alpha,r) of
86 *>
87 *> | (beta A - alpha B) r | / ( ulp max( |beta A|, |alpha B| ) )
88 *>
89 *> Test Matrices
90 *> ---- --------
91 *>
92 *> The sizes of the test matrices are specified by an array
93 *> NN(1:NSIZES); the value of each element NN(j) specifies one size.
94 *> The "types" are specified by a logical array DOTYPE( 1:NTYPES ); if
95 *> DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
96 *> Currently, the list of possible types is:
97 *>
98 *> (1) ( 0, 0 ) (a pair of zero matrices)
99 *>
100 *> (2) ( I, 0 ) (an identity and a zero matrix)
101 *>
102 *> (3) ( 0, I ) (an identity and a zero matrix)
103 *>
104 *> (4) ( I, I ) (a pair of identity matrices)
105 *>
106 *> t t
107 *> (5) ( J , J ) (a pair of transposed Jordan blocks)
108 *>
109 *> t ( I 0 )
110 *> (6) ( X, Y ) where X = ( J 0 ) and Y = ( t )
111 *> ( 0 I ) ( 0 J )
112 *> and I is a k x k identity and J a (k+1)x(k+1)
113 *> Jordan block; k=(N-1)/2
114 *>
115 *> (7) ( D, I ) where D is diag( 0, 1,..., N-1 ) (a diagonal
116 *> matrix with those diagonal entries.)
117 *> (8) ( I, D )
118 *>
119 *> (9) ( big*D, small*I ) where "big" is near overflow and small=1/big
120 *>
121 *> (10) ( small*D, big*I )
122 *>
123 *> (11) ( big*I, small*D )
124 *>
125 *> (12) ( small*I, big*D )
126 *>
127 *> (13) ( big*D, big*I )
128 *>
129 *> (14) ( small*D, small*I )
130 *>
131 *> (15) ( D1, D2 ) where D1 is diag( 0, 0, 1, ..., N-3, 0 ) and
132 *> D2 is diag( 0, N-3, N-4,..., 1, 0, 0 )
133 *> t t
134 *> (16) Q ( J , J ) Z where Q and Z are random unitary matrices.
135 *>
136 *> (17) Q ( T1, T2 ) Z where T1 and T2 are upper triangular matrices
137 *> with random O(1) entries above the diagonal
138 *> and diagonal entries diag(T1) =
139 *> ( 0, 0, 1, ..., N-3, 0 ) and diag(T2) =
140 *> ( 0, N-3, N-4,..., 1, 0, 0 )
141 *>
142 *> (18) Q ( T1, T2 ) Z diag(T1) = ( 0, 0, 1, 1, s, ..., s, 0 )
143 *> diag(T2) = ( 0, 1, 0, 1,..., 1, 0 )
144 *> s = machine precision.
145 *>
146 *> (19) Q ( T1, T2 ) Z diag(T1)=( 0,0,1,1, 1-d, ..., 1-(N-5)*d=s, 0 )
147 *> diag(T2) = ( 0, 1, 0, 1, ..., 1, 0 )
148 *>
149 *> N-5
150 *> (20) Q ( T1, T2 ) Z diag(T1)=( 0, 0, 1, 1, a, ..., a =s, 0 )
151 *> diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )
152 *>
153 *> (21) Q ( T1, T2 ) Z diag(T1)=( 0, 0, 1, r1, r2, ..., r(N-4), 0 )
154 *> diag(T2) = ( 0, 1, 0, 1, ..., 1, 0, 0 )
155 *> where r1,..., r(N-4) are random.
156 *>
157 *> (22) Q ( big*T1, small*T2 ) Z diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
158 *> diag(T2) = ( 0, 1, ..., 1, 0, 0 )
159 *>
160 *> (23) Q ( small*T1, big*T2 ) Z diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
161 *> diag(T2) = ( 0, 1, ..., 1, 0, 0 )
162 *>
163 *> (24) Q ( small*T1, small*T2 ) Z diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
164 *> diag(T2) = ( 0, 1, ..., 1, 0, 0 )
165 *>
166 *> (25) Q ( big*T1, big*T2 ) Z diag(T1) = ( 0, 0, 1, ..., N-3, 0 )
167 *> diag(T2) = ( 0, 1, ..., 1, 0, 0 )
168 *>
169 *> (26) Q ( T1, T2 ) Z where T1 and T2 are random upper-triangular
170 *> matrices.
171 *> \endverbatim
172 *
173 * Arguments:
174 * ==========
175 *
176 *> \param[in] NSIZES
177 *> \verbatim
178 *> NSIZES is INTEGER
179 *> The number of sizes of matrices to use. If it is zero,
180 *> CDRVGG does nothing. It must be at least zero.
181 *> \endverbatim
182 *>
183 *> \param[in] NN
184 *> \verbatim
185 *> NN is INTEGER array, dimension (NSIZES)
186 *> An array containing the sizes to be used for the matrices.
187 *> Zero values will be skipped. The values must be at least
188 *> zero.
189 *> \endverbatim
190 *>
191 *> \param[in] NTYPES
192 *> \verbatim
193 *> NTYPES is INTEGER
194 *> The number of elements in DOTYPE. If it is zero, CDRVGG
195 *> does nothing. It must be at least zero. If it is MAXTYP+1
196 *> and NSIZES is 1, then an additional type, MAXTYP+1 is
197 *> defined, which is to use whatever matrix is in A. This
198 *> is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
199 *> DOTYPE(MAXTYP+1) is .TRUE. .
200 *> \endverbatim
201 *>
202 *> \param[in] DOTYPE
203 *> \verbatim
204 *> DOTYPE is LOGICAL array, dimension (NTYPES)
205 *> If DOTYPE(j) is .TRUE., then for each size in NN a
206 *> matrix of that size and of type j will be generated.
207 *> If NTYPES is smaller than the maximum number of types
208 *> defined (PARAMETER MAXTYP), then types NTYPES+1 through
209 *> MAXTYP will not be generated. If NTYPES is larger
210 *> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
211 *> will be ignored.
212 *> \endverbatim
213 *>
214 *> \param[in,out] ISEED
215 *> \verbatim
216 *> ISEED is INTEGER array, dimension (4)
217 *> On entry ISEED specifies the seed of the random number
218 *> generator. The array elements should be between 0 and 4095;
219 *> if not they will be reduced mod 4096. Also, ISEED(4) must
220 *> be odd. The random number generator uses a linear
221 *> congruential sequence limited to small integers, and so
222 *> should produce machine independent random numbers. The
223 *> values of ISEED are changed on exit, and can be used in the
224 *> next call to CDRVGG to continue the same random number
225 *> sequence.
226 *> \endverbatim
227 *>
228 *> \param[in] THRESH
229 *> \verbatim
230 *> THRESH is REAL
231 *> A test will count as "failed" if the "error", computed as
232 *> described above, exceeds THRESH. Note that the error is
233 *> scaled to be O(1), so THRESH should be a reasonably small
234 *> multiple of 1, e.g., 10 or 100. In particular, it should
235 *> not depend on the precision (single vs. double) or the size
236 *> of the matrix. It must be at least zero.
237 *> \endverbatim
238 *>
239 *> \param[in] THRSHN
240 *> \verbatim
241 *> THRSHN is REAL
242 *> Threshhold for reporting eigenvector normalization error.
243 *> If the normalization of any eigenvector differs from 1 by
244 *> more than THRSHN*ulp, then a special error message will be
245 *> printed. (This is handled separately from the other tests,
246 *> since only a compiler or programming error should cause an
247 *> error message, at least if THRSHN is at least 5--10.)
248 *> \endverbatim
249 *>
250 *> \param[in] NOUNIT
251 *> \verbatim
252 *> NOUNIT is INTEGER
253 *> The FORTRAN unit number for printing out error messages
254 *> (e.g., if a routine returns IINFO not equal to 0.)
255 *> \endverbatim
256 *>
257 *> \param[in,out] A
258 *> \verbatim
259 *> A is COMPLEX array, dimension (LDA, max(NN))
260 *> Used to hold the original A matrix. Used as input only
261 *> if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
262 *> DOTYPE(MAXTYP+1)=.TRUE.
263 *> \endverbatim
264 *>
265 *> \param[in] LDA
266 *> \verbatim
267 *> LDA is INTEGER
268 *> The leading dimension of A, B, S, T, S2, and T2.
269 *> It must be at least 1 and at least max( NN ).
270 *> \endverbatim
271 *>
272 *> \param[in,out] B
273 *> \verbatim
274 *> B is COMPLEX array, dimension (LDA, max(NN))
275 *> Used to hold the original B matrix. Used as input only
276 *> if NTYPES=MAXTYP+1, DOTYPE(1:MAXTYP)=.FALSE., and
277 *> DOTYPE(MAXTYP+1)=.TRUE.
278 *> \endverbatim
279 *>
280 *> \param[out] S
281 *> \verbatim
282 *> S is COMPLEX array, dimension (LDA, max(NN))
283 *> The upper triangular matrix computed from A by CGEGS.
284 *> \endverbatim
285 *>
286 *> \param[out] T
287 *> \verbatim
288 *> T is COMPLEX array, dimension (LDA, max(NN))
289 *> The upper triangular matrix computed from B by CGEGS.
290 *> \endverbatim
291 *>
292 *> \param[out] S2
293 *> \verbatim
294 *> S2 is COMPLEX array, dimension (LDA, max(NN))
295 *> The matrix computed from A by CGEGV. This will be the
296 *> Schur (upper triangular) form of some matrix related to A,
297 *> but will not, in general, be the same as S.
298 *> \endverbatim
299 *>
300 *> \param[out] T2
301 *> \verbatim
302 *> T2 is COMPLEX array, dimension (LDA, max(NN))
303 *> The matrix computed from B by CGEGV. This will be the
304 *> Schur form of some matrix related to B, but will not, in
305 *> general, be the same as T.
306 *> \endverbatim
307 *>
308 *> \param[out] Q
309 *> \verbatim
310 *> Q is COMPLEX array, dimension (LDQ, max(NN))
311 *> The (left) unitary matrix computed by CGEGS.
312 *> \endverbatim
313 *>
314 *> \param[in] LDQ
315 *> \verbatim
316 *> LDQ is INTEGER
317 *> The leading dimension of Q, Z, VL, and VR. It must
318 *> be at least 1 and at least max( NN ).
319 *> \endverbatim
320 *>
321 *> \param[out] Z
322 *> \verbatim
323 *> Z is COMPLEX array, dimension (LDQ, max(NN))
324 *> The (right) unitary matrix computed by CGEGS.
325 *> \endverbatim
326 *>
327 *> \param[out] ALPHA1
328 *> \verbatim
329 *> ALPHA1 is COMPLEX array, dimension (max(NN))
330 *> \endverbatim
331 *>
332 *> \param[out] BETA1
333 *> \verbatim
334 *> BETA1 is COMPLEX array, dimension (max(NN))
335 *>
336 *> The generalized eigenvalues of (A,B) computed by CGEGS.
337 *> ALPHA1(k) / BETA1(k) is the k-th generalized eigenvalue of
338 *> the matrices in A and B.
339 *> \endverbatim
340 *>
341 *> \param[out] ALPHA2
342 *> \verbatim
343 *> ALPHA2 is COMPLEX array, dimension (max(NN))
344 *> \endverbatim
345 *>
346 *> \param[out] BETA2
347 *> \verbatim
348 *> BETA2 is COMPLEX array, dimension (max(NN))
349 *>
350 *> The generalized eigenvalues of (A,B) computed by CGEGV.
351 *> ALPHA2(k) / BETA2(k) is the k-th generalized eigenvalue of
352 *> the matrices in A and B.
353 *> \endverbatim
354 *>
355 *> \param[out] VL
356 *> \verbatim
357 *> VL is COMPLEX array, dimension (LDQ, max(NN))
358 *> The (lower triangular) left eigenvector matrix for the
359 *> matrices in A and B.
360 *> \endverbatim
361 *>
362 *> \param[out] VR
363 *> \verbatim
364 *> VR is COMPLEX array, dimension (LDQ, max(NN))
365 *> The (upper triangular) right eigenvector matrix for the
366 *> matrices in A and B.
367 *> \endverbatim
368 *>
369 *> \param[out] WORK
370 *> \verbatim
371 *> WORK is COMPLEX array, dimension (LWORK)
372 *> \endverbatim
373 *>
374 *> \param[in] LWORK
375 *> \verbatim
376 *> LWORK is INTEGER
377 *> The number of entries in WORK. This must be at least
378 *> MAX( 2*N, N*(NB+1), (k+1)*(2*k+N+1) ), where "k" is the
379 *> sum of the blocksize and number-of-shifts for CHGEQZ, and
380 *> NB is the greatest of the blocksizes for CGEQRF, CUNMQR,
381 *> and CUNGQR. (The blocksizes and the number-of-shifts are
382 *> retrieved through calls to ILAENV.)
383 *> \endverbatim
384 *>
385 *> \param[out] RWORK
386 *> \verbatim
387 *> RWORK is REAL array, dimension (8*N)
388 *> \endverbatim
389 *>
390 *> \param[out] RESULT
391 *> \verbatim
392 *> RESULT is REAL array, dimension (7)
393 *> The values computed by the tests described above.
394 *> The values are currently limited to 1/ulp, to avoid
395 *> overflow.
396 *> \endverbatim
397 *>
398 *> \param[out] INFO
399 *> \verbatim
400 *> INFO is INTEGER
401 *> = 0: successful exit
402 *> < 0: if INFO = -i, the i-th argument had an illegal value.
403 *> > 0: A routine returned an error code. INFO is the
404 *> absolute value of the INFO value returned.
405 *> \endverbatim
406 *
407 * Authors:
408 * ========
409 *
410 *> \author Univ. of Tennessee
411 *> \author Univ. of California Berkeley
412 *> \author Univ. of Colorado Denver
413 *> \author NAG Ltd.
414 *
415 *> \date November 2011
416 *
417 *> \ingroup complex_eig
418 *
419 * =====================================================================
420  SUBROUTINE cdrvgg( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
421  $ thrshn, nounit, a, lda, b, s, t, s2, t2, q,
422  $ ldq, z, alpha1, beta1, alpha2, beta2, vl, vr,
423  $ work, lwork, rwork, result, info )
424 *
425 * -- LAPACK test routine (version 3.4.0) --
426 * -- LAPACK is a software package provided by Univ. of Tennessee, --
427 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
428 * November 2011
429 *
430 * .. Scalar Arguments ..
431  INTEGER info, lda, ldq, lwork, nounit, nsizes, ntypes
432  REAL thresh, thrshn
433 * ..
434 * .. Array Arguments ..
435 *
436 * =====================================================================
437 *
438  LOGICAL dotype( * )
439  INTEGER iseed( 4 ), nn( * )
440  REAL result( * ), rwork( * )
441  COMPLEX a( lda, * ), alpha1( * ), alpha2( * ),
442  $ b( lda, * ), beta1( * ), beta2( * ),
443  $ q( ldq, * ), s( lda, * ), s2( lda, * ),
444  $ t( lda, * ), t2( lda, * ), vl( ldq, * ),
445  $ vr( ldq, * ), work( * ), z( ldq, * )
446 * ..
447 * .. Parameters ..
448  REAL zero, one
449  parameter( zero = 0.0e+0, one = 1.0e+0 )
450  COMPLEX czero, cone
451  parameter( czero = ( 0.0e+0, 0.0e+0 ),
452  $ cone = ( 1.0e+0, 0.0e+0 ) )
453  INTEGER maxtyp
454  parameter( maxtyp = 26 )
455 * ..
456 * .. Local Scalars ..
457  LOGICAL badnn
458  INTEGER i1, iadd, iinfo, in, j, jc, jr, jsize, jtype,
459  $ lwkopt, mtypes, n, n1, nb, nbz, nerrs, nmats,
460  $ nmax, ns, ntest, ntestt
461  REAL safmax, safmin, temp1, temp2, ulp, ulpinv
462  COMPLEX ctemp, x
463 * ..
464 * .. Local Arrays ..
465  LOGICAL lasign( maxtyp ), lbsign( maxtyp )
466  INTEGER ioldsd( 4 ), kadd( 6 ), kamagn( maxtyp ),
467  $ katype( maxtyp ), kazero( maxtyp ),
468  $ kbmagn( maxtyp ), kbtype( maxtyp ),
469  $ kbzero( maxtyp ), kclass( maxtyp ),
470  $ ktrian( maxtyp ), kz1( 6 ), kz2( 6 )
471  REAL dumma( 4 ), rmagn( 0: 3 )
472 * ..
473 * .. External Functions ..
474  INTEGER ilaenv
475  REAL slamch
476  COMPLEX clarnd
477  EXTERNAL ilaenv, slamch, clarnd
478 * ..
479 * .. External Subroutines ..
480  EXTERNAL alasvm, cgegs, cgegv, cget51, cget52, clacpy,
482 * ..
483 * .. Intrinsic Functions ..
484  INTRINSIC abs, aimag, conjg, max, min, REAL, sign
485 * ..
486 * .. Statement Functions ..
487  REAL abs1
488 * ..
489 * .. Statement Function definitions ..
490  abs1( x ) = abs( REAL( X ) ) + abs( aimag( x ) )
491 * ..
492 * .. Data statements ..
493  DATA kclass / 15*1, 10*2, 1*3 /
494  DATA kz1 / 0, 1, 2, 1, 3, 3 /
495  DATA kz2 / 0, 0, 1, 2, 1, 1 /
496  DATA kadd / 0, 0, 0, 0, 3, 2 /
497  DATA katype / 0, 1, 0, 1, 2, 3, 4, 1, 4, 4, 1, 1, 4,
498  $ 4, 4, 2, 4, 5, 8, 7, 9, 4*4, 0 /
499  DATA kbtype / 0, 0, 1, 1, 2, -3, 1, 4, 1, 1, 4, 4,
500  $ 1, 1, -4, 2, -4, 8*8, 0 /
501  DATA kazero / 6*1, 2, 1, 2*2, 2*1, 2*2, 3, 1, 3,
502  $ 4*5, 4*3, 1 /
503  DATA kbzero / 6*1, 1, 2, 2*1, 2*2, 2*1, 4, 1, 4,
504  $ 4*6, 4*4, 1 /
505  DATA kamagn / 8*1, 2, 3, 2, 3, 2, 3, 7*1, 2, 3, 3,
506  $ 2, 1 /
507  DATA kbmagn / 8*1, 3, 2, 3, 2, 2, 3, 7*1, 3, 2, 3,
508  $ 2, 1 /
509  DATA ktrian / 16*0, 10*1 /
510  DATA lasign / 6*.false., .true., .false., 2*.true.,
511  $ 2*.false., 3*.true., .false., .true.,
512  $ 3*.false., 5*.true., .false. /
513  DATA lbsign / 7*.false., .true., 2*.false.,
514  $ 2*.true., 2*.false., .true., .false., .true.,
515  $ 9*.false. /
516 * ..
517 * .. Executable Statements ..
518 *
519 * Check for errors
520 *
521  info = 0
522 *
523  badnn = .false.
524  nmax = 1
525  DO 10 j = 1, nsizes
526  nmax = max( nmax, nn( j ) )
527  IF( nn( j ).LT.0 )
528  $ badnn = .true.
529  10 CONTINUE
530 *
531 * Maximum blocksize and shift -- we assume that blocksize and number
532 * of shifts are monotone increasing functions of N.
533 *
534  nb = max( 1, ilaenv( 1, 'CGEQRF', ' ', nmax, nmax, -1, -1 ),
535  $ ilaenv( 1, 'CUNMQR', 'LC', nmax, nmax, nmax, -1 ),
536  $ ilaenv( 1, 'CUNGQR', ' ', nmax, nmax, nmax, -1 ) )
537  nbz = ilaenv( 1, 'CHGEQZ', 'SII', nmax, 1, nmax, 0 )
538  ns = ilaenv( 4, 'CHGEQZ', 'SII', nmax, 1, nmax, 0 )
539  i1 = nbz + ns
540  lwkopt = max( 2*nmax, nmax*( nb+1 ), ( 2*i1+nmax+1 )*( i1+1 ) )
541 *
542 * Check for errors
543 *
544  IF( nsizes.LT.0 ) THEN
545  info = -1
546  ELSE IF( badnn ) THEN
547  info = -2
548  ELSE IF( ntypes.LT.0 ) THEN
549  info = -3
550  ELSE IF( thresh.LT.zero ) THEN
551  info = -6
552  ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
553  info = -10
554  ELSE IF( ldq.LE.1 .OR. ldq.LT.nmax ) THEN
555  info = -19
556  ELSE IF( lwkopt.GT.lwork ) THEN
557  info = -30
558  END IF
559 *
560  IF( info.NE.0 ) THEN
561  CALL xerbla( 'CDRVGG', -info )
562  RETURN
563  END IF
564 *
565 * Quick return if possible
566 *
567  IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
568  $ RETURN
569 *
570  ulp = slamch( 'Precision' )
571  safmin = slamch( 'Safe minimum' )
572  safmin = safmin / ulp
573  safmax = one / safmin
574  CALL slabad( safmin, safmax )
575  ulpinv = one / ulp
576 *
577 * The values RMAGN(2:3) depend on N, see below.
578 *
579  rmagn( 0 ) = zero
580  rmagn( 1 ) = one
581 *
582 * Loop over sizes, types
583 *
584  ntestt = 0
585  nerrs = 0
586  nmats = 0
587 *
588  DO 160 jsize = 1, nsizes
589  n = nn( jsize )
590  n1 = max( 1, n )
591  rmagn( 2 ) = safmax*ulp / REAL( n1 )
592  rmagn( 3 ) = safmin*ulpinv*n1
593 *
594  IF( nsizes.NE.1 ) THEN
595  mtypes = min( maxtyp, ntypes )
596  ELSE
597  mtypes = min( maxtyp+1, ntypes )
598  END IF
599 *
600  DO 150 jtype = 1, mtypes
601  IF( .NOT.dotype( jtype ) )
602  $ go to 150
603  nmats = nmats + 1
604  ntest = 0
605 *
606 * Save ISEED in case of an error.
607 *
608  DO 20 j = 1, 4
609  ioldsd( j ) = iseed( j )
610  20 CONTINUE
611 *
612 * Initialize RESULT
613 *
614  DO 30 j = 1, 7
615  result( j ) = zero
616  30 CONTINUE
617 *
618 * Compute A and B
619 *
620 * Description of control parameters:
621 *
622 * KCLASS: =1 means w/o rotation, =2 means w/ rotation,
623 * =3 means random.
624 * KATYPE: the "type" to be passed to CLATM4 for computing A.
625 * KAZERO: the pattern of zeros on the diagonal for A:
626 * =1: ( xxx ), =2: (0, xxx ) =3: ( 0, 0, xxx, 0 ),
627 * =4: ( 0, xxx, 0, 0 ), =5: ( 0, 0, 1, xxx, 0 ),
628 * =6: ( 0, 1, 0, xxx, 0 ). (xxx means a string of
629 * non-zero entries.)
630 * KAMAGN: the magnitude of the matrix: =0: zero, =1: O(1),
631 * =2: large, =3: small.
632 * LASIGN: .TRUE. if the diagonal elements of A are to be
633 * multiplied by a random magnitude 1 number.
634 * KBTYPE, KBZERO, KBMAGN, IBSIGN: the same, but for B.
635 * KTRIAN: =0: don't fill in the upper triangle, =1: do.
636 * KZ1, KZ2, KADD: used to implement KAZERO and KBZERO.
637 * RMAGN: used to implement KAMAGN and KBMAGN.
638 *
639  IF( mtypes.GT.maxtyp )
640  $ go to 110
641  iinfo = 0
642  IF( kclass( jtype ).LT.3 ) THEN
643 *
644 * Generate A (w/o rotation)
645 *
646  IF( abs( katype( jtype ) ).EQ.3 ) THEN
647  in = 2*( ( n-1 ) / 2 ) + 1
648  IF( in.NE.n )
649  $ CALL claset( 'Full', n, n, czero, czero, a, lda )
650  ELSE
651  in = n
652  END IF
653  CALL clatm4( katype( jtype ), in, kz1( kazero( jtype ) ),
654  $ kz2( kazero( jtype ) ), lasign( jtype ),
655  $ rmagn( kamagn( jtype ) ), ulp,
656  $ rmagn( ktrian( jtype )*kamagn( jtype ) ), 2,
657  $ iseed, a, lda )
658  iadd = kadd( kazero( jtype ) )
659  IF( iadd.GT.0 .AND. iadd.LE.n )
660  $ a( iadd, iadd ) = rmagn( kamagn( jtype ) )
661 *
662 * Generate B (w/o rotation)
663 *
664  IF( abs( kbtype( jtype ) ).EQ.3 ) THEN
665  in = 2*( ( n-1 ) / 2 ) + 1
666  IF( in.NE.n )
667  $ CALL claset( 'Full', n, n, czero, czero, b, lda )
668  ELSE
669  in = n
670  END IF
671  CALL clatm4( kbtype( jtype ), in, kz1( kbzero( jtype ) ),
672  $ kz2( kbzero( jtype ) ), lbsign( jtype ),
673  $ rmagn( kbmagn( jtype ) ), one,
674  $ rmagn( ktrian( jtype )*kbmagn( jtype ) ), 2,
675  $ iseed, b, lda )
676  iadd = kadd( kbzero( jtype ) )
677  IF( iadd.NE.0 .AND. iadd.LE.n )
678  $ b( iadd, iadd ) = rmagn( kbmagn( jtype ) )
679 *
680  IF( kclass( jtype ).EQ.2 .AND. n.GT.0 ) THEN
681 *
682 * Include rotations
683 *
684 * Generate Q, Z as Householder transformations times
685 * a diagonal matrix.
686 *
687  DO 50 jc = 1, n - 1
688  DO 40 jr = jc, n
689  q( jr, jc ) = clarnd( 3, iseed )
690  z( jr, jc ) = clarnd( 3, iseed )
691  40 CONTINUE
692  CALL clarfg( n+1-jc, q( jc, jc ), q( jc+1, jc ), 1,
693  $ work( jc ) )
694  work( 2*n+jc ) = sign( one, REAL( Q( JC, JC ) ) )
695  q( jc, jc ) = cone
696  CALL clarfg( n+1-jc, z( jc, jc ), z( jc+1, jc ), 1,
697  $ work( n+jc ) )
698  work( 3*n+jc ) = sign( one, REAL( Z( JC, JC ) ) )
699  z( jc, jc ) = cone
700  50 CONTINUE
701  ctemp = clarnd( 3, iseed )
702  q( n, n ) = cone
703  work( n ) = czero
704  work( 3*n ) = ctemp / abs( ctemp )
705  ctemp = clarnd( 3, iseed )
706  z( n, n ) = cone
707  work( 2*n ) = czero
708  work( 4*n ) = ctemp / abs( ctemp )
709 *
710 * Apply the diagonal matrices
711 *
712  DO 70 jc = 1, n
713  DO 60 jr = 1, n
714  a( jr, jc ) = work( 2*n+jr )*
715  $ conjg( work( 3*n+jc ) )*
716  $ a( jr, jc )
717  b( jr, jc ) = work( 2*n+jr )*
718  $ conjg( work( 3*n+jc ) )*
719  $ b( jr, jc )
720  60 CONTINUE
721  70 CONTINUE
722  CALL cunm2r( 'L', 'N', n, n, n-1, q, ldq, work, a,
723  $ lda, work( 2*n+1 ), iinfo )
724  IF( iinfo.NE.0 )
725  $ go to 100
726  CALL cunm2r( 'R', 'C', n, n, n-1, z, ldq, work( n+1 ),
727  $ a, lda, work( 2*n+1 ), iinfo )
728  IF( iinfo.NE.0 )
729  $ go to 100
730  CALL cunm2r( 'L', 'N', n, n, n-1, q, ldq, work, b,
731  $ lda, work( 2*n+1 ), iinfo )
732  IF( iinfo.NE.0 )
733  $ go to 100
734  CALL cunm2r( 'R', 'C', n, n, n-1, z, ldq, work( n+1 ),
735  $ b, lda, work( 2*n+1 ), iinfo )
736  IF( iinfo.NE.0 )
737  $ go to 100
738  END IF
739  ELSE
740 *
741 * Random matrices
742 *
743  DO 90 jc = 1, n
744  DO 80 jr = 1, n
745  a( jr, jc ) = rmagn( kamagn( jtype ) )*
746  $ clarnd( 4, iseed )
747  b( jr, jc ) = rmagn( kbmagn( jtype ) )*
748  $ clarnd( 4, iseed )
749  80 CONTINUE
750  90 CONTINUE
751  END IF
752 *
753  100 CONTINUE
754 *
755  IF( iinfo.NE.0 ) THEN
756  WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
757  $ ioldsd
758  info = abs( iinfo )
759  RETURN
760  END IF
761 *
762  110 CONTINUE
763 *
764 * Call CGEGS to compute H, T, Q, Z, alpha, and beta.
765 *
766  CALL clacpy( ' ', n, n, a, lda, s, lda )
767  CALL clacpy( ' ', n, n, b, lda, t, lda )
768  ntest = 1
769  result( 1 ) = ulpinv
770 *
771  CALL cgegs( 'V', 'V', n, s, lda, t, lda, alpha1, beta1, q,
772  $ ldq, z, ldq, work, lwork, rwork, iinfo )
773  IF( iinfo.NE.0 ) THEN
774  WRITE( nounit, fmt = 9999 )'CGEGS', iinfo, n, jtype,
775  $ ioldsd
776  info = abs( iinfo )
777  go to 130
778  END IF
779 *
780  ntest = 4
781 *
782 * Do tests 1--4
783 *
784  CALL cget51( 1, n, a, lda, s, lda, q, ldq, z, ldq, work,
785  $ rwork, result( 1 ) )
786  CALL cget51( 1, n, b, lda, t, lda, q, ldq, z, ldq, work,
787  $ rwork, result( 2 ) )
788  CALL cget51( 3, n, b, lda, t, lda, q, ldq, q, ldq, work,
789  $ rwork, result( 3 ) )
790  CALL cget51( 3, n, b, lda, t, lda, z, ldq, z, ldq, work,
791  $ rwork, result( 4 ) )
792 *
793 * Do test 5: compare eigenvalues with diagonals.
794 *
795  temp1 = zero
796 *
797  DO 120 j = 1, n
798  temp2 = ( abs1( alpha1( j )-s( j, j ) ) /
799  $ max( safmin, abs1( alpha1( j ) ), abs1( s( j,
800  $ j ) ) )+abs1( beta1( j )-t( j, j ) ) /
801  $ max( safmin, abs1( beta1( j ) ), abs1( t( j,
802  $ j ) ) ) ) / ulp
803  temp1 = max( temp1, temp2 )
804  120 CONTINUE
805  result( 5 ) = temp1
806 *
807 * Call CGEGV to compute S2, T2, VL, and VR, do tests.
808 *
809 * Eigenvalues and Eigenvectors
810 *
811  CALL clacpy( ' ', n, n, a, lda, s2, lda )
812  CALL clacpy( ' ', n, n, b, lda, t2, lda )
813  ntest = 6
814  result( 6 ) = ulpinv
815 *
816  CALL cgegv( 'V', 'V', n, s2, lda, t2, lda, alpha2, beta2,
817  $ vl, ldq, vr, ldq, work, lwork, rwork, iinfo )
818  IF( iinfo.NE.0 ) THEN
819  WRITE( nounit, fmt = 9999 )'CGEGV', iinfo, n, jtype,
820  $ ioldsd
821  info = abs( iinfo )
822  go to 130
823  END IF
824 *
825  ntest = 7
826 *
827 * Do Tests 6 and 7
828 *
829  CALL cget52( .true., n, a, lda, b, lda, vl, ldq, alpha2,
830  $ beta2, work, rwork, dumma( 1 ) )
831  result( 6 ) = dumma( 1 )
832  IF( dumma( 2 ).GT.thrshn ) THEN
833  WRITE( nounit, fmt = 9998 )'Left', 'CGEGV', dumma( 2 ),
834  $ n, jtype, ioldsd
835  END IF
836 *
837  CALL cget52( .false., n, a, lda, b, lda, vr, ldq, alpha2,
838  $ beta2, work, rwork, dumma( 1 ) )
839  result( 7 ) = dumma( 1 )
840  IF( dumma( 2 ).GT.thresh ) THEN
841  WRITE( nounit, fmt = 9998 )'Right', 'CGEGV', dumma( 2 ),
842  $ n, jtype, ioldsd
843  END IF
844 *
845 * End of Loop -- Check for RESULT(j) > THRESH
846 *
847  130 CONTINUE
848 *
849  ntestt = ntestt + ntest
850 *
851 * Print out tests which fail.
852 *
853  DO 140 jr = 1, ntest
854  IF( result( jr ).GE.thresh ) THEN
855 *
856 * If this is the first test to fail,
857 * print a header to the data file.
858 *
859  IF( nerrs.EQ.0 ) THEN
860  WRITE( nounit, fmt = 9997 )'CGG'
861 *
862 * Matrix types
863 *
864  WRITE( nounit, fmt = 9996 )
865  WRITE( nounit, fmt = 9995 )
866  WRITE( nounit, fmt = 9994 )'Unitary'
867 *
868 * Tests performed
869 *
870  WRITE( nounit, fmt = 9993 )'unitary', '*',
871  $ 'conjugate transpose', ( '*', j = 1, 5 )
872 *
873  END IF
874  nerrs = nerrs + 1
875  IF( result( jr ).LT.10000.0 ) THEN
876  WRITE( nounit, fmt = 9992 )n, jtype, ioldsd, jr,
877  $ result( jr )
878  ELSE
879  WRITE( nounit, fmt = 9991 )n, jtype, ioldsd, jr,
880  $ result( jr )
881  END IF
882  END IF
883  140 CONTINUE
884 *
885  150 CONTINUE
886  160 CONTINUE
887 *
888 * Summary
889 *
890  CALL alasvm( 'CGG', nounit, nerrs, ntestt, 0 )
891  RETURN
892 *
893  9999 FORMAT( ' CDRVGG: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
894  $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
895 *
896  9998 FORMAT( ' CDRVGG: ', a, ' Eigenvectors from ', a, ' incorrectly ',
897  $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
898  $ 'N=', i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5,
899  $ ')' )
900 *
901  9997 FORMAT( / 1x, a3,
902  $ ' -- Complex Generalized eigenvalue problem driver' )
903 *
904  9996 FORMAT( ' Matrix types (see CDRVGG for details): ' )
905 *
906  9995 FORMAT( ' Special Matrices:', 23x,
907  $ '(J''=transposed Jordan block)',
908  $ / ' 1=(0,0) 2=(I,0) 3=(0,I) 4=(I,I) 5=(J'',J'') ',
909  $ '6=(diag(J'',I), diag(I,J''))', / ' Diagonal Matrices: ( ',
910  $ 'D=diag(0,1,2,...) )', / ' 7=(D,I) 9=(large*D, small*I',
911  $ ') 11=(large*I, small*D) 13=(large*D, large*I)', /
912  $ ' 8=(I,D) 10=(small*D, large*I) 12=(small*I, large*D) ',
913  $ ' 14=(small*D, small*I)', / ' 15=(D, reversed D)' )
914  9994 FORMAT( ' Matrices Rotated by Random ', a, ' Matrices U, V:',
915  $ / ' 16=Transposed Jordan Blocks 19=geometric ',
916  $ 'alpha, beta=0,1', / ' 17=arithm. alpha&beta ',
917  $ ' 20=arithmetic alpha, beta=0,1', / ' 18=clustered ',
918  $ 'alpha, beta=0,1 21=random alpha, beta=0,1',
919  $ / ' Large & Small Matrices:', / ' 22=(large, small) ',
920  $ '23=(small,large) 24=(small,small) 25=(large,large)',
921  $ / ' 26=random O(1) matrices.' )
922 *
923  9993 FORMAT( / ' Tests performed: (S is Schur, T is triangular, ',
924  $ 'Q and Z are ', a, ',', / 20x,
925  $ 'l and r are the appropriate left and right', / 19x,
926  $ 'eigenvectors, resp., a is alpha, b is beta, and', / 19x, a,
927  $ ' means ', a, '.)', / ' 1 = | A - Q S Z', a,
928  $ ' | / ( |A| n ulp ) 2 = | B - Q T Z', a,
929  $ ' | / ( |B| n ulp )', / ' 3 = | I - QQ', a,
930  $ ' | / ( n ulp ) 4 = | I - ZZ', a,
931  $ ' | / ( n ulp )', /
932  $ ' 5 = difference between (alpha,beta) and diagonals of',
933  $ ' (S,T)', / ' 6 = max | ( b A - a B )', a,
934  $ ' l | / const. 7 = max | ( b A - a B ) r | / const.',
935  $ / 1x )
936  9992 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
937  $ 4( i4, ',' ), ' result ', i3, ' is', 0p, f8.2 )
938  9991 FORMAT( ' Matrix order=', i5, ', type=', i2, ', seed=',
939  $ 4( i4, ',' ), ' result ', i3, ' is', 1p, e10.3 )
940 *
941 * End of CDRVGG
942 *
943  END