1 SUBROUTINE ZDRGVX( NSIZE, THRESH, NIN, NOUT, A, LDA, B, AI, BI,
2 $ ALPHA, BETA, VL, VR, ILO, IHI, LSCALE, RSCALE,
3 $ S, DTRU, DIF, DIFTRU, WORK, LWORK, RWORK,
4 $ IWORK, LIWORK, RESULT, BWORK, INFO )
5 *
6 * -- LAPACK test routine (version 3.1) --
7 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
8 * November 2006
9 *
10 * .. Scalar Arguments ..
11 INTEGER IHI, ILO, INFO, LDA, LIWORK, LWORK, NIN, NOUT,
12 $ NSIZE
13 DOUBLE PRECISION THRESH
14 * ..
15 * .. Array Arguments ..
16 LOGICAL BWORK( * )
17 INTEGER IWORK( * )
18 DOUBLE PRECISION DIF( * ), DIFTRU( * ), DTRU( * ), LSCALE( * ),
19 $ RESULT( 4 ), RSCALE( * ), RWORK( * ), S( * )
20 COMPLEX*16 A( LDA, * ), AI( LDA, * ), ALPHA( * ),
21 $ B( LDA, * ), BETA( * ), BI( LDA, * ),
22 $ VL( LDA, * ), VR( LDA, * ), WORK( * )
23 * ..
24 *
25 * Purpose
26 * =======
27 *
28 * ZDRGVX checks the nonsymmetric generalized eigenvalue problem
29 * expert driver ZGGEVX.
30 *
31 * ZGGEVX computes the generalized eigenvalues, (optionally) the left
32 * and/or right eigenvectors, (optionally) computes a balancing
33 * transformation to improve the conditioning, and (optionally)
34 * reciprocal condition numbers for the eigenvalues and eigenvectors.
35 *
36 * When ZDRGVX is called with NSIZE > 0, two types of test matrix pairs
37 * are generated by the subroutine DLATM6 and test the driver ZGGEVX.
38 * The test matrices have the known exact condition numbers for
39 * eigenvalues. For the condition numbers of the eigenvectors
40 * corresponding the first and last eigenvalues are also know
41 * ``exactly'' (see ZLATM6).
42 * For each matrix pair, the following tests will be performed and
43 * compared with the threshhold THRESH.
44 *
45 * (1) max over all left eigenvalue/-vector pairs (beta/alpha,l) of
46 *
47 * | l**H * (beta A - alpha B) | / ( ulp max( |beta A|, |alpha B| ) )
48 *
49 * where l**H is the conjugate tranpose of l.
50 *
51 * (2) max over all right eigenvalue/-vector pairs (beta/alpha,r) of
52 *
53 * | (beta A - alpha B) r | / ( ulp max( |beta A|, |alpha B| ) )
54 *
55 * (3) The condition number S(i) of eigenvalues computed by ZGGEVX
56 * differs less than a factor THRESH from the exact S(i) (see
57 * ZLATM6).
58 *
59 * (4) DIF(i) computed by ZTGSNA differs less than a factor 10*THRESH
60 * from the exact value (for the 1st and 5th vectors only).
61 *
62 * Test Matrices
63 * =============
64 *
65 * Two kinds of test matrix pairs
66 * (A, B) = inverse(YH) * (Da, Db) * inverse(X)
67 * are used in the tests:
68 *
69 * 1: Da = 1+a 0 0 0 0 Db = 1 0 0 0 0
70 * 0 2+a 0 0 0 0 1 0 0 0
71 * 0 0 3+a 0 0 0 0 1 0 0
72 * 0 0 0 4+a 0 0 0 0 1 0
73 * 0 0 0 0 5+a , 0 0 0 0 1 , and
74 *
75 * 2: Da = 1 -1 0 0 0 Db = 1 0 0 0 0
76 * 1 1 0 0 0 0 1 0 0 0
77 * 0 0 1 0 0 0 0 1 0 0
78 * 0 0 0 1+a 1+b 0 0 0 1 0
79 * 0 0 0 -1-b 1+a , 0 0 0 0 1 .
80 *
81 * In both cases the same inverse(YH) and inverse(X) are used to compute
82 * (A, B), giving the exact eigenvectors to (A,B) as (YH, X):
83 *
84 * YH: = 1 0 -y y -y X = 1 0 -x -x x
85 * 0 1 -y y -y 0 1 x -x -x
86 * 0 0 1 0 0 0 0 1 0 0
87 * 0 0 0 1 0 0 0 0 1 0
88 * 0 0 0 0 1, 0 0 0 0 1 , where
89 *
90 * a, b, x and y will have all values independently of each other from
91 * { sqrt(sqrt(ULP)), 0.1, 1, 10, 1/sqrt(sqrt(ULP)) }.
92 *
93 * Arguments
94 * =========
95 *
96 * NSIZE (input) INTEGER
97 * The number of sizes of matrices to use. NSIZE must be at
98 * least zero. If it is zero, no randomly generated matrices
99 * are tested, but any test matrices read from NIN will be
100 * tested. If it is not zero, then N = 5.
101 *
102 * THRESH (input) DOUBLE PRECISION
103 * A test will count as "failed" if the "error", computed as
104 * described above, exceeds THRESH. Note that the error
105 * is scaled to be O(1), so THRESH should be a reasonably
106 * small multiple of 1, e.g., 10 or 100. In particular,
107 * it should not depend on the precision (single vs. double)
108 * or the size of the matrix. It must be at least zero.
109 *
110 * NIN (input) INTEGER
111 * The FORTRAN unit number for reading in the data file of
112 * problems to solve.
113 *
114 * NOUT (input) INTEGER
115 * The FORTRAN unit number for printing out error messages
116 * (e.g., if a routine returns IINFO not equal to 0.)
117 *
118 * A (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
119 * Used to hold the matrix whose eigenvalues are to be
120 * computed. On exit, A contains the last matrix actually used.
121 *
122 * LDA (input) INTEGER
123 * The leading dimension of A, B, AI, BI, Ao, and Bo.
124 * It must be at least 1 and at least NSIZE.
125 *
126 * B (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
127 * Used to hold the matrix whose eigenvalues are to be
128 * computed. On exit, B contains the last matrix actually used.
129 *
130 * AI (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
131 * Copy of A, modified by ZGGEVX.
132 *
133 * BI (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
134 * Copy of B, modified by ZGGEVX.
135 *
136 * ALPHA (workspace) COMPLEX*16 array, dimension (NSIZE)
137 * BETA (workspace) COMPLEX*16 array, dimension (NSIZE)
138 * On exit, ALPHA/BETA are the eigenvalues.
139 *
140 * VL (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
141 * VL holds the left eigenvectors computed by ZGGEVX.
142 *
143 * VR (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
144 * VR holds the right eigenvectors computed by ZGGEVX.
145 *
146 * ILO (output/workspace) INTEGER
147 *
148 * IHI (output/workspace) INTEGER
149 *
150 * LSCALE (output/workspace) DOUBLE PRECISION array, dimension (N)
151 *
152 * RSCALE (output/workspace) DOUBLE PRECISION array, dimension (N)
153 *
154 * S (output/workspace) DOUBLE PRECISION array, dimension (N)
155 *
156 * DTRU (output/workspace) DOUBLE PRECISION array, dimension (N)
157 *
158 * DIF (output/workspace) DOUBLE PRECISION array, dimension (N)
159 *
160 * DIFTRU (output/workspace) DOUBLE PRECISION array, dimension (N)
161 *
162 * WORK (workspace) COMPLEX*16 array, dimension (LWORK)
163 *
164 * LWORK (input) INTEGER
165 * Leading dimension of WORK. LWORK >= 2*N*N + 2*N
166 *
167 * RWORK (workspace) DOUBLE PRECISION array, dimension (6*N)
168 *
169 * IWORK (workspace) INTEGER array, dimension (LIWORK)
170 *
171 * LIWORK (input) INTEGER
172 * Leading dimension of IWORK. LIWORK >= N+2.
173 *
174 * RESULT (output/workspace) DOUBLE PRECISION array, dimension (4)
175 *
176 * BWORK (workspace) LOGICAL array, dimension (N)
177 *
178 * INFO (output) INTEGER
179 * = 0: successful exit
180 * < 0: if INFO = -i, the i-th argument had an illegal value.
181 * > 0: A routine returned an error code.
182 *
183 * =====================================================================
184 *
185 * .. Parameters ..
186 DOUBLE PRECISION ZERO, ONE, TEN, TNTH
187 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0, TEN = 1.0D+1,
188 $ TNTH = 1.0D-1 )
189 * ..
190 * .. Local Scalars ..
191 INTEGER I, IPTYPE, IWA, IWB, IWX, IWY, J, LINFO,
192 $ MAXWRK, MINWRK, N, NERRS, NMAX, NPTKNT, NTESTT
193 DOUBLE PRECISION ABNORM, ANORM, BNORM, RATIO1, RATIO2, THRSH2,
194 $ ULP, ULPINV
195 * ..
196 * .. Local Arrays ..
197 COMPLEX*16 WEIGHT( 5 )
198 * ..
199 * .. External Functions ..
200 INTEGER ILAENV
201 DOUBLE PRECISION DLAMCH, ZLANGE
202 EXTERNAL ILAENV, DLAMCH, ZLANGE
203 * ..
204 * .. External Subroutines ..
205 EXTERNAL ALASVM, XERBLA, ZGET52, ZGGEVX, ZLACPY, ZLATM6
206 * ..
207 * .. Intrinsic Functions ..
208 INTRINSIC ABS, DCMPLX, MAX, SQRT
209 * ..
210 * .. Executable Statements ..
211 *
212 * Check for errors
213 *
214 INFO = 0
215 *
216 NMAX = 5
217 *
218 IF( NSIZE.LT.0 ) THEN
219 INFO = -1
220 ELSE IF( THRESH.LT.ZERO ) THEN
221 INFO = -2
222 ELSE IF( NIN.LE.0 ) THEN
223 INFO = -3
224 ELSE IF( NOUT.LE.0 ) THEN
225 INFO = -4
226 ELSE IF( LDA.LT.1 .OR. LDA.LT.NMAX ) THEN
227 INFO = -6
228 ELSE IF( LIWORK.LT.NMAX+2 ) THEN
229 INFO = -26
230 END IF
231 *
232 * Compute workspace
233 * (Note: Comments in the code beginning "Workspace:" describe the
234 * minimal amount of workspace needed at that point in the code,
235 * as well as the preferred amount for good performance.
236 * NB refers to the optimal block size for the immediately
237 * following subroutine, as returned by ILAENV.)
238 *
239 MINWRK = 1
240 IF( INFO.EQ.0 .AND. LWORK.GE.1 ) THEN
241 MINWRK = 2*NMAX*( NMAX+1 )
242 MAXWRK = NMAX*( 1+ILAENV( 1, 'ZGEQRF', ' ', NMAX, 1, NMAX,
243 $ 0 ) )
244 MAXWRK = MAX( MAXWRK, 2*NMAX*( NMAX+1 ) )
245 WORK( 1 ) = MAXWRK
246 END IF
247 *
248 IF( LWORK.LT.MINWRK )
249 $ INFO = -23
250 *
251 IF( INFO.NE.0 ) THEN
252 CALL XERBLA( 'ZDRGVX', -INFO )
253 RETURN
254 END IF
255 *
256 N = 5
257 ULP = DLAMCH( 'P' )
258 ULPINV = ONE / ULP
259 THRSH2 = TEN*THRESH
260 NERRS = 0
261 NPTKNT = 0
262 NTESTT = 0
263 *
264 IF( NSIZE.EQ.0 )
265 $ GO TO 90
266 *
267 * Parameters used for generating test matrices.
268 *
269 WEIGHT( 1 ) = DCMPLX( SQRT( SQRT( ULP ) ), ZERO )
270 WEIGHT( 2 ) = DCMPLX( TNTH, ZERO )
271 WEIGHT( 3 ) = ONE
272 WEIGHT( 4 ) = ONE / WEIGHT( 2 )
273 WEIGHT( 5 ) = ONE / WEIGHT( 1 )
274 *
275 DO 80 IPTYPE = 1, 2
276 DO 70 IWA = 1, 5
277 DO 60 IWB = 1, 5
278 DO 50 IWX = 1, 5
279 DO 40 IWY = 1, 5
280 *
281 * generated a pair of test matrix
282 *
283 CALL ZLATM6( IPTYPE, 5, A, LDA, B, VR, LDA, VL,
284 $ LDA, WEIGHT( IWA ), WEIGHT( IWB ),
285 $ WEIGHT( IWX ), WEIGHT( IWY ), DTRU,
286 $ DIFTRU )
287 *
288 * Compute eigenvalues/eigenvectors of (A, B).
289 * Compute eigenvalue/eigenvector condition numbers
290 * using computed eigenvectors.
291 *
292 CALL ZLACPY( 'F', N, N, A, LDA, AI, LDA )
293 CALL ZLACPY( 'F', N, N, B, LDA, BI, LDA )
294 *
295 CALL ZGGEVX( 'N', 'V', 'V', 'B', N, AI, LDA, BI,
296 $ LDA, ALPHA, BETA, VL, LDA, VR, LDA,
297 $ ILO, IHI, LSCALE, RSCALE, ANORM,
298 $ BNORM, S, DIF, WORK, LWORK, RWORK,
299 $ IWORK, BWORK, LINFO )
300 IF( LINFO.NE.0 ) THEN
301 WRITE( NOUT, FMT = 9999 )'ZGGEVX', LINFO, N,
302 $ IPTYPE, IWA, IWB, IWX, IWY
303 GO TO 30
304 END IF
305 *
306 * Compute the norm(A, B)
307 *
308 CALL ZLACPY( 'Full', N, N, AI, LDA, WORK, N )
309 CALL ZLACPY( 'Full', N, N, BI, LDA, WORK( N*N+1 ),
310 $ N )
311 ABNORM = ZLANGE( 'Fro', N, 2*N, WORK, N, RWORK )
312 *
313 * Tests (1) and (2)
314 *
315 RESULT( 1 ) = ZERO
316 CALL ZGET52( .TRUE., N, A, LDA, B, LDA, VL, LDA,
317 $ ALPHA, BETA, WORK, RWORK,
318 $ RESULT( 1 ) )
319 IF( RESULT( 2 ).GT.THRESH ) THEN
320 WRITE( NOUT, FMT = 9998 )'Left', 'ZGGEVX',
321 $ RESULT( 2 ), N, IPTYPE, IWA, IWB, IWX, IWY
322 END IF
323 *
324 RESULT( 2 ) = ZERO
325 CALL ZGET52( .FALSE., N, A, LDA, B, LDA, VR, LDA,
326 $ ALPHA, BETA, WORK, RWORK,
327 $ RESULT( 2 ) )
328 IF( RESULT( 3 ).GT.THRESH ) THEN
329 WRITE( NOUT, FMT = 9998 )'Right', 'ZGGEVX',
330 $ RESULT( 3 ), N, IPTYPE, IWA, IWB, IWX, IWY
331 END IF
332 *
333 * Test (3)
334 *
335 RESULT( 3 ) = ZERO
336 DO 10 I = 1, N
337 IF( S( I ).EQ.ZERO ) THEN
338 IF( DTRU( I ).GT.ABNORM*ULP )
339 $ RESULT( 3 ) = ULPINV
340 ELSE IF( DTRU( I ).EQ.ZERO ) THEN
341 IF( S( I ).GT.ABNORM*ULP )
342 $ RESULT( 3 ) = ULPINV
343 ELSE
344 RWORK( I ) = MAX( ABS( DTRU( I ) / S( I ) ),
345 $ ABS( S( I ) / DTRU( I ) ) )
346 RESULT( 3 ) = MAX( RESULT( 3 ), RWORK( I ) )
347 END IF
348 10 CONTINUE
349 *
350 * Test (4)
351 *
352 RESULT( 4 ) = ZERO
353 IF( DIF( 1 ).EQ.ZERO ) THEN
354 IF( DIFTRU( 1 ).GT.ABNORM*ULP )
355 $ RESULT( 4 ) = ULPINV
356 ELSE IF( DIFTRU( 1 ).EQ.ZERO ) THEN
357 IF( DIF( 1 ).GT.ABNORM*ULP )
358 $ RESULT( 4 ) = ULPINV
359 ELSE IF( DIF( 5 ).EQ.ZERO ) THEN
360 IF( DIFTRU( 5 ).GT.ABNORM*ULP )
361 $ RESULT( 4 ) = ULPINV
362 ELSE IF( DIFTRU( 5 ).EQ.ZERO ) THEN
363 IF( DIF( 5 ).GT.ABNORM*ULP )
364 $ RESULT( 4 ) = ULPINV
365 ELSE
366 RATIO1 = MAX( ABS( DIFTRU( 1 ) / DIF( 1 ) ),
367 $ ABS( DIF( 1 ) / DIFTRU( 1 ) ) )
368 RATIO2 = MAX( ABS( DIFTRU( 5 ) / DIF( 5 ) ),
369 $ ABS( DIF( 5 ) / DIFTRU( 5 ) ) )
370 RESULT( 4 ) = MAX( RATIO1, RATIO2 )
371 END IF
372 *
373 NTESTT = NTESTT + 4
374 *
375 * Print out tests which fail.
376 *
377 DO 20 J = 1, 4
378 IF( ( RESULT( J ).GE.THRSH2 .AND. J.GE.4 ) .OR.
379 $ ( RESULT( J ).GE.THRESH .AND. J.LE.3 ) )
380 $ THEN
381 *
382 * If this is the first test to fail,
383 * print a header to the data file.
384 *
385 IF( NERRS.EQ.0 ) THEN
386 WRITE( NOUT, FMT = 9997 )'ZXV'
387 *
388 * Print out messages for built-in examples
389 *
390 * Matrix types
391 *
392 WRITE( NOUT, FMT = 9995 )
393 WRITE( NOUT, FMT = 9994 )
394 WRITE( NOUT, FMT = 9993 )
395 *
396 * Tests performed
397 *
398 WRITE( NOUT, FMT = 9992 )'''',
399 $ 'transpose', ''''
400 *
401 END IF
402 NERRS = NERRS + 1
403 IF( RESULT( J ).LT.10000.0D0 ) THEN
404 WRITE( NOUT, FMT = 9991 )IPTYPE, IWA,
405 $ IWB, IWX, IWY, J, RESULT( J )
406 ELSE
407 WRITE( NOUT, FMT = 9990 )IPTYPE, IWA,
408 $ IWB, IWX, IWY, J, RESULT( J )
409 END IF
410 END IF
411 20 CONTINUE
412 *
413 30 CONTINUE
414 *
415 40 CONTINUE
416 50 CONTINUE
417 60 CONTINUE
418 70 CONTINUE
419 80 CONTINUE
420 *
421 GO TO 150
422 *
423 90 CONTINUE
424 *
425 * Read in data from file to check accuracy of condition estimation
426 * Read input data until N=0
427 *
428 READ( NIN, FMT = *, END = 150 )N
429 IF( N.EQ.0 )
430 $ GO TO 150
431 DO 100 I = 1, N
432 READ( NIN, FMT = * )( A( I, J ), J = 1, N )
433 100 CONTINUE
434 DO 110 I = 1, N
435 READ( NIN, FMT = * )( B( I, J ), J = 1, N )
436 110 CONTINUE
437 READ( NIN, FMT = * )( DTRU( I ), I = 1, N )
438 READ( NIN, FMT = * )( DIFTRU( I ), I = 1, N )
439 *
440 NPTKNT = NPTKNT + 1
441 *
442 * Compute eigenvalues/eigenvectors of (A, B).
443 * Compute eigenvalue/eigenvector condition numbers
444 * using computed eigenvectors.
445 *
446 CALL ZLACPY( 'F', N, N, A, LDA, AI, LDA )
447 CALL ZLACPY( 'F', N, N, B, LDA, BI, LDA )
448 *
449 CALL ZGGEVX( 'N', 'V', 'V', 'B', N, AI, LDA, BI, LDA, ALPHA, BETA,
450 $ VL, LDA, VR, LDA, ILO, IHI, LSCALE, RSCALE, ANORM,
451 $ BNORM, S, DIF, WORK, LWORK, RWORK, IWORK, BWORK,
452 $ LINFO )
453 *
454 IF( LINFO.NE.0 ) THEN
455 WRITE( NOUT, FMT = 9987 )'ZGGEVX', LINFO, N, NPTKNT
456 GO TO 140
457 END IF
458 *
459 * Compute the norm(A, B)
460 *
461 CALL ZLACPY( 'Full', N, N, AI, LDA, WORK, N )
462 CALL ZLACPY( 'Full', N, N, BI, LDA, WORK( N*N+1 ), N )
463 ABNORM = ZLANGE( 'Fro', N, 2*N, WORK, N, RWORK )
464 *
465 * Tests (1) and (2)
466 *
467 RESULT( 1 ) = ZERO
468 CALL ZGET52( .TRUE., N, A, LDA, B, LDA, VL, LDA, ALPHA, BETA,
469 $ WORK, RWORK, RESULT( 1 ) )
470 IF( RESULT( 2 ).GT.THRESH ) THEN
471 WRITE( NOUT, FMT = 9986 )'Left', 'ZGGEVX', RESULT( 2 ), N,
472 $ NPTKNT
473 END IF
474 *
475 RESULT( 2 ) = ZERO
476 CALL ZGET52( .FALSE., N, A, LDA, B, LDA, VR, LDA, ALPHA, BETA,
477 $ WORK, RWORK, RESULT( 2 ) )
478 IF( RESULT( 3 ).GT.THRESH ) THEN
479 WRITE( NOUT, FMT = 9986 )'Right', 'ZGGEVX', RESULT( 3 ), N,
480 $ NPTKNT
481 END IF
482 *
483 * Test (3)
484 *
485 RESULT( 3 ) = ZERO
486 DO 120 I = 1, N
487 IF( S( I ).EQ.ZERO ) THEN
488 IF( DTRU( I ).GT.ABNORM*ULP )
489 $ RESULT( 3 ) = ULPINV
490 ELSE IF( DTRU( I ).EQ.ZERO ) THEN
491 IF( S( I ).GT.ABNORM*ULP )
492 $ RESULT( 3 ) = ULPINV
493 ELSE
494 RWORK( I ) = MAX( ABS( DTRU( I ) / S( I ) ),
495 $ ABS( S( I ) / DTRU( I ) ) )
496 RESULT( 3 ) = MAX( RESULT( 3 ), RWORK( I ) )
497 END IF
498 120 CONTINUE
499 *
500 * Test (4)
501 *
502 RESULT( 4 ) = ZERO
503 IF( DIF( 1 ).EQ.ZERO ) THEN
504 IF( DIFTRU( 1 ).GT.ABNORM*ULP )
505 $ RESULT( 4 ) = ULPINV
506 ELSE IF( DIFTRU( 1 ).EQ.ZERO ) THEN
507 IF( DIF( 1 ).GT.ABNORM*ULP )
508 $ RESULT( 4 ) = ULPINV
509 ELSE IF( DIF( 5 ).EQ.ZERO ) THEN
510 IF( DIFTRU( 5 ).GT.ABNORM*ULP )
511 $ RESULT( 4 ) = ULPINV
512 ELSE IF( DIFTRU( 5 ).EQ.ZERO ) THEN
513 IF( DIF( 5 ).GT.ABNORM*ULP )
514 $ RESULT( 4 ) = ULPINV
515 ELSE
516 RATIO1 = MAX( ABS( DIFTRU( 1 ) / DIF( 1 ) ),
517 $ ABS( DIF( 1 ) / DIFTRU( 1 ) ) )
518 RATIO2 = MAX( ABS( DIFTRU( 5 ) / DIF( 5 ) ),
519 $ ABS( DIF( 5 ) / DIFTRU( 5 ) ) )
520 RESULT( 4 ) = MAX( RATIO1, RATIO2 )
521 END IF
522 *
523 NTESTT = NTESTT + 4
524 *
525 * Print out tests which fail.
526 *
527 DO 130 J = 1, 4
528 IF( RESULT( J ).GE.THRSH2 ) THEN
529 *
530 * If this is the first test to fail,
531 * print a header to the data file.
532 *
533 IF( NERRS.EQ.0 ) THEN
534 WRITE( NOUT, FMT = 9997 )'ZXV'
535 *
536 * Print out messages for built-in examples
537 *
538 * Matrix types
539 *
540 WRITE( NOUT, FMT = 9996 )
541 *
542 * Tests performed
543 *
544 WRITE( NOUT, FMT = 9992 )'''', 'transpose', ''''
545 *
546 END IF
547 NERRS = NERRS + 1
548 IF( RESULT( J ).LT.10000.0D0 ) THEN
549 WRITE( NOUT, FMT = 9989 )NPTKNT, N, J, RESULT( J )
550 ELSE
551 WRITE( NOUT, FMT = 9988 )NPTKNT, N, J, RESULT( J )
552 END IF
553 END IF
554 130 CONTINUE
555 *
556 140 CONTINUE
557 *
558 GO TO 90
559 150 CONTINUE
560 *
561 * Summary
562 *
563 CALL ALASVM( 'ZXV', NOUT, NERRS, NTESTT, 0 )
564 *
565 WORK( 1 ) = MAXWRK
566 *
567 RETURN
568 *
569 9999 FORMAT( ' ZDRGVX: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
570 $ I6, ', JTYPE=', I6, ')' )
571 *
572 9998 FORMAT( ' ZDRGVX: ', A, ' Eigenvectors from ', A, ' incorrectly ',
573 $ 'normalized.', / ' Bits of error=', 0P, G10.3, ',', 9X,
574 $ 'N=', I6, ', JTYPE=', I6, ', IWA=', I5, ', IWB=', I5,
575 $ ', IWX=', I5, ', IWY=', I5 )
576 *
577 9997 FORMAT( / 1X, A3, ' -- Complex Expert Eigenvalue/vector',
578 $ ' problem driver' )
579 *
580 9996 FORMAT( 'Input Example' )
581 *
582 9995 FORMAT( ' Matrix types: ', / )
583 *
584 9994 FORMAT( ' TYPE 1: Da is diagonal, Db is identity, ',
585 $ / ' A = Y^(-H) Da X^(-1), B = Y^(-H) Db X^(-1) ',
586 $ / ' YH and X are left and right eigenvectors. ', / )
587 *
588 9993 FORMAT( ' TYPE 2: Da is quasi-diagonal, Db is identity, ',
589 $ / ' A = Y^(-H) Da X^(-1), B = Y^(-H) Db X^(-1) ',
590 $ / ' YH and X are left and right eigenvectors. ', / )
591 *
592 9992 FORMAT( / ' Tests performed: ', / 4X,
593 $ ' a is alpha, b is beta, l is a left eigenvector, ', / 4X,
594 $ ' r is a right eigenvector and ', A, ' means ', A, '.',
595 $ / ' 1 = max | ( b A - a B )', A, ' l | / const.',
596 $ / ' 2 = max | ( b A - a B ) r | / const.',
597 $ / ' 3 = max ( Sest/Stru, Stru/Sest ) ',
598 $ ' over all eigenvalues', /
599 $ ' 4 = max( DIFest/DIFtru, DIFtru/DIFest ) ',
600 $ ' over the 1st and 5th eigenvectors', / )
601 *
602 9991 FORMAT( ' Type=', I2, ',', ' IWA=', I2, ', IWB=', I2, ', IWX=',
603 $ I2, ', IWY=', I2, ', result ', I2, ' is', 0P, F8.2 )
604 *
605 9990 FORMAT( ' Type=', I2, ',', ' IWA=', I2, ', IWB=', I2, ', IWX=',
606 $ I2, ', IWY=', I2, ', result ', I2, ' is', 1P, D10.3 )
607 *
608 9989 FORMAT( ' Input example #', I2, ', matrix order=', I4, ',',
609 $ ' result ', I2, ' is', 0P, F8.2 )
610 *
611 9988 FORMAT( ' Input example #', I2, ', matrix order=', I4, ',',
612 $ ' result ', I2, ' is', 1P, D10.3 )
613 *
614 9987 FORMAT( ' ZDRGVX: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
615 $ I6, ', Input example #', I2, ')' )
616 *
617 9986 FORMAT( ' ZDRGVX: ', A, ' Eigenvectors from ', A, ' incorrectly ',
618 $ 'normalized.', / ' Bits of error=', 0P, G10.3, ',', 9X,
619 $ 'N=', I6, ', Input Example #', I2, ')' )
620 *
621 * End of ZDRGVX
622 *
623 END
2 $ ALPHA, BETA, VL, VR, ILO, IHI, LSCALE, RSCALE,
3 $ S, DTRU, DIF, DIFTRU, WORK, LWORK, RWORK,
4 $ IWORK, LIWORK, RESULT, BWORK, INFO )
5 *
6 * -- LAPACK test routine (version 3.1) --
7 * Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
8 * November 2006
9 *
10 * .. Scalar Arguments ..
11 INTEGER IHI, ILO, INFO, LDA, LIWORK, LWORK, NIN, NOUT,
12 $ NSIZE
13 DOUBLE PRECISION THRESH
14 * ..
15 * .. Array Arguments ..
16 LOGICAL BWORK( * )
17 INTEGER IWORK( * )
18 DOUBLE PRECISION DIF( * ), DIFTRU( * ), DTRU( * ), LSCALE( * ),
19 $ RESULT( 4 ), RSCALE( * ), RWORK( * ), S( * )
20 COMPLEX*16 A( LDA, * ), AI( LDA, * ), ALPHA( * ),
21 $ B( LDA, * ), BETA( * ), BI( LDA, * ),
22 $ VL( LDA, * ), VR( LDA, * ), WORK( * )
23 * ..
24 *
25 * Purpose
26 * =======
27 *
28 * ZDRGVX checks the nonsymmetric generalized eigenvalue problem
29 * expert driver ZGGEVX.
30 *
31 * ZGGEVX computes the generalized eigenvalues, (optionally) the left
32 * and/or right eigenvectors, (optionally) computes a balancing
33 * transformation to improve the conditioning, and (optionally)
34 * reciprocal condition numbers for the eigenvalues and eigenvectors.
35 *
36 * When ZDRGVX is called with NSIZE > 0, two types of test matrix pairs
37 * are generated by the subroutine DLATM6 and test the driver ZGGEVX.
38 * The test matrices have the known exact condition numbers for
39 * eigenvalues. For the condition numbers of the eigenvectors
40 * corresponding the first and last eigenvalues are also know
41 * ``exactly'' (see ZLATM6).
42 * For each matrix pair, the following tests will be performed and
43 * compared with the threshhold THRESH.
44 *
45 * (1) max over all left eigenvalue/-vector pairs (beta/alpha,l) of
46 *
47 * | l**H * (beta A - alpha B) | / ( ulp max( |beta A|, |alpha B| ) )
48 *
49 * where l**H is the conjugate tranpose of l.
50 *
51 * (2) max over all right eigenvalue/-vector pairs (beta/alpha,r) of
52 *
53 * | (beta A - alpha B) r | / ( ulp max( |beta A|, |alpha B| ) )
54 *
55 * (3) The condition number S(i) of eigenvalues computed by ZGGEVX
56 * differs less than a factor THRESH from the exact S(i) (see
57 * ZLATM6).
58 *
59 * (4) DIF(i) computed by ZTGSNA differs less than a factor 10*THRESH
60 * from the exact value (for the 1st and 5th vectors only).
61 *
62 * Test Matrices
63 * =============
64 *
65 * Two kinds of test matrix pairs
66 * (A, B) = inverse(YH) * (Da, Db) * inverse(X)
67 * are used in the tests:
68 *
69 * 1: Da = 1+a 0 0 0 0 Db = 1 0 0 0 0
70 * 0 2+a 0 0 0 0 1 0 0 0
71 * 0 0 3+a 0 0 0 0 1 0 0
72 * 0 0 0 4+a 0 0 0 0 1 0
73 * 0 0 0 0 5+a , 0 0 0 0 1 , and
74 *
75 * 2: Da = 1 -1 0 0 0 Db = 1 0 0 0 0
76 * 1 1 0 0 0 0 1 0 0 0
77 * 0 0 1 0 0 0 0 1 0 0
78 * 0 0 0 1+a 1+b 0 0 0 1 0
79 * 0 0 0 -1-b 1+a , 0 0 0 0 1 .
80 *
81 * In both cases the same inverse(YH) and inverse(X) are used to compute
82 * (A, B), giving the exact eigenvectors to (A,B) as (YH, X):
83 *
84 * YH: = 1 0 -y y -y X = 1 0 -x -x x
85 * 0 1 -y y -y 0 1 x -x -x
86 * 0 0 1 0 0 0 0 1 0 0
87 * 0 0 0 1 0 0 0 0 1 0
88 * 0 0 0 0 1, 0 0 0 0 1 , where
89 *
90 * a, b, x and y will have all values independently of each other from
91 * { sqrt(sqrt(ULP)), 0.1, 1, 10, 1/sqrt(sqrt(ULP)) }.
92 *
93 * Arguments
94 * =========
95 *
96 * NSIZE (input) INTEGER
97 * The number of sizes of matrices to use. NSIZE must be at
98 * least zero. If it is zero, no randomly generated matrices
99 * are tested, but any test matrices read from NIN will be
100 * tested. If it is not zero, then N = 5.
101 *
102 * THRESH (input) DOUBLE PRECISION
103 * A test will count as "failed" if the "error", computed as
104 * described above, exceeds THRESH. Note that the error
105 * is scaled to be O(1), so THRESH should be a reasonably
106 * small multiple of 1, e.g., 10 or 100. In particular,
107 * it should not depend on the precision (single vs. double)
108 * or the size of the matrix. It must be at least zero.
109 *
110 * NIN (input) INTEGER
111 * The FORTRAN unit number for reading in the data file of
112 * problems to solve.
113 *
114 * NOUT (input) INTEGER
115 * The FORTRAN unit number for printing out error messages
116 * (e.g., if a routine returns IINFO not equal to 0.)
117 *
118 * A (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
119 * Used to hold the matrix whose eigenvalues are to be
120 * computed. On exit, A contains the last matrix actually used.
121 *
122 * LDA (input) INTEGER
123 * The leading dimension of A, B, AI, BI, Ao, and Bo.
124 * It must be at least 1 and at least NSIZE.
125 *
126 * B (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
127 * Used to hold the matrix whose eigenvalues are to be
128 * computed. On exit, B contains the last matrix actually used.
129 *
130 * AI (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
131 * Copy of A, modified by ZGGEVX.
132 *
133 * BI (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
134 * Copy of B, modified by ZGGEVX.
135 *
136 * ALPHA (workspace) COMPLEX*16 array, dimension (NSIZE)
137 * BETA (workspace) COMPLEX*16 array, dimension (NSIZE)
138 * On exit, ALPHA/BETA are the eigenvalues.
139 *
140 * VL (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
141 * VL holds the left eigenvectors computed by ZGGEVX.
142 *
143 * VR (workspace) COMPLEX*16 array, dimension (LDA, NSIZE)
144 * VR holds the right eigenvectors computed by ZGGEVX.
145 *
146 * ILO (output/workspace) INTEGER
147 *
148 * IHI (output/workspace) INTEGER
149 *
150 * LSCALE (output/workspace) DOUBLE PRECISION array, dimension (N)
151 *
152 * RSCALE (output/workspace) DOUBLE PRECISION array, dimension (N)
153 *
154 * S (output/workspace) DOUBLE PRECISION array, dimension (N)
155 *
156 * DTRU (output/workspace) DOUBLE PRECISION array, dimension (N)
157 *
158 * DIF (output/workspace) DOUBLE PRECISION array, dimension (N)
159 *
160 * DIFTRU (output/workspace) DOUBLE PRECISION array, dimension (N)
161 *
162 * WORK (workspace) COMPLEX*16 array, dimension (LWORK)
163 *
164 * LWORK (input) INTEGER
165 * Leading dimension of WORK. LWORK >= 2*N*N + 2*N
166 *
167 * RWORK (workspace) DOUBLE PRECISION array, dimension (6*N)
168 *
169 * IWORK (workspace) INTEGER array, dimension (LIWORK)
170 *
171 * LIWORK (input) INTEGER
172 * Leading dimension of IWORK. LIWORK >= N+2.
173 *
174 * RESULT (output/workspace) DOUBLE PRECISION array, dimension (4)
175 *
176 * BWORK (workspace) LOGICAL array, dimension (N)
177 *
178 * INFO (output) INTEGER
179 * = 0: successful exit
180 * < 0: if INFO = -i, the i-th argument had an illegal value.
181 * > 0: A routine returned an error code.
182 *
183 * =====================================================================
184 *
185 * .. Parameters ..
186 DOUBLE PRECISION ZERO, ONE, TEN, TNTH
187 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0, TEN = 1.0D+1,
188 $ TNTH = 1.0D-1 )
189 * ..
190 * .. Local Scalars ..
191 INTEGER I, IPTYPE, IWA, IWB, IWX, IWY, J, LINFO,
192 $ MAXWRK, MINWRK, N, NERRS, NMAX, NPTKNT, NTESTT
193 DOUBLE PRECISION ABNORM, ANORM, BNORM, RATIO1, RATIO2, THRSH2,
194 $ ULP, ULPINV
195 * ..
196 * .. Local Arrays ..
197 COMPLEX*16 WEIGHT( 5 )
198 * ..
199 * .. External Functions ..
200 INTEGER ILAENV
201 DOUBLE PRECISION DLAMCH, ZLANGE
202 EXTERNAL ILAENV, DLAMCH, ZLANGE
203 * ..
204 * .. External Subroutines ..
205 EXTERNAL ALASVM, XERBLA, ZGET52, ZGGEVX, ZLACPY, ZLATM6
206 * ..
207 * .. Intrinsic Functions ..
208 INTRINSIC ABS, DCMPLX, MAX, SQRT
209 * ..
210 * .. Executable Statements ..
211 *
212 * Check for errors
213 *
214 INFO = 0
215 *
216 NMAX = 5
217 *
218 IF( NSIZE.LT.0 ) THEN
219 INFO = -1
220 ELSE IF( THRESH.LT.ZERO ) THEN
221 INFO = -2
222 ELSE IF( NIN.LE.0 ) THEN
223 INFO = -3
224 ELSE IF( NOUT.LE.0 ) THEN
225 INFO = -4
226 ELSE IF( LDA.LT.1 .OR. LDA.LT.NMAX ) THEN
227 INFO = -6
228 ELSE IF( LIWORK.LT.NMAX+2 ) THEN
229 INFO = -26
230 END IF
231 *
232 * Compute workspace
233 * (Note: Comments in the code beginning "Workspace:" describe the
234 * minimal amount of workspace needed at that point in the code,
235 * as well as the preferred amount for good performance.
236 * NB refers to the optimal block size for the immediately
237 * following subroutine, as returned by ILAENV.)
238 *
239 MINWRK = 1
240 IF( INFO.EQ.0 .AND. LWORK.GE.1 ) THEN
241 MINWRK = 2*NMAX*( NMAX+1 )
242 MAXWRK = NMAX*( 1+ILAENV( 1, 'ZGEQRF', ' ', NMAX, 1, NMAX,
243 $ 0 ) )
244 MAXWRK = MAX( MAXWRK, 2*NMAX*( NMAX+1 ) )
245 WORK( 1 ) = MAXWRK
246 END IF
247 *
248 IF( LWORK.LT.MINWRK )
249 $ INFO = -23
250 *
251 IF( INFO.NE.0 ) THEN
252 CALL XERBLA( 'ZDRGVX', -INFO )
253 RETURN
254 END IF
255 *
256 N = 5
257 ULP = DLAMCH( 'P' )
258 ULPINV = ONE / ULP
259 THRSH2 = TEN*THRESH
260 NERRS = 0
261 NPTKNT = 0
262 NTESTT = 0
263 *
264 IF( NSIZE.EQ.0 )
265 $ GO TO 90
266 *
267 * Parameters used for generating test matrices.
268 *
269 WEIGHT( 1 ) = DCMPLX( SQRT( SQRT( ULP ) ), ZERO )
270 WEIGHT( 2 ) = DCMPLX( TNTH, ZERO )
271 WEIGHT( 3 ) = ONE
272 WEIGHT( 4 ) = ONE / WEIGHT( 2 )
273 WEIGHT( 5 ) = ONE / WEIGHT( 1 )
274 *
275 DO 80 IPTYPE = 1, 2
276 DO 70 IWA = 1, 5
277 DO 60 IWB = 1, 5
278 DO 50 IWX = 1, 5
279 DO 40 IWY = 1, 5
280 *
281 * generated a pair of test matrix
282 *
283 CALL ZLATM6( IPTYPE, 5, A, LDA, B, VR, LDA, VL,
284 $ LDA, WEIGHT( IWA ), WEIGHT( IWB ),
285 $ WEIGHT( IWX ), WEIGHT( IWY ), DTRU,
286 $ DIFTRU )
287 *
288 * Compute eigenvalues/eigenvectors of (A, B).
289 * Compute eigenvalue/eigenvector condition numbers
290 * using computed eigenvectors.
291 *
292 CALL ZLACPY( 'F', N, N, A, LDA, AI, LDA )
293 CALL ZLACPY( 'F', N, N, B, LDA, BI, LDA )
294 *
295 CALL ZGGEVX( 'N', 'V', 'V', 'B', N, AI, LDA, BI,
296 $ LDA, ALPHA, BETA, VL, LDA, VR, LDA,
297 $ ILO, IHI, LSCALE, RSCALE, ANORM,
298 $ BNORM, S, DIF, WORK, LWORK, RWORK,
299 $ IWORK, BWORK, LINFO )
300 IF( LINFO.NE.0 ) THEN
301 WRITE( NOUT, FMT = 9999 )'ZGGEVX', LINFO, N,
302 $ IPTYPE, IWA, IWB, IWX, IWY
303 GO TO 30
304 END IF
305 *
306 * Compute the norm(A, B)
307 *
308 CALL ZLACPY( 'Full', N, N, AI, LDA, WORK, N )
309 CALL ZLACPY( 'Full', N, N, BI, LDA, WORK( N*N+1 ),
310 $ N )
311 ABNORM = ZLANGE( 'Fro', N, 2*N, WORK, N, RWORK )
312 *
313 * Tests (1) and (2)
314 *
315 RESULT( 1 ) = ZERO
316 CALL ZGET52( .TRUE., N, A, LDA, B, LDA, VL, LDA,
317 $ ALPHA, BETA, WORK, RWORK,
318 $ RESULT( 1 ) )
319 IF( RESULT( 2 ).GT.THRESH ) THEN
320 WRITE( NOUT, FMT = 9998 )'Left', 'ZGGEVX',
321 $ RESULT( 2 ), N, IPTYPE, IWA, IWB, IWX, IWY
322 END IF
323 *
324 RESULT( 2 ) = ZERO
325 CALL ZGET52( .FALSE., N, A, LDA, B, LDA, VR, LDA,
326 $ ALPHA, BETA, WORK, RWORK,
327 $ RESULT( 2 ) )
328 IF( RESULT( 3 ).GT.THRESH ) THEN
329 WRITE( NOUT, FMT = 9998 )'Right', 'ZGGEVX',
330 $ RESULT( 3 ), N, IPTYPE, IWA, IWB, IWX, IWY
331 END IF
332 *
333 * Test (3)
334 *
335 RESULT( 3 ) = ZERO
336 DO 10 I = 1, N
337 IF( S( I ).EQ.ZERO ) THEN
338 IF( DTRU( I ).GT.ABNORM*ULP )
339 $ RESULT( 3 ) = ULPINV
340 ELSE IF( DTRU( I ).EQ.ZERO ) THEN
341 IF( S( I ).GT.ABNORM*ULP )
342 $ RESULT( 3 ) = ULPINV
343 ELSE
344 RWORK( I ) = MAX( ABS( DTRU( I ) / S( I ) ),
345 $ ABS( S( I ) / DTRU( I ) ) )
346 RESULT( 3 ) = MAX( RESULT( 3 ), RWORK( I ) )
347 END IF
348 10 CONTINUE
349 *
350 * Test (4)
351 *
352 RESULT( 4 ) = ZERO
353 IF( DIF( 1 ).EQ.ZERO ) THEN
354 IF( DIFTRU( 1 ).GT.ABNORM*ULP )
355 $ RESULT( 4 ) = ULPINV
356 ELSE IF( DIFTRU( 1 ).EQ.ZERO ) THEN
357 IF( DIF( 1 ).GT.ABNORM*ULP )
358 $ RESULT( 4 ) = ULPINV
359 ELSE IF( DIF( 5 ).EQ.ZERO ) THEN
360 IF( DIFTRU( 5 ).GT.ABNORM*ULP )
361 $ RESULT( 4 ) = ULPINV
362 ELSE IF( DIFTRU( 5 ).EQ.ZERO ) THEN
363 IF( DIF( 5 ).GT.ABNORM*ULP )
364 $ RESULT( 4 ) = ULPINV
365 ELSE
366 RATIO1 = MAX( ABS( DIFTRU( 1 ) / DIF( 1 ) ),
367 $ ABS( DIF( 1 ) / DIFTRU( 1 ) ) )
368 RATIO2 = MAX( ABS( DIFTRU( 5 ) / DIF( 5 ) ),
369 $ ABS( DIF( 5 ) / DIFTRU( 5 ) ) )
370 RESULT( 4 ) = MAX( RATIO1, RATIO2 )
371 END IF
372 *
373 NTESTT = NTESTT + 4
374 *
375 * Print out tests which fail.
376 *
377 DO 20 J = 1, 4
378 IF( ( RESULT( J ).GE.THRSH2 .AND. J.GE.4 ) .OR.
379 $ ( RESULT( J ).GE.THRESH .AND. J.LE.3 ) )
380 $ THEN
381 *
382 * If this is the first test to fail,
383 * print a header to the data file.
384 *
385 IF( NERRS.EQ.0 ) THEN
386 WRITE( NOUT, FMT = 9997 )'ZXV'
387 *
388 * Print out messages for built-in examples
389 *
390 * Matrix types
391 *
392 WRITE( NOUT, FMT = 9995 )
393 WRITE( NOUT, FMT = 9994 )
394 WRITE( NOUT, FMT = 9993 )
395 *
396 * Tests performed
397 *
398 WRITE( NOUT, FMT = 9992 )'''',
399 $ 'transpose', ''''
400 *
401 END IF
402 NERRS = NERRS + 1
403 IF( RESULT( J ).LT.10000.0D0 ) THEN
404 WRITE( NOUT, FMT = 9991 )IPTYPE, IWA,
405 $ IWB, IWX, IWY, J, RESULT( J )
406 ELSE
407 WRITE( NOUT, FMT = 9990 )IPTYPE, IWA,
408 $ IWB, IWX, IWY, J, RESULT( J )
409 END IF
410 END IF
411 20 CONTINUE
412 *
413 30 CONTINUE
414 *
415 40 CONTINUE
416 50 CONTINUE
417 60 CONTINUE
418 70 CONTINUE
419 80 CONTINUE
420 *
421 GO TO 150
422 *
423 90 CONTINUE
424 *
425 * Read in data from file to check accuracy of condition estimation
426 * Read input data until N=0
427 *
428 READ( NIN, FMT = *, END = 150 )N
429 IF( N.EQ.0 )
430 $ GO TO 150
431 DO 100 I = 1, N
432 READ( NIN, FMT = * )( A( I, J ), J = 1, N )
433 100 CONTINUE
434 DO 110 I = 1, N
435 READ( NIN, FMT = * )( B( I, J ), J = 1, N )
436 110 CONTINUE
437 READ( NIN, FMT = * )( DTRU( I ), I = 1, N )
438 READ( NIN, FMT = * )( DIFTRU( I ), I = 1, N )
439 *
440 NPTKNT = NPTKNT + 1
441 *
442 * Compute eigenvalues/eigenvectors of (A, B).
443 * Compute eigenvalue/eigenvector condition numbers
444 * using computed eigenvectors.
445 *
446 CALL ZLACPY( 'F', N, N, A, LDA, AI, LDA )
447 CALL ZLACPY( 'F', N, N, B, LDA, BI, LDA )
448 *
449 CALL ZGGEVX( 'N', 'V', 'V', 'B', N, AI, LDA, BI, LDA, ALPHA, BETA,
450 $ VL, LDA, VR, LDA, ILO, IHI, LSCALE, RSCALE, ANORM,
451 $ BNORM, S, DIF, WORK, LWORK, RWORK, IWORK, BWORK,
452 $ LINFO )
453 *
454 IF( LINFO.NE.0 ) THEN
455 WRITE( NOUT, FMT = 9987 )'ZGGEVX', LINFO, N, NPTKNT
456 GO TO 140
457 END IF
458 *
459 * Compute the norm(A, B)
460 *
461 CALL ZLACPY( 'Full', N, N, AI, LDA, WORK, N )
462 CALL ZLACPY( 'Full', N, N, BI, LDA, WORK( N*N+1 ), N )
463 ABNORM = ZLANGE( 'Fro', N, 2*N, WORK, N, RWORK )
464 *
465 * Tests (1) and (2)
466 *
467 RESULT( 1 ) = ZERO
468 CALL ZGET52( .TRUE., N, A, LDA, B, LDA, VL, LDA, ALPHA, BETA,
469 $ WORK, RWORK, RESULT( 1 ) )
470 IF( RESULT( 2 ).GT.THRESH ) THEN
471 WRITE( NOUT, FMT = 9986 )'Left', 'ZGGEVX', RESULT( 2 ), N,
472 $ NPTKNT
473 END IF
474 *
475 RESULT( 2 ) = ZERO
476 CALL ZGET52( .FALSE., N, A, LDA, B, LDA, VR, LDA, ALPHA, BETA,
477 $ WORK, RWORK, RESULT( 2 ) )
478 IF( RESULT( 3 ).GT.THRESH ) THEN
479 WRITE( NOUT, FMT = 9986 )'Right', 'ZGGEVX', RESULT( 3 ), N,
480 $ NPTKNT
481 END IF
482 *
483 * Test (3)
484 *
485 RESULT( 3 ) = ZERO
486 DO 120 I = 1, N
487 IF( S( I ).EQ.ZERO ) THEN
488 IF( DTRU( I ).GT.ABNORM*ULP )
489 $ RESULT( 3 ) = ULPINV
490 ELSE IF( DTRU( I ).EQ.ZERO ) THEN
491 IF( S( I ).GT.ABNORM*ULP )
492 $ RESULT( 3 ) = ULPINV
493 ELSE
494 RWORK( I ) = MAX( ABS( DTRU( I ) / S( I ) ),
495 $ ABS( S( I ) / DTRU( I ) ) )
496 RESULT( 3 ) = MAX( RESULT( 3 ), RWORK( I ) )
497 END IF
498 120 CONTINUE
499 *
500 * Test (4)
501 *
502 RESULT( 4 ) = ZERO
503 IF( DIF( 1 ).EQ.ZERO ) THEN
504 IF( DIFTRU( 1 ).GT.ABNORM*ULP )
505 $ RESULT( 4 ) = ULPINV
506 ELSE IF( DIFTRU( 1 ).EQ.ZERO ) THEN
507 IF( DIF( 1 ).GT.ABNORM*ULP )
508 $ RESULT( 4 ) = ULPINV
509 ELSE IF( DIF( 5 ).EQ.ZERO ) THEN
510 IF( DIFTRU( 5 ).GT.ABNORM*ULP )
511 $ RESULT( 4 ) = ULPINV
512 ELSE IF( DIFTRU( 5 ).EQ.ZERO ) THEN
513 IF( DIF( 5 ).GT.ABNORM*ULP )
514 $ RESULT( 4 ) = ULPINV
515 ELSE
516 RATIO1 = MAX( ABS( DIFTRU( 1 ) / DIF( 1 ) ),
517 $ ABS( DIF( 1 ) / DIFTRU( 1 ) ) )
518 RATIO2 = MAX( ABS( DIFTRU( 5 ) / DIF( 5 ) ),
519 $ ABS( DIF( 5 ) / DIFTRU( 5 ) ) )
520 RESULT( 4 ) = MAX( RATIO1, RATIO2 )
521 END IF
522 *
523 NTESTT = NTESTT + 4
524 *
525 * Print out tests which fail.
526 *
527 DO 130 J = 1, 4
528 IF( RESULT( J ).GE.THRSH2 ) THEN
529 *
530 * If this is the first test to fail,
531 * print a header to the data file.
532 *
533 IF( NERRS.EQ.0 ) THEN
534 WRITE( NOUT, FMT = 9997 )'ZXV'
535 *
536 * Print out messages for built-in examples
537 *
538 * Matrix types
539 *
540 WRITE( NOUT, FMT = 9996 )
541 *
542 * Tests performed
543 *
544 WRITE( NOUT, FMT = 9992 )'''', 'transpose', ''''
545 *
546 END IF
547 NERRS = NERRS + 1
548 IF( RESULT( J ).LT.10000.0D0 ) THEN
549 WRITE( NOUT, FMT = 9989 )NPTKNT, N, J, RESULT( J )
550 ELSE
551 WRITE( NOUT, FMT = 9988 )NPTKNT, N, J, RESULT( J )
552 END IF
553 END IF
554 130 CONTINUE
555 *
556 140 CONTINUE
557 *
558 GO TO 90
559 150 CONTINUE
560 *
561 * Summary
562 *
563 CALL ALASVM( 'ZXV', NOUT, NERRS, NTESTT, 0 )
564 *
565 WORK( 1 ) = MAXWRK
566 *
567 RETURN
568 *
569 9999 FORMAT( ' ZDRGVX: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
570 $ I6, ', JTYPE=', I6, ')' )
571 *
572 9998 FORMAT( ' ZDRGVX: ', A, ' Eigenvectors from ', A, ' incorrectly ',
573 $ 'normalized.', / ' Bits of error=', 0P, G10.3, ',', 9X,
574 $ 'N=', I6, ', JTYPE=', I6, ', IWA=', I5, ', IWB=', I5,
575 $ ', IWX=', I5, ', IWY=', I5 )
576 *
577 9997 FORMAT( / 1X, A3, ' -- Complex Expert Eigenvalue/vector',
578 $ ' problem driver' )
579 *
580 9996 FORMAT( 'Input Example' )
581 *
582 9995 FORMAT( ' Matrix types: ', / )
583 *
584 9994 FORMAT( ' TYPE 1: Da is diagonal, Db is identity, ',
585 $ / ' A = Y^(-H) Da X^(-1), B = Y^(-H) Db X^(-1) ',
586 $ / ' YH and X are left and right eigenvectors. ', / )
587 *
588 9993 FORMAT( ' TYPE 2: Da is quasi-diagonal, Db is identity, ',
589 $ / ' A = Y^(-H) Da X^(-1), B = Y^(-H) Db X^(-1) ',
590 $ / ' YH and X are left and right eigenvectors. ', / )
591 *
592 9992 FORMAT( / ' Tests performed: ', / 4X,
593 $ ' a is alpha, b is beta, l is a left eigenvector, ', / 4X,
594 $ ' r is a right eigenvector and ', A, ' means ', A, '.',
595 $ / ' 1 = max | ( b A - a B )', A, ' l | / const.',
596 $ / ' 2 = max | ( b A - a B ) r | / const.',
597 $ / ' 3 = max ( Sest/Stru, Stru/Sest ) ',
598 $ ' over all eigenvalues', /
599 $ ' 4 = max( DIFest/DIFtru, DIFtru/DIFest ) ',
600 $ ' over the 1st and 5th eigenvectors', / )
601 *
602 9991 FORMAT( ' Type=', I2, ',', ' IWA=', I2, ', IWB=', I2, ', IWX=',
603 $ I2, ', IWY=', I2, ', result ', I2, ' is', 0P, F8.2 )
604 *
605 9990 FORMAT( ' Type=', I2, ',', ' IWA=', I2, ', IWB=', I2, ', IWX=',
606 $ I2, ', IWY=', I2, ', result ', I2, ' is', 1P, D10.3 )
607 *
608 9989 FORMAT( ' Input example #', I2, ', matrix order=', I4, ',',
609 $ ' result ', I2, ' is', 0P, F8.2 )
610 *
611 9988 FORMAT( ' Input example #', I2, ', matrix order=', I4, ',',
612 $ ' result ', I2, ' is', 1P, D10.3 )
613 *
614 9987 FORMAT( ' ZDRGVX: ', A, ' returned INFO=', I6, '.', / 9X, 'N=',
615 $ I6, ', Input example #', I2, ')' )
616 *
617 9986 FORMAT( ' ZDRGVX: ', A, ' Eigenvectors from ', A, ' incorrectly ',
618 $ 'normalized.', / ' Bits of error=', 0P, G10.3, ',', 9X,
619 $ 'N=', I6, ', Input Example #', I2, ')' )
620 *
621 * End of ZDRGVX
622 *
623 END