LAPACK  3.10.1
LAPACK: Linear Algebra PACKage
dchkhs.f
Go to the documentation of this file.
1 *> \brief \b DCHKHS
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 DCHKHS( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
12 * NOUNIT, A, LDA, H, T1, T2, U, LDU, Z, UZ, WR1,
13 * WI1, WR2, WI2, WR3, WI3, EVECTL, EVECTR, EVECTY,
14 * EVECTX, UU, TAU, WORK, NWORK, IWORK, SELECT,
15 * RESULT, INFO )
16 *
17 * .. Scalar Arguments ..
18 * INTEGER INFO, LDA, LDU, NOUNIT, NSIZES, NTYPES, NWORK
19 * DOUBLE PRECISION THRESH
20 * ..
21 * .. Array Arguments ..
22 * LOGICAL DOTYPE( * ), SELECT( * )
23 * INTEGER ISEED( 4 ), IWORK( * ), NN( * )
24 * DOUBLE PRECISION A( LDA, * ), EVECTL( LDU, * ),
25 * $ EVECTR( LDU, * ), EVECTX( LDU, * ),
26 * $ EVECTY( LDU, * ), H( LDA, * ), RESULT( 14 ),
27 * $ T1( LDA, * ), T2( LDA, * ), TAU( * ),
28 * $ U( LDU, * ), UU( LDU, * ), UZ( LDU, * ),
29 * $ WI1( * ), WI2( * ), WI3( * ), WORK( * ),
30 * $ WR1( * ), WR2( * ), WR3( * ), Z( LDU, * )
31 * ..
32 *
33 *
34 *> \par Purpose:
35 * =============
36 *>
37 *> \verbatim
38 *>
39 *> DCHKHS checks the nonsymmetric eigenvalue problem routines.
40 *>
41 *> DGEHRD factors A as U H U' , where ' means transpose,
42 *> H is hessenberg, and U is an orthogonal matrix.
43 *>
44 *> DORGHR generates the orthogonal matrix U.
45 *>
46 *> DORMHR multiplies a matrix by the orthogonal matrix U.
47 *>
48 *> DHSEQR factors H as Z T Z' , where Z is orthogonal and
49 *> T is "quasi-triangular", and the eigenvalue vector W.
50 *>
51 *> DTREVC computes the left and right eigenvector matrices
52 *> L and R for T.
53 *>
54 *> DHSEIN computes the left and right eigenvector matrices
55 *> Y and X for H, using inverse iteration.
56 *>
57 *> When DCHKHS is called, a number of matrix "sizes" ("n's") and a
58 *> number of matrix "types" are specified. For each size ("n")
59 *> and each type of matrix, one matrix will be generated and used
60 *> to test the nonsymmetric eigenroutines. For each matrix, 14
61 *> tests will be performed:
62 *>
63 *> (1) | A - U H U**T | / ( |A| n ulp )
64 *>
65 *> (2) | I - UU**T | / ( n ulp )
66 *>
67 *> (3) | H - Z T Z**T | / ( |H| n ulp )
68 *>
69 *> (4) | I - ZZ**T | / ( n ulp )
70 *>
71 *> (5) | A - UZ H (UZ)**T | / ( |A| n ulp )
72 *>
73 *> (6) | I - UZ (UZ)**T | / ( n ulp )
74 *>
75 *> (7) | T(Z computed) - T(Z not computed) | / ( |T| ulp )
76 *>
77 *> (8) | W(Z computed) - W(Z not computed) | / ( |W| ulp )
78 *>
79 *> (9) | TR - RW | / ( |T| |R| ulp )
80 *>
81 *> (10) | L**H T - W**H L | / ( |T| |L| ulp )
82 *>
83 *> (11) | HX - XW | / ( |H| |X| ulp )
84 *>
85 *> (12) | Y**H H - W**H Y | / ( |H| |Y| ulp )
86 *>
87 *> (13) | AX - XW | / ( |A| |X| ulp )
88 *>
89 *> (14) | Y**H A - W**H Y | / ( |A| |Y| ulp )
90 *>
91 *> The "sizes" are specified by an array NN(1:NSIZES); the value of
92 *> each element NN(j) specifies one size.
93 *> The "types" are specified by a logical array DOTYPE( 1:NTYPES );
94 *> if DOTYPE(j) is .TRUE., then matrix type "j" will be generated.
95 *> Currently, the list of possible types is:
96 *>
97 *> (1) The zero matrix.
98 *> (2) The identity matrix.
99 *> (3) A (transposed) Jordan block, with 1's on the diagonal.
100 *>
101 *> (4) A diagonal matrix with evenly spaced entries
102 *> 1, ..., ULP and random signs.
103 *> (ULP = (first number larger than 1) - 1 )
104 *> (5) A diagonal matrix with geometrically spaced entries
105 *> 1, ..., ULP and random signs.
106 *> (6) A diagonal matrix with "clustered" entries 1, ULP, ..., ULP
107 *> and random signs.
108 *>
109 *> (7) Same as (4), but multiplied by SQRT( overflow threshold )
110 *> (8) Same as (4), but multiplied by SQRT( 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 SQRT( overflow threshold )
153 *> (18) Same as (16), but multiplied by SQRT( underflow threshold )
154 *>
155 *> (19) Nonsymmetric matrix with random entries chosen from (-1,1).
156 *> (20) Same as (19), but multiplied by SQRT( overflow threshold )
157 *> (21) Same as (19), but multiplied by SQRT( underflow threshold )
158 *> \endverbatim
159 *
160 * Arguments:
161 * ==========
162 *
163 *> \verbatim
164 *> NSIZES - INTEGER
165 *> The number of sizes of matrices to use. If it is zero,
166 *> DCHKHS does nothing. It must be at least zero.
167 *> Not modified.
168 *>
169 *> NN - INTEGER array, dimension (NSIZES)
170 *> An array containing the sizes to be used for the matrices.
171 *> Zero values will be skipped. The values must be at least
172 *> zero.
173 *> Not modified.
174 *>
175 *> NTYPES - INTEGER
176 *> The number of elements in DOTYPE. If it is zero, DCHKHS
177 *> does nothing. It must be at least zero. If it is MAXTYP+1
178 *> and NSIZES is 1, then an additional type, MAXTYP+1 is
179 *> defined, which is to use whatever matrix is in A. This
180 *> is only useful if DOTYPE(1:MAXTYP) is .FALSE. and
181 *> DOTYPE(MAXTYP+1) is .TRUE. .
182 *> Not modified.
183 *>
184 *> DOTYPE - LOGICAL array, dimension (NTYPES)
185 *> If DOTYPE(j) is .TRUE., then for each size in NN a
186 *> matrix of that size and of type j will be generated.
187 *> If NTYPES is smaller than the maximum number of types
188 *> defined (PARAMETER MAXTYP), then types NTYPES+1 through
189 *> MAXTYP will not be generated. If NTYPES is larger
190 *> than MAXTYP, DOTYPE(MAXTYP+1) through DOTYPE(NTYPES)
191 *> will be ignored.
192 *> Not modified.
193 *>
194 *> ISEED - INTEGER array, dimension (4)
195 *> On entry ISEED specifies the seed of the random number
196 *> generator. The array elements should be between 0 and 4095;
197 *> if not they will be reduced mod 4096. Also, ISEED(4) must
198 *> be odd. The random number generator uses a linear
199 *> congruential sequence limited to small integers, and so
200 *> should produce machine independent random numbers. The
201 *> values of ISEED are changed on exit, and can be used in the
202 *> next call to DCHKHS to continue the same random number
203 *> sequence.
204 *> Modified.
205 *>
206 *> THRESH - DOUBLE PRECISION
207 *> A test will count as "failed" if the "error", computed as
208 *> described above, exceeds THRESH. Note that the error
209 *> is scaled to be O(1), so THRESH should be a reasonably
210 *> small multiple of 1, e.g., 10 or 100. In particular,
211 *> it should not depend on the precision (single vs. double)
212 *> or the size of the matrix. It must be at least zero.
213 *> Not modified.
214 *>
215 *> NOUNIT - INTEGER
216 *> The FORTRAN unit number for printing out error messages
217 *> (e.g., if a routine returns IINFO not equal to 0.)
218 *> Not modified.
219 *>
220 *> A - DOUBLE PRECISION array, dimension (LDA,max(NN))
221 *> Used to hold the matrix whose eigenvalues are to be
222 *> computed. On exit, A contains the last matrix actually
223 *> used.
224 *> Modified.
225 *>
226 *> LDA - INTEGER
227 *> The leading dimension of A, H, T1 and T2. It must be at
228 *> least 1 and at least max( NN ).
229 *> Not modified.
230 *>
231 *> H - DOUBLE PRECISION array, dimension (LDA,max(NN))
232 *> The upper hessenberg matrix computed by DGEHRD. On exit,
233 *> H contains the Hessenberg form of the matrix in A.
234 *> Modified.
235 *>
236 *> T1 - DOUBLE PRECISION array, dimension (LDA,max(NN))
237 *> The Schur (="quasi-triangular") matrix computed by DHSEQR
238 *> if Z is computed. On exit, T1 contains the Schur form of
239 *> the matrix in A.
240 *> Modified.
241 *>
242 *> T2 - DOUBLE PRECISION array, dimension (LDA,max(NN))
243 *> The Schur matrix computed by DHSEQR when Z is not computed.
244 *> This should be identical to T1.
245 *> Modified.
246 *>
247 *> LDU - INTEGER
248 *> The leading dimension of U, Z, UZ and UU. It must be at
249 *> least 1 and at least max( NN ).
250 *> Not modified.
251 *>
252 *> U - DOUBLE PRECISION array, dimension (LDU,max(NN))
253 *> The orthogonal matrix computed by DGEHRD.
254 *> Modified.
255 *>
256 *> Z - DOUBLE PRECISION array, dimension (LDU,max(NN))
257 *> The orthogonal matrix computed by DHSEQR.
258 *> Modified.
259 *>
260 *> UZ - DOUBLE PRECISION array, dimension (LDU,max(NN))
261 *> The product of U times Z.
262 *> Modified.
263 *>
264 *> WR1 - DOUBLE PRECISION array, dimension (max(NN))
265 *> WI1 - DOUBLE PRECISION array, dimension (max(NN))
266 *> The real and imaginary parts of the eigenvalues of A,
267 *> as computed when Z is computed.
268 *> On exit, WR1 + WI1*i are the eigenvalues of the matrix in A.
269 *> Modified.
270 *>
271 *> WR2 - DOUBLE PRECISION array, dimension (max(NN))
272 *> WI2 - DOUBLE PRECISION array, dimension (max(NN))
273 *> The real and imaginary parts of the eigenvalues of A,
274 *> as computed when T is computed but not Z.
275 *> On exit, WR2 + WI2*i are the eigenvalues of the matrix in A.
276 *> Modified.
277 *>
278 *> WR3 - DOUBLE PRECISION array, dimension (max(NN))
279 *> WI3 - DOUBLE PRECISION array, dimension (max(NN))
280 *> Like WR1, WI1, these arrays contain the eigenvalues of A,
281 *> but those computed when DHSEQR only computes the
282 *> eigenvalues, i.e., not the Schur vectors and no more of the
283 *> Schur form than is necessary for computing the
284 *> eigenvalues.
285 *> Modified.
286 *>
287 *> EVECTL - DOUBLE PRECISION array, dimension (LDU,max(NN))
288 *> The (upper triangular) left eigenvector matrix for the
289 *> matrix in T1. For complex conjugate pairs, the real part
290 *> is stored in one row and the imaginary part in the next.
291 *> Modified.
292 *>
293 *> EVEZTR - DOUBLE PRECISION array, dimension (LDU,max(NN))
294 *> The (upper triangular) right eigenvector matrix for the
295 *> matrix in T1. For complex conjugate pairs, the real part
296 *> is stored in one column and the imaginary part in the next.
297 *> Modified.
298 *>
299 *> EVECTY - DOUBLE PRECISION array, dimension (LDU,max(NN))
300 *> The left eigenvector matrix for the
301 *> matrix in H. For complex conjugate pairs, the real part
302 *> is stored in one row and the imaginary part in the next.
303 *> Modified.
304 *>
305 *> EVECTX - DOUBLE PRECISION array, dimension (LDU,max(NN))
306 *> The right eigenvector matrix for the
307 *> matrix in H. For complex conjugate pairs, the real part
308 *> is stored in one column and the imaginary part in the next.
309 *> Modified.
310 *>
311 *> UU - DOUBLE PRECISION array, dimension (LDU,max(NN))
312 *> Details of the orthogonal matrix computed by DGEHRD.
313 *> Modified.
314 *>
315 *> TAU - DOUBLE PRECISION array, dimension(max(NN))
316 *> Further details of the orthogonal matrix computed by DGEHRD.
317 *> Modified.
318 *>
319 *> WORK - DOUBLE PRECISION array, dimension (NWORK)
320 *> Workspace.
321 *> Modified.
322 *>
323 *> NWORK - INTEGER
324 *> The number of entries in WORK. NWORK >= 4*NN(j)*NN(j) + 2.
325 *>
326 *> IWORK - INTEGER array, dimension (max(NN))
327 *> Workspace.
328 *> Modified.
329 *>
330 *> SELECT - LOGICAL array, dimension (max(NN))
331 *> Workspace.
332 *> Modified.
333 *>
334 *> RESULT - DOUBLE PRECISION array, dimension (14)
335 *> The values computed by the fourteen tests described above.
336 *> The values are currently limited to 1/ulp, to avoid
337 *> overflow.
338 *> Modified.
339 *>
340 *> INFO - INTEGER
341 *> If 0, then everything ran OK.
342 *> -1: NSIZES < 0
343 *> -2: Some NN(j) < 0
344 *> -3: NTYPES < 0
345 *> -6: THRESH < 0
346 *> -9: LDA < 1 or LDA < NMAX, where NMAX is max( NN(j) ).
347 *> -14: LDU < 1 or LDU < NMAX.
348 *> -28: NWORK too small.
349 *> If DLATMR, SLATMS, or SLATME returns an error code, the
350 *> absolute value of it is returned.
351 *> If 1, then DHSEQR could not find all the shifts.
352 *> If 2, then the EISPACK code (for small blocks) failed.
353 *> If >2, then 30*N iterations were not enough to find an
354 *> eigenvalue or to decompose the problem.
355 *> Modified.
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 *> MTEST The number of tests defined: care must be taken
365 *> that (1) the size of RESULT, (2) the number of
366 *> tests actually performed, and (3) MTEST agree.
367 *> NTEST The number of tests performed on this matrix
368 *> so far. This should be less than MTEST, and
369 *> equal to it by the last test. It will be less
370 *> if any of the routines being tested indicates
371 *> that it could not compute the matrices that
372 *> would be tested.
373 *> NMAX Largest value in NN.
374 *> NMATS The number of matrices generated so far.
375 *> NERRS The number of tests which have exceeded THRESH
376 *> so far (computed by DLAFTS).
377 *> COND, CONDS,
378 *> IMODE Values to be passed to the matrix generators.
379 *> ANORM Norm of A; passed to matrix generators.
380 *>
381 *> OVFL, UNFL Overflow and underflow thresholds.
382 *> ULP, ULPINV Finest relative precision and its inverse.
383 *> RTOVFL, RTUNFL,
384 *> RTULP, RTULPI Square roots of the previous 4 values.
385 *>
386 *> The following four arrays decode JTYPE:
387 *> KTYPE(j) The general type (1-10) for type "j".
388 *> KMODE(j) The MODE value to be passed to the matrix
389 *> generator for type "j".
390 *> KMAGN(j) The order of magnitude ( O(1),
391 *> O(overflow^(1/2) ), O(underflow^(1/2) )
392 *> KCONDS(j) Selects whether CONDS is to be 1 or
393 *> 1/sqrt(ulp). (0 means irrelevant.)
394 *> \endverbatim
395 *
396 * Authors:
397 * ========
398 *
399 *> \author Univ. of Tennessee
400 *> \author Univ. of California Berkeley
401 *> \author Univ. of Colorado Denver
402 *> \author NAG Ltd.
403 *
404 *> \ingroup double_eig
405 *
406 * =====================================================================
407  SUBROUTINE dchkhs( NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH,
408  $ NOUNIT, A, LDA, H, T1, T2, U, LDU, Z, UZ, WR1,
409  $ WI1, WR2, WI2, WR3, WI3, EVECTL, EVECTR,
410  $ EVECTY, EVECTX, UU, TAU, WORK, NWORK, IWORK,
411  $ SELECT, RESULT, INFO )
412 *
413 * -- LAPACK test routine --
414 * -- LAPACK is a software package provided by Univ. of Tennessee, --
415 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
416 *
417 * .. Scalar Arguments ..
418  INTEGER INFO, LDA, LDU, NOUNIT, NSIZES, NTYPES, NWORK
419  DOUBLE PRECISION THRESH
420 * ..
421 * .. Array Arguments ..
422  LOGICAL DOTYPE( * ), SELECT( * )
423  INTEGER ISEED( 4 ), IWORK( * ), NN( * )
424  DOUBLE PRECISION A( LDA, * ), EVECTL( LDU, * ),
425  $ EVECTR( LDU, * ), EVECTX( LDU, * ),
426  $ evecty( ldu, * ), h( lda, * ), result( 14 ),
427  $ t1( lda, * ), t2( lda, * ), tau( * ),
428  $ u( ldu, * ), uu( ldu, * ), uz( ldu, * ),
429  $ wi1( * ), wi2( * ), wi3( * ), work( * ),
430  $ wr1( * ), wr2( * ), wr3( * ), z( ldu, * )
431 * ..
432 *
433 * =====================================================================
434 *
435 * .. Parameters ..
436  DOUBLE PRECISION ZERO, ONE
437  PARAMETER ( ZERO = 0.0d0, one = 1.0d0 )
438  INTEGER MAXTYP
439  PARAMETER ( MAXTYP = 21 )
440 * ..
441 * .. Local Scalars ..
442  LOGICAL BADNN, MATCH
443  INTEGER I, IHI, IINFO, ILO, IMODE, IN, ITYPE, J, JCOL,
444  $ JJ, JSIZE, JTYPE, K, MTYPES, N, N1, NERRS,
445  $ NMATS, NMAX, NSELC, NSELR, NTEST, NTESTT
446  DOUBLE PRECISION ANINV, ANORM, COND, CONDS, OVFL, RTOVFL, RTULP,
447  $ rtulpi, rtunfl, temp1, temp2, ulp, ulpinv, unfl
448 * ..
449 * .. Local Arrays ..
450  CHARACTER ADUMMA( 1 )
451  INTEGER IDUMMA( 1 ), IOLDSD( 4 ), KCONDS( MAXTYP ),
452  $ KMAGN( MAXTYP ), KMODE( MAXTYP ),
453  $ KTYPE( MAXTYP )
454  DOUBLE PRECISION DUMMA( 6 )
455 * ..
456 * .. External Functions ..
457  DOUBLE PRECISION DLAMCH
458  EXTERNAL dlamch
459 * ..
460 * .. External Subroutines ..
461  EXTERNAL dcopy, dgehrd, dgemm, dget10, dget22, dhsein,
464  $ dtrevc, xerbla
465 * ..
466 * .. Intrinsic Functions ..
467  INTRINSIC abs, dble, max, min, sqrt
468 * ..
469 * .. Data statements ..
470  DATA ktype / 1, 2, 3, 5*4, 4*6, 6*6, 3*9 /
471  DATA kmagn / 3*1, 1, 1, 1, 2, 3, 4*1, 1, 1, 1, 1, 2,
472  $ 3, 1, 2, 3 /
473  DATA kmode / 3*0, 4, 3, 1, 4, 4, 4, 3, 1, 5, 4, 3,
474  $ 1, 5, 5, 5, 4, 3, 1 /
475  DATA kconds / 3*0, 5*0, 4*1, 6*2, 3*0 /
476 * ..
477 * .. Executable Statements ..
478 *
479 * Check for errors
480 *
481  ntestt = 0
482  info = 0
483 *
484  badnn = .false.
485  nmax = 0
486  DO 10 j = 1, nsizes
487  nmax = max( nmax, nn( j ) )
488  IF( nn( j ).LT.0 )
489  $ badnn = .true.
490  10 CONTINUE
491 *
492 * Check for errors
493 *
494  IF( nsizes.LT.0 ) THEN
495  info = -1
496  ELSE IF( badnn ) THEN
497  info = -2
498  ELSE IF( ntypes.LT.0 ) THEN
499  info = -3
500  ELSE IF( thresh.LT.zero ) THEN
501  info = -6
502  ELSE IF( lda.LE.1 .OR. lda.LT.nmax ) THEN
503  info = -9
504  ELSE IF( ldu.LE.1 .OR. ldu.LT.nmax ) THEN
505  info = -14
506  ELSE IF( 4*nmax*nmax+2.GT.nwork ) THEN
507  info = -28
508  END IF
509 *
510  IF( info.NE.0 ) THEN
511  CALL xerbla( 'DCHKHS', -info )
512  RETURN
513  END IF
514 *
515 * Quick return if possible
516 *
517  IF( nsizes.EQ.0 .OR. ntypes.EQ.0 )
518  $ RETURN
519 *
520 * More important constants
521 *
522  unfl = dlamch( 'Safe minimum' )
523  ovfl = dlamch( 'Overflow' )
524  CALL dlabad( unfl, ovfl )
525  ulp = dlamch( 'Epsilon' )*dlamch( 'Base' )
526  ulpinv = one / ulp
527  rtunfl = sqrt( unfl )
528  rtovfl = sqrt( ovfl )
529  rtulp = sqrt( ulp )
530  rtulpi = one / rtulp
531 *
532 * Loop over sizes, types
533 *
534  nerrs = 0
535  nmats = 0
536 *
537  DO 270 jsize = 1, nsizes
538  n = nn( jsize )
539  IF( n.EQ.0 )
540  $ GO TO 270
541  n1 = max( 1, n )
542  aninv = one / dble( n1 )
543 *
544  IF( nsizes.NE.1 ) THEN
545  mtypes = min( maxtyp, ntypes )
546  ELSE
547  mtypes = min( maxtyp+1, ntypes )
548  END IF
549 *
550  DO 260 jtype = 1, mtypes
551  IF( .NOT.dotype( jtype ) )
552  $ GO TO 260
553  nmats = nmats + 1
554  ntest = 0
555 *
556 * Save ISEED in case of an error.
557 *
558  DO 20 j = 1, 4
559  ioldsd( j ) = iseed( j )
560  20 CONTINUE
561 *
562 * Initialize RESULT
563 *
564  DO 30 j = 1, 14
565  result( j ) = zero
566  30 CONTINUE
567 *
568 * Compute "A"
569 *
570 * Control parameters:
571 *
572 * KMAGN KCONDS KMODE KTYPE
573 * =1 O(1) 1 clustered 1 zero
574 * =2 large large clustered 2 identity
575 * =3 small exponential Jordan
576 * =4 arithmetic diagonal, (w/ eigenvalues)
577 * =5 random log symmetric, w/ eigenvalues
578 * =6 random general, w/ eigenvalues
579 * =7 random diagonal
580 * =8 random symmetric
581 * =9 random general
582 * =10 random triangular
583 *
584  IF( mtypes.GT.maxtyp )
585  $ GO TO 100
586 *
587  itype = ktype( jtype )
588  imode = kmode( jtype )
589 *
590 * Compute norm
591 *
592  GO TO ( 40, 50, 60 )kmagn( jtype )
593 *
594  40 CONTINUE
595  anorm = one
596  GO TO 70
597 *
598  50 CONTINUE
599  anorm = ( rtovfl*ulp )*aninv
600  GO TO 70
601 *
602  60 CONTINUE
603  anorm = rtunfl*n*ulpinv
604  GO TO 70
605 *
606  70 CONTINUE
607 *
608  CALL dlaset( 'Full', lda, n, zero, zero, a, lda )
609  iinfo = 0
610  cond = ulpinv
611 *
612 * Special Matrices
613 *
614  IF( itype.EQ.1 ) THEN
615 *
616 * Zero
617 *
618  iinfo = 0
619 *
620  ELSE IF( itype.EQ.2 ) THEN
621 *
622 * Identity
623 *
624  DO 80 jcol = 1, n
625  a( jcol, jcol ) = anorm
626  80 CONTINUE
627 *
628  ELSE IF( itype.EQ.3 ) THEN
629 *
630 * Jordan Block
631 *
632  DO 90 jcol = 1, n
633  a( jcol, jcol ) = anorm
634  IF( jcol.GT.1 )
635  $ a( jcol, jcol-1 ) = one
636  90 CONTINUE
637 *
638  ELSE IF( itype.EQ.4 ) THEN
639 *
640 * Diagonal Matrix, [Eigen]values Specified
641 *
642  CALL dlatms( n, n, 'S', iseed, 'S', work, imode, cond,
643  $ anorm, 0, 0, 'N', a, lda, work( n+1 ),
644  $ iinfo )
645 *
646  ELSE IF( itype.EQ.5 ) THEN
647 *
648 * Symmetric, eigenvalues specified
649 *
650  CALL dlatms( n, n, 'S', iseed, 'S', work, imode, cond,
651  $ anorm, n, n, 'N', a, lda, work( n+1 ),
652  $ iinfo )
653 *
654  ELSE IF( itype.EQ.6 ) THEN
655 *
656 * General, eigenvalues specified
657 *
658  IF( kconds( jtype ).EQ.1 ) THEN
659  conds = one
660  ELSE IF( kconds( jtype ).EQ.2 ) THEN
661  conds = rtulpi
662  ELSE
663  conds = zero
664  END IF
665 *
666  adumma( 1 ) = ' '
667  CALL dlatme( n, 'S', iseed, work, imode, cond, one,
668  $ adumma, 'T', 'T', 'T', work( n+1 ), 4,
669  $ conds, n, n, anorm, a, lda, work( 2*n+1 ),
670  $ iinfo )
671 *
672  ELSE IF( itype.EQ.7 ) THEN
673 *
674 * Diagonal, random eigenvalues
675 *
676  CALL dlatmr( n, n, 'S', iseed, 'S', work, 6, one, one,
677  $ 'T', 'N', work( n+1 ), 1, one,
678  $ work( 2*n+1 ), 1, one, 'N', idumma, 0, 0,
679  $ zero, anorm, 'NO', a, lda, iwork, iinfo )
680 *
681  ELSE IF( itype.EQ.8 ) THEN
682 *
683 * Symmetric, random eigenvalues
684 *
685  CALL dlatmr( n, n, 'S', iseed, 'S', work, 6, one, one,
686  $ 'T', 'N', work( n+1 ), 1, one,
687  $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
688  $ zero, anorm, 'NO', a, lda, iwork, iinfo )
689 *
690  ELSE IF( itype.EQ.9 ) THEN
691 *
692 * General, random eigenvalues
693 *
694  CALL dlatmr( n, n, 'S', iseed, 'N', work, 6, one, one,
695  $ 'T', 'N', work( n+1 ), 1, one,
696  $ work( 2*n+1 ), 1, one, 'N', idumma, n, n,
697  $ zero, anorm, 'NO', a, lda, iwork, iinfo )
698 *
699  ELSE IF( itype.EQ.10 ) THEN
700 *
701 * Triangular, random eigenvalues
702 *
703  CALL dlatmr( n, n, 'S', iseed, 'N', work, 6, one, one,
704  $ 'T', 'N', work( n+1 ), 1, one,
705  $ work( 2*n+1 ), 1, one, 'N', idumma, n, 0,
706  $ zero, anorm, 'NO', a, lda, iwork, iinfo )
707 *
708  ELSE
709 *
710  iinfo = 1
711  END IF
712 *
713  IF( iinfo.NE.0 ) THEN
714  WRITE( nounit, fmt = 9999 )'Generator', iinfo, n, jtype,
715  $ ioldsd
716  info = abs( iinfo )
717  RETURN
718  END IF
719 *
720  100 CONTINUE
721 *
722 * Call DGEHRD to compute H and U, do tests.
723 *
724  CALL dlacpy( ' ', n, n, a, lda, h, lda )
725 *
726  ntest = 1
727 *
728  ilo = 1
729  ihi = n
730 *
731  CALL dgehrd( n, ilo, ihi, h, lda, work, work( n+1 ),
732  $ nwork-n, iinfo )
733 *
734  IF( iinfo.NE.0 ) THEN
735  result( 1 ) = ulpinv
736  WRITE( nounit, fmt = 9999 )'DGEHRD', iinfo, n, jtype,
737  $ ioldsd
738  info = abs( iinfo )
739  GO TO 250
740  END IF
741 *
742  DO 120 j = 1, n - 1
743  uu( j+1, j ) = zero
744  DO 110 i = j + 2, n
745  u( i, j ) = h( i, j )
746  uu( i, j ) = h( i, j )
747  h( i, j ) = zero
748  110 CONTINUE
749  120 CONTINUE
750  CALL dcopy( n-1, work, 1, tau, 1 )
751  CALL dorghr( n, ilo, ihi, u, ldu, work, work( n+1 ),
752  $ nwork-n, iinfo )
753  ntest = 2
754 *
755  CALL dhst01( n, ilo, ihi, a, lda, h, lda, u, ldu, work,
756  $ nwork, result( 1 ) )
757 *
758 * Call DHSEQR to compute T1, T2 and Z, do tests.
759 *
760 * Eigenvalues only (WR3,WI3)
761 *
762  CALL dlacpy( ' ', n, n, h, lda, t2, lda )
763  ntest = 3
764  result( 3 ) = ulpinv
765 *
766  CALL dhseqr( 'E', 'N', n, ilo, ihi, t2, lda, wr3, wi3, uz,
767  $ ldu, work, nwork, iinfo )
768  IF( iinfo.NE.0 ) THEN
769  WRITE( nounit, fmt = 9999 )'DHSEQR(E)', iinfo, n, jtype,
770  $ ioldsd
771  IF( iinfo.LE.n+2 ) THEN
772  info = abs( iinfo )
773  GO TO 250
774  END IF
775  END IF
776 *
777 * Eigenvalues (WR2,WI2) and Full Schur Form (T2)
778 *
779  CALL dlacpy( ' ', n, n, h, lda, t2, lda )
780 *
781  CALL dhseqr( 'S', 'N', n, ilo, ihi, t2, lda, wr2, wi2, uz,
782  $ ldu, work, nwork, iinfo )
783  IF( iinfo.NE.0 .AND. iinfo.LE.n+2 ) THEN
784  WRITE( nounit, fmt = 9999 )'DHSEQR(S)', iinfo, n, jtype,
785  $ ioldsd
786  info = abs( iinfo )
787  GO TO 250
788  END IF
789 *
790 * Eigenvalues (WR1,WI1), Schur Form (T1), and Schur vectors
791 * (UZ)
792 *
793  CALL dlacpy( ' ', n, n, h, lda, t1, lda )
794  CALL dlacpy( ' ', n, n, u, ldu, uz, ldu )
795 *
796  CALL dhseqr( 'S', 'V', n, ilo, ihi, t1, lda, wr1, wi1, uz,
797  $ ldu, work, nwork, iinfo )
798  IF( iinfo.NE.0 .AND. iinfo.LE.n+2 ) THEN
799  WRITE( nounit, fmt = 9999 )'DHSEQR(V)', iinfo, n, jtype,
800  $ ioldsd
801  info = abs( iinfo )
802  GO TO 250
803  END IF
804 *
805 * Compute Z = U' UZ
806 *
807  CALL dgemm( 'T', 'N', n, n, n, one, u, ldu, uz, ldu, zero,
808  $ z, ldu )
809  ntest = 8
810 *
811 * Do Tests 3: | H - Z T Z' | / ( |H| n ulp )
812 * and 4: | I - Z Z' | / ( n ulp )
813 *
814  CALL dhst01( n, ilo, ihi, h, lda, t1, lda, z, ldu, work,
815  $ nwork, result( 3 ) )
816 *
817 * Do Tests 5: | A - UZ T (UZ)' | / ( |A| n ulp )
818 * and 6: | I - UZ (UZ)' | / ( n ulp )
819 *
820  CALL dhst01( n, ilo, ihi, a, lda, t1, lda, uz, ldu, work,
821  $ nwork, result( 5 ) )
822 *
823 * Do Test 7: | T2 - T1 | / ( |T| n ulp )
824 *
825  CALL dget10( n, n, t2, lda, t1, lda, work, result( 7 ) )
826 *
827 * Do Test 8: | W2 - W1 | / ( max(|W1|,|W2|) ulp )
828 *
829  temp1 = zero
830  temp2 = zero
831  DO 130 j = 1, n
832  temp1 = max( temp1, abs( wr1( j ) )+abs( wi1( j ) ),
833  $ abs( wr2( j ) )+abs( wi2( j ) ) )
834  temp2 = max( temp2, abs( wr1( j )-wr2( j ) )+
835  & abs( wi1( j )-wi2( j ) ) )
836  130 CONTINUE
837 *
838  result( 8 ) = temp2 / max( unfl, ulp*max( temp1, temp2 ) )
839 *
840 * Compute the Left and Right Eigenvectors of T
841 *
842 * Compute the Right eigenvector Matrix:
843 *
844  ntest = 9
845  result( 9 ) = ulpinv
846 *
847 * Select last max(N/4,1) real, max(N/4,1) complex eigenvectors
848 *
849  nselc = 0
850  nselr = 0
851  j = n
852  140 CONTINUE
853  IF( wi1( j ).EQ.zero ) THEN
854  IF( nselr.LT.max( n / 4, 1 ) ) THEN
855  nselr = nselr + 1
856  SELECT( j ) = .true.
857  ELSE
858  SELECT( j ) = .false.
859  END IF
860  j = j - 1
861  ELSE
862  IF( nselc.LT.max( n / 4, 1 ) ) THEN
863  nselc = nselc + 1
864  SELECT( j ) = .true.
865  SELECT( j-1 ) = .false.
866  ELSE
867  SELECT( j ) = .false.
868  SELECT( j-1 ) = .false.
869  END IF
870  j = j - 2
871  END IF
872  IF( j.GT.0 )
873  $ GO TO 140
874 *
875  CALL dtrevc( 'Right', 'All', SELECT, n, t1, lda, dumma, ldu,
876  $ evectr, ldu, n, in, work, iinfo )
877  IF( iinfo.NE.0 ) THEN
878  WRITE( nounit, fmt = 9999 )'DTREVC(R,A)', iinfo, n,
879  $ jtype, ioldsd
880  info = abs( iinfo )
881  GO TO 250
882  END IF
883 *
884 * Test 9: | TR - RW | / ( |T| |R| ulp )
885 *
886  CALL dget22( 'N', 'N', 'N', n, t1, lda, evectr, ldu, wr1,
887  $ wi1, work, dumma( 1 ) )
888  result( 9 ) = dumma( 1 )
889  IF( dumma( 2 ).GT.thresh ) THEN
890  WRITE( nounit, fmt = 9998 )'Right', 'DTREVC',
891  $ dumma( 2 ), n, jtype, ioldsd
892  END IF
893 *
894 * Compute selected right eigenvectors and confirm that
895 * they agree with previous right eigenvectors
896 *
897  CALL dtrevc( 'Right', 'Some', SELECT, n, t1, lda, dumma,
898  $ ldu, evectl, ldu, n, in, work, iinfo )
899  IF( iinfo.NE.0 ) THEN
900  WRITE( nounit, fmt = 9999 )'DTREVC(R,S)', iinfo, n,
901  $ jtype, ioldsd
902  info = abs( iinfo )
903  GO TO 250
904  END IF
905 *
906  k = 1
907  match = .true.
908  DO 170 j = 1, n
909  IF( SELECT( j ) .AND. wi1( j ).EQ.zero ) THEN
910  DO 150 jj = 1, n
911  IF( evectr( jj, j ).NE.evectl( jj, k ) ) THEN
912  match = .false.
913  GO TO 180
914  END IF
915  150 CONTINUE
916  k = k + 1
917  ELSE IF( SELECT( j ) .AND. wi1( j ).NE.zero ) THEN
918  DO 160 jj = 1, n
919  IF( evectr( jj, j ).NE.evectl( jj, k ) .OR.
920  $ evectr( jj, j+1 ).NE.evectl( jj, k+1 ) ) THEN
921  match = .false.
922  GO TO 180
923  END IF
924  160 CONTINUE
925  k = k + 2
926  END IF
927  170 CONTINUE
928  180 CONTINUE
929  IF( .NOT.match )
930  $ WRITE( nounit, fmt = 9997 )'Right', 'DTREVC', n, jtype,
931  $ ioldsd
932 *
933 * Compute the Left eigenvector Matrix:
934 *
935  ntest = 10
936  result( 10 ) = ulpinv
937  CALL dtrevc( 'Left', 'All', SELECT, n, t1, lda, evectl, ldu,
938  $ dumma, ldu, n, in, work, iinfo )
939  IF( iinfo.NE.0 ) THEN
940  WRITE( nounit, fmt = 9999 )'DTREVC(L,A)', iinfo, n,
941  $ jtype, ioldsd
942  info = abs( iinfo )
943  GO TO 250
944  END IF
945 *
946 * Test 10: | LT - WL | / ( |T| |L| ulp )
947 *
948  CALL dget22( 'Trans', 'N', 'Conj', n, t1, lda, evectl, ldu,
949  $ wr1, wi1, work, dumma( 3 ) )
950  result( 10 ) = dumma( 3 )
951  IF( dumma( 4 ).GT.thresh ) THEN
952  WRITE( nounit, fmt = 9998 )'Left', 'DTREVC', dumma( 4 ),
953  $ n, jtype, ioldsd
954  END IF
955 *
956 * Compute selected left eigenvectors and confirm that
957 * they agree with previous left eigenvectors
958 *
959  CALL dtrevc( 'Left', 'Some', SELECT, n, t1, lda, evectr,
960  $ ldu, dumma, ldu, n, in, work, iinfo )
961  IF( iinfo.NE.0 ) THEN
962  WRITE( nounit, fmt = 9999 )'DTREVC(L,S)', iinfo, n,
963  $ jtype, ioldsd
964  info = abs( iinfo )
965  GO TO 250
966  END IF
967 *
968  k = 1
969  match = .true.
970  DO 210 j = 1, n
971  IF( SELECT( j ) .AND. wi1( j ).EQ.zero ) THEN
972  DO 190 jj = 1, n
973  IF( evectl( jj, j ).NE.evectr( jj, k ) ) THEN
974  match = .false.
975  GO TO 220
976  END IF
977  190 CONTINUE
978  k = k + 1
979  ELSE IF( SELECT( j ) .AND. wi1( j ).NE.zero ) THEN
980  DO 200 jj = 1, n
981  IF( evectl( jj, j ).NE.evectr( jj, k ) .OR.
982  $ evectl( jj, j+1 ).NE.evectr( jj, k+1 ) ) THEN
983  match = .false.
984  GO TO 220
985  END IF
986  200 CONTINUE
987  k = k + 2
988  END IF
989  210 CONTINUE
990  220 CONTINUE
991  IF( .NOT.match )
992  $ WRITE( nounit, fmt = 9997 )'Left', 'DTREVC', n, jtype,
993  $ ioldsd
994 *
995 * Call DHSEIN for Right eigenvectors of H, do test 11
996 *
997  ntest = 11
998  result( 11 ) = ulpinv
999  DO 230 j = 1, n
1000  SELECT( j ) = .true.
1001  230 CONTINUE
1002 *
1003  CALL dhsein( 'Right', 'Qr', 'Ninitv', SELECT, n, h, lda,
1004  $ wr3, wi3, dumma, ldu, evectx, ldu, n1, in,
1005  $ work, iwork, iwork, iinfo )
1006  IF( iinfo.NE.0 ) THEN
1007  WRITE( nounit, fmt = 9999 )'DHSEIN(R)', iinfo, n, jtype,
1008  $ ioldsd
1009  info = abs( iinfo )
1010  IF( iinfo.LT.0 )
1011  $ GO TO 250
1012  ELSE
1013 *
1014 * Test 11: | HX - XW | / ( |H| |X| ulp )
1015 *
1016 * (from inverse iteration)
1017 *
1018  CALL dget22( 'N', 'N', 'N', n, h, lda, evectx, ldu, wr3,
1019  $ wi3, work, dumma( 1 ) )
1020  IF( dumma( 1 ).LT.ulpinv )
1021  $ result( 11 ) = dumma( 1 )*aninv
1022  IF( dumma( 2 ).GT.thresh ) THEN
1023  WRITE( nounit, fmt = 9998 )'Right', 'DHSEIN',
1024  $ dumma( 2 ), n, jtype, ioldsd
1025  END IF
1026  END IF
1027 *
1028 * Call DHSEIN for Left eigenvectors of H, do test 12
1029 *
1030  ntest = 12
1031  result( 12 ) = ulpinv
1032  DO 240 j = 1, n
1033  SELECT( j ) = .true.
1034  240 CONTINUE
1035 *
1036  CALL dhsein( 'Left', 'Qr', 'Ninitv', SELECT, n, h, lda, wr3,
1037  $ wi3, evecty, ldu, dumma, ldu, n1, in, work,
1038  $ iwork, iwork, iinfo )
1039  IF( iinfo.NE.0 ) THEN
1040  WRITE( nounit, fmt = 9999 )'DHSEIN(L)', iinfo, n, jtype,
1041  $ ioldsd
1042  info = abs( iinfo )
1043  IF( iinfo.LT.0 )
1044  $ GO TO 250
1045  ELSE
1046 *
1047 * Test 12: | YH - WY | / ( |H| |Y| ulp )
1048 *
1049 * (from inverse iteration)
1050 *
1051  CALL dget22( 'C', 'N', 'C', n, h, lda, evecty, ldu, wr3,
1052  $ wi3, work, dumma( 3 ) )
1053  IF( dumma( 3 ).LT.ulpinv )
1054  $ result( 12 ) = dumma( 3 )*aninv
1055  IF( dumma( 4 ).GT.thresh ) THEN
1056  WRITE( nounit, fmt = 9998 )'Left', 'DHSEIN',
1057  $ dumma( 4 ), n, jtype, ioldsd
1058  END IF
1059  END IF
1060 *
1061 * Call DORMHR for Right eigenvectors of A, do test 13
1062 *
1063  ntest = 13
1064  result( 13 ) = ulpinv
1065 *
1066  CALL dormhr( 'Left', 'No transpose', n, n, ilo, ihi, uu,
1067  $ ldu, tau, evectx, ldu, work, nwork, iinfo )
1068  IF( iinfo.NE.0 ) THEN
1069  WRITE( nounit, fmt = 9999 )'DORMHR(R)', iinfo, n, jtype,
1070  $ ioldsd
1071  info = abs( iinfo )
1072  IF( iinfo.LT.0 )
1073  $ GO TO 250
1074  ELSE
1075 *
1076 * Test 13: | AX - XW | / ( |A| |X| ulp )
1077 *
1078 * (from inverse iteration)
1079 *
1080  CALL dget22( 'N', 'N', 'N', n, a, lda, evectx, ldu, wr3,
1081  $ wi3, work, dumma( 1 ) )
1082  IF( dumma( 1 ).LT.ulpinv )
1083  $ result( 13 ) = dumma( 1 )*aninv
1084  END IF
1085 *
1086 * Call DORMHR for Left eigenvectors of A, do test 14
1087 *
1088  ntest = 14
1089  result( 14 ) = ulpinv
1090 *
1091  CALL dormhr( 'Left', 'No transpose', n, n, ilo, ihi, uu,
1092  $ ldu, tau, evecty, ldu, work, nwork, iinfo )
1093  IF( iinfo.NE.0 ) THEN
1094  WRITE( nounit, fmt = 9999 )'DORMHR(L)', iinfo, n, jtype,
1095  $ ioldsd
1096  info = abs( iinfo )
1097  IF( iinfo.LT.0 )
1098  $ GO TO 250
1099  ELSE
1100 *
1101 * Test 14: | YA - WY | / ( |A| |Y| ulp )
1102 *
1103 * (from inverse iteration)
1104 *
1105  CALL dget22( 'C', 'N', 'C', n, a, lda, evecty, ldu, wr3,
1106  $ wi3, work, dumma( 3 ) )
1107  IF( dumma( 3 ).LT.ulpinv )
1108  $ result( 14 ) = dumma( 3 )*aninv
1109  END IF
1110 *
1111 * End of Loop -- Check for RESULT(j) > THRESH
1112 *
1113  250 CONTINUE
1114 *
1115  ntestt = ntestt + ntest
1116  CALL dlafts( 'DHS', n, n, jtype, ntest, result, ioldsd,
1117  $ thresh, nounit, nerrs )
1118 *
1119  260 CONTINUE
1120  270 CONTINUE
1121 *
1122 * Summary
1123 *
1124  CALL dlasum( 'DHS', nounit, nerrs, ntestt )
1125 *
1126  RETURN
1127 *
1128  9999 FORMAT( ' DCHKHS: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
1129  $ i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1130  9998 FORMAT( ' DCHKHS: ', a, ' Eigenvectors from ', a, ' incorrectly ',
1131  $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
1132  $ 'N=', i6, ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5,
1133  $ ')' )
1134  9997 FORMAT( ' DCHKHS: Selected ', a, ' Eigenvectors from ', a,
1135  $ ' do not match other eigenvectors ', 9x, 'N=', i6,
1136  $ ', JTYPE=', i6, ', ISEED=(', 3( i5, ',' ), i5, ')' )
1137 *
1138 * End of DCHKHS
1139 *
1140  END
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine dlacpy(UPLO, M, N, A, LDA, B, LDB)
DLACPY copies all or part of one two-dimensional array to another.
Definition: dlacpy.f:103
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:110
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine dcopy(N, DX, INCX, DY, INCY)
DCOPY
Definition: dcopy.f:82
subroutine dgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
DGEMM
Definition: dgemm.f:187
subroutine dhst01(N, ILO, IHI, A, LDA, H, LDH, Q, LDQ, WORK, LWORK, RESULT)
DHST01
Definition: dhst01.f:134
subroutine dlasum(TYPE, IOUNIT, IE, NRUN)
DLASUM
Definition: dlasum.f:43
subroutine dget22(TRANSA, TRANSE, TRANSW, N, A, LDA, E, LDE, WR, WI, WORK, RESULT)
DGET22
Definition: dget22.f:168
subroutine dlafts(TYPE, M, N, IMAT, NTESTS, RESULT, ISEED, THRESH, IOUNIT, IE)
DLAFTS
Definition: dlafts.f:99
subroutine dchkhs(NSIZES, NN, NTYPES, DOTYPE, ISEED, THRESH, NOUNIT, A, LDA, H, T1, T2, U, LDU, Z, UZ, WR1, WI1, WR2, WI2, WR3, WI3, EVECTL, EVECTR, EVECTY, EVECTX, UU, TAU, WORK, NWORK, IWORK, SELECT, RESULT, INFO)
DCHKHS
Definition: dchkhs.f:412
subroutine dget10(M, N, A, LDA, B, LDB, WORK, RESULT)
DGET10
Definition: dget10.f:93
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:471
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:332
subroutine dlatms(M, N, DIST, ISEED, SYM, D, MODE, COND, DMAX, KL, KU, PACK, A, LDA, WORK, INFO)
DLATMS
Definition: dlatms.f:321
subroutine dgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
DGEHRD
Definition: dgehrd.f:167
subroutine dhsein(SIDE, EIGSRC, INITV, SELECT, N, H, LDH, WR, WI, VL, LDVL, VR, LDVR, MM, M, WORK, IFAILL, IFAILR, INFO)
DHSEIN
Definition: dhsein.f:263
subroutine dtrevc(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, INFO)
DTREVC
Definition: dtrevc.f:222
subroutine dhseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, WR, WI, Z, LDZ, WORK, LWORK, INFO)
DHSEQR
Definition: dhseqr.f:316
subroutine dormhr(SIDE, TRANS, M, N, ILO, IHI, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
DORMHR
Definition: dormhr.f:178
subroutine dorghr(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
DORGHR
Definition: dorghr.f:126