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