LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
ddrgvx.f
Go to the documentation of this file.
1 *> \brief \b DDRGVX
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 DDRGVX( NSIZE, THRESH, NIN, NOUT, A, LDA, B, AI, BI,
12 * ALPHAR, ALPHAI, BETA, VL, VR, ILO, IHI, LSCALE,
13 * RSCALE, S, DTRU, DIF, DIFTRU, WORK, LWORK,
14 * IWORK, LIWORK, RESULT, BWORK, INFO )
15 *
16 * .. Scalar Arguments ..
17 * INTEGER IHI, ILO, INFO, LDA, LIWORK, LWORK, NIN, NOUT,
18 * $ NSIZE
19 * DOUBLE PRECISION THRESH
20 * ..
21 * .. Array Arguments ..
22 * LOGICAL BWORK( * )
23 * INTEGER IWORK( * )
24 * DOUBLE PRECISION A( LDA, * ), AI( LDA, * ), ALPHAI( * ),
25 * $ ALPHAR( * ), B( LDA, * ), BETA( * ),
26 * $ BI( LDA, * ), DIF( * ), DIFTRU( * ), DTRU( * ),
27 * $ LSCALE( * ), RESULT( 4 ), RSCALE( * ), S( * ),
28 * $ VL( LDA, * ), VR( LDA, * ), WORK( * )
29 * ..
30 *
31 *
32 *> \par Purpose:
33 * =============
34 *>
35 *> \verbatim
36 *>
37 *> DDRGVX checks the nonsymmetric generalized eigenvalue problem
38 *> expert driver DGGEVX.
39 *>
40 *> DGGEVX computes the generalized eigenvalues, (optionally) the left
41 *> and/or right eigenvectors, (optionally) computes a balancing
42 *> transformation to improve the conditioning, and (optionally)
43 *> reciprocal condition numbers for the eigenvalues and eigenvectors.
44 *>
45 *> When DDRGVX is called with NSIZE > 0, two types of test matrix pairs
46 *> are generated by the subroutine DLATM6 and test the driver DGGEVX.
47 *> The test matrices have the known exact condition numbers for
48 *> eigenvalues. For the condition numbers of the eigenvectors
49 *> corresponding the first and last eigenvalues are also know
50 *> ``exactly'' (see DLATM6).
51 *>
52 *> For each matrix pair, the following tests will be performed and
53 *> compared with the threshhold THRESH.
54 *>
55 *> (1) max over all left eigenvalue/-vector pairs (beta/alpha,l) of
56 *>
57 *> | l**H * (beta A - alpha B) | / ( ulp max( |beta A|, |alpha B| ) )
58 *>
59 *> where l**H is the conjugate tranpose of l.
60 *>
61 *> (2) max over all right eigenvalue/-vector pairs (beta/alpha,r) of
62 *>
63 *> | (beta A - alpha B) r | / ( ulp max( |beta A|, |alpha B| ) )
64 *>
65 *> (3) The condition number S(i) of eigenvalues computed by DGGEVX
66 *> differs less than a factor THRESH from the exact S(i) (see
67 *> DLATM6).
68 *>
69 *> (4) DIF(i) computed by DTGSNA differs less than a factor 10*THRESH
70 *> from the exact value (for the 1st and 5th vectors only).
71 *>
72 *> Test Matrices
73 *> =============
74 *>
75 *> Two kinds of test matrix pairs
76 *>
77 *> (A, B) = inverse(YH) * (Da, Db) * inverse(X)
78 *>
79 *> are used in the tests:
80 *>
81 *> 1: Da = 1+a 0 0 0 0 Db = 1 0 0 0 0
82 *> 0 2+a 0 0 0 0 1 0 0 0
83 *> 0 0 3+a 0 0 0 0 1 0 0
84 *> 0 0 0 4+a 0 0 0 0 1 0
85 *> 0 0 0 0 5+a , 0 0 0 0 1 , and
86 *>
87 *> 2: Da = 1 -1 0 0 0 Db = 1 0 0 0 0
88 *> 1 1 0 0 0 0 1 0 0 0
89 *> 0 0 1 0 0 0 0 1 0 0
90 *> 0 0 0 1+a 1+b 0 0 0 1 0
91 *> 0 0 0 -1-b 1+a , 0 0 0 0 1 .
92 *>
93 *> In both cases the same inverse(YH) and inverse(X) are used to compute
94 *> (A, B), giving the exact eigenvectors to (A,B) as (YH, X):
95 *>
96 *> YH: = 1 0 -y y -y X = 1 0 -x -x x
97 *> 0 1 -y y -y 0 1 x -x -x
98 *> 0 0 1 0 0 0 0 1 0 0
99 *> 0 0 0 1 0 0 0 0 1 0
100 *> 0 0 0 0 1, 0 0 0 0 1 , where
101 *>
102 *> a, b, x and y will have all values independently of each other from
103 *> { sqrt(sqrt(ULP)), 0.1, 1, 10, 1/sqrt(sqrt(ULP)) }.
104 *> \endverbatim
105 *
106 * Arguments:
107 * ==========
108 *
109 *> \param[in] NSIZE
110 *> \verbatim
111 *> NSIZE is INTEGER
112 *> The number of sizes of matrices to use. NSIZE must be at
113 *> least zero. If it is zero, no randomly generated matrices
114 *> are tested, but any test matrices read from NIN will be
115 *> tested.
116 *> \endverbatim
117 *>
118 *> \param[in] THRESH
119 *> \verbatim
120 *> THRESH is DOUBLE PRECISION
121 *> A test will count as "failed" if the "error", computed as
122 *> described above, exceeds THRESH. Note that the error
123 *> is scaled to be O(1), so THRESH should be a reasonably
124 *> small multiple of 1, e.g., 10 or 100. In particular,
125 *> it should not depend on the precision (single vs. double)
126 *> or the size of the matrix. It must be at least zero.
127 *> \endverbatim
128 *>
129 *> \param[in] NIN
130 *> \verbatim
131 *> NIN is INTEGER
132 *> The FORTRAN unit number for reading in the data file of
133 *> problems to solve.
134 *> \endverbatim
135 *>
136 *> \param[in] NOUT
137 *> \verbatim
138 *> NOUT is INTEGER
139 *> The FORTRAN unit number for printing out error messages
140 *> (e.g., if a routine returns IINFO not equal to 0.)
141 *> \endverbatim
142 *>
143 *> \param[out] A
144 *> \verbatim
145 *> A is DOUBLE PRECISION array, dimension (LDA, NSIZE)
146 *> Used to hold the matrix whose eigenvalues are to be
147 *> computed. On exit, A contains the last matrix actually used.
148 *> \endverbatim
149 *>
150 *> \param[in] LDA
151 *> \verbatim
152 *> LDA is INTEGER
153 *> The leading dimension of A, B, AI, BI, Ao, and Bo.
154 *> It must be at least 1 and at least NSIZE.
155 *> \endverbatim
156 *>
157 *> \param[out] B
158 *> \verbatim
159 *> B is DOUBLE PRECISION array, dimension (LDA, NSIZE)
160 *> Used to hold the matrix whose eigenvalues are to be
161 *> computed. On exit, B contains the last matrix actually used.
162 *> \endverbatim
163 *>
164 *> \param[out] AI
165 *> \verbatim
166 *> AI is DOUBLE PRECISION array, dimension (LDA, NSIZE)
167 *> Copy of A, modified by DGGEVX.
168 *> \endverbatim
169 *>
170 *> \param[out] BI
171 *> \verbatim
172 *> BI is DOUBLE PRECISION array, dimension (LDA, NSIZE)
173 *> Copy of B, modified by DGGEVX.
174 *> \endverbatim
175 *>
176 *> \param[out] ALPHAR
177 *> \verbatim
178 *> ALPHAR is DOUBLE PRECISION array, dimension (NSIZE)
179 *> \endverbatim
180 *>
181 *> \param[out] ALPHAI
182 *> \verbatim
183 *> ALPHAI is DOUBLE PRECISION array, dimension (NSIZE)
184 *> \endverbatim
185 *>
186 *> \param[out] BETA
187 *> \verbatim
188 *> BETA is DOUBLE PRECISION array, dimension (NSIZE)
189 *>
190 *> On exit, (ALPHAR + ALPHAI*i)/BETA are the eigenvalues.
191 *> \endverbatim
192 *>
193 *> \param[out] VL
194 *> \verbatim
195 *> VL is DOUBLE PRECISION array, dimension (LDA, NSIZE)
196 *> VL holds the left eigenvectors computed by DGGEVX.
197 *> \endverbatim
198 *>
199 *> \param[out] VR
200 *> \verbatim
201 *> VR is DOUBLE PRECISION array, dimension (LDA, NSIZE)
202 *> VR holds the right eigenvectors computed by DGGEVX.
203 *> \endverbatim
204 *>
205 *> \param[out] ILO
206 *> \verbatim
207 *> ILO is INTEGER
208 *> \endverbatim
209 *>
210 *> \param[out] IHI
211 *> \verbatim
212 *> IHI is INTEGER
213 *> \endverbatim
214 *>
215 *> \param[out] LSCALE
216 *> \verbatim
217 *> LSCALE is DOUBLE PRECISION array, dimension (N)
218 *> \endverbatim
219 *>
220 *> \param[out] RSCALE
221 *> \verbatim
222 *> RSCALE is DOUBLE PRECISION array, dimension (N)
223 *> \endverbatim
224 *>
225 *> \param[out] S
226 *> \verbatim
227 *> S is DOUBLE PRECISION array, dimension (N)
228 *> \endverbatim
229 *>
230 *> \param[out] DTRU
231 *> \verbatim
232 *> DTRU is DOUBLE PRECISION array, dimension (N)
233 *> \endverbatim
234 *>
235 *> \param[out] DIF
236 *> \verbatim
237 *> DIF is DOUBLE PRECISION array, dimension (N)
238 *> \endverbatim
239 *>
240 *> \param[out] DIFTRU
241 *> \verbatim
242 *> DIFTRU is DOUBLE PRECISION array, dimension (N)
243 *> \endverbatim
244 *>
245 *> \param[out] WORK
246 *> \verbatim
247 *> WORK is DOUBLE PRECISION array, dimension (LWORK)
248 *> \endverbatim
249 *>
250 *> \param[in] LWORK
251 *> \verbatim
252 *> LWORK is INTEGER
253 *> Leading dimension of WORK. LWORK >= 2*N*N+12*N+16.
254 *> \endverbatim
255 *>
256 *> \param[out] IWORK
257 *> \verbatim
258 *> IWORK is INTEGER array, dimension (LIWORK)
259 *> \endverbatim
260 *>
261 *> \param[in] LIWORK
262 *> \verbatim
263 *> LIWORK is INTEGER
264 *> Leading dimension of IWORK. Must be at least N+6.
265 *> \endverbatim
266 *>
267 *> \param[out] RESULT
268 *> \verbatim
269 *> RESULT is DOUBLE PRECISION array, dimension (4)
270 *> \endverbatim
271 *>
272 *> \param[out] BWORK
273 *> \verbatim
274 *> BWORK is LOGICAL array, dimension (N)
275 *> \endverbatim
276 *>
277 *> \param[out] INFO
278 *> \verbatim
279 *> INFO is INTEGER
280 *> = 0: successful exit
281 *> < 0: if INFO = -i, the i-th argument had an illegal value.
282 *> > 0: A routine returned an error code.
283 *> \endverbatim
284 *
285 * Authors:
286 * ========
287 *
288 *> \author Univ. of Tennessee
289 *> \author Univ. of California Berkeley
290 *> \author Univ. of Colorado Denver
291 *> \author NAG Ltd.
292 *
293 *> \date April 2012
294 *
295 *> \ingroup double_eig
296 *
297 * =====================================================================
298  SUBROUTINE ddrgvx( NSIZE, THRESH, NIN, NOUT, A, LDA, B, AI, BI,
299  $ alphar, alphai, beta, vl, vr, ilo, ihi, lscale,
300  $ rscale, s, dtru, dif, diftru, work, lwork,
301  $ iwork, liwork, result, bwork, info )
302 *
303 * -- LAPACK test routine (version 3.4.1) --
304 * -- LAPACK is a software package provided by Univ. of Tennessee, --
305 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
306 * April 2012
307 *
308 * .. Scalar Arguments ..
309  INTEGER ihi, ilo, info, lda, liwork, lwork, nin, nout,
310  $ nsize
311  DOUBLE PRECISION thresh
312 * ..
313 * .. Array Arguments ..
314  LOGICAL bwork( * )
315  INTEGER iwork( * )
316  DOUBLE PRECISION a( lda, * ), ai( lda, * ), alphai( * ),
317  $ alphar( * ), b( lda, * ), beta( * ),
318  $ bi( lda, * ), dif( * ), diftru( * ), dtru( * ),
319  $ lscale( * ), result( 4 ), rscale( * ), s( * ),
320  $ vl( lda, * ), vr( lda, * ), work( * )
321 * ..
322 *
323 * =====================================================================
324 *
325 * .. Parameters ..
326  DOUBLE PRECISION zero, one, ten, tnth, half
327  parameter( zero = 0.0d+0, one = 1.0d+0, ten = 1.0d+1,
328  $ tnth = 1.0d-1, half = 0.5d+0 )
329 * ..
330 * .. Local Scalars ..
331  INTEGER i, iptype, iwa, iwb, iwx, iwy, j, linfo,
332  $ maxwrk, minwrk, n, nerrs, nmax, nptknt, ntestt
333  DOUBLE PRECISION abnorm, anorm, bnorm, ratio1, ratio2, thrsh2,
334  $ ulp, ulpinv
335 * ..
336 * .. Local Arrays ..
337  DOUBLE PRECISION weight( 5 )
338 * ..
339 * .. External Functions ..
340  INTEGER ilaenv
341  DOUBLE PRECISION dlamch, dlange
342  EXTERNAL ilaenv, dlamch, dlange
343 * ..
344 * .. External Subroutines ..
345  EXTERNAL alasvm, dget52, dggevx, dlacpy, dlatm6, xerbla
346 * ..
347 * .. Intrinsic Functions ..
348  INTRINSIC abs, max, sqrt
349 * ..
350 * .. Executable Statements ..
351 *
352 * Check for errors
353 *
354  info = 0
355 *
356  nmax = 5
357 *
358  IF( nsize.LT.0 ) THEN
359  info = -1
360  ELSE IF( thresh.LT.zero ) THEN
361  info = -2
362  ELSE IF( nin.LE.0 ) THEN
363  info = -3
364  ELSE IF( nout.LE.0 ) THEN
365  info = -4
366  ELSE IF( lda.LT.1 .OR. lda.LT.nmax ) THEN
367  info = -6
368  ELSE IF( liwork.LT.nmax+6 ) THEN
369  info = -26
370  END IF
371 *
372 * Compute workspace
373 * (Note: Comments in the code beginning "Workspace:" describe the
374 * minimal amount of workspace needed at that point in the code,
375 * as well as the preferred amount for good performance.
376 * NB refers to the optimal block size for the immediately
377 * following subroutine, as returned by ILAENV.)
378 *
379  minwrk = 1
380  IF( info.EQ.0 .AND. lwork.GE.1 ) THEN
381  minwrk = 2*nmax*nmax + 12*nmax + 16
382  maxwrk = 6*nmax + nmax*ilaenv( 1, 'DGEQRF', ' ', nmax, 1, nmax,
383  $ 0 )
384  maxwrk = max( maxwrk, 2*nmax*nmax+12*nmax+16 )
385  work( 1 ) = maxwrk
386  END IF
387 *
388  IF( lwork.LT.minwrk )
389  $ info = -24
390 *
391  IF( info.NE.0 ) THEN
392  CALL xerbla( 'DDRGVX', -info )
393  return
394  END IF
395 *
396  n = 5
397  ulp = dlamch( 'P' )
398  ulpinv = one / ulp
399  thrsh2 = ten*thresh
400  nerrs = 0
401  nptknt = 0
402  ntestt = 0
403 *
404  IF( nsize.EQ.0 )
405  $ go to 90
406 *
407 * Parameters used for generating test matrices.
408 *
409  weight( 1 ) = tnth
410  weight( 2 ) = half
411  weight( 3 ) = one
412  weight( 4 ) = one / weight( 2 )
413  weight( 5 ) = one / weight( 1 )
414 *
415  DO 80 iptype = 1, 2
416  DO 70 iwa = 1, 5
417  DO 60 iwb = 1, 5
418  DO 50 iwx = 1, 5
419  DO 40 iwy = 1, 5
420 *
421 * generated a test matrix pair
422 *
423  CALL dlatm6( iptype, 5, a, lda, b, vr, lda, vl,
424  $ lda, weight( iwa ), weight( iwb ),
425  $ weight( iwx ), weight( iwy ), dtru,
426  $ diftru )
427 *
428 * Compute eigenvalues/eigenvectors of (A, B).
429 * Compute eigenvalue/eigenvector condition numbers
430 * using computed eigenvectors.
431 *
432  CALL dlacpy( 'F', n, n, a, lda, ai, lda )
433  CALL dlacpy( 'F', n, n, b, lda, bi, lda )
434 *
435  CALL dggevx( 'N', 'V', 'V', 'B', n, ai, lda, bi,
436  $ lda, alphar, alphai, beta, vl, lda,
437  $ vr, lda, ilo, ihi, lscale, rscale,
438  $ anorm, bnorm, s, dif, work, lwork,
439  $ iwork, bwork, linfo )
440  IF( linfo.NE.0 ) THEN
441  result( 1 ) = ulpinv
442  WRITE( nout, fmt = 9999 )'DGGEVX', linfo, n,
443  $ iptype
444  go to 30
445  END IF
446 *
447 * Compute the norm(A, B)
448 *
449  CALL dlacpy( 'Full', n, n, ai, lda, work, n )
450  CALL dlacpy( 'Full', n, n, bi, lda, work( n*n+1 ),
451  $ n )
452  abnorm = dlange( 'Fro', n, 2*n, work, n, work )
453 *
454 * Tests (1) and (2)
455 *
456  result( 1 ) = zero
457  CALL dget52( .true., n, a, lda, b, lda, vl, lda,
458  $ alphar, alphai, beta, work,
459  $ result( 1 ) )
460  IF( result( 2 ).GT.thresh ) THEN
461  WRITE( nout, fmt = 9998 )'Left', 'DGGEVX',
462  $ result( 2 ), n, iptype, iwa, iwb, iwx, iwy
463  END IF
464 *
465  result( 2 ) = zero
466  CALL dget52( .false., n, a, lda, b, lda, vr, lda,
467  $ alphar, alphai, beta, work,
468  $ result( 2 ) )
469  IF( result( 3 ).GT.thresh ) THEN
470  WRITE( nout, fmt = 9998 )'Right', 'DGGEVX',
471  $ result( 3 ), n, iptype, iwa, iwb, iwx, iwy
472  END IF
473 *
474 * Test (3)
475 *
476  result( 3 ) = zero
477  DO 10 i = 1, n
478  IF( s( i ).EQ.zero ) THEN
479  IF( dtru( i ).GT.abnorm*ulp )
480  $ result( 3 ) = ulpinv
481  ELSE IF( dtru( i ).EQ.zero ) THEN
482  IF( s( i ).GT.abnorm*ulp )
483  $ result( 3 ) = ulpinv
484  ELSE
485  work( i ) = max( abs( dtru( i ) / s( i ) ),
486  $ abs( s( i ) / dtru( i ) ) )
487  result( 3 ) = max( result( 3 ), work( i ) )
488  END IF
489  10 continue
490 *
491 * Test (4)
492 *
493  result( 4 ) = zero
494  IF( dif( 1 ).EQ.zero ) THEN
495  IF( diftru( 1 ).GT.abnorm*ulp )
496  $ result( 4 ) = ulpinv
497  ELSE IF( diftru( 1 ).EQ.zero ) THEN
498  IF( dif( 1 ).GT.abnorm*ulp )
499  $ result( 4 ) = ulpinv
500  ELSE IF( dif( 5 ).EQ.zero ) THEN
501  IF( diftru( 5 ).GT.abnorm*ulp )
502  $ result( 4 ) = ulpinv
503  ELSE IF( diftru( 5 ).EQ.zero ) THEN
504  IF( dif( 5 ).GT.abnorm*ulp )
505  $ result( 4 ) = ulpinv
506  ELSE
507  ratio1 = max( abs( diftru( 1 ) / dif( 1 ) ),
508  $ abs( dif( 1 ) / diftru( 1 ) ) )
509  ratio2 = max( abs( diftru( 5 ) / dif( 5 ) ),
510  $ abs( dif( 5 ) / diftru( 5 ) ) )
511  result( 4 ) = max( ratio1, ratio2 )
512  END IF
513 *
514  ntestt = ntestt + 4
515 *
516 * Print out tests which fail.
517 *
518  DO 20 j = 1, 4
519  IF( ( result( j ).GE.thrsh2 .AND. j.GE.4 ) .OR.
520  $ ( result( j ).GE.thresh .AND. j.LE.3 ) )
521  $ THEN
522 *
523 * If this is the first test to fail,
524 * print a header to the data file.
525 *
526  IF( nerrs.EQ.0 ) THEN
527  WRITE( nout, fmt = 9997 )'DXV'
528 *
529 * Print out messages for built-in examples
530 *
531 * Matrix types
532 *
533  WRITE( nout, fmt = 9995 )
534  WRITE( nout, fmt = 9994 )
535  WRITE( nout, fmt = 9993 )
536 *
537 * Tests performed
538 *
539  WRITE( nout, fmt = 9992 )'''',
540  $ 'transpose', ''''
541 *
542  END IF
543  nerrs = nerrs + 1
544  IF( result( j ).LT.10000.0d0 ) THEN
545  WRITE( nout, fmt = 9991 )iptype, iwa,
546  $ iwb, iwx, iwy, j, result( j )
547  ELSE
548  WRITE( nout, fmt = 9990 )iptype, iwa,
549  $ iwb, iwx, iwy, j, result( j )
550  END IF
551  END IF
552  20 continue
553 *
554  30 continue
555 *
556  40 continue
557  50 continue
558  60 continue
559  70 continue
560  80 continue
561 *
562  go to 150
563 *
564  90 continue
565 *
566 * Read in data from file to check accuracy of condition estimation
567 * Read input data until N=0
568 *
569  READ( nin, fmt = *, END = 150 )n
570  IF( n.EQ.0 )
571  $ go to 150
572  DO 100 i = 1, n
573  READ( nin, fmt = * )( a( i, j ), j = 1, n )
574  100 continue
575  DO 110 i = 1, n
576  READ( nin, fmt = * )( b( i, j ), j = 1, n )
577  110 continue
578  READ( nin, fmt = * )( dtru( i ), i = 1, n )
579  READ( nin, fmt = * )( diftru( i ), i = 1, n )
580 *
581  nptknt = nptknt + 1
582 *
583 * Compute eigenvalues/eigenvectors of (A, B).
584 * Compute eigenvalue/eigenvector condition numbers
585 * using computed eigenvectors.
586 *
587  CALL dlacpy( 'F', n, n, a, lda, ai, lda )
588  CALL dlacpy( 'F', n, n, b, lda, bi, lda )
589 *
590  CALL dggevx( 'N', 'V', 'V', 'B', n, ai, lda, bi, lda, alphar,
591  $ alphai, beta, vl, lda, vr, lda, ilo, ihi, lscale,
592  $ rscale, anorm, bnorm, s, dif, work, lwork, iwork,
593  $ bwork, linfo )
594 *
595  IF( linfo.NE.0 ) THEN
596  result( 1 ) = ulpinv
597  WRITE( nout, fmt = 9987 )'DGGEVX', linfo, n, nptknt
598  go to 140
599  END IF
600 *
601 * Compute the norm(A, B)
602 *
603  CALL dlacpy( 'Full', n, n, ai, lda, work, n )
604  CALL dlacpy( 'Full', n, n, bi, lda, work( n*n+1 ), n )
605  abnorm = dlange( 'Fro', n, 2*n, work, n, work )
606 *
607 * Tests (1) and (2)
608 *
609  result( 1 ) = zero
610  CALL dget52( .true., n, a, lda, b, lda, vl, lda, alphar, alphai,
611  $ beta, work, result( 1 ) )
612  IF( result( 2 ).GT.thresh ) THEN
613  WRITE( nout, fmt = 9986 )'Left', 'DGGEVX', result( 2 ), n,
614  $ nptknt
615  END IF
616 *
617  result( 2 ) = zero
618  CALL dget52( .false., n, a, lda, b, lda, vr, lda, alphar, alphai,
619  $ beta, work, result( 2 ) )
620  IF( result( 3 ).GT.thresh ) THEN
621  WRITE( nout, fmt = 9986 )'Right', 'DGGEVX', result( 3 ), n,
622  $ nptknt
623  END IF
624 *
625 * Test (3)
626 *
627  result( 3 ) = zero
628  DO 120 i = 1, n
629  IF( s( i ).EQ.zero ) THEN
630  IF( dtru( i ).GT.abnorm*ulp )
631  $ result( 3 ) = ulpinv
632  ELSE IF( dtru( i ).EQ.zero ) THEN
633  IF( s( i ).GT.abnorm*ulp )
634  $ result( 3 ) = ulpinv
635  ELSE
636  work( i ) = max( abs( dtru( i ) / s( i ) ),
637  $ abs( s( i ) / dtru( i ) ) )
638  result( 3 ) = max( result( 3 ), work( i ) )
639  END IF
640  120 continue
641 *
642 * Test (4)
643 *
644  result( 4 ) = zero
645  IF( dif( 1 ).EQ.zero ) THEN
646  IF( diftru( 1 ).GT.abnorm*ulp )
647  $ result( 4 ) = ulpinv
648  ELSE IF( diftru( 1 ).EQ.zero ) THEN
649  IF( dif( 1 ).GT.abnorm*ulp )
650  $ result( 4 ) = ulpinv
651  ELSE IF( dif( 5 ).EQ.zero ) THEN
652  IF( diftru( 5 ).GT.abnorm*ulp )
653  $ result( 4 ) = ulpinv
654  ELSE IF( diftru( 5 ).EQ.zero ) THEN
655  IF( dif( 5 ).GT.abnorm*ulp )
656  $ result( 4 ) = ulpinv
657  ELSE
658  ratio1 = max( abs( diftru( 1 ) / dif( 1 ) ),
659  $ abs( dif( 1 ) / diftru( 1 ) ) )
660  ratio2 = max( abs( diftru( 5 ) / dif( 5 ) ),
661  $ abs( dif( 5 ) / diftru( 5 ) ) )
662  result( 4 ) = max( ratio1, ratio2 )
663  END IF
664 *
665  ntestt = ntestt + 4
666 *
667 * Print out tests which fail.
668 *
669  DO 130 j = 1, 4
670  IF( result( j ).GE.thrsh2 ) THEN
671 *
672 * If this is the first test to fail,
673 * print a header to the data file.
674 *
675  IF( nerrs.EQ.0 ) THEN
676  WRITE( nout, fmt = 9997 )'DXV'
677 *
678 * Print out messages for built-in examples
679 *
680 * Matrix types
681 *
682  WRITE( nout, fmt = 9996 )
683 *
684 * Tests performed
685 *
686  WRITE( nout, fmt = 9992 )'''', 'transpose', ''''
687 *
688  END IF
689  nerrs = nerrs + 1
690  IF( result( j ).LT.10000.0d0 ) THEN
691  WRITE( nout, fmt = 9989 )nptknt, n, j, result( j )
692  ELSE
693  WRITE( nout, fmt = 9988 )nptknt, n, j, result( j )
694  END IF
695  END IF
696  130 continue
697 *
698  140 continue
699 *
700  go to 90
701  150 continue
702 *
703 * Summary
704 *
705  CALL alasvm( 'DXV', nout, nerrs, ntestt, 0 )
706 *
707  work( 1 ) = maxwrk
708 *
709  return
710 *
711  9999 format( ' DDRGVX: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
712  $ i6, ', JTYPE=', i6, ')' )
713 *
714  9998 format( ' DDRGVX: ', a, ' Eigenvectors from ', a, ' incorrectly ',
715  $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
716  $ 'N=', i6, ', JTYPE=', i6, ', IWA=', i5, ', IWB=', i5,
717  $ ', IWX=', i5, ', IWY=', i5 )
718 *
719  9997 format( / 1x, a3, ' -- Real Expert Eigenvalue/vector',
720  $ ' problem driver' )
721 *
722  9996 format( ' Input Example' )
723 *
724  9995 format( ' Matrix types: ', / )
725 *
726  9994 format( ' TYPE 1: Da is diagonal, Db is identity, ',
727  $ / ' A = Y^(-H) Da X^(-1), B = Y^(-H) Db X^(-1) ',
728  $ / ' YH and X are left and right eigenvectors. ', / )
729 *
730  9993 format( ' TYPE 2: Da is quasi-diagonal, Db is identity, ',
731  $ / ' A = Y^(-H) Da X^(-1), B = Y^(-H) Db X^(-1) ',
732  $ / ' YH and X are left and right eigenvectors. ', / )
733 *
734  9992 format( / ' Tests performed: ', / 4x,
735  $ ' a is alpha, b is beta, l is a left eigenvector, ', / 4x,
736  $ ' r is a right eigenvector and ', a, ' means ', a, '.',
737  $ / ' 1 = max | ( b A - a B )', a, ' l | / const.',
738  $ / ' 2 = max | ( b A - a B ) r | / const.',
739  $ / ' 3 = max ( Sest/Stru, Stru/Sest ) ',
740  $ ' over all eigenvalues', /
741  $ ' 4 = max( DIFest/DIFtru, DIFtru/DIFest ) ',
742  $ ' over the 1st and 5th eigenvectors', / )
743 *
744  9991 format( ' Type=', i2, ',', ' IWA=', i2, ', IWB=', i2, ', IWX=',
745  $ i2, ', IWY=', i2, ', result ', i2, ' is', 0p, f8.2 )
746  9990 format( ' Type=', i2, ',', ' IWA=', i2, ', IWB=', i2, ', IWX=',
747  $ i2, ', IWY=', i2, ', result ', i2, ' is', 1p, d10.3 )
748  9989 format( ' Input example #', i2, ', matrix order=', i4, ',',
749  $ ' result ', i2, ' is', 0p, f8.2 )
750  9988 format( ' Input example #', i2, ', matrix order=', i4, ',',
751  $ ' result ', i2, ' is', 1p, d10.3 )
752  9987 format( ' DDRGVX: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
753  $ i6, ', Input example #', i2, ')' )
754 *
755  9986 format( ' DDRGVX: ', a, ' Eigenvectors from ', a, ' incorrectly ',
756  $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
757  $ 'N=', i6, ', Input Example #', i2, ')' )
758 *
759 *
760 * End of DDRGVX
761 *
762  END