LAPACK  3.10.1
LAPACK: Linear Algebra PACKage
cget23.f
Go to the documentation of this file.
1 *> \brief \b CGET23
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 CGET23( COMP, ISRT, BALANC, JTYPE, THRESH, ISEED,
12 * NOUNIT, N, A, LDA, H, W, W1, VL, LDVL, VR,
13 * LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN,
14 * RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT,
15 * WORK, LWORK, RWORK, INFO )
16 *
17 * .. Scalar Arguments ..
18 * LOGICAL COMP
19 * CHARACTER BALANC
20 * INTEGER INFO, ISRT, JTYPE, LDA, LDLRE, LDVL, LDVR,
21 * $ LWORK, N, NOUNIT
22 * REAL THRESH
23 * ..
24 * .. Array Arguments ..
25 * INTEGER ISEED( 4 )
26 * REAL RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
27 * $ RCNDV1( * ), RCONDE( * ), RCONDV( * ),
28 * $ RESULT( 11 ), RWORK( * ), SCALE( * ),
29 * $ SCALE1( * )
30 * COMPLEX A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
31 * $ VL( LDVL, * ), VR( LDVR, * ), W( * ), W1( * ),
32 * $ WORK( * )
33 * ..
34 *
35 *
36 *> \par Purpose:
37 * =============
38 *>
39 *> \verbatim
40 *>
41 *> CGET23 checks the nonsymmetric eigenvalue problem driver CGEEVX.
42 *> If COMP = .FALSE., the first 8 of the following tests will be
43 *> performed on the input matrix A, and also test 9 if LWORK is
44 *> sufficiently large.
45 *> if COMP is .TRUE. all 11 tests will be performed.
46 *>
47 *> (1) | A * VR - VR * W | / ( n |A| ulp )
48 *>
49 *> Here VR is the matrix of unit right eigenvectors.
50 *> W is a diagonal matrix with diagonal entries W(j).
51 *>
52 *> (2) | A**H * VL - VL * W**H | / ( n |A| ulp )
53 *>
54 *> Here VL is the matrix of unit left eigenvectors, A**H is the
55 *> conjugate transpose of A, and W is as above.
56 *>
57 *> (3) | |VR(i)| - 1 | / ulp and largest component real
58 *>
59 *> VR(i) denotes the i-th column of VR.
60 *>
61 *> (4) | |VL(i)| - 1 | / ulp and largest component real
62 *>
63 *> VL(i) denotes the i-th column of VL.
64 *>
65 *> (5) 0 if W(full) = W(partial), 1/ulp otherwise
66 *>
67 *> W(full) denotes the eigenvalues computed when VR, VL, RCONDV
68 *> and RCONDE are also computed, and W(partial) denotes the
69 *> eigenvalues computed when only some of VR, VL, RCONDV, and
70 *> RCONDE are computed.
71 *>
72 *> (6) 0 if VR(full) = VR(partial), 1/ulp otherwise
73 *>
74 *> VR(full) denotes the right eigenvectors computed when VL, RCONDV
75 *> and RCONDE are computed, and VR(partial) denotes the result
76 *> when only some of VL and RCONDV are computed.
77 *>
78 *> (7) 0 if VL(full) = VL(partial), 1/ulp otherwise
79 *>
80 *> VL(full) denotes the left eigenvectors computed when VR, RCONDV
81 *> and RCONDE are computed, and VL(partial) denotes the result
82 *> when only some of VR and RCONDV are computed.
83 *>
84 *> (8) 0 if SCALE, ILO, IHI, ABNRM (full) =
85 *> SCALE, ILO, IHI, ABNRM (partial)
86 *> 1/ulp otherwise
87 *>
88 *> SCALE, ILO, IHI and ABNRM describe how the matrix is balanced.
89 *> (full) is when VR, VL, RCONDE and RCONDV are also computed, and
90 *> (partial) is when some are not computed.
91 *>
92 *> (9) 0 if RCONDV(full) = RCONDV(partial), 1/ulp otherwise
93 *>
94 *> RCONDV(full) denotes the reciprocal condition numbers of the
95 *> right eigenvectors computed when VR, VL and RCONDE are also
96 *> computed. RCONDV(partial) denotes the reciprocal condition
97 *> numbers when only some of VR, VL and RCONDE are computed.
98 *>
99 *> (10) |RCONDV - RCDVIN| / cond(RCONDV)
100 *>
101 *> RCONDV is the reciprocal right eigenvector condition number
102 *> computed by CGEEVX and RCDVIN (the precomputed true value)
103 *> is supplied as input. cond(RCONDV) is the condition number of
104 *> RCONDV, and takes errors in computing RCONDV into account, so
105 *> that the resulting quantity should be O(ULP). cond(RCONDV) is
106 *> essentially given by norm(A)/RCONDE.
107 *>
108 *> (11) |RCONDE - RCDEIN| / cond(RCONDE)
109 *>
110 *> RCONDE is the reciprocal eigenvalue condition number
111 *> computed by CGEEVX and RCDEIN (the precomputed true value)
112 *> is supplied as input. cond(RCONDE) is the condition number
113 *> of RCONDE, and takes errors in computing RCONDE into account,
114 *> so that the resulting quantity should be O(ULP). cond(RCONDE)
115 *> is essentially given by norm(A)/RCONDV.
116 *> \endverbatim
117 *
118 * Arguments:
119 * ==========
120 *
121 *> \param[in] COMP
122 *> \verbatim
123 *> COMP is LOGICAL
124 *> COMP describes which input tests to perform:
125 *> = .FALSE. if the computed condition numbers are not to
126 *> be tested against RCDVIN and RCDEIN
127 *> = .TRUE. if they are to be compared
128 *> \endverbatim
129 *>
130 *> \param[in] ISRT
131 *> \verbatim
132 *> ISRT is INTEGER
133 *> If COMP = .TRUE., ISRT indicates in how the eigenvalues
134 *> corresponding to values in RCDVIN and RCDEIN are ordered:
135 *> = 0 means the eigenvalues are sorted by
136 *> increasing real part
137 *> = 1 means the eigenvalues are sorted by
138 *> increasing imaginary part
139 *> If COMP = .FALSE., ISRT is not referenced.
140 *> \endverbatim
141 *>
142 *> \param[in] BALANC
143 *> \verbatim
144 *> BALANC is CHARACTER
145 *> Describes the balancing option to be tested.
146 *> = 'N' for no permuting or diagonal scaling
147 *> = 'P' for permuting but no diagonal scaling
148 *> = 'S' for no permuting but diagonal scaling
149 *> = 'B' for permuting and diagonal scaling
150 *> \endverbatim
151 *>
152 *> \param[in] JTYPE
153 *> \verbatim
154 *> JTYPE is INTEGER
155 *> Type of input matrix. Used to label output if error occurs.
156 *> \endverbatim
157 *>
158 *> \param[in] THRESH
159 *> \verbatim
160 *> THRESH is REAL
161 *> A test will count as "failed" if the "error", computed as
162 *> described above, exceeds THRESH. Note that the error
163 *> is scaled to be O(1), so THRESH should be a reasonably
164 *> small multiple of 1, e.g., 10 or 100. In particular,
165 *> it should not depend on the precision (single vs. double)
166 *> or the size of the matrix. It must be at least zero.
167 *> \endverbatim
168 *>
169 *> \param[in] ISEED
170 *> \verbatim
171 *> ISEED is INTEGER array, dimension (4)
172 *> If COMP = .FALSE., the random number generator seed
173 *> used to produce matrix.
174 *> If COMP = .TRUE., ISEED(1) = the number of the example.
175 *> Used to label output if error occurs.
176 *> \endverbatim
177 *>
178 *> \param[in] NOUNIT
179 *> \verbatim
180 *> NOUNIT is INTEGER
181 *> The FORTRAN unit number for printing out error messages
182 *> (e.g., if a routine returns INFO not equal to 0.)
183 *> \endverbatim
184 *>
185 *> \param[in] N
186 *> \verbatim
187 *> N is INTEGER
188 *> The dimension of A. N must be at least 0.
189 *> \endverbatim
190 *>
191 *> \param[in,out] A
192 *> \verbatim
193 *> A is COMPLEX array, dimension (LDA,N)
194 *> Used to hold the matrix whose eigenvalues are to be
195 *> computed.
196 *> \endverbatim
197 *>
198 *> \param[in] LDA
199 *> \verbatim
200 *> LDA is INTEGER
201 *> The leading dimension of A, and H. LDA must be at
202 *> least 1 and at least N.
203 *> \endverbatim
204 *>
205 *> \param[out] H
206 *> \verbatim
207 *> H is COMPLEX array, dimension (LDA,N)
208 *> Another copy of the test matrix A, modified by CGEEVX.
209 *> \endverbatim
210 *>
211 *> \param[out] W
212 *> \verbatim
213 *> W is COMPLEX array, dimension (N)
214 *> Contains the eigenvalues of A.
215 *> \endverbatim
216 *>
217 *> \param[out] W1
218 *> \verbatim
219 *> W1 is COMPLEX array, dimension (N)
220 *> Like W, this array contains the eigenvalues of A,
221 *> but those computed when CGEEVX only computes a partial
222 *> eigendecomposition, i.e. not the eigenvalues and left
223 *> and right eigenvectors.
224 *> \endverbatim
225 *>
226 *> \param[out] VL
227 *> \verbatim
228 *> VL is COMPLEX array, dimension (LDVL,N)
229 *> VL holds the computed left eigenvectors.
230 *> \endverbatim
231 *>
232 *> \param[in] LDVL
233 *> \verbatim
234 *> LDVL is INTEGER
235 *> Leading dimension of VL. Must be at least max(1,N).
236 *> \endverbatim
237 *>
238 *> \param[out] VR
239 *> \verbatim
240 *> VR is COMPLEX array, dimension (LDVR,N)
241 *> VR holds the computed right eigenvectors.
242 *> \endverbatim
243 *>
244 *> \param[in] LDVR
245 *> \verbatim
246 *> LDVR is INTEGER
247 *> Leading dimension of VR. Must be at least max(1,N).
248 *> \endverbatim
249 *>
250 *> \param[out] LRE
251 *> \verbatim
252 *> LRE is COMPLEX array, dimension (LDLRE,N)
253 *> LRE holds the computed right or left eigenvectors.
254 *> \endverbatim
255 *>
256 *> \param[in] LDLRE
257 *> \verbatim
258 *> LDLRE is INTEGER
259 *> Leading dimension of LRE. Must be at least max(1,N).
260 *> \endverbatim
261 *>
262 *> \param[out] RCONDV
263 *> \verbatim
264 *> RCONDV is REAL array, dimension (N)
265 *> RCONDV holds the computed reciprocal condition numbers
266 *> for eigenvectors.
267 *> \endverbatim
268 *>
269 *> \param[out] RCNDV1
270 *> \verbatim
271 *> RCNDV1 is REAL array, dimension (N)
272 *> RCNDV1 holds more computed reciprocal condition numbers
273 *> for eigenvectors.
274 *> \endverbatim
275 *>
276 *> \param[in] RCDVIN
277 *> \verbatim
278 *> RCDVIN is REAL array, dimension (N)
279 *> When COMP = .TRUE. RCDVIN holds the precomputed reciprocal
280 *> condition numbers for eigenvectors to be compared with
281 *> RCONDV.
282 *> \endverbatim
283 *>
284 *> \param[out] RCONDE
285 *> \verbatim
286 *> RCONDE is REAL array, dimension (N)
287 *> RCONDE holds the computed reciprocal condition numbers
288 *> for eigenvalues.
289 *> \endverbatim
290 *>
291 *> \param[out] RCNDE1
292 *> \verbatim
293 *> RCNDE1 is REAL array, dimension (N)
294 *> RCNDE1 holds more computed reciprocal condition numbers
295 *> for eigenvalues.
296 *> \endverbatim
297 *>
298 *> \param[in] RCDEIN
299 *> \verbatim
300 *> RCDEIN is REAL array, dimension (N)
301 *> When COMP = .TRUE. RCDEIN holds the precomputed reciprocal
302 *> condition numbers for eigenvalues to be compared with
303 *> RCONDE.
304 *> \endverbatim
305 *>
306 *> \param[out] SCALE
307 *> \verbatim
308 *> SCALE is REAL array, dimension (N)
309 *> Holds information describing balancing of matrix.
310 *> \endverbatim
311 *>
312 *> \param[out] SCALE1
313 *> \verbatim
314 *> SCALE1 is REAL array, dimension (N)
315 *> Holds information describing balancing of matrix.
316 *> \endverbatim
317 *>
318 *> \param[out] RESULT
319 *> \verbatim
320 *> RESULT is REAL array, dimension (11)
321 *> The values computed by the 11 tests described above.
322 *> The values are currently limited to 1/ulp, to avoid
323 *> overflow.
324 *> \endverbatim
325 *>
326 *> \param[out] WORK
327 *> \verbatim
328 *> WORK is COMPLEX array, dimension (LWORK)
329 *> \endverbatim
330 *>
331 *> \param[in] LWORK
332 *> \verbatim
333 *> LWORK is INTEGER
334 *> The number of entries in WORK. This must be at least
335 *> 2*N, and 2*N+N**2 if tests 9, 10 or 11 are to be performed.
336 *> \endverbatim
337 *>
338 *> \param[out] RWORK
339 *> \verbatim
340 *> RWORK is REAL array, dimension (2*N)
341 *> \endverbatim
342 *>
343 *> \param[out] INFO
344 *> \verbatim
345 *> INFO is INTEGER
346 *> If 0, successful exit.
347 *> If <0, input parameter -INFO had an incorrect value.
348 *> If >0, CGEEVX returned an error code, the absolute
349 *> value of which is returned.
350 *> \endverbatim
351 *
352 * Authors:
353 * ========
354 *
355 *> \author Univ. of Tennessee
356 *> \author Univ. of California Berkeley
357 *> \author Univ. of Colorado Denver
358 *> \author NAG Ltd.
359 *
360 *> \ingroup complex_eig
361 *
362 * =====================================================================
363  SUBROUTINE cget23( COMP, ISRT, BALANC, JTYPE, THRESH, ISEED,
364  $ NOUNIT, N, A, LDA, H, W, W1, VL, LDVL, VR,
365  $ LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN,
366  $ RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT,
367  $ WORK, LWORK, RWORK, INFO )
368 *
369 * -- LAPACK test routine --
370 * -- LAPACK is a software package provided by Univ. of Tennessee, --
371 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
372 *
373 * .. Scalar Arguments ..
374  LOGICAL COMP
375  CHARACTER BALANC
376  INTEGER INFO, ISRT, JTYPE, LDA, LDLRE, LDVL, LDVR,
377  $ LWORK, N, NOUNIT
378  REAL THRESH
379 * ..
380 * .. Array Arguments ..
381  INTEGER ISEED( 4 )
382  REAL RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
383  $ RCNDV1( * ), RCONDE( * ), RCONDV( * ),
384  $ RESULT( 11 ), RWORK( * ), SCALE( * ),
385  $ scale1( * )
386  COMPLEX A( LDA, * ), H( LDA, * ), LRE( LDLRE, * ),
387  $ VL( LDVL, * ), VR( LDVR, * ), W( * ), W1( * ),
388  $ WORK( * )
389 * ..
390 *
391 * =====================================================================
392 *
393 * .. Parameters ..
394  REAL ZERO, ONE, TWO
395  PARAMETER ( ZERO = 0.0e0, one = 1.0e0, two = 2.0e0 )
396  REAL EPSIN
397  PARAMETER ( EPSIN = 5.9605e-8 )
398 * ..
399 * .. Local Scalars ..
400  LOGICAL BALOK, NOBAL
401  CHARACTER SENSE
402  INTEGER I, IHI, IHI1, IINFO, ILO, ILO1, ISENS, ISENSM,
403  $ J, JJ, KMIN
404  REAL ABNRM, ABNRM1, EPS, SMLNUM, TNRM, TOL, TOLIN,
405  $ ulp, ulpinv, v, vmax, vmx, vricmp, vrimin,
406  $ vrmx, vtst
407  COMPLEX CTMP
408 * ..
409 * .. Local Arrays ..
410  CHARACTER SENS( 2 )
411  REAL RES( 2 )
412  COMPLEX CDUM( 1 )
413 * ..
414 * .. External Functions ..
415  LOGICAL LSAME
416  REAL SCNRM2, SLAMCH
417  EXTERNAL LSAME, SCNRM2, SLAMCH
418 * ..
419 * .. External Subroutines ..
420  EXTERNAL cgeevx, cget22, clacpy, xerbla
421 * ..
422 * .. Intrinsic Functions ..
423  INTRINSIC abs, aimag, max, min, real
424 * ..
425 * .. Data statements ..
426  DATA sens / 'N', 'V' /
427 * ..
428 * .. Executable Statements ..
429 *
430 * Check for errors
431 *
432  nobal = lsame( balanc, 'N' )
433  balok = nobal .OR. lsame( balanc, 'P' ) .OR.
434  $ lsame( balanc, 'S' ) .OR. lsame( balanc, 'B' )
435  info = 0
436  IF( isrt.NE.0 .AND. isrt.NE.1 ) THEN
437  info = -2
438  ELSE IF( .NOT.balok ) THEN
439  info = -3
440  ELSE IF( thresh.LT.zero ) THEN
441  info = -5
442  ELSE IF( nounit.LE.0 ) THEN
443  info = -7
444  ELSE IF( n.LT.0 ) THEN
445  info = -8
446  ELSE IF( lda.LT.1 .OR. lda.LT.n ) THEN
447  info = -10
448  ELSE IF( ldvl.LT.1 .OR. ldvl.LT.n ) THEN
449  info = -15
450  ELSE IF( ldvr.LT.1 .OR. ldvr.LT.n ) THEN
451  info = -17
452  ELSE IF( ldlre.LT.1 .OR. ldlre.LT.n ) THEN
453  info = -19
454  ELSE IF( lwork.LT.2*n .OR. ( comp .AND. lwork.LT.2*n+n*n ) ) THEN
455  info = -30
456  END IF
457 *
458  IF( info.NE.0 ) THEN
459  CALL xerbla( 'CGET23', -info )
460  RETURN
461  END IF
462 *
463 * Quick return if nothing to do
464 *
465  DO 10 i = 1, 11
466  result( i ) = -one
467  10 CONTINUE
468 *
469  IF( n.EQ.0 )
470  $ RETURN
471 *
472 * More Important constants
473 *
474  ulp = slamch( 'Precision' )
475  smlnum = slamch( 'S' )
476  ulpinv = one / ulp
477 *
478 * Compute eigenvalues and eigenvectors, and test them
479 *
480  IF( lwork.GE.2*n+n*n ) THEN
481  sense = 'B'
482  isensm = 2
483  ELSE
484  sense = 'E'
485  isensm = 1
486  END IF
487  CALL clacpy( 'F', n, n, a, lda, h, lda )
488  CALL cgeevx( balanc, 'V', 'V', sense, n, h, lda, w, vl, ldvl, vr,
489  $ ldvr, ilo, ihi, scale, abnrm, rconde, rcondv, work,
490  $ lwork, rwork, iinfo )
491  IF( iinfo.NE.0 ) THEN
492  result( 1 ) = ulpinv
493  IF( jtype.NE.22 ) THEN
494  WRITE( nounit, fmt = 9998 )'CGEEVX1', iinfo, n, jtype,
495  $ balanc, iseed
496  ELSE
497  WRITE( nounit, fmt = 9999 )'CGEEVX1', iinfo, n, iseed( 1 )
498  END IF
499  info = abs( iinfo )
500  RETURN
501  END IF
502 *
503 * Do Test (1)
504 *
505  CALL cget22( 'N', 'N', 'N', n, a, lda, vr, ldvr, w, work, rwork,
506  $ res )
507  result( 1 ) = res( 1 )
508 *
509 * Do Test (2)
510 *
511  CALL cget22( 'C', 'N', 'C', n, a, lda, vl, ldvl, w, work, rwork,
512  $ res )
513  result( 2 ) = res( 1 )
514 *
515 * Do Test (3)
516 *
517  DO 30 j = 1, n
518  tnrm = scnrm2( n, vr( 1, j ), 1 )
519  result( 3 ) = max( result( 3 ),
520  $ min( ulpinv, abs( tnrm-one ) / ulp ) )
521  vmx = zero
522  vrmx = zero
523  DO 20 jj = 1, n
524  vtst = abs( vr( jj, j ) )
525  IF( vtst.GT.vmx )
526  $ vmx = vtst
527  IF( aimag( vr( jj, j ) ).EQ.zero .AND.
528  $ abs( real( vr( jj, j ) ) ).GT.vrmx )
529  $ vrmx = abs( real( vr( jj, j ) ) )
530  20 CONTINUE
531  IF( vrmx / vmx.LT.one-two*ulp )
532  $ result( 3 ) = ulpinv
533  30 CONTINUE
534 *
535 * Do Test (4)
536 *
537  DO 50 j = 1, n
538  tnrm = scnrm2( n, vl( 1, j ), 1 )
539  result( 4 ) = max( result( 4 ),
540  $ min( ulpinv, abs( tnrm-one ) / ulp ) )
541  vmx = zero
542  vrmx = zero
543  DO 40 jj = 1, n
544  vtst = abs( vl( jj, j ) )
545  IF( vtst.GT.vmx )
546  $ vmx = vtst
547  IF( aimag( vl( jj, j ) ).EQ.zero .AND.
548  $ abs( real( vl( jj, j ) ) ).GT.vrmx )
549  $ vrmx = abs( real( vl( jj, j ) ) )
550  40 CONTINUE
551  IF( vrmx / vmx.LT.one-two*ulp )
552  $ result( 4 ) = ulpinv
553  50 CONTINUE
554 *
555 * Test for all options of computing condition numbers
556 *
557  DO 200 isens = 1, isensm
558 *
559  sense = sens( isens )
560 *
561 * Compute eigenvalues only, and test them
562 *
563  CALL clacpy( 'F', n, n, a, lda, h, lda )
564  CALL cgeevx( balanc, 'N', 'N', sense, n, h, lda, w1, cdum, 1,
565  $ cdum, 1, ilo1, ihi1, scale1, abnrm1, rcnde1,
566  $ rcndv1, work, lwork, rwork, iinfo )
567  IF( iinfo.NE.0 ) THEN
568  result( 1 ) = ulpinv
569  IF( jtype.NE.22 ) THEN
570  WRITE( nounit, fmt = 9998 )'CGEEVX2', iinfo, n, jtype,
571  $ balanc, iseed
572  ELSE
573  WRITE( nounit, fmt = 9999 )'CGEEVX2', iinfo, n,
574  $ iseed( 1 )
575  END IF
576  info = abs( iinfo )
577  GO TO 190
578  END IF
579 *
580 * Do Test (5)
581 *
582  DO 60 j = 1, n
583  IF( w( j ).NE.w1( j ) )
584  $ result( 5 ) = ulpinv
585  60 CONTINUE
586 *
587 * Do Test (8)
588 *
589  IF( .NOT.nobal ) THEN
590  DO 70 j = 1, n
591  IF( scale( j ).NE.scale1( j ) )
592  $ result( 8 ) = ulpinv
593  70 CONTINUE
594  IF( ilo.NE.ilo1 )
595  $ result( 8 ) = ulpinv
596  IF( ihi.NE.ihi1 )
597  $ result( 8 ) = ulpinv
598  IF( abnrm.NE.abnrm1 )
599  $ result( 8 ) = ulpinv
600  END IF
601 *
602 * Do Test (9)
603 *
604  IF( isens.EQ.2 .AND. n.GT.1 ) THEN
605  DO 80 j = 1, n
606  IF( rcondv( j ).NE.rcndv1( j ) )
607  $ result( 9 ) = ulpinv
608  80 CONTINUE
609  END IF
610 *
611 * Compute eigenvalues and right eigenvectors, and test them
612 *
613  CALL clacpy( 'F', n, n, a, lda, h, lda )
614  CALL cgeevx( balanc, 'N', 'V', sense, n, h, lda, w1, cdum, 1,
615  $ lre, ldlre, ilo1, ihi1, scale1, abnrm1, rcnde1,
616  $ rcndv1, work, lwork, rwork, iinfo )
617  IF( iinfo.NE.0 ) THEN
618  result( 1 ) = ulpinv
619  IF( jtype.NE.22 ) THEN
620  WRITE( nounit, fmt = 9998 )'CGEEVX3', iinfo, n, jtype,
621  $ balanc, iseed
622  ELSE
623  WRITE( nounit, fmt = 9999 )'CGEEVX3', iinfo, n,
624  $ iseed( 1 )
625  END IF
626  info = abs( iinfo )
627  GO TO 190
628  END IF
629 *
630 * Do Test (5) again
631 *
632  DO 90 j = 1, n
633  IF( w( j ).NE.w1( j ) )
634  $ result( 5 ) = ulpinv
635  90 CONTINUE
636 *
637 * Do Test (6)
638 *
639  DO 110 j = 1, n
640  DO 100 jj = 1, n
641  IF( vr( j, jj ).NE.lre( j, jj ) )
642  $ result( 6 ) = ulpinv
643  100 CONTINUE
644  110 CONTINUE
645 *
646 * Do Test (8) again
647 *
648  IF( .NOT.nobal ) THEN
649  DO 120 j = 1, n
650  IF( scale( j ).NE.scale1( j ) )
651  $ result( 8 ) = ulpinv
652  120 CONTINUE
653  IF( ilo.NE.ilo1 )
654  $ result( 8 ) = ulpinv
655  IF( ihi.NE.ihi1 )
656  $ result( 8 ) = ulpinv
657  IF( abnrm.NE.abnrm1 )
658  $ result( 8 ) = ulpinv
659  END IF
660 *
661 * Do Test (9) again
662 *
663  IF( isens.EQ.2 .AND. n.GT.1 ) THEN
664  DO 130 j = 1, n
665  IF( rcondv( j ).NE.rcndv1( j ) )
666  $ result( 9 ) = ulpinv
667  130 CONTINUE
668  END IF
669 *
670 * Compute eigenvalues and left eigenvectors, and test them
671 *
672  CALL clacpy( 'F', n, n, a, lda, h, lda )
673  CALL cgeevx( balanc, 'V', 'N', sense, n, h, lda, w1, lre,
674  $ ldlre, cdum, 1, ilo1, ihi1, scale1, abnrm1,
675  $ rcnde1, rcndv1, work, lwork, rwork, iinfo )
676  IF( iinfo.NE.0 ) THEN
677  result( 1 ) = ulpinv
678  IF( jtype.NE.22 ) THEN
679  WRITE( nounit, fmt = 9998 )'CGEEVX4', iinfo, n, jtype,
680  $ balanc, iseed
681  ELSE
682  WRITE( nounit, fmt = 9999 )'CGEEVX4', iinfo, n,
683  $ iseed( 1 )
684  END IF
685  info = abs( iinfo )
686  GO TO 190
687  END IF
688 *
689 * Do Test (5) again
690 *
691  DO 140 j = 1, n
692  IF( w( j ).NE.w1( j ) )
693  $ result( 5 ) = ulpinv
694  140 CONTINUE
695 *
696 * Do Test (7)
697 *
698  DO 160 j = 1, n
699  DO 150 jj = 1, n
700  IF( vl( j, jj ).NE.lre( j, jj ) )
701  $ result( 7 ) = ulpinv
702  150 CONTINUE
703  160 CONTINUE
704 *
705 * Do Test (8) again
706 *
707  IF( .NOT.nobal ) THEN
708  DO 170 j = 1, n
709  IF( scale( j ).NE.scale1( j ) )
710  $ result( 8 ) = ulpinv
711  170 CONTINUE
712  IF( ilo.NE.ilo1 )
713  $ result( 8 ) = ulpinv
714  IF( ihi.NE.ihi1 )
715  $ result( 8 ) = ulpinv
716  IF( abnrm.NE.abnrm1 )
717  $ result( 8 ) = ulpinv
718  END IF
719 *
720 * Do Test (9) again
721 *
722  IF( isens.EQ.2 .AND. n.GT.1 ) THEN
723  DO 180 j = 1, n
724  IF( rcondv( j ).NE.rcndv1( j ) )
725  $ result( 9 ) = ulpinv
726  180 CONTINUE
727  END IF
728 *
729  190 CONTINUE
730 *
731  200 CONTINUE
732 *
733 * If COMP, compare condition numbers to precomputed ones
734 *
735  IF( comp ) THEN
736  CALL clacpy( 'F', n, n, a, lda, h, lda )
737  CALL cgeevx( 'N', 'V', 'V', 'B', n, h, lda, w, vl, ldvl, vr,
738  $ ldvr, ilo, ihi, scale, abnrm, rconde, rcondv,
739  $ work, lwork, rwork, iinfo )
740  IF( iinfo.NE.0 ) THEN
741  result( 1 ) = ulpinv
742  WRITE( nounit, fmt = 9999 )'CGEEVX5', iinfo, n, iseed( 1 )
743  info = abs( iinfo )
744  GO TO 250
745  END IF
746 *
747 * Sort eigenvalues and condition numbers lexicographically
748 * to compare with inputs
749 *
750  DO 220 i = 1, n - 1
751  kmin = i
752  IF( isrt.EQ.0 ) THEN
753  vrimin = real( w( i ) )
754  ELSE
755  vrimin = aimag( w( i ) )
756  END IF
757  DO 210 j = i + 1, n
758  IF( isrt.EQ.0 ) THEN
759  vricmp = real( w( j ) )
760  ELSE
761  vricmp = aimag( w( j ) )
762  END IF
763  IF( vricmp.LT.vrimin ) THEN
764  kmin = j
765  vrimin = vricmp
766  END IF
767  210 CONTINUE
768  ctmp = w( kmin )
769  w( kmin ) = w( i )
770  w( i ) = ctmp
771  vrimin = rconde( kmin )
772  rconde( kmin ) = rconde( i )
773  rconde( i ) = vrimin
774  vrimin = rcondv( kmin )
775  rcondv( kmin ) = rcondv( i )
776  rcondv( i ) = vrimin
777  220 CONTINUE
778 *
779 * Compare condition numbers for eigenvectors
780 * taking their condition numbers into account
781 *
782  result( 10 ) = zero
783  eps = max( epsin, ulp )
784  v = max( real( n )*eps*abnrm, smlnum )
785  IF( abnrm.EQ.zero )
786  $ v = one
787  DO 230 i = 1, n
788  IF( v.GT.rcondv( i )*rconde( i ) ) THEN
789  tol = rcondv( i )
790  ELSE
791  tol = v / rconde( i )
792  END IF
793  IF( v.GT.rcdvin( i )*rcdein( i ) ) THEN
794  tolin = rcdvin( i )
795  ELSE
796  tolin = v / rcdein( i )
797  END IF
798  tol = max( tol, smlnum / eps )
799  tolin = max( tolin, smlnum / eps )
800  IF( eps*( rcdvin( i )-tolin ).GT.rcondv( i )+tol ) THEN
801  vmax = one / eps
802  ELSE IF( rcdvin( i )-tolin.GT.rcondv( i )+tol ) THEN
803  vmax = ( rcdvin( i )-tolin ) / ( rcondv( i )+tol )
804  ELSE IF( rcdvin( i )+tolin.LT.eps*( rcondv( i )-tol ) ) THEN
805  vmax = one / eps
806  ELSE IF( rcdvin( i )+tolin.LT.rcondv( i )-tol ) THEN
807  vmax = ( rcondv( i )-tol ) / ( rcdvin( i )+tolin )
808  ELSE
809  vmax = one
810  END IF
811  result( 10 ) = max( result( 10 ), vmax )
812  230 CONTINUE
813 *
814 * Compare condition numbers for eigenvalues
815 * taking their condition numbers into account
816 *
817  result( 11 ) = zero
818  DO 240 i = 1, n
819  IF( v.GT.rcondv( i ) ) THEN
820  tol = one
821  ELSE
822  tol = v / rcondv( i )
823  END IF
824  IF( v.GT.rcdvin( i ) ) THEN
825  tolin = one
826  ELSE
827  tolin = v / rcdvin( i )
828  END IF
829  tol = max( tol, smlnum / eps )
830  tolin = max( tolin, smlnum / eps )
831  IF( eps*( rcdein( i )-tolin ).GT.rconde( i )+tol ) THEN
832  vmax = one / eps
833  ELSE IF( rcdein( i )-tolin.GT.rconde( i )+tol ) THEN
834  vmax = ( rcdein( i )-tolin ) / ( rconde( i )+tol )
835  ELSE IF( rcdein( i )+tolin.LT.eps*( rconde( i )-tol ) ) THEN
836  vmax = one / eps
837  ELSE IF( rcdein( i )+tolin.LT.rconde( i )-tol ) THEN
838  vmax = ( rconde( i )-tol ) / ( rcdein( i )+tolin )
839  ELSE
840  vmax = one
841  END IF
842  result( 11 ) = max( result( 11 ), vmax )
843  240 CONTINUE
844  250 CONTINUE
845 *
846  END IF
847 *
848  9999 FORMAT( ' CGET23: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
849  $ i6, ', INPUT EXAMPLE NUMBER = ', i4 )
850  9998 FORMAT( ' CGET23: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
851  $ i6, ', JTYPE=', i6, ', BALANC = ', a, ', ISEED=(',
852  $ 3( i5, ',' ), i5, ')' )
853 *
854  RETURN
855 *
856 * End of CGET23
857 *
858  END
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine cget23(COMP, ISRT, BALANC, JTYPE, THRESH, ISEED, NOUNIT, N, A, LDA, H, W, W1, VL, LDVL, VR, LDVR, LRE, LDLRE, RCONDV, RCNDV1, RCDVIN, RCONDE, RCNDE1, RCDEIN, SCALE, SCALE1, RESULT, WORK, LWORK, RWORK, INFO)
CGET23
Definition: cget23.f:368
subroutine cget22(TRANSA, TRANSE, TRANSW, N, A, LDA, E, LDE, W, WORK, RWORK, RESULT)
CGET22
Definition: cget22.f:144
subroutine cgeevx(BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, W, VL, LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV, WORK, LWORK, RWORK, INFO)
CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices
Definition: cgeevx.f:288
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