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