OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
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*> \ingroup complex_eig
292*
293* =====================================================================
294 SUBROUTINE cdrgvx( NSIZE, THRESH, NIN, NOUT, A, LDA, B, AI, BI,
295 $ ALPHA, BETA, VL, VR, ILO, IHI, LSCALE, RSCALE,
296 $ S, STRU, DIF, DIFTRU, WORK, LWORK, RWORK,
297 $ IWORK, LIWORK, RESULT, BWORK, INFO )
298*
299* -- LAPACK test routine --
300* -- LAPACK is a software package provided by Univ. of Tennessee, --
301* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
302*
303* .. Scalar Arguments ..
304 INTEGER IHI, ILO, INFO, LDA, LIWORK, LWORK, NIN, NOUT,
305 $ NSIZE
306 REAL THRESH
307* ..
308* .. Array Arguments ..
309 LOGICAL BWORK( * )
310 INTEGER IWORK( * )
311 REAL DIF( * ), DIFTRU( * ), LSCALE( * ),
312 $ result( 4 ), rscale( * ), rwork( * ), s( * ),
313 $ stru( * )
314 COMPLEX A( LDA, * ), AI( LDA, * ), ALPHA( * ),
315 $ B( LDA, * ), BETA( * ), BI( LDA, * ),
316 $ VL( LDA, * ), VR( LDA, * ), WORK( * )
317* ..
318*
319* =====================================================================
320*
321* .. Parameters ..
322 REAL ZERO, ONE, TEN, TNTH, HALF
323 PARAMETER ( ZERO = 0.0e+0, one = 1.0e+0, ten = 1.0e+1,
324 $ tnth = 1.0e-1, half = 0.5e+0 )
325* ..
326* .. Local Scalars ..
327 INTEGER I, IPTYPE, IWA, IWB, IWX, IWY, J, LINFO,
328 $ MAXWRK, MINWRK, N, NERRS, NMAX, NPTKNT, NTESTT
329 REAL ABNORM, ANORM, BNORM, RATIO1, RATIO2, THRSH2,
330 $ ulp, ulpinv
331* ..
332* .. Local Arrays ..
333 COMPLEX WEIGHT( 5 )
334* ..
335* .. External Functions ..
336 INTEGER ILAENV
337 REAL CLANGE, SLAMCH
338 EXTERNAL ILAENV, CLANGE, SLAMCH
339* ..
340* .. External Subroutines ..
341 EXTERNAL alasvm, cget52, cggevx, clacpy, clatm6, xerbla
342* ..
343* .. Intrinsic Functions ..
344 INTRINSIC abs, cmplx, max, sqrt
345* ..
346* .. Executable Statements ..
347*
348* Check for errors
349*
350 info = 0
351*
352 nmax = 5
353*
354 IF( nsize.LT.0 ) THEN
355 info = -1
356 ELSE IF( thresh.LT.zero ) THEN
357 info = -2
358 ELSE IF( nin.LE.0 ) THEN
359 info = -3
360 ELSE IF( nout.LE.0 ) THEN
361 info = -4
362 ELSE IF( lda.LT.1 .OR. lda.LT.nmax ) THEN
363 info = -6
364 ELSE IF( liwork.LT.nmax+2 ) THEN
365 info = -26
366 END IF
367*
368* Compute workspace
369* (Note: Comments in the code beginning "Workspace:" describe the
370* minimal amount of workspace needed at that point in the code,
371* as well as the preferred amount for good performance.
372* NB refers to the optimal block size for the immediately
373* following subroutine, as returned by ILAENV.)
374*
375 minwrk = 1
376 IF( info.EQ.0 .AND. lwork.GE.1 ) THEN
377 minwrk = 2*nmax*( nmax+1 )
378 maxwrk = nmax*( 1+ilaenv( 1, 'CGEQRF', ' ', nmax, 1, nmax,
379 $ 0 ) )
380 maxwrk = max( maxwrk, 2*nmax*( nmax+1 ) )
381 work( 1 ) = maxwrk
382 END IF
383*
384 IF( lwork.LT.minwrk )
385 $ info = -23
386*
387 IF( info.NE.0 ) THEN
388 CALL xerbla( 'CDRGVX', -info )
389 RETURN
390 END IF
391*
392 n = 5
393 ulp = slamch( 'P' )
394 ulpinv = one / ulp
395 thrsh2 = ten*thresh
396 nerrs = 0
397 nptknt = 0
398 ntestt = 0
399*
400 IF( nsize.EQ.0 )
401 $ GO TO 90
402*
403* Parameters used for generating test matrices.
404*
405 weight( 1 ) = cmplx( tnth, zero )
406 weight( 2 ) = cmplx( half, zero )
407 weight( 3 ) = one
408 weight( 4 ) = one / weight( 2 )
409 weight( 5 ) = one / weight( 1 )
410*
411 DO 80 iptype = 1, 2
412 DO 70 iwa = 1, 5
413 DO 60 iwb = 1, 5
414 DO 50 iwx = 1, 5
415 DO 40 iwy = 1, 5
416*
417* generated a pair of test matrix
418*
419 CALL clatm6( iptype, 5, a, lda, b, vr, lda, vl,
420 $ lda, weight( iwa ), weight( iwb ),
421 $ weight( iwx ), weight( iwy ), stru,
422 $ diftru )
423*
424* Compute eigenvalues/eigenvectors of (A, B).
425* Compute eigenvalue/eigenvector condition numbers
426* using computed eigenvectors.
427*
428 CALL clacpy( 'F', n, n, a, lda, ai, lda )
429 CALL clacpy( 'F', n, n, b, lda, bi, lda )
430*
431 CALL cggevx( 'N', 'V', 'V', 'B', n, ai, lda, bi,
432 $ lda, alpha, beta, vl, lda, vr, lda,
433 $ ilo, ihi, lscale, rscale, anorm,
434 $ bnorm, s, dif, work, lwork, rwork,
435 $ iwork, bwork, linfo )
436 IF( linfo.NE.0 ) THEN
437 WRITE( nout, fmt = 9999 )'CGGEVX', linfo, n,
438 $ iptype, iwa, iwb, iwx, iwy
439 GO TO 30
440 END IF
441*
442* Compute the norm(A, B)
443*
444 CALL clacpy( 'Full', n, n, ai, lda, work, n )
445 CALL clacpy( 'Full', n, n, bi, lda, work( n*n+1 ),
446 $ n )
447 abnorm = clange( 'Fro', n, 2*n, work, n, rwork )
448*
449* Tests (1) and (2)
450*
451 result( 1 ) = zero
452 CALL cget52( .true., n, a, lda, b, lda, vl, lda,
453 $ alpha, beta, work, rwork,
454 $ result( 1 ) )
455 IF( result( 2 ).GT.thresh ) THEN
456 WRITE( nout, fmt = 9998 )'Left', 'CGGEVX',
457 $ result( 2 ), n, iptype, iwa, iwb, iwx, iwy
458 END IF
459*
460 result( 2 ) = zero
461 CALL cget52( .false., n, a, lda, b, lda, vr, lda,
462 $ alpha, beta, work, rwork,
463 $ result( 2 ) )
464 IF( result( 3 ).GT.thresh ) THEN
465 WRITE( nout, fmt = 9998 )'Right', 'CGGEVX',
466 $ result( 3 ), n, iptype, iwa, iwb, iwx, iwy
467 END IF
468*
469* Test (3)
470*
471 result( 3 ) = zero
472 DO 10 i = 1, n
473 IF( s( i ).EQ.zero ) THEN
474 IF( stru( i ).GT.abnorm*ulp )
475 $ result( 3 ) = ulpinv
476 ELSE IF( stru( i ).EQ.zero ) THEN
477 IF( s( i ).GT.abnorm*ulp )
478 $ result( 3 ) = ulpinv
479 ELSE
480 rwork( i ) = max( abs( stru( i ) / s( i ) ),
481 $ abs( s( i ) / stru( i ) ) )
482 result( 3 ) = max( result( 3 ), rwork( i ) )
483 END IF
484 10 CONTINUE
485*
486* Test (4)
487*
488 result( 4 ) = zero
489 IF( dif( 1 ).EQ.zero ) THEN
490 IF( diftru( 1 ).GT.abnorm*ulp )
491 $ result( 4 ) = ulpinv
492 ELSE IF( diftru( 1 ).EQ.zero ) THEN
493 IF( dif( 1 ).GT.abnorm*ulp )
494 $ result( 4 ) = ulpinv
495 ELSE IF( dif( 5 ).EQ.zero ) THEN
496 IF( diftru( 5 ).GT.abnorm*ulp )
497 $ result( 4 ) = ulpinv
498 ELSE IF( diftru( 5 ).EQ.zero ) THEN
499 IF( dif( 5 ).GT.abnorm*ulp )
500 $ result( 4 ) = ulpinv
501 ELSE
502 ratio1 = max( abs( diftru( 1 ) / dif( 1 ) ),
503 $ abs( dif( 1 ) / diftru( 1 ) ) )
504 ratio2 = max( abs( diftru( 5 ) / dif( 5 ) ),
505 $ abs( dif( 5 ) / diftru( 5 ) ) )
506 result( 4 ) = max( ratio1, ratio2 )
507 END IF
508*
509 ntestt = ntestt + 4
510*
511* Print out tests which fail.
512*
513 DO 20 j = 1, 4
514 IF( ( result( j ).GE.thrsh2 .AND. j.GE.4 ) .OR.
515 $ ( result( j ).GE.thresh .AND. j.LE.3 ) )
516 $ THEN
517*
518* If this is the first test to fail,
519* print a header to the data file.
520*
521 IF( nerrs.EQ.0 ) THEN
522 WRITE( nout, fmt = 9997 )'CXV'
523*
524* Print out messages for built-in examples
525*
526* Matrix types
527*
528 WRITE( nout, fmt = 9995 )
529 WRITE( nout, fmt = 9994 )
530 WRITE( nout, fmt = 9993 )
531*
532* Tests performed
533*
534 WRITE( nout, fmt = 9992 )'''',
535 $ 'transpose', ''''
536*
537 END IF
538 NERRS = NERRS + 1
539.LT. IF( RESULT( J )10000.0 ) THEN
540 WRITE( NOUT, FMT = 9991 )IPTYPE, IWA,
541 $ IWB, IWX, IWY, J, RESULT( J )
542 ELSE
543 WRITE( NOUT, FMT = 9990 )IPTYPE, IWA,
544 $ IWB, IWX, IWY, J, RESULT( J )
545 END IF
546 END IF
547 20 CONTINUE
548*
549 30 CONTINUE
550*
551 40 CONTINUE
552 50 CONTINUE
553 60 CONTINUE
554 70 CONTINUE
555 80 CONTINUE
556*
557 GO TO 150
558*
559 90 CONTINUE
560*
561* Read in data from file to check accuracy of condition estimation
562* Read input data until N=0
563*
564 READ( NIN, FMT = *, END = 150 )N
565.EQ. IF( N0 )
566 $ GO TO 150
567 DO 100 I = 1, N
568 READ( NIN, FMT = * )( A( I, J ), J = 1, N )
569 100 CONTINUE
570 DO 110 I = 1, N
571 READ( NIN, FMT = * )( B( I, J ), J = 1, N )
572 110 CONTINUE
573 READ( NIN, FMT = * )( STRU( I ), I = 1, N )
574 READ( NIN, FMT = * )( DIFTRU( I ), I = 1, N )
575*
576 NPTKNT = NPTKNT + 1
577*
578* Compute eigenvalues/eigenvectors of (A, B).
579* Compute eigenvalue/eigenvector condition numbers
580* using computed eigenvectors.
581*
582 CALL CLACPY( 'f', N, N, A, LDA, AI, LDA )
583 CALL CLACPY( 'f', N, N, B, LDA, BI, LDA )
584*
585 CALL CGGEVX( 'n', 'v', 'v', 'b', N, AI, LDA, BI, LDA, ALPHA, BETA,
586 $ VL, LDA, VR, LDA, ILO, IHI, LSCALE, RSCALE, ANORM,
587 $ BNORM, S, DIF, WORK, LWORK, RWORK, IWORK, BWORK,
588 $ LINFO )
589*
590.NE. IF( LINFO0 ) THEN
591 WRITE( NOUT, FMT = 9987 )'cggevx', LINFO, N, NPTKNT
592 GO TO 140
593 END IF
594*
595* Compute the norm(A, B)
596*
597 CALL CLACPY( 'full', N, N, AI, LDA, WORK, N )
598 CALL CLACPY( 'full', N, N, BI, LDA, WORK( N*N+1 ), N )
599 ABNORM = CLANGE( 'fro', N, 2*N, WORK, N, RWORK )
600*
601* Tests (1) and (2)
602*
603 RESULT( 1 ) = ZERO
604 CALL CGET52( .TRUE., N, A, LDA, B, LDA, VL, LDA, ALPHA, BETA,
605 $ WORK, RWORK, RESULT( 1 ) )
606.GT. IF( RESULT( 2 )THRESH ) THEN
607 WRITE( NOUT, FMT = 9986 )'left', 'cggevx', RESULT( 2 ), N,
608 $ NPTKNT
609 END IF
610*
611 RESULT( 2 ) = ZERO
612 CALL CGET52( .FALSE., N, A, LDA, B, LDA, VR, LDA, ALPHA, BETA,
613 $ WORK, RWORK, RESULT( 2 ) )
614.GT. IF( RESULT( 3 )THRESH ) THEN
615 WRITE( NOUT, FMT = 9986 )'right', 'cggevx', RESULT( 3 ), N,
616 $ NPTKNT
617 END IF
618*
619* Test (3)
620*
621 RESULT( 3 ) = ZERO
622 DO 120 I = 1, N
623.EQ. IF( S( I )ZERO ) THEN
624.GT. IF( STRU( I )ABNORM*ULP )
625 $ RESULT( 3 ) = ULPINV
626.EQ. ELSE IF( STRU( I )ZERO ) THEN
627.GT. IF( S( I )ABNORM*ULP )
628 $ RESULT( 3 ) = ULPINV
629 ELSE
630 RWORK( I ) = MAX( ABS( STRU( I ) / S( I ) ),
631 $ ABS( S( I ) / STRU( I ) ) )
632 RESULT( 3 ) = MAX( RESULT( 3 ), RWORK( I ) )
633 END IF
634 120 CONTINUE
635*
636* Test (4)
637*
638 RESULT( 4 ) = ZERO
639.EQ. IF( DIF( 1 )ZERO ) THEN
640.GT. IF( DIFTRU( 1 )ABNORM*ULP )
641 $ RESULT( 4 ) = ULPINV
642.EQ. ELSE IF( DIFTRU( 1 )ZERO ) THEN
643.GT. IF( DIF( 1 )ABNORM*ULP )
644 $ RESULT( 4 ) = ULPINV
645.EQ. ELSE IF( DIF( 5 )ZERO ) THEN
646.GT. IF( DIFTRU( 5 )ABNORM*ULP )
647 $ RESULT( 4 ) = ULPINV
648.EQ. ELSE IF( DIFTRU( 5 )ZERO ) THEN
649.GT. IF( DIF( 5 )ABNORM*ULP )
650 $ RESULT( 4 ) = ULPINV
651 ELSE
652 RATIO1 = MAX( ABS( DIFTRU( 1 ) / DIF( 1 ) ),
653 $ ABS( DIF( 1 ) / DIFTRU( 1 ) ) )
654 RATIO2 = MAX( ABS( DIFTRU( 5 ) / DIF( 5 ) ),
655 $ ABS( DIF( 5 ) / DIFTRU( 5 ) ) )
656 RESULT( 4 ) = MAX( RATIO1, RATIO2 )
657 END IF
658*
659 NTESTT = NTESTT + 4
660*
661* Print out tests which fail.
662*
663 DO 130 J = 1, 4
664.GE. IF( RESULT( J )THRSH2 ) THEN
665*
666* If this is the first test to fail,
667* print a header to the data file.
668*
669.EQ. IF( NERRS0 ) THEN
670 WRITE( NOUT, FMT = 9997 )'cxv'
671*
672* Print out messages for built-in examples
673*
674* Matrix types
675*
676 WRITE( NOUT, FMT = 9996 )
677*
678* Tests performed
679*
680 WRITE( NOUT, FMT = 9992 )'''', 'transpose', ''''
681*
682 END IF
683 NERRS = NERRS + 1
684.LT. IF( RESULT( J )10000.0 ) THEN
685 WRITE( NOUT, FMT = 9989 )NPTKNT, N, J, RESULT( J )
686 ELSE
687 WRITE( NOUT, FMT = 9988 )NPTKNT, N, J, RESULT( J )
688 END IF
689 END IF
690 130 CONTINUE
691*
692 140 CONTINUE
693*
694 GO TO 90
695 150 CONTINUE
696*
697* Summary
698*
699 CALL ALASVM( 'cxv', NOUT, NERRS, NTESTT, 0 )
700*
701 WORK( 1 ) = MAXWRK
702*
703 RETURN
704*
705 9999 FORMAT( ' cdrgvx: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
706 $ i6, ', JTYPE=', i6, ')' )
707*
708 9998 FORMAT( ' CDRGVX: ', a, ' Eigenvectors from ', a, ' incorrectly ',
709 $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
710 $ 'N=', i6, ', JTYPE=', i6, ', IWA=', i5, ', IWB=', i5,
711 $ ', IWX=', i5, ', IWY=', i5 )
712*
713 9997 FORMAT( / 1x, a3, ' -- Complex Expert Eigenvalue/vector',
714 $ ' problem driver' )
715*
716 9996 FORMAT( 'Input Example' )
717*
718 9995 FORMAT( ' Matrix types: ', / )
719*
720 9994 FORMAT( ' TYPE 1: Da is diagonal, Db is identity, ',
721 $ / ' a = y^(-h) da x^(-1), b = y^(-h) db x^(-1) ',
722 $ / ' yh and x are left and right eigenvectors. ', / )
723*
724 9993 FORMAT( ' TYPE 2: da is quasi-diagonal, db is identity, ',
725 $ / ' a = y^(-h) da x^(-1), b = y^(-h) db x^(-1) ',
726 $ / ' yh and x are left and right eigenvectors. ', / )
727*
728 9992 FORMAT( / ' tests performed: ', / 4X,
729 $ ' a is alpha, b is beta, l is a left eigenvector, ', / 4X,
730 $ ' r is a right eigenvector and ', A, ' means ', A, '.',
731 $ / ' 1 = max | ( b a - a b )', A, ' l | / const.',
732 $ / ' 2 = max | ( b a - a b ) r | / const.',
733 $ / ' 3 = max( sest/stru, stru/sest ) ',
734 $ ' over all eigenvalues', /
735 $ ' 4 = max( difest/diftru, diftru/difest ) ',
736 $ ' over the 1st and 5th eigenvectors', / )
737*
738 9991 FORMAT( ' type=', I2, ',', ' iwa=', I2, ', iwb=', I2, ', iwx=',
739 $ I2, ', iwy=', I2, ', result ', I2, ' is', 0P, F8.2 )
740*
741 9990 FORMAT( ' type=', I2, ',', ' iwa=', I2, ', iwb=', I2, ', iwx=',
742 $ I2, ', iwy=', I2, ', result ', I2, ' is', 1P, E10.3 )
743*
744 9989 FORMAT( ' input example #', I2, ', matrix order=', I4, ',',
745 $ ' result ', i2, ' is', 0p, f8.2 )
746*
747 9988 FORMAT( ' Input example #', i2, ', matrix order=', i4, ',',
748 $ ' result ', i2, ' is', 1p, e10.3 )
749*
750 9987 FORMAT( ' CDRGVX: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
751 $ i6, ', Input example #', i2, ')' )
752*
753 9986 FORMAT( ' CDRGVX: ', a, ' Eigenvectors from ', a, ' incorrectly ',
754 $ 'normalized.', / ' Bits of error=', 0p, g10.3, ',', 9x,
755 $ 'N=', i6, ', Input Example #', i2, ')' )
756*
757* End of CDRGVX
758*
759 END
cmplx
float cmplx[2]
Definition pblas.h:136
the
end diagonal values have been computed in the(sparse) matrix id.SOL
xerbla
subroutine xerbla(srname, info)
XERBLA
Definition xerbla.f:60
alasvm
subroutine alasvm(type, nout, nfail, nrun, nerrs)
ALASVM
Definition alasvm.f:73
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:374
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:103
cget52
subroutine cget52(left, n, a, lda, b, ldb, e, lde, alpha, beta, work, rwork, result)
CGET52
Definition cget52.f:161
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:298
clatm6
subroutine clatm6(type, n, a, lda, b, x, ldx, y, ldy, alpha, beta, wx, wy, s, dif)
CLATM6
Definition clatm6.f:174
max
#define max(a, b)
Definition macros.h:21