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