LAPACK  3.6.0
LAPACK: Linear Algebra PACKage
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
subroutine alasvm(TYPE, NOUT, NFAIL, NRUN, NERRS)
ALASVM
Definition: alasvm.f:75
subroutine sdrgvx(NSIZE, THRESH, NIN, NOUT, A, LDA, B, AI, BI, ALPHAR, ALPHAI, BETA, VL, VR, ILO, IHI, LSCALE, RSCALE, S, STRU, DIF, DIFTRU, WORK, LWORK, IWORK, LIWORK, RESULT, BWORK, INFO)
SDRGVX
Definition: sdrgvx.f:303
subroutine sget52(LEFT, N, A, LDA, B, LDB, E, LDE, ALPHAR, ALPHAI, BETA, WORK, RESULT)
SGET52
Definition: sget52.f:201
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine slatm6(TYPE, N, A, LDA, B, X, LDX, Y, LDY, ALPHA, BETA, WX, WY, S, DIF)
SLATM6
Definition: slatm6.f:178
subroutine sggevx(BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, B, LDB, ALPHAR, ALPHAI, BETA, VL, LDVL, VR, LDVR, ILO, IHI, LSCALE, RSCALE, ABNRM, BBNRM, RCONDE, RCONDV, WORK, LWORK, IWORK, BWORK, INFO)
SGGEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices ...
Definition: sggevx.f:393
subroutine slacpy(UPLO, M, N, A, LDA, B, LDB)
SLACPY copies all or part of one two-dimensional array to another.
Definition: slacpy.f:105