LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
zget23.f
Go to the documentation of this file.
1 *> \brief \b ZGET23
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 ZGET23( 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 * DOUBLE PRECISION THRESH
23 * ..
24 * .. Array Arguments ..
25 * INTEGER ISEED( 4 )
26 * DOUBLE PRECISION RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
27 * $ RCNDV1( * ), RCONDE( * ), RCONDV( * ),
28 * $ RESULT( 11 ), RWORK( * ), SCALE( * ),
29 * $ SCALE1( * )
30 * COMPLEX*16 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 *> ZGET23 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 ZGEEVX 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 ZGEEVX 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 DOUBLE PRECISION
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*16 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*16 array, dimension (LDA,N)
208 *> Another copy of the test matrix A, modified by ZGEEVX.
209 *> \endverbatim
210 *>
211 *> \param[out] W
212 *> \verbatim
213 *> W is COMPLEX*16 array, dimension (N)
214 *> Contains the eigenvalues of A.
215 *> \endverbatim
216 *>
217 *> \param[out] W1
218 *> \verbatim
219 *> W1 is COMPLEX*16 array, dimension (N)
220 *> Like W, this array contains the eigenvalues of A,
221 *> but those computed when ZGEEVX 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*16 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*16 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*16 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION 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 DOUBLE PRECISION array, dimension (N)
309 *> Holds information describing balancing of matrix.
310 *> \endverbatim
311 *>
312 *> \param[out] SCALE1
313 *> \verbatim
314 *> SCALE1 is DOUBLE PRECISION array, dimension (N)
315 *> Holds information describing balancing of matrix.
316 *> \endverbatim
317 *>
318 *> \param[out] RESULT
319 *> \verbatim
320 *> RESULT is DOUBLE PRECISION 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*16 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 DOUBLE PRECISION 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, ZGEEVX 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 *> \date November 2011
361 *
362 *> \ingroup complex16_eig
363 *
364 * =====================================================================
365  SUBROUTINE zget23( COMP, ISRT, BALANC, JTYPE, THRESH, ISEED,
366  $ nounit, n, a, lda, h, w, w1, vl, ldvl, vr,
367  $ ldvr, lre, ldlre, rcondv, rcndv1, rcdvin,
368  $ rconde, rcnde1, rcdein, scale, scale1, result,
369  $ work, lwork, rwork, info )
370 *
371 * -- LAPACK test routine (version 3.4.0) --
372 * -- LAPACK is a software package provided by Univ. of Tennessee, --
373 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
374 * November 2011
375 *
376 * .. Scalar Arguments ..
377  LOGICAL COMP
378  CHARACTER BALANC
379  INTEGER INFO, ISRT, JTYPE, LDA, LDLRE, LDVL, LDVR,
380  $ lwork, n, nounit
381  DOUBLE PRECISION THRESH
382 * ..
383 * .. Array Arguments ..
384  INTEGER ISEED( 4 )
385  DOUBLE PRECISION RCDEIN( * ), RCDVIN( * ), RCNDE1( * ),
386  $ rcndv1( * ), rconde( * ), rcondv( * ),
387  $ result( 11 ), rwork( * ), scale( * ),
388  $ scale1( * )
389  COMPLEX*16 A( lda, * ), H( lda, * ), LRE( ldlre, * ),
390  $ vl( ldvl, * ), vr( ldvr, * ), w( * ), w1( * ),
391  $ work( * )
392 * ..
393 *
394 * =====================================================================
395 *
396 * .. Parameters ..
397  DOUBLE PRECISION ZERO, ONE, TWO
398  parameter ( zero = 0.0d0, one = 1.0d0, two = 2.0d0 )
399  DOUBLE PRECISION EPSIN
400  parameter ( epsin = 5.9605d-8 )
401 * ..
402 * .. Local Scalars ..
403  LOGICAL BALOK, NOBAL
404  CHARACTER SENSE
405  INTEGER I, IHI, IHI1, IINFO, ILO, ILO1, ISENS, ISENSM,
406  $ j, jj, kmin
407  DOUBLE PRECISION ABNRM, ABNRM1, EPS, SMLNUM, TNRM, TOL, TOLIN,
408  $ ulp, ulpinv, v, vmax, vmx, vricmp, vrimin,
409  $ vrmx, vtst
410  COMPLEX*16 CTMP
411 * ..
412 * .. Local Arrays ..
413  CHARACTER SENS( 2 )
414  DOUBLE PRECISION RES( 2 )
415  COMPLEX*16 CDUM( 1 )
416 * ..
417 * .. External Functions ..
418  LOGICAL LSAME
419  DOUBLE PRECISION DLAMCH, DZNRM2
420  EXTERNAL lsame, dlamch, dznrm2
421 * ..
422 * .. External Subroutines ..
423  EXTERNAL xerbla, zgeevx, zget22, zlacpy
424 * ..
425 * .. Intrinsic Functions ..
426  INTRINSIC abs, dble, dimag, max, min
427 * ..
428 * .. Data statements ..
429  DATA sens / 'N', 'V' /
430 * ..
431 * .. Executable Statements ..
432 *
433 * Check for errors
434 *
435  nobal = lsame( balanc, 'N' )
436  balok = nobal .OR. lsame( balanc, 'P' ) .OR.
437  $ lsame( balanc, 'S' ) .OR. lsame( balanc, 'B' )
438  info = 0
439  IF( isrt.NE.0 .AND. isrt.NE.1 ) THEN
440  info = -2
441  ELSE IF( .NOT.balok ) THEN
442  info = -3
443  ELSE IF( thresh.LT.zero ) THEN
444  info = -5
445  ELSE IF( nounit.LE.0 ) THEN
446  info = -7
447  ELSE IF( n.LT.0 ) THEN
448  info = -8
449  ELSE IF( lda.LT.1 .OR. lda.LT.n ) THEN
450  info = -10
451  ELSE IF( ldvl.LT.1 .OR. ldvl.LT.n ) THEN
452  info = -15
453  ELSE IF( ldvr.LT.1 .OR. ldvr.LT.n ) THEN
454  info = -17
455  ELSE IF( ldlre.LT.1 .OR. ldlre.LT.n ) THEN
456  info = -19
457  ELSE IF( lwork.LT.2*n .OR. ( comp .AND. lwork.LT.2*n+n*n ) ) THEN
458  info = -30
459  END IF
460 *
461  IF( info.NE.0 ) THEN
462  CALL xerbla( 'ZGET23', -info )
463  RETURN
464  END IF
465 *
466 * Quick return if nothing to do
467 *
468  DO 10 i = 1, 11
469  result( i ) = -one
470  10 CONTINUE
471 *
472  IF( n.EQ.0 )
473  $ RETURN
474 *
475 * More Important constants
476 *
477  ulp = dlamch( 'Precision' )
478  smlnum = dlamch( 'S' )
479  ulpinv = one / ulp
480 *
481 * Compute eigenvalues and eigenvectors, and test them
482 *
483  IF( lwork.GE.2*n+n*n ) THEN
484  sense = 'B'
485  isensm = 2
486  ELSE
487  sense = 'E'
488  isensm = 1
489  END IF
490  CALL zlacpy( 'F', n, n, a, lda, h, lda )
491  CALL zgeevx( balanc, 'V', 'V', sense, n, h, lda, w, vl, ldvl, vr,
492  $ ldvr, ilo, ihi, scale, abnrm, rconde, rcondv, work,
493  $ lwork, rwork, iinfo )
494  IF( iinfo.NE.0 ) THEN
495  result( 1 ) = ulpinv
496  IF( jtype.NE.22 ) THEN
497  WRITE( nounit, fmt = 9998 )'ZGEEVX1', iinfo, n, jtype,
498  $ balanc, iseed
499  ELSE
500  WRITE( nounit, fmt = 9999 )'ZGEEVX1', iinfo, n, iseed( 1 )
501  END IF
502  info = abs( iinfo )
503  RETURN
504  END IF
505 *
506 * Do Test (1)
507 *
508  CALL zget22( 'N', 'N', 'N', n, a, lda, vr, ldvr, w, work, rwork,
509  $ res )
510  result( 1 ) = res( 1 )
511 *
512 * Do Test (2)
513 *
514  CALL zget22( 'C', 'N', 'C', n, a, lda, vl, ldvl, w, work, rwork,
515  $ res )
516  result( 2 ) = res( 1 )
517 *
518 * Do Test (3)
519 *
520  DO 30 j = 1, n
521  tnrm = dznrm2( n, vr( 1, j ), 1 )
522  result( 3 ) = max( result( 3 ),
523  $ min( ulpinv, abs( tnrm-one ) / ulp ) )
524  vmx = zero
525  vrmx = zero
526  DO 20 jj = 1, n
527  vtst = abs( vr( jj, j ) )
528  IF( vtst.GT.vmx )
529  $ vmx = vtst
530  IF( dimag( vr( jj, j ) ).EQ.zero .AND.
531  $ abs( dble( vr( jj, j ) ) ).GT.vrmx )
532  $ vrmx = abs( dble( vr( jj, j ) ) )
533  20 CONTINUE
534  IF( vrmx / vmx.LT.one-two*ulp )
535  $ result( 3 ) = ulpinv
536  30 CONTINUE
537 *
538 * Do Test (4)
539 *
540  DO 50 j = 1, n
541  tnrm = dznrm2( n, vl( 1, j ), 1 )
542  result( 4 ) = max( result( 4 ),
543  $ min( ulpinv, abs( tnrm-one ) / ulp ) )
544  vmx = zero
545  vrmx = zero
546  DO 40 jj = 1, n
547  vtst = abs( vl( jj, j ) )
548  IF( vtst.GT.vmx )
549  $ vmx = vtst
550  IF( dimag( vl( jj, j ) ).EQ.zero .AND.
551  $ abs( dble( vl( jj, j ) ) ).GT.vrmx )
552  $ vrmx = abs( dble( vl( jj, j ) ) )
553  40 CONTINUE
554  IF( vrmx / vmx.LT.one-two*ulp )
555  $ result( 4 ) = ulpinv
556  50 CONTINUE
557 *
558 * Test for all options of computing condition numbers
559 *
560  DO 200 isens = 1, isensm
561 *
562  sense = sens( isens )
563 *
564 * Compute eigenvalues only, and test them
565 *
566  CALL zlacpy( 'F', n, n, a, lda, h, lda )
567  CALL zgeevx( balanc, 'N', 'N', sense, n, h, lda, w1, cdum, 1,
568  $ cdum, 1, ilo1, ihi1, scale1, abnrm1, rcnde1,
569  $ rcndv1, work, lwork, rwork, iinfo )
570  IF( iinfo.NE.0 ) THEN
571  result( 1 ) = ulpinv
572  IF( jtype.NE.22 ) THEN
573  WRITE( nounit, fmt = 9998 )'ZGEEVX2', iinfo, n, jtype,
574  $ balanc, iseed
575  ELSE
576  WRITE( nounit, fmt = 9999 )'ZGEEVX2', iinfo, n,
577  $ iseed( 1 )
578  END IF
579  info = abs( iinfo )
580  GO TO 190
581  END IF
582 *
583 * Do Test (5)
584 *
585  DO 60 j = 1, n
586  IF( w( j ).NE.w1( j ) )
587  $ result( 5 ) = ulpinv
588  60 CONTINUE
589 *
590 * Do Test (8)
591 *
592  IF( .NOT.nobal ) THEN
593  DO 70 j = 1, n
594  IF( scale( j ).NE.scale1( j ) )
595  $ result( 8 ) = ulpinv
596  70 CONTINUE
597  IF( ilo.NE.ilo1 )
598  $ result( 8 ) = ulpinv
599  IF( ihi.NE.ihi1 )
600  $ result( 8 ) = ulpinv
601  IF( abnrm.NE.abnrm1 )
602  $ result( 8 ) = ulpinv
603  END IF
604 *
605 * Do Test (9)
606 *
607  IF( isens.EQ.2 .AND. n.GT.1 ) THEN
608  DO 80 j = 1, n
609  IF( rcondv( j ).NE.rcndv1( j ) )
610  $ result( 9 ) = ulpinv
611  80 CONTINUE
612  END IF
613 *
614 * Compute eigenvalues and right eigenvectors, and test them
615 *
616  CALL zlacpy( 'F', n, n, a, lda, h, lda )
617  CALL zgeevx( balanc, 'N', 'V', sense, n, h, lda, w1, cdum, 1,
618  $ lre, ldlre, ilo1, ihi1, scale1, abnrm1, rcnde1,
619  $ rcndv1, work, lwork, rwork, iinfo )
620  IF( iinfo.NE.0 ) THEN
621  result( 1 ) = ulpinv
622  IF( jtype.NE.22 ) THEN
623  WRITE( nounit, fmt = 9998 )'ZGEEVX3', iinfo, n, jtype,
624  $ balanc, iseed
625  ELSE
626  WRITE( nounit, fmt = 9999 )'ZGEEVX3', iinfo, n,
627  $ iseed( 1 )
628  END IF
629  info = abs( iinfo )
630  GO TO 190
631  END IF
632 *
633 * Do Test (5) again
634 *
635  DO 90 j = 1, n
636  IF( w( j ).NE.w1( j ) )
637  $ result( 5 ) = ulpinv
638  90 CONTINUE
639 *
640 * Do Test (6)
641 *
642  DO 110 j = 1, n
643  DO 100 jj = 1, n
644  IF( vr( j, jj ).NE.lre( j, jj ) )
645  $ result( 6 ) = ulpinv
646  100 CONTINUE
647  110 CONTINUE
648 *
649 * Do Test (8) again
650 *
651  IF( .NOT.nobal ) THEN
652  DO 120 j = 1, n
653  IF( scale( j ).NE.scale1( j ) )
654  $ result( 8 ) = ulpinv
655  120 CONTINUE
656  IF( ilo.NE.ilo1 )
657  $ result( 8 ) = ulpinv
658  IF( ihi.NE.ihi1 )
659  $ result( 8 ) = ulpinv
660  IF( abnrm.NE.abnrm1 )
661  $ result( 8 ) = ulpinv
662  END IF
663 *
664 * Do Test (9) again
665 *
666  IF( isens.EQ.2 .AND. n.GT.1 ) THEN
667  DO 130 j = 1, n
668  IF( rcondv( j ).NE.rcndv1( j ) )
669  $ result( 9 ) = ulpinv
670  130 CONTINUE
671  END IF
672 *
673 * Compute eigenvalues and left eigenvectors, and test them
674 *
675  CALL zlacpy( 'F', n, n, a, lda, h, lda )
676  CALL zgeevx( balanc, 'V', 'N', sense, n, h, lda, w1, lre,
677  $ ldlre, cdum, 1, ilo1, ihi1, scale1, abnrm1,
678  $ rcnde1, rcndv1, work, lwork, rwork, iinfo )
679  IF( iinfo.NE.0 ) THEN
680  result( 1 ) = ulpinv
681  IF( jtype.NE.22 ) THEN
682  WRITE( nounit, fmt = 9998 )'ZGEEVX4', iinfo, n, jtype,
683  $ balanc, iseed
684  ELSE
685  WRITE( nounit, fmt = 9999 )'ZGEEVX4', iinfo, n,
686  $ iseed( 1 )
687  END IF
688  info = abs( iinfo )
689  GO TO 190
690  END IF
691 *
692 * Do Test (5) again
693 *
694  DO 140 j = 1, n
695  IF( w( j ).NE.w1( j ) )
696  $ result( 5 ) = ulpinv
697  140 CONTINUE
698 *
699 * Do Test (7)
700 *
701  DO 160 j = 1, n
702  DO 150 jj = 1, n
703  IF( vl( j, jj ).NE.lre( j, jj ) )
704  $ result( 7 ) = ulpinv
705  150 CONTINUE
706  160 CONTINUE
707 *
708 * Do Test (8) again
709 *
710  IF( .NOT.nobal ) THEN
711  DO 170 j = 1, n
712  IF( scale( j ).NE.scale1( j ) )
713  $ result( 8 ) = ulpinv
714  170 CONTINUE
715  IF( ilo.NE.ilo1 )
716  $ result( 8 ) = ulpinv
717  IF( ihi.NE.ihi1 )
718  $ result( 8 ) = ulpinv
719  IF( abnrm.NE.abnrm1 )
720  $ result( 8 ) = ulpinv
721  END IF
722 *
723 * Do Test (9) again
724 *
725  IF( isens.EQ.2 .AND. n.GT.1 ) THEN
726  DO 180 j = 1, n
727  IF( rcondv( j ).NE.rcndv1( j ) )
728  $ result( 9 ) = ulpinv
729  180 CONTINUE
730  END IF
731 *
732  190 CONTINUE
733 *
734  200 CONTINUE
735 *
736 * If COMP, compare condition numbers to precomputed ones
737 *
738  IF( comp ) THEN
739  CALL zlacpy( 'F', n, n, a, lda, h, lda )
740  CALL zgeevx( 'N', 'V', 'V', 'B', n, h, lda, w, vl, ldvl, vr,
741  $ ldvr, ilo, ihi, scale, abnrm, rconde, rcondv,
742  $ work, lwork, rwork, iinfo )
743  IF( iinfo.NE.0 ) THEN
744  result( 1 ) = ulpinv
745  WRITE( nounit, fmt = 9999 )'ZGEEVX5', iinfo, n, iseed( 1 )
746  info = abs( iinfo )
747  GO TO 250
748  END IF
749 *
750 * Sort eigenvalues and condition numbers lexicographically
751 * to compare with inputs
752 *
753  DO 220 i = 1, n - 1
754  kmin = i
755  IF( isrt.EQ.0 ) THEN
756  vrimin = dble( w( i ) )
757  ELSE
758  vrimin = dimag( w( i ) )
759  END IF
760  DO 210 j = i + 1, n
761  IF( isrt.EQ.0 ) THEN
762  vricmp = dble( w( j ) )
763  ELSE
764  vricmp = dimag( w( j ) )
765  END IF
766  IF( vricmp.LT.vrimin ) THEN
767  kmin = j
768  vrimin = vricmp
769  END IF
770  210 CONTINUE
771  ctmp = w( kmin )
772  w( kmin ) = w( i )
773  w( i ) = ctmp
774  vrimin = rconde( kmin )
775  rconde( kmin ) = rconde( i )
776  rconde( i ) = vrimin
777  vrimin = rcondv( kmin )
778  rcondv( kmin ) = rcondv( i )
779  rcondv( i ) = vrimin
780  220 CONTINUE
781 *
782 * Compare condition numbers for eigenvectors
783 * taking their condition numbers into account
784 *
785  result( 10 ) = zero
786  eps = max( epsin, ulp )
787  v = max( dble( n )*eps*abnrm, smlnum )
788  IF( abnrm.EQ.zero )
789  $ v = one
790  DO 230 i = 1, n
791  IF( v.GT.rcondv( i )*rconde( i ) ) THEN
792  tol = rcondv( i )
793  ELSE
794  tol = v / rconde( i )
795  END IF
796  IF( v.GT.rcdvin( i )*rcdein( i ) ) THEN
797  tolin = rcdvin( i )
798  ELSE
799  tolin = v / rcdein( i )
800  END IF
801  tol = max( tol, smlnum / eps )
802  tolin = max( tolin, smlnum / eps )
803  IF( eps*( rcdvin( i )-tolin ).GT.rcondv( i )+tol ) THEN
804  vmax = one / eps
805  ELSE IF( rcdvin( i )-tolin.GT.rcondv( i )+tol ) THEN
806  vmax = ( rcdvin( i )-tolin ) / ( rcondv( i )+tol )
807  ELSE IF( rcdvin( i )+tolin.LT.eps*( rcondv( i )-tol ) ) THEN
808  vmax = one / eps
809  ELSE IF( rcdvin( i )+tolin.LT.rcondv( i )-tol ) THEN
810  vmax = ( rcondv( i )-tol ) / ( rcdvin( i )+tolin )
811  ELSE
812  vmax = one
813  END IF
814  result( 10 ) = max( result( 10 ), vmax )
815  230 CONTINUE
816 *
817 * Compare condition numbers for eigenvalues
818 * taking their condition numbers into account
819 *
820  result( 11 ) = zero
821  DO 240 i = 1, n
822  IF( v.GT.rcondv( i ) ) THEN
823  tol = one
824  ELSE
825  tol = v / rcondv( i )
826  END IF
827  IF( v.GT.rcdvin( i ) ) THEN
828  tolin = one
829  ELSE
830  tolin = v / rcdvin( i )
831  END IF
832  tol = max( tol, smlnum / eps )
833  tolin = max( tolin, smlnum / eps )
834  IF( eps*( rcdein( i )-tolin ).GT.rconde( i )+tol ) THEN
835  vmax = one / eps
836  ELSE IF( rcdein( i )-tolin.GT.rconde( i )+tol ) THEN
837  vmax = ( rcdein( i )-tolin ) / ( rconde( i )+tol )
838  ELSE IF( rcdein( i )+tolin.LT.eps*( rconde( i )-tol ) ) THEN
839  vmax = one / eps
840  ELSE IF( rcdein( i )+tolin.LT.rconde( i )-tol ) THEN
841  vmax = ( rconde( i )-tol ) / ( rcdein( i )+tolin )
842  ELSE
843  vmax = one
844  END IF
845  result( 11 ) = max( result( 11 ), vmax )
846  240 CONTINUE
847  250 CONTINUE
848 *
849  END IF
850 *
851  9999 FORMAT( ' ZGET23: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
852  $ i6, ', INPUT EXAMPLE NUMBER = ', i4 )
853  9998 FORMAT( ' ZGET23: ', a, ' returned INFO=', i6, '.', / 9x, 'N=',
854  $ i6, ', JTYPE=', i6, ', BALANC = ', a, ', ISEED=(',
855  $ 3( i5, ',' ), i5, ')' )
856 *
857  RETURN
858 *
859 * End of ZGET23
860 *
861  END
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:105
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zgeevx(BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, W, VL, LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV, WORK, LWORK, RWORK, INFO)
ZGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices ...
Definition: zgeevx.f:289
subroutine zget22(TRANSA, TRANSE, TRANSW, N, A, LDA, E, LDE, W, WORK, RWORK, RESULT)
ZGET22
Definition: zget22.f:145
subroutine zget23(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)
ZGET23
Definition: zget23.f:370