LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Classes Files Functions Variables Typedefs Macros
cdrvvx.f
Go to the documentation of this file.
1 *> \brief \b CDRVVX
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 CDRVVX( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
12 * NIUNIT, NOUNIT, A, LDA, H, W, W1, VL, LDVL, VR,
13 * LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN,
14 * RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT,
15 * WORK, NWORK, RWORK, INFO )
16 *
17 * .. Scalar Arguments ..
18 * INTEGER INFO, LDA, LDLRE, LDVL, LDVR, NIUNIT, NOUNIT,
19 * $ NSIZES, NTYPES, NWORK
20 * REAL THRESH
21 * ..
22 * .. Array Arguments ..
23 * LOGICAL DOTYPE( * )
24 * INTEGER ISEED( 4 ), NN( * )
25 * REAL RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
26 * $ RCNDV1( * ), RCONDE( * ), RCONDV( * ),
27 * $ RESULT( 11 ), RWORK( * ), SCALE( * ),
28 * $ SCALE1( * )
29 * COMPLEX A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
30 * $ VL( LDVL, * ), VR( LDVR, * ), W( * ), W1( * ),
31 * $ WORK( * )
32 * ..
33 *
34 *
35 *> \par Purpose:
36 * =============
37 *>
38 *> \verbatim
39 *>
40 *> CDRVVX checks the nonsymmetric eigenvalue problem expert driver
41 *> CGEEVX.
42 *>
43 *> CDRVVX uses both test matrices generated randomly depending on
44 *> data supplied in the calling sequence, as well as on data
45 *> read from an input file and including precomputed condition
46 *> numbers to which it compares the ones it computes.
47 *>
48 *> When CDRVVX is called, a number of matrix "sizes" ("n's") and a
49 *> number of matrix "types" are specified in the calling sequence.
50 *> For each size ("n") and each type of matrix, one matrix will be
51 *> generated and used to test the nonsymmetric eigenroutines. For
52 *> each matrix, 9 tests will be performed:
53 *>
54 *> (1) | A * VR - VR * W | / ( n |A| ulp )
55 *>
56 *> Here VR is the matrix of unit right eigenvectors.
57 *> W is a diagonal matrix with diagonal entries W(j).
58 *>
59 *> (2) | A**H * VL - VL * W**H | / ( n |A| ulp )
60 *>
61 *> Here VL is the matrix of unit left eigenvectors, A**H is the
62 *> conjugate transpose of A, and W is as above.
63 *>
64 *> (3) | |VR(i)| - 1 | / ulp and largest component real
65 *>
66 *> VR(i) denotes the i-th column of VR.
67 *>
68 *> (4) | |VL(i)| - 1 | / ulp and largest component real
69 *>
70 *> VL(i) denotes the i-th column of VL.
71 *>
72 *> (5) W(full) = W(partial)
73 *>
74 *> W(full) denotes the eigenvalues computed when VR, VL, RCONDV
75 *> and RCONDE are also computed, and W(partial) denotes the
76 *> eigenvalues computed when only some of VR, VL, RCONDV, and
77 *> RCONDE are computed.
78 *>
79 *> (6) VR(full) = VR(partial)
80 *>
81 *> VR(full) denotes the right eigenvectors computed when VL, RCONDV
82 *> and RCONDE are computed, and VR(partial) denotes the result
83 *> when only some of VL and RCONDV are computed.
84 *>
85 *> (7) VL(full) = VL(partial)
86 *>
87 *> VL(full) denotes the left eigenvectors computed when VR, RCONDV
88 *> and RCONDE are computed, and VL(partial) denotes the result
89 *> when only some of VR and RCONDV are computed.
90 *>
91 *> (8) 0 if SCALE, ILO, IHI, ABNRM (full) =
92 *> SCALE, ILO, IHI, ABNRM (partial)
93 *> 1/ulp otherwise
94 *>
95 *> SCALE, ILO, IHI and ABNRM describe how the matrix is balanced.
96 *> (full) is when VR, VL, RCONDE and RCONDV are also computed, and
97 *> (partial) is when some are not computed.
98 *>
99 *> (9) RCONDV(full) = RCONDV(partial)
100 *>
101 *> RCONDV(full) denotes the reciprocal condition numbers of the
102 *> right eigenvectors computed when VR, VL and RCONDE are also
103 *> computed. RCONDV(partial) denotes the reciprocal condition
104 *> numbers when only some of VR, VL and RCONDE are computed.
105 *>
106 *> The "sizes" are specified by an array NN(1:NSIZES); the value of
107 *> each element NN(j) specifies one size.
108 *> The "types" are specified by a logical array DOTYPE( 1:NTYPES );
109 *> if DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
110 *> Currently, the list of possible types is:
111 *>
112 *> (1) The zero matrix.
113 *> (2) The identity matrix.
114 *> (3) A (transposed) Jordan block, with 1's on the diagonal.
115 *>
116 *> (4) A diagonal matrix with evenly spaced entries
117 *> 1, ..., ULP and random complex angles.
118 *> (ULP = (first number larger than 1) - 1 )
119 *> (5) A diagonal matrix with geometrically spaced entries
120 *> 1, ..., ULP and random complex angles.
121 *> (6) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
122 *> and random complex angles.
123 *>
124 *> (7) Same as (4), but multiplied by a constant near
125 *> the overflow threshold
126 *> (8) Same as (4), but multiplied by a constant near
127 *> the underflow threshold
128 *>
129 *> (9) A matrix of the form U' T U, where U is unitary and
130 *> T has evenly spaced entries 1, ..., ULP with random complex
131 *> angles on the diagonal and random O(1) entries in the upper
132 *> triangle.
133 *>
134 *> (10) A matrix of the form U' T U, where U is unitary and
135 *> T has geometrically spaced entries 1, ..., ULP with random
136 *> complex angles on the diagonal and random O(1) entries in
137 *> the upper triangle.
138 *>
139 *> (11) A matrix of the form U' T U, where U is unitary and
140 *> T has "clustered" entries 1, ULP,..., ULP with random
141 *> complex angles on the diagonal and random O(1) entries in
142 *> the upper triangle.
143 *>
144 *> (12) A matrix of the form U' T U, where U is unitary and
145 *> T has complex eigenvalues randomly chosen from
146 *> ULP < |z| < 1 and random O(1) entries in the upper
147 *> triangle.
148 *>
149 *> (13) A matrix of the form X' T X, where X has condition
150 *> SQRT( ULP ) and T has evenly spaced entries 1, ..., ULP
151 *> with random complex angles on the diagonal and random O(1)
152 *> entries in the upper triangle.
153 *>
154 *> (14) A matrix of the form X' T X, where X has condition
155 *> SQRT( ULP ) and T has geometrically spaced entries
156 *> 1, ..., ULP with random complex angles on the diagonal
157 *> and random O(1) entries in the upper triangle.
158 *>
159 *> (15) A matrix of the form X' T X, where X has condition
160 *> SQRT( ULP ) and T has "clustered" entries 1, ULP,..., ULP
161 *> with random complex angles on the diagonal and random O(1)
162 *> entries in the upper triangle.
163 *>
164 *> (16) A matrix of the form X' T X, where X has condition
165 *> SQRT( ULP ) and T has complex eigenvalues randomly chosen
166 *> from ULP < |z| < 1 and random O(1) entries in the upper
167 *> triangle.
168 *>
169 *> (17) Same as (16), but multiplied by a constant
170 *> near the overflow threshold
171 *> (18) Same as (16), but multiplied by a constant
172 *> near the underflow threshold
173 *>
174 *> (19) Nonsymmetric matrix with random entries chosen from |z| < 1
175 *> If N is at least 4, all entries in first two rows and last
176 *> row, and first column and last two columns are zero.
177 *> (20) Same as (19), but multiplied by a constant
178 *> near the overflow threshold
179 *> (21) Same as (19), but multiplied by a constant
180 *> near the underflow threshold
181 *>
182 *> In addition, an input file will be read from logical unit number
183 *> NIUNIT. The file contains matrices along with precomputed
184 *> eigenvalues and reciprocal condition numbers for the eigenvalues
185 *> and right eigenvectors. For these matrices, in addition to tests
186 *> (1) to (9) we will compute the following two tests:
187 *>
188 *> (10) |RCONDV - RCDVIN| / cond(RCONDV)
189 *>
190 *> RCONDV is the reciprocal right eigenvector condition number
191 *> computed by CGEEVX and RCDVIN (the precomputed true value)
192 *> is supplied as input. cond(RCONDV) is the condition number of
193 *> RCONDV, and takes errors in computing RCONDV into account, so
194 *> that the resulting quantity should be O(ULP). cond(RCONDV) is
195 *> essentially given by norm(A)/RCONDE.
196 *>
197 *> (11) |RCONDE - RCDEIN| / cond(RCONDE)
198 *>
199 *> RCONDE is the reciprocal eigenvalue condition number
200 *> computed by CGEEVX and RCDEIN (the precomputed true value)
201 *> is supplied as input. cond(RCONDE) is the condition number
202 *> of RCONDE, and takes errors in computing RCONDE into account,
203 *> so that the resulting quantity should be O(ULP). cond(RCONDE)
204 *> is essentially given by norm(A)/RCONDV.
205 *> \endverbatim
206 *
207 * Arguments:
208 * ==========
209 *
210 *> \param[in] NSIZES
211 *> \verbatim
212 *> NSIZES is INTEGER
213 *> The number of sizes of matrices to use. NSIZES must be at
214 *> least zero. If it is zero, no randomly generated matrices
215 *> are tested, but any test matrices read from NIUNIT will be
216 *> tested.
217 *> \endverbatim
218 *>
219 *> \param[in] NN
220 *> \verbatim
221 *> NN is INTEGER array, dimension (NSIZES)
222 *> An array containing the sizes to be used for the matrices.
223 *> Zero values will be skipped. The values must be at least
224 *> zero.
225 *> \endverbatim
226 *>
227 *> \param[in] NTYPES
228 *> \verbatim
229 *> NTYPES is INTEGER
230 *> The number of elements in DOTYPE. NTYPES must be at least
231 *> zero. If it is zero, no randomly generated test matrices
232 *> are tested, but and test matrices read from NIUNIT will be
233 *> tested. If it is MAXTYP+1 and NSIZES is 1, then an
234 *> additional type, MAXTYP+1 is defined, which is to use
235 *> whatever matrix is in A. This is only useful if
236 *> DOTYPE(1:MAXTYP) is .FALSE. and DOTYPE(MAXTYP+1) is .TRUE. .
237 *> \endverbatim
238 *>
239 *> \param[in] DOTYPE
240 *> \verbatim
241 *> DOTYPE is LOGICAL array, dimension (NTYPES)
242 *> If DOTYPE(j) is .TRUE., then for each size in NN a
243 *> matrix of that size and of type j will be generated.
244 *> If NTYPES is smaller than the maximum number of types
245 *> defined (PARAMETER MAXTYP), then types NTYPES+1 through
246 *> MAXTYP will not be generated. If NTYPES is larger
247 *> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
248 *> will be ignored.
249 *> \endverbatim
250 *>
251 *> \param[in,out] ISEED
252 *> \verbatim
253 *> ISEED is INTEGER array, dimension (4)
254 *> On entry ISEED specifies the seed of the random number
255 *> generator. The array elements should be between 0 and 4095;
256 *> if not they will be reduced mod 4096. Also, ISEED(4) must
257 *> be odd. The random number generator uses a linear
258 *> congruential sequence limited to small integers, and so
259 *> should produce machine independent random numbers. The
260 *> values of ISEED are changed on exit, and can be used in the
261 *> next call to CDRVVX to continue the same random number
262 *> sequence.
263 *> \endverbatim
264 *>
265 *> \param[in] THRESH
266 *> \verbatim
267 *> THRESH is REAL
268 *> A test will count as "failed" if the "error", computed as
269 *> described above, exceeds THRESH. Note that the error
270 *> is scaled to be O(1), so THRESH should be a reasonably
271 *> small multiple of 1, e.g., 10 or 100. In particular,
272 *> it should not depend on the precision (single vs. double)
273 *> or the size of the matrix. It must be at least zero.
274 *> \endverbatim
275 *>
276 *> \param[in] NIUNIT
277 *> \verbatim
278 *> NIUNIT is INTEGER
279 *> The FORTRAN unit number for reading in the data file of
280 *> problems to solve.
281 *> \endverbatim
282 *>
283 *> \param[in] NOUNIT
284 *> \verbatim
285 *> NOUNIT is INTEGER
286 *> The FORTRAN unit number for printing out error messages
287 *> (e.g., if a routine returns INFO not equal to 0.)
288 *> \endverbatim
289 *>
290 *> \param[out] A
291 *> \verbatim
292 *> A is COMPLEX array, dimension (LDA, max(NN,12))
293 *> Used to hold the matrix whose eigenvalues are to be
294 *> computed. On exit, A contains the last matrix actually used.
295 *> \endverbatim
296 *>
297 *> \param[in] LDA
298 *> \verbatim
299 *> LDA is INTEGER
300 *> The leading dimension of A, and H. LDA must be at
301 *> least 1 and at least max( NN, 12 ). (12 is the
302 *> dimension of the largest matrix on the precomputed
303 *> input file.)
304 *> \endverbatim
305 *>
306 *> \param[out] H
307 *> \verbatim
308 *> H is COMPLEX array, dimension (LDA, max(NN,12))
309 *> Another copy of the test matrix A, modified by CGEEVX.
310 *> \endverbatim
311 *>
312 *> \param[out] W
313 *> \verbatim
314 *> W is COMPLEX array, dimension (max(NN,12))
315 *> Contains the eigenvalues of A.
316 *> \endverbatim
317 *>
318 *> \param[out] W1
319 *> \verbatim
320 *> W1 is COMPLEX array, dimension (max(NN,12))
321 *> Like W, this array contains the eigenvalues of A,
322 *> but those computed when CGEEVX only computes a partial
323 *> eigendecomposition, i.e. not the eigenvalues and left
324 *> and right eigenvectors.
325 *> \endverbatim
326 *>
327 *> \param[out] VL
328 *> \verbatim
329 *> VL is COMPLEX array, dimension (LDVL, max(NN,12))
330 *> VL holds the computed left eigenvectors.
331 *> \endverbatim
332 *>
333 *> \param[in] LDVL
334 *> \verbatim
335 *> LDVL is INTEGER
336 *> Leading dimension of VL. Must be at least max(1,max(NN,12)).
337 *> \endverbatim
338 *>
339 *> \param[out] VR
340 *> \verbatim
341 *> VR is COMPLEX array, dimension (LDVR, max(NN,12))
342 *> VR holds the computed right eigenvectors.
343 *> \endverbatim
344 *>
345 *> \param[in] LDVR
346 *> \verbatim
347 *> LDVR is INTEGER
348 *> Leading dimension of VR. Must be at least max(1,max(NN,12)).
349 *> \endverbatim
350 *>
351 *> \param[out] LRE
352 *> \verbatim
353 *> LRE is COMPLEX array, dimension (LDLRE, max(NN,12))
354 *> LRE holds the computed right or left eigenvectors.
355 *> \endverbatim
356 *>
357 *> \param[in] LDLRE
358 *> \verbatim
359 *> LDLRE is INTEGER
360 *> Leading dimension of LRE. Must be at least max(1,max(NN,12))
361 *> \endverbatim
362 *>
363 *> \param[out] RCONDV
364 *> \verbatim
365 *> RCONDV is REAL array, dimension (N)
366 *> RCONDV holds the computed reciprocal condition numbers
367 *> for eigenvectors.
368 *> \endverbatim
369 *>
370 *> \param[out] RCNDV1
371 *> \verbatim
372 *> RCNDV1 is REAL array, dimension (N)
373 *> RCNDV1 holds more computed reciprocal condition numbers
374 *> for eigenvectors.
375 *> \endverbatim
376 *>
377 *> \param[in] RCDVIN
378 *> \verbatim
379 *> RCDVIN is REAL array, dimension (N)
380 *> When COMP = .TRUE. RCDVIN holds the precomputed reciprocal
381 *> condition numbers for eigenvectors to be compared with
382 *> RCONDV.
383 *> \endverbatim
384 *>
385 *> \param[out] RCONDE
386 *> \verbatim
387 *> RCONDE is REAL array, dimension (N)
388 *> RCONDE holds the computed reciprocal condition numbers
389 *> for eigenvalues.
390 *> \endverbatim
391 *>
392 *> \param[out] RCNDE1
393 *> \verbatim
394 *> RCNDE1 is REAL array, dimension (N)
395 *> RCNDE1 holds more computed reciprocal condition numbers
396 *> for eigenvalues.
397 *> \endverbatim
398 *>
399 *> \param[in] RCDEIN
400 *> \verbatim
401 *> RCDEIN is REAL array, dimension (N)
402 *> When COMP = .TRUE. RCDEIN holds the precomputed reciprocal
403 *> condition numbers for eigenvalues to be compared with
404 *> RCONDE.
405 *> \endverbatim
406 *>
407 *> \param[out] SCALE
408 *> \verbatim
409 *> SCALE is REAL array, dimension (N)
410 *> Holds information describing balancing of matrix.
411 *> \endverbatim
412 *>
413 *> \param[out] SCALE1
414 *> \verbatim
415 *> SCALE1 is REAL array, dimension (N)
416 *> Holds information describing balancing of matrix.
417 *> \endverbatim
418 *>
419 *> \param[out] RESULT
420 *> \verbatim
421 *> RESULT is REAL array, dimension (11)
422 *> The values computed by the seven tests described above.
423 *> The values are currently limited to 1/ulp, to avoid
424 *> overflow.
425 *> \endverbatim
426 *>
427 *> \param[out] WORK
428 *> \verbatim
429 *> WORK is COMPLEX array, dimension (NWORK)
430 *> \endverbatim
431 *>
432 *> \param[in] NWORK
433 *> \verbatim
434 *> NWORK is INTEGER
435 *> The number of entries in WORK. This must be at least
436 *> max(6*12+2*12**2,6*NN(j)+2*NN(j)**2) =
437 *> max( 360 ,6*NN(j)+2*NN(j)**2) for all j.
438 *> \endverbatim
439 *>
440 *> \param[out] RWORK
441 *> \verbatim
442 *> RWORK is REAL array, dimension (2*max(NN,12))
443 *> \endverbatim
444 *>
445 *> \param[out] INFO
446 *> \verbatim
447 *> INFO is INTEGER
448 *> If 0, then successful exit.
449 *> If <0, then input paramter -INFO is incorrect.
450 *> If >0, CLATMR, CLATMS, CLATME or CGET23 returned an error
451 *> code, and INFO is its absolute value.
452 *>
453 *>-----------------------------------------------------------------------
454 *>
455 *> Some Local Variables and Parameters:
456 *> ---- ----- --------- --- ----------
457 *>
458 *> ZERO, ONE Real 0 and 1.
459 *> MAXTYP The number of types defined.
460 *> NMAX Largest value in NN or 12.
461 *> NERRS The number of tests which have exceeded THRESH
462 *> COND, CONDS,
463 *> IMODE Values to be passed to the matrix generators.
464 *> ANORM Norm of A; passed to matrix generators.
465 *>
466 *> OVFL, UNFL Overflow and underflow thresholds.
467 *> ULP, ULPINV Finest relative precision and its inverse.
468 *> RTULP, RTULPI Square roots of the previous 4 values.
469 *>
470 *> The following four arrays decode JTYPE:
471 *> KTYPE(j) The general type (1-10) for type "j".
472 *> KMODE(j) The MODE value to be passed to the matrix
473 *> generator for type "j".
474 *> KMAGN(j) The order of magnitude ( O(1),
475 *> O(overflow^(1/2) ), O(underflow^(1/2) )
476 *> KCONDS(j) Selectw whether CONDS is to be 1 or
477 *> 1/sqrt(ulp). (0 means irrelevant.)
478 *> \endverbatim
479 *
480 * Authors:
481 * ========
482 *
483 *> \author Univ. of Tennessee
484 *> \author Univ. of California Berkeley
485 *> \author Univ. of Colorado Denver
486 *> \author NAG Ltd.
487 *
488 *> \date November 2011
489 *
490 *> \ingroup complex_eig
491 *
492 * =====================================================================
493  SUBROUTINE cdrvvx( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
494  $ niunit, nounit, a, lda, h, w, w1, vl, ldvl, vr,
495  $ ldvr, lre, ldlre, rcondv, rcndv1, rcdvin,
496  $ rconde, rcnde1, rcdein, scale, scale1, result,
497  $ work, nwork, rwork, info )
498 *
499 * -- LAPACK test routine (version 3.4.0) --
500 * -- LAPACK is a software package provided by Univ. of Tennessee, --
501 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
502 * November 2011
503 *
504 * .. Scalar Arguments ..
505  INTEGER info, lda, ldlre, ldvl, ldvr, niunit, nounit,
506  $ nsizes, ntypes, nwork
507  REAL thresh
508 * ..
509 * .. Array Arguments ..
510  LOGICAL dotype( * )
511  INTEGER iseed( 4 ), nn( * )
512  REAL rcdein( * ), rcdvin( * ), rcnde1( * ),
513  $ rcndv1( * ), rconde( * ), rcondv( * ),
514  $ result( 11 ), rwork( * ), scale( * ),
515  $ scale1( * )
516  COMPLEX a( lda, * ), h( lda, * ), lre( ldlre, * ),
517  $ vl( ldvl, * ), vr( ldvr, * ), w( * ), w1( * ),
518  $ work( * )
519 * ..
520 *
521 * =====================================================================
522 *
523 * .. Parameters ..
524  COMPLEX czero
525  parameter( czero = ( 0.0e+0, 0.0e+0 ) )
526  COMPLEX cone
527  parameter( cone = ( 1.0e+0, 0.0e+0 ) )
528  REAL zero, one
529  parameter( zero = 0.0e+0, one = 1.0e+0 )
530  INTEGER maxtyp
531  parameter( maxtyp = 21 )
532 * ..
533 * .. Local Scalars ..
534  LOGICAL badnn
535  CHARACTER balanc
536  CHARACTER*3 path
537  INTEGER i, ibal, iinfo, imode, isrt, itype, iwk, j,
538  $ jcol, jsize, jtype, mtypes, n, nerrs,
539  $ nfail, nmax, nnwork, ntest, ntestf, ntestt
540  REAL anorm, cond, conds, ovfl, rtulp, rtulpi, ulp,
541  $ ulpinv, unfl, wi, wr
542 * ..
543 * .. Local Arrays ..
544  CHARACTER bal( 4 )
545  INTEGER idumma( 1 ), ioldsd( 4 ), kconds( maxtyp ),
546  $ kmagn( maxtyp ), kmode( maxtyp ),
547  $ ktype( maxtyp )
548 * ..
549 * .. External Functions ..
550  REAL slamch
551  EXTERNAL slamch
552 * ..
553 * .. External Subroutines ..
554  EXTERNAL cget23, clatme, clatmr, clatms, claset, slabad,
555  $ slasum, xerbla
556 * ..
557 * .. Intrinsic Functions ..
558  INTRINSIC abs, cmplx, max, min, sqrt
559 * ..
560 * .. Data statements ..
561  DATA ktype / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
562  DATA kmagn / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
563  $ 3, 1, 2, 3 /
564  DATA kmode / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
565  $ 1, 5, 5, 5, 4, 3, 1 /
566  DATA kconds / 3*0, 5*0, 4*1, 6*2, 3*0 /
567  DATA bal / 'N', 'P', 'S', 'B' /
568 * ..
569 * .. Executable Statements ..
570 *
571  path( 1: 1 ) = 'Complex precision'
572  path( 2: 3 ) = 'VX'
573 *
574 * Check for errors
575 *
576  ntestt = 0
577  ntestf = 0
578  info = 0
579 *
580 * Important constants
581 *
582  badnn = .false.
583 *
584 * 7 is the largest dimension in the input file of precomputed
585 * problems
586 *
587  nmax = 7
588  DO 10 j = 1, nsizes
589  nmax = max( nmax, nn( j ) )
590  IF( nn( j ).LT.0 )
591  $ badnn = .true.
592  10 CONTINUE
593 *
594 * Check for errors
595 *
596  IF( nsizes.LT.0 ) THEN
597  info = -1
598  ELSE IF( badnn ) THEN
599  info = -2
600  ELSE IF( ntypes.LT.0 ) THEN
601  info = -3
602  ELSE IF( thresh.LT.zero ) THEN
603  info = -6
604  ELSE IF( lda.LT.1 .OR. lda.LT.nmax ) THEN
605  info = -10
606  ELSE IF( ldvl.LT.1 .OR. ldvl.LT.nmax ) THEN
607  info = -15
608  ELSE IF( ldvr.LT.1 .OR. ldvr.LT.nmax ) THEN
609  info = -17
610  ELSE IF( ldlre.LT.1 .OR. ldlre.LT.nmax ) THEN
611  info = -19
612  ELSE IF( 6*nmax+2*nmax**2.GT.nwork ) THEN
613  info = -30
614  END IF
615 *
616  IF( info.NE.0 ) THEN
617  CALL xerbla( 'CDRVVX', -info )
618  RETURN
619  END IF
620 *
621 * If nothing to do check on NIUNIT
622 *
623  IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
624  $ go to 160
625 *
626 * More Important constants
627 *
628  unfl = slamch( 'Safe minimum' )
629  ovfl = one / unfl
630  CALL slabad( unfl, ovfl )
631  ulp = slamch( 'Precision' )
632  ulpinv = one / ulp
633  rtulp = sqrt( ulp )
634  rtulpi = one / rtulp
635 *
636 * Loop over sizes, types
637 *
638  nerrs = 0
639 *
640  DO 150 jsize = 1, nsizes
641  n = nn( jsize )
642  IF( nsizes.NE.1 ) THEN
643  mtypes = min( maxtyp, ntypes )
644  ELSE
645  mtypes = min( maxtyp+1, ntypes )
646  END IF
647 *
648  DO 140 jtype = 1, mtypes
649  IF( .NOT.dotype( jtype ) )
650  $ go to 140
651 *
652 * Save ISEED in case of an error.
653 *
654  DO 20 j = 1, 4
655  ioldsd( j ) = iseed( j )
656  20 CONTINUE
657 *
658 * Compute "A"
659 *
660 * Control parameters:
661 *
662 * KMAGN KCONDS KMODE KTYPE
663 * =1 O(1) 1 clustered 1 zero
664 * =2 large large clustered 2 identity
665 * =3 small exponential Jordan
666 * =4 arithmetic diagonal, (w/ eigenvalues)
667 * =5 random log symmetric, w/ eigenvalues
668 * =6 random general, w/ eigenvalues
669 * =7 random diagonal
670 * =8 random symmetric
671 * =9 random general
672 * =10 random triangular
673 *
674  IF( mtypes.GT.maxtyp )
675  $ go to 90
676 *
677  itype = ktype( jtype )
678  imode = kmode( jtype )
679 *
680 * Compute norm
681 *
682  go to( 30, 40, 50 )kmagn( jtype )
683 *
684  30 CONTINUE
685  anorm = one
686  go to 60
687 *
688  40 CONTINUE
689  anorm = ovfl*ulp
690  go to 60
691 *
692  50 CONTINUE
693  anorm = unfl*ulpinv
694  go to 60
695 *
696  60 CONTINUE
697 *
698  CALL claset( 'Full', lda, n, czero, czero, a, lda )
699  iinfo = 0
700  cond = ulpinv
701 *
702 * Special Matrices -- Identity & Jordan block
703 *
704 * Zero
705 *
706  IF( itype.EQ.1 ) THEN
707  iinfo = 0
708 *
709  ELSE IF( itype.EQ.2 ) THEN
710 *
711 * Identity
712 *
713  DO 70 jcol = 1, n
714  a( jcol, jcol ) = anorm
715  70 CONTINUE
716 *
717  ELSE IF( itype.EQ.3 ) THEN
718 *
719 * Jordan Block
720 *
721  DO 80 jcol = 1, n
722  a( jcol, jcol ) = anorm
723  IF( jcol.GT.1 )
724  $ a( jcol, jcol-1 ) = one
725  80 CONTINUE
726 *
727  ELSE IF( itype.EQ.4 ) THEN
728 *
729 * Diagonal Matrix, [Eigen]values Specified
730 *
731  CALL clatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
732  $ anorm, 0, 0, 'N', a, lda, work( n+1 ),
733  $ iinfo )
734 *
735  ELSE IF( itype.EQ.5 ) THEN
736 *
737 * Symmetric, eigenvalues specified
738 *
739  CALL clatms( n, n, 'S', iseed, 'H', rwork, imode, cond,
740  $ anorm, n, n, 'N', a, lda, work( n+1 ),
741  $ iinfo )
742 *
743  ELSE IF( itype.EQ.6 ) THEN
744 *
745 * General, eigenvalues specified
746 *
747  IF( kconds( jtype ).EQ.1 ) THEN
748  conds = one
749  ELSE IF( kconds( jtype ).EQ.2 ) THEN
750  conds = rtulpi
751  ELSE
752  conds = zero
753  END IF
754 *
755  CALL clatme( n, 'D', iseed, work, imode, cond, cone,
756  $ 'T', 'T', 'T', rwork, 4, conds, n, n, anorm,
757  $ a, lda, work( 2*n+1 ), iinfo )
758 *
759  ELSE IF( itype.EQ.7 ) THEN
760 *
761 * Diagonal, random eigenvalues
762 *
763  CALL clatmr( n, n, 'D', iseed, 'S', work, 6, one, cone,
764  $ 'T', 'N', work( n+1 ), 1, one,
765  $ work( 2*n+1 ), 1, one, 'N', idumma, 0, 0,
766  $ zero, anorm, 'NO', a, lda, idumma, iinfo )
767 *
768  ELSE IF( itype.EQ.8 ) THEN
769 *
770 * Symmetric, random eigenvalues
771 *
772  CALL clatmr( n, n, 'D', iseed, 'H', work, 6, one, cone,
773  $ 'T', 'N', work( n+1 ), 1, one,
774  $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
775  $ zero, anorm, 'NO', a, lda, idumma, iinfo )
776 *
777  ELSE IF( itype.EQ.9 ) THEN
778 *
779 * General, random eigenvalues
780 *
781  CALL clatmr( n, n, 'D', iseed, 'N', work, 6, one, cone,
782  $ 'T', 'N', work( n+1 ), 1, one,
783  $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
784  $ zero, anorm, 'NO', a, lda, idumma, iinfo )
785  IF( n.GE.4 ) THEN
786  CALL claset( 'Full', 2, n, czero, czero, a, lda )
787  CALL claset( 'Full', n-3, 1, czero, czero, a( 3, 1 ),
788  $ lda )
789  CALL claset( 'Full', n-3, 2, czero, czero,
790  $ a( 3, n-1 ), lda )
791  CALL claset( 'Full', 1, n, czero, czero, a( n, 1 ),
792  $ lda )
793  END IF
794 *
795  ELSE IF( itype.EQ.10 ) THEN
796 *
797 * Triangular, random eigenvalues
798 *
799  CALL clatmr( n, n, 'D', iseed, 'N', work, 6, one, cone,
800  $ 'T', 'N', work( n+1 ), 1, one,
801  $ work( 2*n+1 ), 1, one, 'N', idumma, n, 0,
802  $ zero, anorm, 'NO', a, lda, idumma, iinfo )
803 *
804  ELSE
805 *
806  iinfo = 1
807  END IF
808 *
809  IF( iinfo.NE.0 ) THEN
810  WRITE( nounit, fmt = 9992 )'Generator', iinfo, n, jtype,
811  $ ioldsd
812  info = abs( iinfo )
813  RETURN
814  END IF
815 *
816  90 CONTINUE
817 *
818 * Test for minimal and generous workspace
819 *
820  DO 130 iwk = 1, 3
821  IF( iwk.EQ.1 ) THEN
822  nnwork = 2*n
823  ELSE IF( iwk.EQ.2 ) THEN
824  nnwork = 2*n + n**2
825  ELSE
826  nnwork = 6*n + 2*n**2
827  END IF
828  nnwork = max( nnwork, 1 )
829 *
830 * Test for all balancing options
831 *
832  DO 120 ibal = 1, 4
833  balanc = bal( ibal )
834 *
835 * Perform tests
836 *
837  CALL cget23( .false., 0, balanc, jtype, thresh,
838  $ ioldsd, nounit, n, a, lda, h, w, w1, vl,
839  $ ldvl, vr, ldvr, lre, ldlre, rcondv,
840  $ rcndv1, rcdvin, rconde, rcnde1, rcdein,
841  $ scale, scale1, result, work, nnwork,
842  $ rwork, info )
843 *
844 * Check for RESULT(j) > THRESH
845 *
846  ntest = 0
847  nfail = 0
848  DO 100 j = 1, 9
849  IF( result( j ).GE.zero )
850  $ ntest = ntest + 1
851  IF( result( j ).GE.thresh )
852  $ nfail = nfail + 1
853  100 CONTINUE
854 *
855  IF( nfail.GT.0 )
856  $ ntestf = ntestf + 1
857  IF( ntestf.EQ.1 ) THEN
858  WRITE( nounit, fmt = 9999 )path
859  WRITE( nounit, fmt = 9998 )
860  WRITE( nounit, fmt = 9997 )
861  WRITE( nounit, fmt = 9996 )
862  WRITE( nounit, fmt = 9995 )thresh
863  ntestf = 2
864  END IF
865 *
866  DO 110 j = 1, 9
867  IF( result( j ).GE.thresh ) THEN
868  WRITE( nounit, fmt = 9994 )balanc, n, iwk,
869  $ ioldsd, jtype, j, result( j )
870  END IF
871  110 CONTINUE
872 *
873  nerrs = nerrs + nfail
874  ntestt = ntestt + ntest
875 *
876  120 CONTINUE
877  130 CONTINUE
878  140 CONTINUE
879  150 CONTINUE
880 *
881  160 CONTINUE
882 *
883 * Read in data from file to check accuracy of condition estimation.
884 * Assume input eigenvalues are sorted lexicographically (increasing
885 * by real part, then decreasing by imaginary part)
886 *
887  jtype = 0
888  170 CONTINUE
889  READ( niunit, fmt = *, END = 220 )n, isrt
890 *
891 * Read input data until N=0
892 *
893  IF( n.EQ.0 )
894  $ go to 220
895  jtype = jtype + 1
896  iseed( 1 ) = jtype
897  DO 180 i = 1, n
898  READ( niunit, fmt = * )( a( i, j ), j = 1, n )
899  180 CONTINUE
900  DO 190 i = 1, n
901  READ( niunit, fmt = * )wr, wi, rcdein( i ), rcdvin( i )
902  w1( i ) = cmplx( wr, wi )
903  190 CONTINUE
904  CALL cget23( .true., isrt, 'N', 22, thresh, iseed, nounit, n, a,
905  $ lda, h, w, w1, vl, ldvl, vr, ldvr, lre, ldlre,
906  $ rcondv, rcndv1, rcdvin, rconde, rcnde1, rcdein,
907  $ scale, scale1, result, work, 6*n+2*n**2, rwork,
908  $ info )
909 *
910 * Check for RESULT(j) > THRESH
911 *
912  ntest = 0
913  nfail = 0
914  DO 200 j = 1, 11
915  IF( result( j ).GE.zero )
916  $ ntest = ntest + 1
917  IF( result( j ).GE.thresh )
918  $ nfail = nfail + 1
919  200 CONTINUE
920 *
921  IF( nfail.GT.0 )
922  $ ntestf = ntestf + 1
923  IF( ntestf.EQ.1 ) THEN
924  WRITE( nounit, fmt = 9999 )path
925  WRITE( nounit, fmt = 9998 )
926  WRITE( nounit, fmt = 9997 )
927  WRITE( nounit, fmt = 9996 )
928  WRITE( nounit, fmt = 9995 )thresh
929  ntestf = 2
930  END IF
931 *
932  DO 210 j = 1, 11
933  IF( result( j ).GE.thresh ) THEN
934  WRITE( nounit, fmt = 9993 )n, jtype, j, result( j )
935  END IF
936  210 CONTINUE
937 *
938  nerrs = nerrs + nfail
939  ntestt = ntestt + ntest
940  go to 170
941  220 CONTINUE
942 *
943 * Summary
944 *
945  CALL slasum( path, nounit, nerrs, ntestt )
946 *
947  9999 FORMAT( / 1x, a3, ' -- Complex Eigenvalue-Eigenvector ',
948  $ 'Decomposition Expert Driver',
949  $ / ' Matrix types (see CDRVVX for details): ' )
950 *
951  9998 FORMAT( / ' Special Matrices:', / ' 1=Zero matrix. ',
952  $ ' ', ' 5=Diagonal: geometr. spaced entries.',
953  $ / ' 2=Identity matrix. ', ' 6=Diagona',
954  $ 'l: clustered entries.', / ' 3=Transposed Jordan block. ',
955  $ ' ', ' 7=Diagonal: large, evenly spaced.', / ' ',
956  $ '4=Diagonal: evenly spaced entries. ', ' 8=Diagonal: s',
957  $ 'mall, evenly spaced.' )
958  9997 FORMAT( ' Dense, Non-Symmetric Matrices:', / ' 9=Well-cond., ev',
959  $ 'enly spaced eigenvals.', ' 14=Ill-cond., geomet. spaced e',
960  $ 'igenals.', / ' 10=Well-cond., geom. spaced eigenvals. ',
961  $ ' 15=Ill-conditioned, clustered e.vals.', / ' 11=Well-cond',
962  $ 'itioned, clustered e.vals. ', ' 16=Ill-cond., random comp',
963  $ 'lex ', / ' 12=Well-cond., random complex ', ' ',
964  $ ' 17=Ill-cond., large rand. complx ', / ' 13=Ill-condi',
965  $ 'tioned, evenly spaced. ', ' 18=Ill-cond., small rand.',
966  $ ' complx ' )
967  9996 FORMAT( ' 19=Matrix with random O(1) entries. ', ' 21=Matrix ',
968  $ 'with small random entries.', / ' 20=Matrix with large ran',
969  $ 'dom entries. ', ' 22=Matrix read from input file', / )
970  9995 FORMAT( ' Tests performed with test threshold =', f8.2,
971  $ / / ' 1 = | A VR - VR W | / ( n |A| ulp ) ',
972  $ / ' 2 = | transpose(A) VL - VL W | / ( n |A| ulp ) ',
973  $ / ' 3 = | |VR(i)| - 1 | / ulp ',
974  $ / ' 4 = | |VL(i)| - 1 | / ulp ',
975  $ / ' 5 = 0 if W same no matter if VR or VL computed,',
976  $ ' 1/ulp otherwise', /
977  $ ' 6 = 0 if VR same no matter what else computed,',
978  $ ' 1/ulp otherwise', /
979  $ ' 7 = 0 if VL same no matter what else computed,',
980  $ ' 1/ulp otherwise', /
981  $ ' 8 = 0 if RCONDV same no matter what else computed,',
982  $ ' 1/ulp otherwise', /
983  $ ' 9 = 0 if SCALE, ILO, IHI, ABNRM same no matter what else',
984  $ ' computed, 1/ulp otherwise',
985  $ / ' 10 = | RCONDV - RCONDV(precomputed) | / cond(RCONDV),',
986  $ / ' 11 = | RCONDE - RCONDE(precomputed) | / cond(RCONDE),' )
987  9994 FORMAT( ' BALANC=''', a1, ''',N=', i4, ',IWK=', i1, ', seed=',
988  $ 4( i4, ',' ), ' type ', i2, ', test(', i2, ')=', g10.3 )
989  9993 FORMAT( ' N=', i5, ', input example =', i3, ', test(', i2, ')=',
990  $ g10.3 )
991  9992 FORMAT( ' CDRVVX: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
992  $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
993 *
994  RETURN
995 *
996 * End of CDRVVX
997 *
998  END