LAPACK  3.9.0
LAPACK: Linear Algebra PACKage
ztrevc3.f
Go to the documentation of this file.
1 *> \brief \b ZTREVC3
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download ZTREVC3 + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ztrevc3.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ztrevc3.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ztrevc3.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE ZTREVC3( SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR,
22 * $ LDVR, MM, M, WORK, LWORK, RWORK, LRWORK, INFO)
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER HOWMNY, SIDE
26 * INTEGER INFO, LDT, LDVL, LDVR, LWORK, M, MM, N
27 * ..
28 * .. Array Arguments ..
29 * LOGICAL SELECT( * )
30 * DOUBLE PRECISION RWORK( * )
31 * COMPLEX*16 T( LDT, * ), VL( LDVL, * ), VR( LDVR, * ),
32 * $ WORK( * )
33 * ..
34 *
35 *
36 *> \par Purpose:
37 * =============
38 *>
39 *> \verbatim
40 *>
41 *> ZTREVC3 computes some or all of the right and/or left eigenvectors of
42 *> a complex upper triangular matrix T.
43 *> Matrices of this type are produced by the Schur factorization of
44 *> a complex general matrix: A = Q*T*Q**H, as computed by ZHSEQR.
45 *>
46 *> The right eigenvector x and the left eigenvector y of T corresponding
47 *> to an eigenvalue w are defined by:
48 *>
49 *> T*x = w*x, (y**H)*T = w*(y**H)
50 *>
51 *> where y**H denotes the conjugate transpose of the vector y.
52 *> The eigenvalues are not input to this routine, but are read directly
53 *> from the diagonal of T.
54 *>
55 *> This routine returns the matrices X and/or Y of right and left
56 *> eigenvectors of T, or the products Q*X and/or Q*Y, where Q is an
57 *> input matrix. If Q is the unitary factor that reduces a matrix A to
58 *> Schur form T, then Q*X and Q*Y are the matrices of right and left
59 *> eigenvectors of A.
60 *>
61 *> This uses a Level 3 BLAS version of the back transformation.
62 *> \endverbatim
63 *
64 * Arguments:
65 * ==========
66 *
67 *> \param[in] SIDE
68 *> \verbatim
69 *> SIDE is CHARACTER*1
70 *> = 'R': compute right eigenvectors only;
71 *> = 'L': compute left eigenvectors only;
72 *> = 'B': compute both right and left eigenvectors.
73 *> \endverbatim
74 *>
75 *> \param[in] HOWMNY
76 *> \verbatim
77 *> HOWMNY is CHARACTER*1
78 *> = 'A': compute all right and/or left eigenvectors;
79 *> = 'B': compute all right and/or left eigenvectors,
80 *> backtransformed using the matrices supplied in
81 *> VR and/or VL;
82 *> = 'S': compute selected right and/or left eigenvectors,
83 *> as indicated by the logical array SELECT.
84 *> \endverbatim
85 *>
86 *> \param[in] SELECT
87 *> \verbatim
88 *> SELECT is LOGICAL array, dimension (N)
89 *> If HOWMNY = 'S', SELECT specifies the eigenvectors to be
90 *> computed.
91 *> The eigenvector corresponding to the j-th eigenvalue is
92 *> computed if SELECT(j) = .TRUE..
93 *> Not referenced if HOWMNY = 'A' or 'B'.
94 *> \endverbatim
95 *>
96 *> \param[in] N
97 *> \verbatim
98 *> N is INTEGER
99 *> The order of the matrix T. N >= 0.
100 *> \endverbatim
101 *>
102 *> \param[in,out] T
103 *> \verbatim
104 *> T is COMPLEX*16 array, dimension (LDT,N)
105 *> The upper triangular matrix T. T is modified, but restored
106 *> on exit.
107 *> \endverbatim
108 *>
109 *> \param[in] LDT
110 *> \verbatim
111 *> LDT is INTEGER
112 *> The leading dimension of the array T. LDT >= max(1,N).
113 *> \endverbatim
114 *>
115 *> \param[in,out] VL
116 *> \verbatim
117 *> VL is COMPLEX*16 array, dimension (LDVL,MM)
118 *> On entry, if SIDE = 'L' or 'B' and HOWMNY = 'B', VL must
119 *> contain an N-by-N matrix Q (usually the unitary matrix Q of
120 *> Schur vectors returned by ZHSEQR).
121 *> On exit, if SIDE = 'L' or 'B', VL contains:
122 *> if HOWMNY = 'A', the matrix Y of left eigenvectors of T;
123 *> if HOWMNY = 'B', the matrix Q*Y;
124 *> if HOWMNY = 'S', the left eigenvectors of T specified by
125 *> SELECT, stored consecutively in the columns
126 *> of VL, in the same order as their
127 *> eigenvalues.
128 *> Not referenced if SIDE = 'R'.
129 *> \endverbatim
130 *>
131 *> \param[in] LDVL
132 *> \verbatim
133 *> LDVL is INTEGER
134 *> The leading dimension of the array VL.
135 *> LDVL >= 1, and if SIDE = 'L' or 'B', LDVL >= N.
136 *> \endverbatim
137 *>
138 *> \param[in,out] VR
139 *> \verbatim
140 *> VR is COMPLEX*16 array, dimension (LDVR,MM)
141 *> On entry, if SIDE = 'R' or 'B' and HOWMNY = 'B', VR must
142 *> contain an N-by-N matrix Q (usually the unitary matrix Q of
143 *> Schur vectors returned by ZHSEQR).
144 *> On exit, if SIDE = 'R' or 'B', VR contains:
145 *> if HOWMNY = 'A', the matrix X of right eigenvectors of T;
146 *> if HOWMNY = 'B', the matrix Q*X;
147 *> if HOWMNY = 'S', the right eigenvectors of T specified by
148 *> SELECT, stored consecutively in the columns
149 *> of VR, in the same order as their
150 *> eigenvalues.
151 *> Not referenced if SIDE = 'L'.
152 *> \endverbatim
153 *>
154 *> \param[in] LDVR
155 *> \verbatim
156 *> LDVR is INTEGER
157 *> The leading dimension of the array VR.
158 *> LDVR >= 1, and if SIDE = 'R' or 'B', LDVR >= N.
159 *> \endverbatim
160 *>
161 *> \param[in] MM
162 *> \verbatim
163 *> MM is INTEGER
164 *> The number of columns in the arrays VL and/or VR. MM >= M.
165 *> \endverbatim
166 *>
167 *> \param[out] M
168 *> \verbatim
169 *> M is INTEGER
170 *> The number of columns in the arrays VL and/or VR actually
171 *> used to store the eigenvectors.
172 *> If HOWMNY = 'A' or 'B', M is set to N.
173 *> Each selected eigenvector occupies one column.
174 *> \endverbatim
175 *>
176 *> \param[out] WORK
177 *> \verbatim
178 *> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
179 *> \endverbatim
180 *>
181 *> \param[in] LWORK
182 *> \verbatim
183 *> LWORK is INTEGER
184 *> The dimension of array WORK. LWORK >= max(1,2*N).
185 *> For optimum performance, LWORK >= N + 2*N*NB, where NB is
186 *> the optimal blocksize.
187 *>
188 *> If LWORK = -1, then a workspace query is assumed; the routine
189 *> only calculates the optimal size of the WORK array, returns
190 *> this value as the first entry of the WORK array, and no error
191 *> message related to LWORK is issued by XERBLA.
192 *> \endverbatim
193 *>
194 *> \param[out] RWORK
195 *> \verbatim
196 *> RWORK is DOUBLE PRECISION array, dimension (LRWORK)
197 *> \endverbatim
198 *>
199 *> \param[in] LRWORK
200 *> \verbatim
201 *> LRWORK is INTEGER
202 *> The dimension of array RWORK. LRWORK >= max(1,N).
203 *>
204 *> If LRWORK = -1, then a workspace query is assumed; the routine
205 *> only calculates the optimal size of the RWORK array, returns
206 *> this value as the first entry of the RWORK array, and no error
207 *> message related to LRWORK is issued by XERBLA.
208 *> \endverbatim
209 *>
210 *> \param[out] INFO
211 *> \verbatim
212 *> INFO is INTEGER
213 *> = 0: successful exit
214 *> < 0: if INFO = -i, the i-th argument had an illegal value
215 *> \endverbatim
216 *
217 * Authors:
218 * ========
219 *
220 *> \author Univ. of Tennessee
221 *> \author Univ. of California Berkeley
222 *> \author Univ. of Colorado Denver
223 *> \author NAG Ltd.
224 *
225 *> \date November 2017
226 *
227 * @precisions fortran z -> c
228 *
229 *> \ingroup complex16OTHERcomputational
230 *
231 *> \par Further Details:
232 * =====================
233 *>
234 *> \verbatim
235 *>
236 *> The algorithm used in this program is basically backward (forward)
237 *> substitution, with scaling to make the the code robust against
238 *> possible overflow.
239 *>
240 *> Each eigenvector is normalized so that the element of largest
241 *> magnitude has magnitude 1; here the magnitude of a complex number
242 *> (x,y) is taken to be |x| + |y|.
243 *> \endverbatim
244 *>
245 * =====================================================================
246  SUBROUTINE ztrevc3( SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR,
247  $ LDVR, MM, M, WORK, LWORK, RWORK, LRWORK, INFO)
248  IMPLICIT NONE
249 *
250 * -- LAPACK computational routine (version 3.8.0) --
251 * -- LAPACK is a software package provided by Univ. of Tennessee, --
252 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
253 * November 2017
254 *
255 * .. Scalar Arguments ..
256  CHARACTER HOWMNY, SIDE
257  INTEGER INFO, LDT, LDVL, LDVR, LWORK, LRWORK, M, MM, N
258 * ..
259 * .. Array Arguments ..
260  LOGICAL SELECT( * )
261  DOUBLE PRECISION RWORK( * )
262  COMPLEX*16 T( LDT, * ), VL( LDVL, * ), VR( LDVR, * ),
263  $ work( * )
264 * ..
265 *
266 * =====================================================================
267 *
268 * .. Parameters ..
269  DOUBLE PRECISION ZERO, ONE
270  parameter( zero = 0.0d+0, one = 1.0d+0 )
271  COMPLEX*16 CZERO, CONE
272  parameter( czero = ( 0.0d+0, 0.0d+0 ),
273  $ cone = ( 1.0d+0, 0.0d+0 ) )
274  INTEGER NBMIN, NBMAX
275  parameter( nbmin = 8, nbmax = 128 )
276 * ..
277 * .. Local Scalars ..
278  LOGICAL ALLV, BOTHV, LEFTV, LQUERY, OVER, RIGHTV, SOMEV
279  INTEGER I, II, IS, J, K, KI, IV, MAXWRK, NB
280  DOUBLE PRECISION OVFL, REMAX, SCALE, SMIN, SMLNUM, ULP, UNFL
281  COMPLEX*16 CDUM
282 * ..
283 * .. External Functions ..
284  LOGICAL LSAME
285  INTEGER ILAENV, IZAMAX
286  DOUBLE PRECISION DLAMCH, DZASUM
287  EXTERNAL lsame, ilaenv, izamax, dlamch, dzasum
288 * ..
289 * .. External Subroutines ..
290  EXTERNAL xerbla, zcopy, zdscal, zgemv, zlatrs,
291  $ zgemm, dlabad, zlaset, zlacpy
292 * ..
293 * .. Intrinsic Functions ..
294  INTRINSIC abs, dble, dcmplx, conjg, aimag, max
295 * ..
296 * .. Statement Functions ..
297  DOUBLE PRECISION CABS1
298 * ..
299 * .. Statement Function definitions ..
300  cabs1( cdum ) = abs( dble( cdum ) ) + abs( aimag( cdum ) )
301 * ..
302 * .. Executable Statements ..
303 *
304 * Decode and test the input parameters
305 *
306  bothv = lsame( side, 'B' )
307  rightv = lsame( side, 'R' ) .OR. bothv
308  leftv = lsame( side, 'L' ) .OR. bothv
309 *
310  allv = lsame( howmny, 'A' )
311  over = lsame( howmny, 'B' )
312  somev = lsame( howmny, 'S' )
313 *
314 * Set M to the number of columns required to store the selected
315 * eigenvectors.
316 *
317  IF( somev ) THEN
318  m = 0
319  DO 10 j = 1, n
320  IF( SELECT( j ) )
321  $ m = m + 1
322  10 CONTINUE
323  ELSE
324  m = n
325  END IF
326 *
327  info = 0
328  nb = ilaenv( 1, 'ZTREVC', side // howmny, n, -1, -1, -1 )
329  maxwrk = n + 2*n*nb
330  work(1) = maxwrk
331  rwork(1) = n
332  lquery = ( lwork.EQ.-1 .OR. lrwork.EQ.-1 )
333  IF( .NOT.rightv .AND. .NOT.leftv ) THEN
334  info = -1
335  ELSE IF( .NOT.allv .AND. .NOT.over .AND. .NOT.somev ) THEN
336  info = -2
337  ELSE IF( n.LT.0 ) THEN
338  info = -4
339  ELSE IF( ldt.LT.max( 1, n ) ) THEN
340  info = -6
341  ELSE IF( ldvl.LT.1 .OR. ( leftv .AND. ldvl.LT.n ) ) THEN
342  info = -8
343  ELSE IF( ldvr.LT.1 .OR. ( rightv .AND. ldvr.LT.n ) ) THEN
344  info = -10
345  ELSE IF( mm.LT.m ) THEN
346  info = -11
347  ELSE IF( lwork.LT.max( 1, 2*n ) .AND. .NOT.lquery ) THEN
348  info = -14
349  ELSE IF ( lrwork.LT.max( 1, n ) .AND. .NOT.lquery ) THEN
350  info = -16
351  END IF
352  IF( info.NE.0 ) THEN
353  CALL xerbla( 'ZTREVC3', -info )
354  RETURN
355  ELSE IF( lquery ) THEN
356  RETURN
357  END IF
358 *
359 * Quick return if possible.
360 *
361  IF( n.EQ.0 )
362  $ RETURN
363 *
364 * Use blocked version of back-transformation if sufficient workspace.
365 * Zero-out the workspace to avoid potential NaN propagation.
366 *
367  IF( over .AND. lwork .GE. n + 2*n*nbmin ) THEN
368  nb = (lwork - n) / (2*n)
369  nb = min( nb, nbmax )
370  CALL zlaset( 'F', n, 1+2*nb, czero, czero, work, n )
371  ELSE
372  nb = 1
373  END IF
374 *
375 * Set the constants to control overflow.
376 *
377  unfl = dlamch( 'Safe minimum' )
378  ovfl = one / unfl
379  CALL dlabad( unfl, ovfl )
380  ulp = dlamch( 'Precision' )
381  smlnum = unfl*( n / ulp )
382 *
383 * Store the diagonal elements of T in working array WORK.
384 *
385  DO 20 i = 1, n
386  work( i ) = t( i, i )
387  20 CONTINUE
388 *
389 * Compute 1-norm of each column of strictly upper triangular
390 * part of T to control overflow in triangular solver.
391 *
392  rwork( 1 ) = zero
393  DO 30 j = 2, n
394  rwork( j ) = dzasum( j-1, t( 1, j ), 1 )
395  30 CONTINUE
396 *
397  IF( rightv ) THEN
398 *
399 * ============================================================
400 * Compute right eigenvectors.
401 *
402 * IV is index of column in current block.
403 * Non-blocked version always uses IV=NB=1;
404 * blocked version starts with IV=NB, goes down to 1.
405 * (Note the "0-th" column is used to store the original diagonal.)
406  iv = nb
407  is = m
408  DO 80 ki = n, 1, -1
409  IF( somev ) THEN
410  IF( .NOT.SELECT( ki ) )
411  $ GO TO 80
412  END IF
413  smin = max( ulp*( cabs1( t( ki, ki ) ) ), smlnum )
414 *
415 * --------------------------------------------------------
416 * Complex right eigenvector
417 *
418  work( ki + iv*n ) = cone
419 *
420 * Form right-hand side.
421 *
422  DO 40 k = 1, ki - 1
423  work( k + iv*n ) = -t( k, ki )
424  40 CONTINUE
425 *
426 * Solve upper triangular system:
427 * [ T(1:KI-1,1:KI-1) - T(KI,KI) ]*X = SCALE*WORK.
428 *
429  DO 50 k = 1, ki - 1
430  t( k, k ) = t( k, k ) - t( ki, ki )
431  IF( cabs1( t( k, k ) ).LT.smin )
432  $ t( k, k ) = smin
433  50 CONTINUE
434 *
435  IF( ki.GT.1 ) THEN
436  CALL zlatrs( 'Upper', 'No transpose', 'Non-unit', 'Y',
437  $ ki-1, t, ldt, work( 1 + iv*n ), scale,
438  $ rwork, info )
439  work( ki + iv*n ) = scale
440  END IF
441 *
442 * Copy the vector x or Q*x to VR and normalize.
443 *
444  IF( .NOT.over ) THEN
445 * ------------------------------
446 * no back-transform: copy x to VR and normalize.
447  CALL zcopy( ki, work( 1 + iv*n ), 1, vr( 1, is ), 1 )
448 *
449  ii = izamax( ki, vr( 1, is ), 1 )
450  remax = one / cabs1( vr( ii, is ) )
451  CALL zdscal( ki, remax, vr( 1, is ), 1 )
452 *
453  DO 60 k = ki + 1, n
454  vr( k, is ) = czero
455  60 CONTINUE
456 *
457  ELSE IF( nb.EQ.1 ) THEN
458 * ------------------------------
459 * version 1: back-transform each vector with GEMV, Q*x.
460  IF( ki.GT.1 )
461  $ CALL zgemv( 'N', n, ki-1, cone, vr, ldvr,
462  $ work( 1 + iv*n ), 1, dcmplx( scale ),
463  $ vr( 1, ki ), 1 )
464 *
465  ii = izamax( n, vr( 1, ki ), 1 )
466  remax = one / cabs1( vr( ii, ki ) )
467  CALL zdscal( n, remax, vr( 1, ki ), 1 )
468 *
469  ELSE
470 * ------------------------------
471 * version 2: back-transform block of vectors with GEMM
472 * zero out below vector
473  DO k = ki + 1, n
474  work( k + iv*n ) = czero
475  END DO
476 *
477 * Columns IV:NB of work are valid vectors.
478 * When the number of vectors stored reaches NB,
479 * or if this was last vector, do the GEMM
480  IF( (iv.EQ.1) .OR. (ki.EQ.1) ) THEN
481  CALL zgemm( 'N', 'N', n, nb-iv+1, ki+nb-iv, cone,
482  $ vr, ldvr,
483  $ work( 1 + (iv)*n ), n,
484  $ czero,
485  $ work( 1 + (nb+iv)*n ), n )
486 * normalize vectors
487  DO k = iv, nb
488  ii = izamax( n, work( 1 + (nb+k)*n ), 1 )
489  remax = one / cabs1( work( ii + (nb+k)*n ) )
490  CALL zdscal( n, remax, work( 1 + (nb+k)*n ), 1 )
491  END DO
492  CALL zlacpy( 'F', n, nb-iv+1,
493  $ work( 1 + (nb+iv)*n ), n,
494  $ vr( 1, ki ), ldvr )
495  iv = nb
496  ELSE
497  iv = iv - 1
498  END IF
499  END IF
500 *
501 * Restore the original diagonal elements of T.
502 *
503  DO 70 k = 1, ki - 1
504  t( k, k ) = work( k )
505  70 CONTINUE
506 *
507  is = is - 1
508  80 CONTINUE
509  END IF
510 *
511  IF( leftv ) THEN
512 *
513 * ============================================================
514 * Compute left eigenvectors.
515 *
516 * IV is index of column in current block.
517 * Non-blocked version always uses IV=1;
518 * blocked version starts with IV=1, goes up to NB.
519 * (Note the "0-th" column is used to store the original diagonal.)
520  iv = 1
521  is = 1
522  DO 130 ki = 1, n
523 *
524  IF( somev ) THEN
525  IF( .NOT.SELECT( ki ) )
526  $ GO TO 130
527  END IF
528  smin = max( ulp*( cabs1( t( ki, ki ) ) ), smlnum )
529 *
530 * --------------------------------------------------------
531 * Complex left eigenvector
532 *
533  work( ki + iv*n ) = cone
534 *
535 * Form right-hand side.
536 *
537  DO 90 k = ki + 1, n
538  work( k + iv*n ) = -conjg( t( ki, k ) )
539  90 CONTINUE
540 *
541 * Solve conjugate-transposed triangular system:
542 * [ T(KI+1:N,KI+1:N) - T(KI,KI) ]**H * X = SCALE*WORK.
543 *
544  DO 100 k = ki + 1, n
545  t( k, k ) = t( k, k ) - t( ki, ki )
546  IF( cabs1( t( k, k ) ).LT.smin )
547  $ t( k, k ) = smin
548  100 CONTINUE
549 *
550  IF( ki.LT.n ) THEN
551  CALL zlatrs( 'Upper', 'Conjugate transpose', 'Non-unit',
552  $ 'Y', n-ki, t( ki+1, ki+1 ), ldt,
553  $ work( ki+1 + iv*n ), scale, rwork, info )
554  work( ki + iv*n ) = scale
555  END IF
556 *
557 * Copy the vector x or Q*x to VL and normalize.
558 *
559  IF( .NOT.over ) THEN
560 * ------------------------------
561 * no back-transform: copy x to VL and normalize.
562  CALL zcopy( n-ki+1, work( ki + iv*n ), 1, vl(ki,is), 1 )
563 *
564  ii = izamax( n-ki+1, vl( ki, is ), 1 ) + ki - 1
565  remax = one / cabs1( vl( ii, is ) )
566  CALL zdscal( n-ki+1, remax, vl( ki, is ), 1 )
567 *
568  DO 110 k = 1, ki - 1
569  vl( k, is ) = czero
570  110 CONTINUE
571 *
572  ELSE IF( nb.EQ.1 ) THEN
573 * ------------------------------
574 * version 1: back-transform each vector with GEMV, Q*x.
575  IF( ki.LT.n )
576  $ CALL zgemv( 'N', n, n-ki, cone, vl( 1, ki+1 ), ldvl,
577  $ work( ki+1 + iv*n ), 1, dcmplx( scale ),
578  $ vl( 1, ki ), 1 )
579 *
580  ii = izamax( n, vl( 1, ki ), 1 )
581  remax = one / cabs1( vl( ii, ki ) )
582  CALL zdscal( n, remax, vl( 1, ki ), 1 )
583 *
584  ELSE
585 * ------------------------------
586 * version 2: back-transform block of vectors with GEMM
587 * zero out above vector
588 * could go from KI-NV+1 to KI-1
589  DO k = 1, ki - 1
590  work( k + iv*n ) = czero
591  END DO
592 *
593 * Columns 1:IV of work are valid vectors.
594 * When the number of vectors stored reaches NB,
595 * or if this was last vector, do the GEMM
596  IF( (iv.EQ.nb) .OR. (ki.EQ.n) ) THEN
597  CALL zgemm( 'N', 'N', n, iv, n-ki+iv, cone,
598  $ vl( 1, ki-iv+1 ), ldvl,
599  $ work( ki-iv+1 + (1)*n ), n,
600  $ czero,
601  $ work( 1 + (nb+1)*n ), n )
602 * normalize vectors
603  DO k = 1, iv
604  ii = izamax( n, work( 1 + (nb+k)*n ), 1 )
605  remax = one / cabs1( work( ii + (nb+k)*n ) )
606  CALL zdscal( n, remax, work( 1 + (nb+k)*n ), 1 )
607  END DO
608  CALL zlacpy( 'F', n, iv,
609  $ work( 1 + (nb+1)*n ), n,
610  $ vl( 1, ki-iv+1 ), ldvl )
611  iv = 1
612  ELSE
613  iv = iv + 1
614  END IF
615  END IF
616 *
617 * Restore the original diagonal elements of T.
618 *
619  DO 120 k = ki + 1, n
620  t( k, k ) = work( k )
621  120 CONTINUE
622 *
623  is = is + 1
624  130 CONTINUE
625  END IF
626 *
627  RETURN
628 *
629 * End of ZTREVC3
630 *
631  END
zgemv
subroutine zgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
ZGEMV
Definition: zgemv.f:160
zdscal
subroutine zdscal(N, DA, ZX, INCX)
ZDSCAL
Definition: zdscal.f:80
zcopy
subroutine zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:83
zgemm
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:189
zlatrs
subroutine zlatrs(UPLO, TRANS, DIAG, NORMIN, N, A, LDA, X, SCALE, CNORM, INFO)
ZLATRS solves a triangular system of equations with the scale factor set to prevent overflow.
Definition: zlatrs.f:241
zlaset
subroutine zlaset(UPLO, M, N, ALPHA, BETA, A, LDA)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition: zlaset.f:108
zlacpy
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
xerbla
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
dlabad
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:76
ztrevc3
subroutine ztrevc3(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, LWORK, RWORK, LRWORK, INFO)
ZTREVC3
Definition: ztrevc3.f:248