1 SUBROUTINE DDRVVX( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
2 $ NIUNIT, NOUNIT, A, LDA, H, WR, WI, WR1, WI1,
3 $ VL, LDVL, VR, LDVR, LRE, LDLRE, RCONDV, RCNDV1,
4 $ RCDVIN, RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1,
5 $ RESULT, WORK, NWORK, IWORK, INFO )
6 *
7 * -- LAPACK test routine (version 3.1) --
8 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
9 * November 2006
10 *
11 * .. Scalar Arguments ..
12 INTEGER INFO, LDA, LDLRE, LDVL, LDVR, NIUNIT, NOUNIT,
13 $ NSIZES, NTYPES, NWORK
14 DOUBLE PRECISION THRESH
15 * ..
16 * .. Array Arguments ..
17 LOGICAL DOTYPE( * )
18 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
19 DOUBLE PRECISION A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
20 $ RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
21 $ RCNDV1( * ), RCONDE( * ), RCONDV( * ),
22 $ RESULT( 11 ), SCALE( * ), SCALE1( * ),
23 $ VL( LDVL, * ), VR( LDVR, * ), WI( * ),
24 $ WI1( * ), WORK( * ), WR( * ), WR1( * )
25 * ..
26 *
27 * Purpose
28 * =======
29 *
30 * DDRVVX checks the nonsymmetric eigenvalue problem expert driver
31 * DGEEVX.
32 *
33 * DDRVVX uses both test matrices generated randomly depending on
34 * data supplied in the calling sequence, as well as on data
35 * read from an input file and including precomputed condition
36 * numbers to which it compares the ones it computes.
37 *
38 * When DDRVVX is called, a number of matrix "sizes" ("n's") and a
39 * number of matrix "types" are specified in the calling sequence.
40 * For each size ("n") and each type of matrix, one matrix will be
41 * generated and used to test the nonsymmetric eigenroutines. For
42 * each matrix, 9 tests will be performed:
43 *
44 * (1) | A * VR - VR * W | / ( n |A| ulp )
45 *
46 * Here VR is the matrix of unit right eigenvectors.
47 * W is a block diagonal matrix, with a 1x1 block for each
48 * real eigenvalue and a 2x2 block for each complex conjugate
49 * pair. If eigenvalues j and j+1 are a complex conjugate pair,
50 * so WR(j) = WR(j+1) = wr and WI(j) = - WI(j+1) = wi, then the
51 * 2 x 2 block corresponding to the pair will be:
52 *
53 * ( wr wi )
54 * ( -wi wr )
55 *
56 * Such a block multiplying an n x 2 matrix ( ur ui ) on the
57 * right will be the same as multiplying ur + i*ui by wr + i*wi.
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 signs.
118 * (ULP = (first number larger than 1) - 1 )
119 * (5) A diagonal matrix with geometrically spaced entries
120 * 1, ..., ULP and random signs.
121 * (6) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
122 * and random signs.
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 orthogonal and
130 * T has evenly spaced entries 1, ..., ULP with random signs
131 * 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 orthogonal and
135 * T has geometrically spaced entries 1, ..., ULP with random
136 * signs on the diagonal and random O(1) entries in the upper
137 * triangle.
138 *
139 * (11) A matrix of the form U' T U, where U is orthogonal and
140 * T has "clustered" entries 1, ULP,..., ULP with random
141 * signs on the diagonal and random O(1) entries in the upper
142 * triangle.
143 *
144 * (12) A matrix of the form U' T U, where U is orthogonal and
145 * T has real or complex conjugate paired eigenvalues randomly
146 * chosen from ( ULP, 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 signs on the diagonal and random O(1) entries
152 * 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 signs on the diagonal and random
157 * 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 signs on the diagonal and random O(1) entries
162 * 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 real or complex conjugate paired
166 * eigenvalues randomly chosen from ( ULP, 1 ) and random
167 * O(1) entries in the upper 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 (-1,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 DGEEVX 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 DGEEVX 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 *
206 * Arguments
207 * ==========
208 *
209 * NSIZES (input) INTEGER
210 * The number of sizes of matrices to use. NSIZES must be at
211 * least zero. If it is zero, no randomly generated matrices
212 * are tested, but any test matrices read from NIUNIT will be
213 * tested.
214 *
215 * NN (input) INTEGER array, dimension (NSIZES)
216 * An array containing the sizes to be used for the matrices.
217 * Zero values will be skipped. The values must be at least
218 * zero.
219 *
220 * NTYPES (input) INTEGER
221 * The number of elements in DOTYPE. NTYPES must be at least
222 * zero. If it is zero, no randomly generated test matrices
223 * are tested, but and test matrices read from NIUNIT will be
224 * tested. If it is MAXTYP+1 and NSIZES is 1, then an
225 * additional type, MAXTYP+1 is defined, which is to use
226 * whatever matrix is in A. This is only useful if
227 * DOTYPE(1:MAXTYP) is .FALSE. and DOTYPE(MAXTYP+1) is .TRUE. .
228 *
229 * DOTYPE (input) LOGICAL array, dimension (NTYPES)
230 * If DOTYPE(j) is .TRUE., then for each size in NN a
231 * matrix of that size and of type j will be generated.
232 * If NTYPES is smaller than the maximum number of types
233 * defined (PARAMETER MAXTYP), then types NTYPES+1 through
234 * MAXTYP will not be generated. If NTYPES is larger
235 * than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
236 * will be ignored.
237 *
238 * ISEED (input/output) INTEGER array, dimension (4)
239 * On entry ISEED specifies the seed of the random number
240 * generator. The array elements should be between 0 and 4095;
241 * if not they will be reduced mod 4096. Also, ISEED(4) must
242 * be odd. The random number generator uses a linear
243 * congruential sequence limited to small integers, and so
244 * should produce machine independent random numbers. The
245 * values of ISEED are changed on exit, and can be used in the
246 * next call to DDRVVX to continue the same random number
247 * sequence.
248 *
249 * THRESH (input) DOUBLE PRECISION
250 * A test will count as "failed" if the "error", computed as
251 * described above, exceeds THRESH. Note that the error
252 * is scaled to be O(1), so THRESH should be a reasonably
253 * small multiple of 1, e.g., 10 or 100. In particular,
254 * it should not depend on the precision (single vs. double)
255 * or the size of the matrix. It must be at least zero.
256 *
257 * NIUNIT (input) INTEGER
258 * The FORTRAN unit number for reading in the data file of
259 * problems to solve.
260 *
261 * NOUNIT (input) INTEGER
262 * The FORTRAN unit number for printing out error messages
263 * (e.g., if a routine returns INFO not equal to 0.)
264 *
265 * A (workspace) DOUBLE PRECISION array, dimension
266 * (LDA, max(NN,12))
267 * Used to hold the matrix whose eigenvalues are to be
268 * computed. On exit, A contains the last matrix actually used.
269 *
270 * LDA (input) INTEGER
271 * The leading dimension of the arrays A and H.
272 * LDA >= max(NN,12), since 12 is the dimension of the largest
273 * matrix in the precomputed input file.
274 *
275 * H (workspace) DOUBLE PRECISION array, dimension
276 * (LDA, max(NN,12))
277 * Another copy of the test matrix A, modified by DGEEVX.
278 *
279 * WR (workspace) DOUBLE PRECISION array, dimension (max(NN))
280 * WI (workspace) DOUBLE PRECISION array, dimension (max(NN))
281 * The real and imaginary parts of the eigenvalues of A.
282 * On exit, WR + WI*i are the eigenvalues of the matrix in A.
283 *
284 * WR1 (workspace) DOUBLE PRECISION array, dimension (max(NN,12))
285 * WI1 (workspace) DOUBLE PRECISION array, dimension (max(NN,12))
286 * Like WR, WI, these arrays contain the eigenvalues of A,
287 * but those computed when DGEEVX only computes a partial
288 * eigendecomposition, i.e. not the eigenvalues and left
289 * and right eigenvectors.
290 *
291 * VL (workspace) DOUBLE PRECISION array, dimension
292 * (LDVL, max(NN,12))
293 * VL holds the computed left eigenvectors.
294 *
295 * LDVL (input) INTEGER
296 * Leading dimension of VL. Must be at least max(1,max(NN,12)).
297 *
298 * VR (workspace) DOUBLE PRECISION array, dimension
299 * (LDVR, max(NN,12))
300 * VR holds the computed right eigenvectors.
301 *
302 * LDVR (input) INTEGER
303 * Leading dimension of VR. Must be at least max(1,max(NN,12)).
304 *
305 * LRE (workspace) DOUBLE PRECISION array, dimension
306 * (LDLRE, max(NN,12))
307 * LRE holds the computed right or left eigenvectors.
308 *
309 * LDLRE (input) INTEGER
310 * Leading dimension of LRE. Must be at least max(1,max(NN,12))
311 *
312 * RCONDV (workspace) DOUBLE PRECISION array, dimension (N)
313 * RCONDV holds the computed reciprocal condition numbers
314 * for eigenvectors.
315 *
316 * RCNDV1 (workspace) DOUBLE PRECISION array, dimension (N)
317 * RCNDV1 holds more computed reciprocal condition numbers
318 * for eigenvectors.
319 *
320 * RCDVIN (workspace) DOUBLE PRECISION array, dimension (N)
321 * When COMP = .TRUE. RCDVIN holds the precomputed reciprocal
322 * condition numbers for eigenvectors to be compared with
323 * RCONDV.
324 *
325 * RCONDE (workspace) DOUBLE PRECISION array, dimension (N)
326 * RCONDE holds the computed reciprocal condition numbers
327 * for eigenvalues.
328 *
329 * RCNDE1 (workspace) DOUBLE PRECISION array, dimension (N)
330 * RCNDE1 holds more computed reciprocal condition numbers
331 * for eigenvalues.
332 *
333 * RCDEIN (workspace) DOUBLE PRECISION array, dimension (N)
334 * When COMP = .TRUE. RCDEIN holds the precomputed reciprocal
335 * condition numbers for eigenvalues to be compared with
336 * RCONDE.
337 *
338 * RESULT (output) DOUBLE PRECISION array, dimension (11)
339 * The values computed by the seven tests described above.
340 * The values are currently limited to 1/ulp, to avoid overflow.
341 *
342 * WORK (workspace) DOUBLE PRECISION array, dimension (NWORK)
343 *
344 * NWORK (input) INTEGER
345 * The number of entries in WORK. This must be at least
346 * max(6*12+2*12**2,6*NN(j)+2*NN(j)**2) =
347 * max( 360 ,6*NN(j)+2*NN(j)**2) for all j.
348 *
349 * IWORK (workspace) INTEGER array, dimension (2*max(NN,12))
350 *
351 * INFO (output) INTEGER
352 * If 0, then successful exit.
353 * If <0, then input paramter -INFO is incorrect.
354 * If >0, DLATMR, SLATMS, SLATME or DGET23 returned an error
355 * code, and INFO is its absolute value.
356 *
357 *-----------------------------------------------------------------------
358 *
359 * Some Local Variables and Parameters:
360 * ---- ----- --------- --- ----------
361 *
362 * ZERO, ONE Real 0 and 1.
363 * MAXTYP The number of types defined.
364 * NMAX Largest value in NN or 12.
365 * NERRS The number of tests which have exceeded THRESH
366 * COND, CONDS,
367 * IMODE Values to be passed to the matrix generators.
368 * ANORM Norm of A; passed to matrix generators.
369 *
370 * OVFL, UNFL Overflow and underflow thresholds.
371 * ULP, ULPINV Finest relative precision and its inverse.
372 * RTULP, RTULPI Square roots of the previous 4 values.
373 *
374 * The following four arrays decode JTYPE:
375 * KTYPE(j) The general type (1-10) for type "j".
376 * KMODE(j) The MODE value to be passed to the matrix
377 * generator for type "j".
378 * KMAGN(j) The order of magnitude ( O(1),
379 * O(overflow^(1/2) ), O(underflow^(1/2) )
380 * KCONDS(j) Selectw whether CONDS is to be 1 or
381 * 1/sqrt(ulp). (0 means irrelevant.)
382 *
383 * =====================================================================
384 *
385 * .. Parameters ..
386 DOUBLE PRECISION ZERO, ONE
387 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
388 INTEGER MAXTYP
389 PARAMETER ( MAXTYP = 21 )
390 * ..
391 * .. Local Scalars ..
392 LOGICAL BADNN
393 CHARACTER BALANC
394 CHARACTER*3 PATH
395 INTEGER I, IBAL, IINFO, IMODE, ITYPE, IWK, J, JCOL,
396 $ JSIZE, JTYPE, MTYPES, N, NERRS, NFAIL, NMAX,
397 $ NNWORK, NTEST, NTESTF, NTESTT
398 DOUBLE PRECISION ANORM, COND, CONDS, OVFL, RTULP, RTULPI, ULP,
399 $ ULPINV, UNFL
400 * ..
401 * .. Local Arrays ..
402 CHARACTER ADUMMA( 1 ), BAL( 4 )
403 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KCONDS( MAXTYP ),
404 $ KMAGN( MAXTYP ), KMODE( MAXTYP ),
405 $ KTYPE( MAXTYP )
406 * ..
407 * .. External Functions ..
408 DOUBLE PRECISION DLAMCH
409 EXTERNAL DLAMCH
410 * ..
411 * .. External Subroutines ..
412 EXTERNAL DGET23, DLABAD, DLASET, DLASUM, DLATME, DLATMR,
413 $ DLATMS, XERBLA
414 * ..
415 * .. Intrinsic Functions ..
416 INTRINSIC ABS, MAX, MIN, SQRT
417 * ..
418 * .. Data statements ..
419 DATA KTYPE / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
420 DATA KMAGN / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
421 $ 3, 1, 2, 3 /
422 DATA KMODE / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
423 $ 1, 5, 5, 5, 4, 3, 1 /
424 DATA KCONDS / 3*0, 5*0, 4*1, 6*2, 3*0 /
425 DATA BAL / 'N', 'P', 'S', 'B' /
426 * ..
427 * .. Executable Statements ..
428 *
429 PATH( 1: 1 ) = 'Double precision'
430 PATH( 2: 3 ) = 'VX'
431 *
432 * Check for errors
433 *
434 NTESTT = 0
435 NTESTF = 0
436 INFO = 0
437 *
438 * Important constants
439 *
440 BADNN = .FALSE.
441 *
442 * 12 is the largest dimension in the input file of precomputed
443 * problems
444 *
445 NMAX = 12
446 DO 10 J = 1, NSIZES
447 NMAX = MAX( NMAX, NN( J ) )
448 IF( NN( J ).LT.0 )
449 $ BADNN = .TRUE.
450 10 CONTINUE
451 *
452 * Check for errors
453 *
454 IF( NSIZES.LT.0 ) THEN
455 INFO = -1
456 ELSE IF( BADNN ) THEN
457 INFO = -2
458 ELSE IF( NTYPES.LT.0 ) THEN
459 INFO = -3
460 ELSE IF( THRESH.LT.ZERO ) THEN
461 INFO = -6
462 ELSE IF( LDA.LT.1 .OR. LDA.LT.NMAX ) THEN
463 INFO = -10
464 ELSE IF( LDVL.LT.1 .OR. LDVL.LT.NMAX ) THEN
465 INFO = -17
466 ELSE IF( LDVR.LT.1 .OR. LDVR.LT.NMAX ) THEN
467 INFO = -19
468 ELSE IF( LDLRE.LT.1 .OR. LDLRE.LT.NMAX ) THEN
469 INFO = -21
470 ELSE IF( 6*NMAX+2*NMAX**2.GT.NWORK ) THEN
471 INFO = -32
472 END IF
473 *
474 IF( INFO.NE.0 ) THEN
475 CALL XERBLA( 'DDRVVX', -INFO )
476 RETURN
477 END IF
478 *
479 * If nothing to do check on NIUNIT
480 *
481 IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 )
482 $ GO TO 160
483 *
484 * More Important constants
485 *
486 UNFL = DLAMCH( 'Safe minimum' )
487 OVFL = ONE / UNFL
488 CALL DLABAD( UNFL, OVFL )
489 ULP = DLAMCH( 'Precision' )
490 ULPINV = ONE / ULP
491 RTULP = SQRT( ULP )
492 RTULPI = ONE / RTULP
493 *
494 * Loop over sizes, types
495 *
496 NERRS = 0
497 *
498 DO 150 JSIZE = 1, NSIZES
499 N = NN( JSIZE )
500 IF( NSIZES.NE.1 ) THEN
501 MTYPES = MIN( MAXTYP, NTYPES )
502 ELSE
503 MTYPES = MIN( MAXTYP+1, NTYPES )
504 END IF
505 *
506 DO 140 JTYPE = 1, MTYPES
507 IF( .NOT.DOTYPE( JTYPE ) )
508 $ GO TO 140
509 *
510 * Save ISEED in case of an error.
511 *
512 DO 20 J = 1, 4
513 IOLDSD( J ) = ISEED( J )
514 20 CONTINUE
515 *
516 * Compute "A"
517 *
518 * Control parameters:
519 *
520 * KMAGN KCONDS KMODE KTYPE
521 * =1 O(1) 1 clustered 1 zero
522 * =2 large large clustered 2 identity
523 * =3 small exponential Jordan
524 * =4 arithmetic diagonal, (w/ eigenvalues)
525 * =5 random log symmetric, w/ eigenvalues
526 * =6 random general, w/ eigenvalues
527 * =7 random diagonal
528 * =8 random symmetric
529 * =9 random general
530 * =10 random triangular
531 *
532 IF( MTYPES.GT.MAXTYP )
533 $ GO TO 90
534 *
535 ITYPE = KTYPE( JTYPE )
536 IMODE = KMODE( JTYPE )
537 *
538 * Compute norm
539 *
540 GO TO ( 30, 40, 50 )KMAGN( JTYPE )
541 *
542 30 CONTINUE
543 ANORM = ONE
544 GO TO 60
545 *
546 40 CONTINUE
547 ANORM = OVFL*ULP
548 GO TO 60
549 *
550 50 CONTINUE
551 ANORM = UNFL*ULPINV
552 GO TO 60
553 *
554 60 CONTINUE
555 *
556 CALL DLASET( 'Full', LDA, N, ZERO, ZERO, A, LDA )
557 IINFO = 0
558 COND = ULPINV
559 *
560 * Special Matrices -- Identity & Jordan block
561 *
562 * Zero
563 *
564 IF( ITYPE.EQ.1 ) THEN
565 IINFO = 0
566 *
567 ELSE IF( ITYPE.EQ.2 ) THEN
568 *
569 * Identity
570 *
571 DO 70 JCOL = 1, N
572 A( JCOL, JCOL ) = ANORM
573 70 CONTINUE
574 *
575 ELSE IF( ITYPE.EQ.3 ) THEN
576 *
577 * Jordan Block
578 *
579 DO 80 JCOL = 1, N
580 A( JCOL, JCOL ) = ANORM
581 IF( JCOL.GT.1 )
582 $ A( JCOL, JCOL-1 ) = ONE
583 80 CONTINUE
584 *
585 ELSE IF( ITYPE.EQ.4 ) THEN
586 *
587 * Diagonal Matrix, [Eigen]values Specified
588 *
589 CALL DLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND,
590 $ ANORM, 0, 0, 'N', A, LDA, WORK( N+1 ),
591 $ IINFO )
592 *
593 ELSE IF( ITYPE.EQ.5 ) THEN
594 *
595 * Symmetric, eigenvalues specified
596 *
597 CALL DLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND,
598 $ ANORM, N, N, 'N', A, LDA, WORK( N+1 ),
599 $ IINFO )
600 *
601 ELSE IF( ITYPE.EQ.6 ) THEN
602 *
603 * General, eigenvalues specified
604 *
605 IF( KCONDS( JTYPE ).EQ.1 ) THEN
606 CONDS = ONE
607 ELSE IF( KCONDS( JTYPE ).EQ.2 ) THEN
608 CONDS = RTULPI
609 ELSE
610 CONDS = ZERO
611 END IF
612 *
613 ADUMMA( 1 ) = ' '
614 CALL DLATME( N, 'S', ISEED, WORK, IMODE, COND, ONE,
615 $ ADUMMA, 'T', 'T', 'T', WORK( N+1 ), 4,
616 $ CONDS, N, N, ANORM, A, LDA, WORK( 2*N+1 ),
617 $ IINFO )
618 *
619 ELSE IF( ITYPE.EQ.7 ) THEN
620 *
621 * Diagonal, random eigenvalues
622 *
623 CALL DLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE,
624 $ 'T', 'N', WORK( N+1 ), 1, ONE,
625 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0,
626 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
627 *
628 ELSE IF( ITYPE.EQ.8 ) THEN
629 *
630 * Symmetric, random eigenvalues
631 *
632 CALL DLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE,
633 $ 'T', 'N', WORK( N+1 ), 1, ONE,
634 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N,
635 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
636 *
637 ELSE IF( ITYPE.EQ.9 ) THEN
638 *
639 * General, random eigenvalues
640 *
641 CALL DLATMR( N, N, 'S', ISEED, 'N', WORK, 6, ONE, ONE,
642 $ 'T', 'N', WORK( N+1 ), 1, ONE,
643 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N,
644 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
645 IF( N.GE.4 ) THEN
646 CALL DLASET( 'Full', 2, N, ZERO, ZERO, A, LDA )
647 CALL DLASET( 'Full', N-3, 1, ZERO, ZERO, A( 3, 1 ),
648 $ LDA )
649 CALL DLASET( 'Full', N-3, 2, ZERO, ZERO, A( 3, N-1 ),
650 $ LDA )
651 CALL DLASET( 'Full', 1, N, ZERO, ZERO, A( N, 1 ),
652 $ LDA )
653 END IF
654 *
655 ELSE IF( ITYPE.EQ.10 ) THEN
656 *
657 * Triangular, random eigenvalues
658 *
659 CALL DLATMR( N, N, 'S', ISEED, 'N', WORK, 6, ONE, ONE,
660 $ 'T', 'N', WORK( N+1 ), 1, ONE,
661 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, 0,
662 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
663 *
664 ELSE
665 *
666 IINFO = 1
667 END IF
668 *
669 IF( IINFO.NE.0 ) THEN
670 WRITE( NOUNIT, FMT = 9992 )'Generator', IINFO, N, JTYPE,
671 $ IOLDSD
672 INFO = ABS( IINFO )
673 RETURN
674 END IF
675 *
676 90 CONTINUE
677 *
678 * Test for minimal and generous workspace
679 *
680 DO 130 IWK = 1, 3
681 IF( IWK.EQ.1 ) THEN
682 NNWORK = 3*N
683 ELSE IF( IWK.EQ.2 ) THEN
684 NNWORK = 6*N + N**2
685 ELSE
686 NNWORK = 6*N + 2*N**2
687 END IF
688 NNWORK = MAX( NNWORK, 1 )
689 *
690 * Test for all balancing options
691 *
692 DO 120 IBAL = 1, 4
693 BALANC = BAL( IBAL )
694 *
695 * Perform tests
696 *
697 CALL DGET23( .FALSE., BALANC, JTYPE, THRESH, IOLDSD,
698 $ NOUNIT, N, A, LDA, H, WR, WI, WR1, WI1,
699 $ VL, LDVL, VR, LDVR, LRE, LDLRE, RCONDV,
700 $ RCNDV1, RCDVIN, RCONDE, RCNDE1, RCDEIN,
701 $ SCALE, SCALE1, RESULT, WORK, NNWORK,
702 $ IWORK, INFO )
703 *
704 * Check for RESULT(j) > THRESH
705 *
706 NTEST = 0
707 NFAIL = 0
708 DO 100 J = 1, 9
709 IF( RESULT( J ).GE.ZERO )
710 $ NTEST = NTEST + 1
711 IF( RESULT( J ).GE.THRESH )
712 $ NFAIL = NFAIL + 1
713 100 CONTINUE
714 *
715 IF( NFAIL.GT.0 )
716 $ NTESTF = NTESTF + 1
717 IF( NTESTF.EQ.1 ) THEN
718 WRITE( NOUNIT, FMT = 9999 )PATH
719 WRITE( NOUNIT, FMT = 9998 )
720 WRITE( NOUNIT, FMT = 9997 )
721 WRITE( NOUNIT, FMT = 9996 )
722 WRITE( NOUNIT, FMT = 9995 )THRESH
723 NTESTF = 2
724 END IF
725 *
726 DO 110 J = 1, 9
727 IF( RESULT( J ).GE.THRESH ) THEN
728 WRITE( NOUNIT, FMT = 9994 )BALANC, N, IWK,
729 $ IOLDSD, JTYPE, J, RESULT( J )
730 END IF
731 110 CONTINUE
732 *
733 NERRS = NERRS + NFAIL
734 NTESTT = NTESTT + NTEST
735 *
736 120 CONTINUE
737 130 CONTINUE
738 140 CONTINUE
739 150 CONTINUE
740 *
741 160 CONTINUE
742 *
743 * Read in data from file to check accuracy of condition estimation.
744 * Assume input eigenvalues are sorted lexicographically (increasing
745 * by real part, then decreasing by imaginary part)
746 *
747 JTYPE = 0
748 170 CONTINUE
749 READ( NIUNIT, FMT = *, END = 220 )N
750 *
751 * Read input data until N=0
752 *
753 IF( N.EQ.0 )
754 $ GO TO 220
755 JTYPE = JTYPE + 1
756 ISEED( 1 ) = JTYPE
757 DO 180 I = 1, N
758 READ( NIUNIT, FMT = * )( A( I, J ), J = 1, N )
759 180 CONTINUE
760 DO 190 I = 1, N
761 READ( NIUNIT, FMT = * )WR1( I ), WI1( I ), RCDEIN( I ),
762 $ RCDVIN( I )
763 190 CONTINUE
764 CALL DGET23( .TRUE., 'N', 22, THRESH, ISEED, NOUNIT, N, A, LDA, H,
765 $ WR, WI, WR1, WI1, VL, LDVL, VR, LDVR, LRE, LDLRE,
766 $ RCONDV, RCNDV1, RCDVIN, RCONDE, RCNDE1, RCDEIN,
767 $ SCALE, SCALE1, RESULT, WORK, 6*N+2*N**2, IWORK,
768 $ INFO )
769 *
770 * Check for RESULT(j) > THRESH
771 *
772 NTEST = 0
773 NFAIL = 0
774 DO 200 J = 1, 11
775 IF( RESULT( J ).GE.ZERO )
776 $ NTEST = NTEST + 1
777 IF( RESULT( J ).GE.THRESH )
778 $ NFAIL = NFAIL + 1
779 200 CONTINUE
780 *
781 IF( NFAIL.GT.0 )
782 $ NTESTF = NTESTF + 1
783 IF( NTESTF.EQ.1 ) THEN
784 WRITE( NOUNIT, FMT = 9999 )PATH
785 WRITE( NOUNIT, FMT = 9998 )
786 WRITE( NOUNIT, FMT = 9997 )
787 WRITE( NOUNIT, FMT = 9996 )
788 WRITE( NOUNIT, FMT = 9995 )THRESH
789 NTESTF = 2
790 END IF
791 *
792 DO 210 J = 1, 11
793 IF( RESULT( J ).GE.THRESH ) THEN
794 WRITE( NOUNIT, FMT = 9993 )N, JTYPE, J, RESULT( J )
795 END IF
796 210 CONTINUE
797 *
798 NERRS = NERRS + NFAIL
799 NTESTT = NTESTT + NTEST
800 GO TO 170
801 220 CONTINUE
802 *
803 * Summary
804 *
805 CALL DLASUM( PATH, NOUNIT, NERRS, NTESTT )
806 *
807 9999 FORMAT( / 1X, A3, ' -- Real Eigenvalue-Eigenvector Decomposition',
808 $ ' Expert Driver', /
809 $ ' Matrix types (see DDRVVX for details): ' )
810 *
811 9998 FORMAT( / ' Special Matrices:', / ' 1=Zero matrix. ',
812 $ ' ', ' 5=Diagonal: geometr. spaced entries.',
813 $ / ' 2=Identity matrix. ', ' 6=Diagona',
814 $ 'l: clustered entries.', / ' 3=Transposed Jordan block. ',
815 $ ' ', ' 7=Diagonal: large, evenly spaced.', / ' ',
816 $ '4=Diagonal: evenly spaced entries. ', ' 8=Diagonal: s',
817 $ 'mall, evenly spaced.' )
818 9997 FORMAT( ' Dense, Non-Symmetric Matrices:', / ' 9=Well-cond., ev',
819 $ 'enly spaced eigenvals.', ' 14=Ill-cond., geomet. spaced e',
820 $ 'igenals.', / ' 10=Well-cond., geom. spaced eigenvals. ',
821 $ ' 15=Ill-conditioned, clustered e.vals.', / ' 11=Well-cond',
822 $ 'itioned, clustered e.vals. ', ' 16=Ill-cond., random comp',
823 $ 'lex ', / ' 12=Well-cond., random complex ', ' ',
824 $ ' 17=Ill-cond., large rand. complx ', / ' 13=Ill-condi',
825 $ 'tioned, evenly spaced. ', ' 18=Ill-cond., small rand.',
826 $ ' complx ' )
827 9996 FORMAT( ' 19=Matrix with random O(1) entries. ', ' 21=Matrix ',
828 $ 'with small random entries.', / ' 20=Matrix with large ran',
829 $ 'dom entries. ', ' 22=Matrix read from input file', / )
830 9995 FORMAT( ' Tests performed with test threshold =', F8.2,
831 $ / / ' 1 = | A VR - VR W | / ( n |A| ulp ) ',
832 $ / ' 2 = | transpose(A) VL - VL W | / ( n |A| ulp ) ',
833 $ / ' 3 = | |VR(i)| - 1 | / ulp ',
834 $ / ' 4 = | |VL(i)| - 1 | / ulp ',
835 $ / ' 5 = 0 if W same no matter if VR or VL computed,',
836 $ ' 1/ulp otherwise', /
837 $ ' 6 = 0 if VR same no matter what else computed,',
838 $ ' 1/ulp otherwise', /
839 $ ' 7 = 0 if VL same no matter what else computed,',
840 $ ' 1/ulp otherwise', /
841 $ ' 8 = 0 if RCONDV same no matter what else computed,',
842 $ ' 1/ulp otherwise', /
843 $ ' 9 = 0 if SCALE, ILO, IHI, ABNRM same no matter what else',
844 $ ' computed, 1/ulp otherwise',
845 $ / ' 10 = | RCONDV - RCONDV(precomputed) | / cond(RCONDV),',
846 $ / ' 11 = | RCONDE - RCONDE(precomputed) | / cond(RCONDE),' )
847 9994 FORMAT( ' BALANC=''', A1, ''',N=', I4, ',IWK=', I1, ', seed=',
848 $ 4( I4, ',' ), ' type ', I2, ', test(', I2, ')=', G10.3 )
849 9993 FORMAT( ' N=', I5, ', input example =', I3, ', test(', I2, ')=',
850 $ G10.3 )
851 9992 FORMAT( ' DDRVVX: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
852 $ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
853 *
854 RETURN
855 *
856 * End of DDRVVX
857 *
858 END
2 $ NIUNIT, NOUNIT, A, LDA, H, WR, WI, WR1, WI1,
3 $ VL, LDVL, VR, LDVR, LRE, LDLRE, RCONDV, RCNDV1,
4 $ RCDVIN, RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1,
5 $ RESULT, WORK, NWORK, IWORK, INFO )
6 *
7 * -- LAPACK test routine (version 3.1) --
8 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
9 * November 2006
10 *
11 * .. Scalar Arguments ..
12 INTEGER INFO, LDA, LDLRE, LDVL, LDVR, NIUNIT, NOUNIT,
13 $ NSIZES, NTYPES, NWORK
14 DOUBLE PRECISION THRESH
15 * ..
16 * .. Array Arguments ..
17 LOGICAL DOTYPE( * )
18 INTEGER ISEED( 4 ), IWORK( * ), NN( * )
19 DOUBLE PRECISION A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
20 $ RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
21 $ RCNDV1( * ), RCONDE( * ), RCONDV( * ),
22 $ RESULT( 11 ), SCALE( * ), SCALE1( * ),
23 $ VL( LDVL, * ), VR( LDVR, * ), WI( * ),
24 $ WI1( * ), WORK( * ), WR( * ), WR1( * )
25 * ..
26 *
27 * Purpose
28 * =======
29 *
30 * DDRVVX checks the nonsymmetric eigenvalue problem expert driver
31 * DGEEVX.
32 *
33 * DDRVVX uses both test matrices generated randomly depending on
34 * data supplied in the calling sequence, as well as on data
35 * read from an input file and including precomputed condition
36 * numbers to which it compares the ones it computes.
37 *
38 * When DDRVVX is called, a number of matrix "sizes" ("n's") and a
39 * number of matrix "types" are specified in the calling sequence.
40 * For each size ("n") and each type of matrix, one matrix will be
41 * generated and used to test the nonsymmetric eigenroutines. For
42 * each matrix, 9 tests will be performed:
43 *
44 * (1) | A * VR - VR * W | / ( n |A| ulp )
45 *
46 * Here VR is the matrix of unit right eigenvectors.
47 * W is a block diagonal matrix, with a 1x1 block for each
48 * real eigenvalue and a 2x2 block for each complex conjugate
49 * pair. If eigenvalues j and j+1 are a complex conjugate pair,
50 * so WR(j) = WR(j+1) = wr and WI(j) = - WI(j+1) = wi, then the
51 * 2 x 2 block corresponding to the pair will be:
52 *
53 * ( wr wi )
54 * ( -wi wr )
55 *
56 * Such a block multiplying an n x 2 matrix ( ur ui ) on the
57 * right will be the same as multiplying ur + i*ui by wr + i*wi.
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 signs.
118 * (ULP = (first number larger than 1) - 1 )
119 * (5) A diagonal matrix with geometrically spaced entries
120 * 1, ..., ULP and random signs.
121 * (6) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
122 * and random signs.
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 orthogonal and
130 * T has evenly spaced entries 1, ..., ULP with random signs
131 * 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 orthogonal and
135 * T has geometrically spaced entries 1, ..., ULP with random
136 * signs on the diagonal and random O(1) entries in the upper
137 * triangle.
138 *
139 * (11) A matrix of the form U' T U, where U is orthogonal and
140 * T has "clustered" entries 1, ULP,..., ULP with random
141 * signs on the diagonal and random O(1) entries in the upper
142 * triangle.
143 *
144 * (12) A matrix of the form U' T U, where U is orthogonal and
145 * T has real or complex conjugate paired eigenvalues randomly
146 * chosen from ( ULP, 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 signs on the diagonal and random O(1) entries
152 * 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 signs on the diagonal and random
157 * 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 signs on the diagonal and random O(1) entries
162 * 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 real or complex conjugate paired
166 * eigenvalues randomly chosen from ( ULP, 1 ) and random
167 * O(1) entries in the upper 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 (-1,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 DGEEVX 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 DGEEVX 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 *
206 * Arguments
207 * ==========
208 *
209 * NSIZES (input) INTEGER
210 * The number of sizes of matrices to use. NSIZES must be at
211 * least zero. If it is zero, no randomly generated matrices
212 * are tested, but any test matrices read from NIUNIT will be
213 * tested.
214 *
215 * NN (input) INTEGER array, dimension (NSIZES)
216 * An array containing the sizes to be used for the matrices.
217 * Zero values will be skipped. The values must be at least
218 * zero.
219 *
220 * NTYPES (input) INTEGER
221 * The number of elements in DOTYPE. NTYPES must be at least
222 * zero. If it is zero, no randomly generated test matrices
223 * are tested, but and test matrices read from NIUNIT will be
224 * tested. If it is MAXTYP+1 and NSIZES is 1, then an
225 * additional type, MAXTYP+1 is defined, which is to use
226 * whatever matrix is in A. This is only useful if
227 * DOTYPE(1:MAXTYP) is .FALSE. and DOTYPE(MAXTYP+1) is .TRUE. .
228 *
229 * DOTYPE (input) LOGICAL array, dimension (NTYPES)
230 * If DOTYPE(j) is .TRUE., then for each size in NN a
231 * matrix of that size and of type j will be generated.
232 * If NTYPES is smaller than the maximum number of types
233 * defined (PARAMETER MAXTYP), then types NTYPES+1 through
234 * MAXTYP will not be generated. If NTYPES is larger
235 * than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
236 * will be ignored.
237 *
238 * ISEED (input/output) INTEGER array, dimension (4)
239 * On entry ISEED specifies the seed of the random number
240 * generator. The array elements should be between 0 and 4095;
241 * if not they will be reduced mod 4096. Also, ISEED(4) must
242 * be odd. The random number generator uses a linear
243 * congruential sequence limited to small integers, and so
244 * should produce machine independent random numbers. The
245 * values of ISEED are changed on exit, and can be used in the
246 * next call to DDRVVX to continue the same random number
247 * sequence.
248 *
249 * THRESH (input) DOUBLE PRECISION
250 * A test will count as "failed" if the "error", computed as
251 * described above, exceeds THRESH. Note that the error
252 * is scaled to be O(1), so THRESH should be a reasonably
253 * small multiple of 1, e.g., 10 or 100. In particular,
254 * it should not depend on the precision (single vs. double)
255 * or the size of the matrix. It must be at least zero.
256 *
257 * NIUNIT (input) INTEGER
258 * The FORTRAN unit number for reading in the data file of
259 * problems to solve.
260 *
261 * NOUNIT (input) INTEGER
262 * The FORTRAN unit number for printing out error messages
263 * (e.g., if a routine returns INFO not equal to 0.)
264 *
265 * A (workspace) DOUBLE PRECISION array, dimension
266 * (LDA, max(NN,12))
267 * Used to hold the matrix whose eigenvalues are to be
268 * computed. On exit, A contains the last matrix actually used.
269 *
270 * LDA (input) INTEGER
271 * The leading dimension of the arrays A and H.
272 * LDA >= max(NN,12), since 12 is the dimension of the largest
273 * matrix in the precomputed input file.
274 *
275 * H (workspace) DOUBLE PRECISION array, dimension
276 * (LDA, max(NN,12))
277 * Another copy of the test matrix A, modified by DGEEVX.
278 *
279 * WR (workspace) DOUBLE PRECISION array, dimension (max(NN))
280 * WI (workspace) DOUBLE PRECISION array, dimension (max(NN))
281 * The real and imaginary parts of the eigenvalues of A.
282 * On exit, WR + WI*i are the eigenvalues of the matrix in A.
283 *
284 * WR1 (workspace) DOUBLE PRECISION array, dimension (max(NN,12))
285 * WI1 (workspace) DOUBLE PRECISION array, dimension (max(NN,12))
286 * Like WR, WI, these arrays contain the eigenvalues of A,
287 * but those computed when DGEEVX only computes a partial
288 * eigendecomposition, i.e. not the eigenvalues and left
289 * and right eigenvectors.
290 *
291 * VL (workspace) DOUBLE PRECISION array, dimension
292 * (LDVL, max(NN,12))
293 * VL holds the computed left eigenvectors.
294 *
295 * LDVL (input) INTEGER
296 * Leading dimension of VL. Must be at least max(1,max(NN,12)).
297 *
298 * VR (workspace) DOUBLE PRECISION array, dimension
299 * (LDVR, max(NN,12))
300 * VR holds the computed right eigenvectors.
301 *
302 * LDVR (input) INTEGER
303 * Leading dimension of VR. Must be at least max(1,max(NN,12)).
304 *
305 * LRE (workspace) DOUBLE PRECISION array, dimension
306 * (LDLRE, max(NN,12))
307 * LRE holds the computed right or left eigenvectors.
308 *
309 * LDLRE (input) INTEGER
310 * Leading dimension of LRE. Must be at least max(1,max(NN,12))
311 *
312 * RCONDV (workspace) DOUBLE PRECISION array, dimension (N)
313 * RCONDV holds the computed reciprocal condition numbers
314 * for eigenvectors.
315 *
316 * RCNDV1 (workspace) DOUBLE PRECISION array, dimension (N)
317 * RCNDV1 holds more computed reciprocal condition numbers
318 * for eigenvectors.
319 *
320 * RCDVIN (workspace) DOUBLE PRECISION array, dimension (N)
321 * When COMP = .TRUE. RCDVIN holds the precomputed reciprocal
322 * condition numbers for eigenvectors to be compared with
323 * RCONDV.
324 *
325 * RCONDE (workspace) DOUBLE PRECISION array, dimension (N)
326 * RCONDE holds the computed reciprocal condition numbers
327 * for eigenvalues.
328 *
329 * RCNDE1 (workspace) DOUBLE PRECISION array, dimension (N)
330 * RCNDE1 holds more computed reciprocal condition numbers
331 * for eigenvalues.
332 *
333 * RCDEIN (workspace) DOUBLE PRECISION array, dimension (N)
334 * When COMP = .TRUE. RCDEIN holds the precomputed reciprocal
335 * condition numbers for eigenvalues to be compared with
336 * RCONDE.
337 *
338 * RESULT (output) DOUBLE PRECISION array, dimension (11)
339 * The values computed by the seven tests described above.
340 * The values are currently limited to 1/ulp, to avoid overflow.
341 *
342 * WORK (workspace) DOUBLE PRECISION array, dimension (NWORK)
343 *
344 * NWORK (input) INTEGER
345 * The number of entries in WORK. This must be at least
346 * max(6*12+2*12**2,6*NN(j)+2*NN(j)**2) =
347 * max( 360 ,6*NN(j)+2*NN(j)**2) for all j.
348 *
349 * IWORK (workspace) INTEGER array, dimension (2*max(NN,12))
350 *
351 * INFO (output) INTEGER
352 * If 0, then successful exit.
353 * If <0, then input paramter -INFO is incorrect.
354 * If >0, DLATMR, SLATMS, SLATME or DGET23 returned an error
355 * code, and INFO is its absolute value.
356 *
357 *-----------------------------------------------------------------------
358 *
359 * Some Local Variables and Parameters:
360 * ---- ----- --------- --- ----------
361 *
362 * ZERO, ONE Real 0 and 1.
363 * MAXTYP The number of types defined.
364 * NMAX Largest value in NN or 12.
365 * NERRS The number of tests which have exceeded THRESH
366 * COND, CONDS,
367 * IMODE Values to be passed to the matrix generators.
368 * ANORM Norm of A; passed to matrix generators.
369 *
370 * OVFL, UNFL Overflow and underflow thresholds.
371 * ULP, ULPINV Finest relative precision and its inverse.
372 * RTULP, RTULPI Square roots of the previous 4 values.
373 *
374 * The following four arrays decode JTYPE:
375 * KTYPE(j) The general type (1-10) for type "j".
376 * KMODE(j) The MODE value to be passed to the matrix
377 * generator for type "j".
378 * KMAGN(j) The order of magnitude ( O(1),
379 * O(overflow^(1/2) ), O(underflow^(1/2) )
380 * KCONDS(j) Selectw whether CONDS is to be 1 or
381 * 1/sqrt(ulp). (0 means irrelevant.)
382 *
383 * =====================================================================
384 *
385 * .. Parameters ..
386 DOUBLE PRECISION ZERO, ONE
387 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 )
388 INTEGER MAXTYP
389 PARAMETER ( MAXTYP = 21 )
390 * ..
391 * .. Local Scalars ..
392 LOGICAL BADNN
393 CHARACTER BALANC
394 CHARACTER*3 PATH
395 INTEGER I, IBAL, IINFO, IMODE, ITYPE, IWK, J, JCOL,
396 $ JSIZE, JTYPE, MTYPES, N, NERRS, NFAIL, NMAX,
397 $ NNWORK, NTEST, NTESTF, NTESTT
398 DOUBLE PRECISION ANORM, COND, CONDS, OVFL, RTULP, RTULPI, ULP,
399 $ ULPINV, UNFL
400 * ..
401 * .. Local Arrays ..
402 CHARACTER ADUMMA( 1 ), BAL( 4 )
403 INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KCONDS( MAXTYP ),
404 $ KMAGN( MAXTYP ), KMODE( MAXTYP ),
405 $ KTYPE( MAXTYP )
406 * ..
407 * .. External Functions ..
408 DOUBLE PRECISION DLAMCH
409 EXTERNAL DLAMCH
410 * ..
411 * .. External Subroutines ..
412 EXTERNAL DGET23, DLABAD, DLASET, DLASUM, DLATME, DLATMR,
413 $ DLATMS, XERBLA
414 * ..
415 * .. Intrinsic Functions ..
416 INTRINSIC ABS, MAX, MIN, SQRT
417 * ..
418 * .. Data statements ..
419 DATA KTYPE / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
420 DATA KMAGN / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
421 $ 3, 1, 2, 3 /
422 DATA KMODE / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
423 $ 1, 5, 5, 5, 4, 3, 1 /
424 DATA KCONDS / 3*0, 5*0, 4*1, 6*2, 3*0 /
425 DATA BAL / 'N', 'P', 'S', 'B' /
426 * ..
427 * .. Executable Statements ..
428 *
429 PATH( 1: 1 ) = 'Double precision'
430 PATH( 2: 3 ) = 'VX'
431 *
432 * Check for errors
433 *
434 NTESTT = 0
435 NTESTF = 0
436 INFO = 0
437 *
438 * Important constants
439 *
440 BADNN = .FALSE.
441 *
442 * 12 is the largest dimension in the input file of precomputed
443 * problems
444 *
445 NMAX = 12
446 DO 10 J = 1, NSIZES
447 NMAX = MAX( NMAX, NN( J ) )
448 IF( NN( J ).LT.0 )
449 $ BADNN = .TRUE.
450 10 CONTINUE
451 *
452 * Check for errors
453 *
454 IF( NSIZES.LT.0 ) THEN
455 INFO = -1
456 ELSE IF( BADNN ) THEN
457 INFO = -2
458 ELSE IF( NTYPES.LT.0 ) THEN
459 INFO = -3
460 ELSE IF( THRESH.LT.ZERO ) THEN
461 INFO = -6
462 ELSE IF( LDA.LT.1 .OR. LDA.LT.NMAX ) THEN
463 INFO = -10
464 ELSE IF( LDVL.LT.1 .OR. LDVL.LT.NMAX ) THEN
465 INFO = -17
466 ELSE IF( LDVR.LT.1 .OR. LDVR.LT.NMAX ) THEN
467 INFO = -19
468 ELSE IF( LDLRE.LT.1 .OR. LDLRE.LT.NMAX ) THEN
469 INFO = -21
470 ELSE IF( 6*NMAX+2*NMAX**2.GT.NWORK ) THEN
471 INFO = -32
472 END IF
473 *
474 IF( INFO.NE.0 ) THEN
475 CALL XERBLA( 'DDRVVX', -INFO )
476 RETURN
477 END IF
478 *
479 * If nothing to do check on NIUNIT
480 *
481 IF( NSIZES.EQ.0 .OR. NTYPES.EQ.0 )
482 $ GO TO 160
483 *
484 * More Important constants
485 *
486 UNFL = DLAMCH( 'Safe minimum' )
487 OVFL = ONE / UNFL
488 CALL DLABAD( UNFL, OVFL )
489 ULP = DLAMCH( 'Precision' )
490 ULPINV = ONE / ULP
491 RTULP = SQRT( ULP )
492 RTULPI = ONE / RTULP
493 *
494 * Loop over sizes, types
495 *
496 NERRS = 0
497 *
498 DO 150 JSIZE = 1, NSIZES
499 N = NN( JSIZE )
500 IF( NSIZES.NE.1 ) THEN
501 MTYPES = MIN( MAXTYP, NTYPES )
502 ELSE
503 MTYPES = MIN( MAXTYP+1, NTYPES )
504 END IF
505 *
506 DO 140 JTYPE = 1, MTYPES
507 IF( .NOT.DOTYPE( JTYPE ) )
508 $ GO TO 140
509 *
510 * Save ISEED in case of an error.
511 *
512 DO 20 J = 1, 4
513 IOLDSD( J ) = ISEED( J )
514 20 CONTINUE
515 *
516 * Compute "A"
517 *
518 * Control parameters:
519 *
520 * KMAGN KCONDS KMODE KTYPE
521 * =1 O(1) 1 clustered 1 zero
522 * =2 large large clustered 2 identity
523 * =3 small exponential Jordan
524 * =4 arithmetic diagonal, (w/ eigenvalues)
525 * =5 random log symmetric, w/ eigenvalues
526 * =6 random general, w/ eigenvalues
527 * =7 random diagonal
528 * =8 random symmetric
529 * =9 random general
530 * =10 random triangular
531 *
532 IF( MTYPES.GT.MAXTYP )
533 $ GO TO 90
534 *
535 ITYPE = KTYPE( JTYPE )
536 IMODE = KMODE( JTYPE )
537 *
538 * Compute norm
539 *
540 GO TO ( 30, 40, 50 )KMAGN( JTYPE )
541 *
542 30 CONTINUE
543 ANORM = ONE
544 GO TO 60
545 *
546 40 CONTINUE
547 ANORM = OVFL*ULP
548 GO TO 60
549 *
550 50 CONTINUE
551 ANORM = UNFL*ULPINV
552 GO TO 60
553 *
554 60 CONTINUE
555 *
556 CALL DLASET( 'Full', LDA, N, ZERO, ZERO, A, LDA )
557 IINFO = 0
558 COND = ULPINV
559 *
560 * Special Matrices -- Identity & Jordan block
561 *
562 * Zero
563 *
564 IF( ITYPE.EQ.1 ) THEN
565 IINFO = 0
566 *
567 ELSE IF( ITYPE.EQ.2 ) THEN
568 *
569 * Identity
570 *
571 DO 70 JCOL = 1, N
572 A( JCOL, JCOL ) = ANORM
573 70 CONTINUE
574 *
575 ELSE IF( ITYPE.EQ.3 ) THEN
576 *
577 * Jordan Block
578 *
579 DO 80 JCOL = 1, N
580 A( JCOL, JCOL ) = ANORM
581 IF( JCOL.GT.1 )
582 $ A( JCOL, JCOL-1 ) = ONE
583 80 CONTINUE
584 *
585 ELSE IF( ITYPE.EQ.4 ) THEN
586 *
587 * Diagonal Matrix, [Eigen]values Specified
588 *
589 CALL DLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND,
590 $ ANORM, 0, 0, 'N', A, LDA, WORK( N+1 ),
591 $ IINFO )
592 *
593 ELSE IF( ITYPE.EQ.5 ) THEN
594 *
595 * Symmetric, eigenvalues specified
596 *
597 CALL DLATMS( N, N, 'S', ISEED, 'S', WORK, IMODE, COND,
598 $ ANORM, N, N, 'N', A, LDA, WORK( N+1 ),
599 $ IINFO )
600 *
601 ELSE IF( ITYPE.EQ.6 ) THEN
602 *
603 * General, eigenvalues specified
604 *
605 IF( KCONDS( JTYPE ).EQ.1 ) THEN
606 CONDS = ONE
607 ELSE IF( KCONDS( JTYPE ).EQ.2 ) THEN
608 CONDS = RTULPI
609 ELSE
610 CONDS = ZERO
611 END IF
612 *
613 ADUMMA( 1 ) = ' '
614 CALL DLATME( N, 'S', ISEED, WORK, IMODE, COND, ONE,
615 $ ADUMMA, 'T', 'T', 'T', WORK( N+1 ), 4,
616 $ CONDS, N, N, ANORM, A, LDA, WORK( 2*N+1 ),
617 $ IINFO )
618 *
619 ELSE IF( ITYPE.EQ.7 ) THEN
620 *
621 * Diagonal, random eigenvalues
622 *
623 CALL DLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE,
624 $ 'T', 'N', WORK( N+1 ), 1, ONE,
625 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, 0, 0,
626 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
627 *
628 ELSE IF( ITYPE.EQ.8 ) THEN
629 *
630 * Symmetric, random eigenvalues
631 *
632 CALL DLATMR( N, N, 'S', ISEED, 'S', WORK, 6, ONE, ONE,
633 $ 'T', 'N', WORK( N+1 ), 1, ONE,
634 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N,
635 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
636 *
637 ELSE IF( ITYPE.EQ.9 ) THEN
638 *
639 * General, random eigenvalues
640 *
641 CALL DLATMR( N, N, 'S', ISEED, 'N', WORK, 6, ONE, ONE,
642 $ 'T', 'N', WORK( N+1 ), 1, ONE,
643 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, N,
644 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
645 IF( N.GE.4 ) THEN
646 CALL DLASET( 'Full', 2, N, ZERO, ZERO, A, LDA )
647 CALL DLASET( 'Full', N-3, 1, ZERO, ZERO, A( 3, 1 ),
648 $ LDA )
649 CALL DLASET( 'Full', N-3, 2, ZERO, ZERO, A( 3, N-1 ),
650 $ LDA )
651 CALL DLASET( 'Full', 1, N, ZERO, ZERO, A( N, 1 ),
652 $ LDA )
653 END IF
654 *
655 ELSE IF( ITYPE.EQ.10 ) THEN
656 *
657 * Triangular, random eigenvalues
658 *
659 CALL DLATMR( N, N, 'S', ISEED, 'N', WORK, 6, ONE, ONE,
660 $ 'T', 'N', WORK( N+1 ), 1, ONE,
661 $ WORK( 2*N+1 ), 1, ONE, 'N', IDUMMA, N, 0,
662 $ ZERO, ANORM, 'NO', A, LDA, IWORK, IINFO )
663 *
664 ELSE
665 *
666 IINFO = 1
667 END IF
668 *
669 IF( IINFO.NE.0 ) THEN
670 WRITE( NOUNIT, FMT = 9992 )'Generator', IINFO, N, JTYPE,
671 $ IOLDSD
672 INFO = ABS( IINFO )
673 RETURN
674 END IF
675 *
676 90 CONTINUE
677 *
678 * Test for minimal and generous workspace
679 *
680 DO 130 IWK = 1, 3
681 IF( IWK.EQ.1 ) THEN
682 NNWORK = 3*N
683 ELSE IF( IWK.EQ.2 ) THEN
684 NNWORK = 6*N + N**2
685 ELSE
686 NNWORK = 6*N + 2*N**2
687 END IF
688 NNWORK = MAX( NNWORK, 1 )
689 *
690 * Test for all balancing options
691 *
692 DO 120 IBAL = 1, 4
693 BALANC = BAL( IBAL )
694 *
695 * Perform tests
696 *
697 CALL DGET23( .FALSE., BALANC, JTYPE, THRESH, IOLDSD,
698 $ NOUNIT, N, A, LDA, H, WR, WI, WR1, WI1,
699 $ VL, LDVL, VR, LDVR, LRE, LDLRE, RCONDV,
700 $ RCNDV1, RCDVIN, RCONDE, RCNDE1, RCDEIN,
701 $ SCALE, SCALE1, RESULT, WORK, NNWORK,
702 $ IWORK, INFO )
703 *
704 * Check for RESULT(j) > THRESH
705 *
706 NTEST = 0
707 NFAIL = 0
708 DO 100 J = 1, 9
709 IF( RESULT( J ).GE.ZERO )
710 $ NTEST = NTEST + 1
711 IF( RESULT( J ).GE.THRESH )
712 $ NFAIL = NFAIL + 1
713 100 CONTINUE
714 *
715 IF( NFAIL.GT.0 )
716 $ NTESTF = NTESTF + 1
717 IF( NTESTF.EQ.1 ) THEN
718 WRITE( NOUNIT, FMT = 9999 )PATH
719 WRITE( NOUNIT, FMT = 9998 )
720 WRITE( NOUNIT, FMT = 9997 )
721 WRITE( NOUNIT, FMT = 9996 )
722 WRITE( NOUNIT, FMT = 9995 )THRESH
723 NTESTF = 2
724 END IF
725 *
726 DO 110 J = 1, 9
727 IF( RESULT( J ).GE.THRESH ) THEN
728 WRITE( NOUNIT, FMT = 9994 )BALANC, N, IWK,
729 $ IOLDSD, JTYPE, J, RESULT( J )
730 END IF
731 110 CONTINUE
732 *
733 NERRS = NERRS + NFAIL
734 NTESTT = NTESTT + NTEST
735 *
736 120 CONTINUE
737 130 CONTINUE
738 140 CONTINUE
739 150 CONTINUE
740 *
741 160 CONTINUE
742 *
743 * Read in data from file to check accuracy of condition estimation.
744 * Assume input eigenvalues are sorted lexicographically (increasing
745 * by real part, then decreasing by imaginary part)
746 *
747 JTYPE = 0
748 170 CONTINUE
749 READ( NIUNIT, FMT = *, END = 220 )N
750 *
751 * Read input data until N=0
752 *
753 IF( N.EQ.0 )
754 $ GO TO 220
755 JTYPE = JTYPE + 1
756 ISEED( 1 ) = JTYPE
757 DO 180 I = 1, N
758 READ( NIUNIT, FMT = * )( A( I, J ), J = 1, N )
759 180 CONTINUE
760 DO 190 I = 1, N
761 READ( NIUNIT, FMT = * )WR1( I ), WI1( I ), RCDEIN( I ),
762 $ RCDVIN( I )
763 190 CONTINUE
764 CALL DGET23( .TRUE., 'N', 22, THRESH, ISEED, NOUNIT, N, A, LDA, H,
765 $ WR, WI, WR1, WI1, VL, LDVL, VR, LDVR, LRE, LDLRE,
766 $ RCONDV, RCNDV1, RCDVIN, RCONDE, RCNDE1, RCDEIN,
767 $ SCALE, SCALE1, RESULT, WORK, 6*N+2*N**2, IWORK,
768 $ INFO )
769 *
770 * Check for RESULT(j) > THRESH
771 *
772 NTEST = 0
773 NFAIL = 0
774 DO 200 J = 1, 11
775 IF( RESULT( J ).GE.ZERO )
776 $ NTEST = NTEST + 1
777 IF( RESULT( J ).GE.THRESH )
778 $ NFAIL = NFAIL + 1
779 200 CONTINUE
780 *
781 IF( NFAIL.GT.0 )
782 $ NTESTF = NTESTF + 1
783 IF( NTESTF.EQ.1 ) THEN
784 WRITE( NOUNIT, FMT = 9999 )PATH
785 WRITE( NOUNIT, FMT = 9998 )
786 WRITE( NOUNIT, FMT = 9997 )
787 WRITE( NOUNIT, FMT = 9996 )
788 WRITE( NOUNIT, FMT = 9995 )THRESH
789 NTESTF = 2
790 END IF
791 *
792 DO 210 J = 1, 11
793 IF( RESULT( J ).GE.THRESH ) THEN
794 WRITE( NOUNIT, FMT = 9993 )N, JTYPE, J, RESULT( J )
795 END IF
796 210 CONTINUE
797 *
798 NERRS = NERRS + NFAIL
799 NTESTT = NTESTT + NTEST
800 GO TO 170
801 220 CONTINUE
802 *
803 * Summary
804 *
805 CALL DLASUM( PATH, NOUNIT, NERRS, NTESTT )
806 *
807 9999 FORMAT( / 1X, A3, ' -- Real Eigenvalue-Eigenvector Decomposition',
808 $ ' Expert Driver', /
809 $ ' Matrix types (see DDRVVX for details): ' )
810 *
811 9998 FORMAT( / ' Special Matrices:', / ' 1=Zero matrix. ',
812 $ ' ', ' 5=Diagonal: geometr. spaced entries.',
813 $ / ' 2=Identity matrix. ', ' 6=Diagona',
814 $ 'l: clustered entries.', / ' 3=Transposed Jordan block. ',
815 $ ' ', ' 7=Diagonal: large, evenly spaced.', / ' ',
816 $ '4=Diagonal: evenly spaced entries. ', ' 8=Diagonal: s',
817 $ 'mall, evenly spaced.' )
818 9997 FORMAT( ' Dense, Non-Symmetric Matrices:', / ' 9=Well-cond., ev',
819 $ 'enly spaced eigenvals.', ' 14=Ill-cond., geomet. spaced e',
820 $ 'igenals.', / ' 10=Well-cond., geom. spaced eigenvals. ',
821 $ ' 15=Ill-conditioned, clustered e.vals.', / ' 11=Well-cond',
822 $ 'itioned, clustered e.vals. ', ' 16=Ill-cond., random comp',
823 $ 'lex ', / ' 12=Well-cond., random complex ', ' ',
824 $ ' 17=Ill-cond., large rand. complx ', / ' 13=Ill-condi',
825 $ 'tioned, evenly spaced. ', ' 18=Ill-cond., small rand.',
826 $ ' complx ' )
827 9996 FORMAT( ' 19=Matrix with random O(1) entries. ', ' 21=Matrix ',
828 $ 'with small random entries.', / ' 20=Matrix with large ran',
829 $ 'dom entries. ', ' 22=Matrix read from input file', / )
830 9995 FORMAT( ' Tests performed with test threshold =', F8.2,
831 $ / / ' 1 = | A VR - VR W | / ( n |A| ulp ) ',
832 $ / ' 2 = | transpose(A) VL - VL W | / ( n |A| ulp ) ',
833 $ / ' 3 = | |VR(i)| - 1 | / ulp ',
834 $ / ' 4 = | |VL(i)| - 1 | / ulp ',
835 $ / ' 5 = 0 if W same no matter if VR or VL computed,',
836 $ ' 1/ulp otherwise', /
837 $ ' 6 = 0 if VR same no matter what else computed,',
838 $ ' 1/ulp otherwise', /
839 $ ' 7 = 0 if VL same no matter what else computed,',
840 $ ' 1/ulp otherwise', /
841 $ ' 8 = 0 if RCONDV same no matter what else computed,',
842 $ ' 1/ulp otherwise', /
843 $ ' 9 = 0 if SCALE, ILO, IHI, ABNRM same no matter what else',
844 $ ' computed, 1/ulp otherwise',
845 $ / ' 10 = | RCONDV - RCONDV(precomputed) | / cond(RCONDV),',
846 $ / ' 11 = | RCONDE - RCONDE(precomputed) | / cond(RCONDE),' )
847 9994 FORMAT( ' BALANC=''', A1, ''',N=', I4, ',IWK=', I1, ', seed=',
848 $ 4( I4, ',' ), ' type ', I2, ', test(', I2, ')=', G10.3 )
849 9993 FORMAT( ' N=', I5, ', input example =', I3, ', test(', I2, ')=',
850 $ G10.3 )
851 9992 FORMAT( ' DDRVVX: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
852 $ I6, ', JTYPE=', I6, ', ISEED=(', 3( I5, ',' ), I5, ')' )
853 *
854 RETURN
855 *
856 * End of DDRVVX
857 *
858 END