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