LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
zlaqr3.f
Go to the documentation of this file.
1 *> \brief \b ZLAQR3 performs the unitary similarity transformation of a Hessenberg matrix to detect and deflate fully converged eigenvalues from a trailing principal submatrix (aggressive early deflation).
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download ZLAQR3 + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlaqr3.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlaqr3.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlaqr3.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE ZLAQR3( WANTT, WANTZ, N, KTOP, KBOT, NW, H, LDH, ILOZ,
22 * IHIZ, Z, LDZ, NS, ND, SH, V, LDV, NH, T, LDT,
23 * NV, WV, LDWV, WORK, LWORK )
24 *
25 * .. Scalar Arguments ..
26 * INTEGER IHIZ, ILOZ, KBOT, KTOP, LDH, LDT, LDV, LDWV,
27 * $ LDZ, LWORK, N, ND, NH, NS, NV, NW
28 * LOGICAL WANTT, WANTZ
29 * ..
30 * .. Array Arguments ..
31 * COMPLEX*16 H( LDH, * ), SH( * ), T( LDT, * ), V( LDV, * ),
32 * $ WORK( * ), WV( LDWV, * ), Z( LDZ, * )
33 * ..
34 *
35 *
36 *> \par Purpose:
37 * =============
38 *>
39 *> \verbatim
40 *>
41 *> Aggressive early deflation:
42 *>
43 *> ZLAQR3 accepts as input an upper Hessenberg matrix
44 *> H and performs an unitary similarity transformation
45 *> designed to detect and deflate fully converged eigenvalues from
46 *> a trailing principal submatrix. On output H has been over-
47 *> written by a new Hessenberg matrix that is a perturbation of
48 *> an unitary similarity transformation of H. It is to be
49 *> hoped that the final version of H has many zero subdiagonal
50 *> entries.
51 *>
52 *> \endverbatim
53 *
54 * Arguments:
55 * ==========
56 *
57 *> \param[in] WANTT
58 *> \verbatim
59 *> WANTT is LOGICAL
60 *> If .TRUE., then the Hessenberg matrix H is fully updated
61 *> so that the triangular Schur factor may be
62 *> computed (in cooperation with the calling subroutine).
63 *> If .FALSE., then only enough of H is updated to preserve
64 *> the eigenvalues.
65 *> \endverbatim
66 *>
67 *> \param[in] WANTZ
68 *> \verbatim
69 *> WANTZ is LOGICAL
70 *> If .TRUE., then the unitary matrix Z is updated so
71 *> so that the unitary Schur factor may be computed
72 *> (in cooperation with the calling subroutine).
73 *> If .FALSE., then Z is not referenced.
74 *> \endverbatim
75 *>
76 *> \param[in] N
77 *> \verbatim
78 *> N is INTEGER
79 *> The order of the matrix H and (if WANTZ is .TRUE.) the
80 *> order of the unitary matrix Z.
81 *> \endverbatim
82 *>
83 *> \param[in] KTOP
84 *> \verbatim
85 *> KTOP is INTEGER
86 *> It is assumed that either KTOP = 1 or H(KTOP,KTOP-1)=0.
87 *> KBOT and KTOP together determine an isolated block
88 *> along the diagonal of the Hessenberg matrix.
89 *> \endverbatim
90 *>
91 *> \param[in] KBOT
92 *> \verbatim
93 *> KBOT is INTEGER
94 *> It is assumed without a check that either
95 *> KBOT = N or H(KBOT+1,KBOT)=0. KBOT and KTOP together
96 *> determine an isolated block along the diagonal of the
97 *> Hessenberg matrix.
98 *> \endverbatim
99 *>
100 *> \param[in] NW
101 *> \verbatim
102 *> NW is INTEGER
103 *> Deflation window size. 1 .LE. NW .LE. (KBOT-KTOP+1).
104 *> \endverbatim
105 *>
106 *> \param[in,out] H
107 *> \verbatim
108 *> H is COMPLEX*16 array, dimension (LDH,N)
109 *> On input the initial N-by-N section of H stores the
110 *> Hessenberg matrix undergoing aggressive early deflation.
111 *> On output H has been transformed by a unitary
112 *> similarity transformation, perturbed, and the returned
113 *> to Hessenberg form that (it is to be hoped) has some
114 *> zero subdiagonal entries.
115 *> \endverbatim
116 *>
117 *> \param[in] LDH
118 *> \verbatim
119 *> LDH is integer
120 *> Leading dimension of H just as declared in the calling
121 *> subroutine. N .LE. LDH
122 *> \endverbatim
123 *>
124 *> \param[in] ILOZ
125 *> \verbatim
126 *> ILOZ is INTEGER
127 *> \endverbatim
128 *>
129 *> \param[in] IHIZ
130 *> \verbatim
131 *> IHIZ is INTEGER
132 *> Specify the rows of Z to which transformations must be
133 *> applied if WANTZ is .TRUE.. 1 .LE. ILOZ .LE. IHIZ .LE. N.
134 *> \endverbatim
135 *>
136 *> \param[in,out] Z
137 *> \verbatim
138 *> Z is COMPLEX*16 array, dimension (LDZ,N)
139 *> IF WANTZ is .TRUE., then on output, the unitary
140 *> similarity transformation mentioned above has been
141 *> accumulated into Z(ILOZ:IHIZ,ILOZ:IHIZ) from the right.
142 *> If WANTZ is .FALSE., then Z is unreferenced.
143 *> \endverbatim
144 *>
145 *> \param[in] LDZ
146 *> \verbatim
147 *> LDZ is integer
148 *> The leading dimension of Z just as declared in the
149 *> calling subroutine. 1 .LE. LDZ.
150 *> \endverbatim
151 *>
152 *> \param[out] NS
153 *> \verbatim
154 *> NS is integer
155 *> The number of unconverged (ie approximate) eigenvalues
156 *> returned in SR and SI that may be used as shifts by the
157 *> calling subroutine.
158 *> \endverbatim
159 *>
160 *> \param[out] ND
161 *> \verbatim
162 *> ND is integer
163 *> The number of converged eigenvalues uncovered by this
164 *> subroutine.
165 *> \endverbatim
166 *>
167 *> \param[out] SH
168 *> \verbatim
169 *> SH is COMPLEX*16 array, dimension KBOT
170 *> On output, approximate eigenvalues that may
171 *> be used for shifts are stored in SH(KBOT-ND-NS+1)
172 *> through SR(KBOT-ND). Converged eigenvalues are
173 *> stored in SH(KBOT-ND+1) through SH(KBOT).
174 *> \endverbatim
175 *>
176 *> \param[out] V
177 *> \verbatim
178 *> V is COMPLEX*16 array, dimension (LDV,NW)
179 *> An NW-by-NW work array.
180 *> \endverbatim
181 *>
182 *> \param[in] LDV
183 *> \verbatim
184 *> LDV is integer scalar
185 *> The leading dimension of V just as declared in the
186 *> calling subroutine. NW .LE. LDV
187 *> \endverbatim
188 *>
189 *> \param[in] NH
190 *> \verbatim
191 *> NH is integer scalar
192 *> The number of columns of T. NH.GE.NW.
193 *> \endverbatim
194 *>
195 *> \param[out] T
196 *> \verbatim
197 *> T is COMPLEX*16 array, dimension (LDT,NW)
198 *> \endverbatim
199 *>
200 *> \param[in] LDT
201 *> \verbatim
202 *> LDT is integer
203 *> The leading dimension of T just as declared in the
204 *> calling subroutine. NW .LE. LDT
205 *> \endverbatim
206 *>
207 *> \param[in] NV
208 *> \verbatim
209 *> NV is integer
210 *> The number of rows of work array WV available for
211 *> workspace. NV.GE.NW.
212 *> \endverbatim
213 *>
214 *> \param[out] WV
215 *> \verbatim
216 *> WV is COMPLEX*16 array, dimension (LDWV,NW)
217 *> \endverbatim
218 *>
219 *> \param[in] LDWV
220 *> \verbatim
221 *> LDWV is integer
222 *> The leading dimension of W just as declared in the
223 *> calling subroutine. NW .LE. LDV
224 *> \endverbatim
225 *>
226 *> \param[out] WORK
227 *> \verbatim
228 *> WORK is COMPLEX*16 array, dimension LWORK.
229 *> On exit, WORK(1) is set to an estimate of the optimal value
230 *> of LWORK for the given values of N, NW, KTOP and KBOT.
231 *> \endverbatim
232 *>
233 *> \param[in] LWORK
234 *> \verbatim
235 *> LWORK is integer
236 *> The dimension of the work array WORK. LWORK = 2*NW
237 *> suffices, but greater efficiency may result from larger
238 *> values of LWORK.
239 *>
240 *> If LWORK = -1, then a workspace query is assumed; ZLAQR3
241 *> only estimates the optimal workspace size for the given
242 *> values of N, NW, KTOP and KBOT. The estimate is returned
243 *> in WORK(1). No error message related to LWORK is issued
244 *> by XERBLA. Neither H nor Z are accessed.
245 *> \endverbatim
246 *
247 * Authors:
248 * ========
249 *
250 *> \author Univ. of Tennessee
251 *> \author Univ. of California Berkeley
252 *> \author Univ. of Colorado Denver
253 *> \author NAG Ltd.
254 *
255 *> \date June 2016
256 *
257 *> \ingroup complex16OTHERauxiliary
258 *
259 *> \par Contributors:
260 * ==================
261 *>
262 *> Karen Braman and Ralph Byers, Department of Mathematics,
263 *> University of Kansas, USA
264 *>
265 * =====================================================================
266  SUBROUTINE zlaqr3( WANTT, WANTZ, N, KTOP, KBOT, NW, H, LDH, ILOZ,
267  $ ihiz, z, ldz, ns, nd, sh, v, ldv, nh, t, ldt,
268  $ nv, wv, ldwv, work, lwork )
269 *
270 * -- LAPACK auxiliary routine (version 3.6.1) --
271 * -- LAPACK is a software package provided by Univ. of Tennessee, --
272 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
273 * June 2016
274 *
275 * .. Scalar Arguments ..
276  INTEGER IHIZ, ILOZ, KBOT, KTOP, LDH, LDT, LDV, LDWV,
277  $ ldz, lwork, n, nd, nh, ns, nv, nw
278  LOGICAL WANTT, WANTZ
279 * ..
280 * .. Array Arguments ..
281  COMPLEX*16 H( ldh, * ), SH( * ), T( ldt, * ), V( ldv, * ),
282  $ work( * ), wv( ldwv, * ), z( ldz, * )
283 * ..
284 *
285 * ================================================================
286 *
287 * .. Parameters ..
288  COMPLEX*16 ZERO, ONE
289  parameter ( zero = ( 0.0d0, 0.0d0 ),
290  $ one = ( 1.0d0, 0.0d0 ) )
291  DOUBLE PRECISION RZERO, RONE
292  parameter ( rzero = 0.0d0, rone = 1.0d0 )
293 * ..
294 * .. Local Scalars ..
295  COMPLEX*16 BETA, CDUM, S, TAU
296  DOUBLE PRECISION FOO, SAFMAX, SAFMIN, SMLNUM, ULP
297  INTEGER I, IFST, ILST, INFO, INFQR, J, JW, KCOL, KLN,
298  $ knt, krow, kwtop, ltop, lwk1, lwk2, lwk3,
299  $ lwkopt, nmin
300 * ..
301 * .. External Functions ..
302  DOUBLE PRECISION DLAMCH
303  INTEGER ILAENV
304  EXTERNAL dlamch, ilaenv
305 * ..
306 * .. External Subroutines ..
307  EXTERNAL dlabad, zcopy, zgehrd, zgemm, zlacpy, zlahqr,
309 * ..
310 * .. Intrinsic Functions ..
311  INTRINSIC abs, dble, dcmplx, dconjg, dimag, int, max, min
312 * ..
313 * .. Statement Functions ..
314  DOUBLE PRECISION CABS1
315 * ..
316 * .. Statement Function definitions ..
317  cabs1( cdum ) = abs( dble( cdum ) ) + abs( dimag( cdum ) )
318 * ..
319 * .. Executable Statements ..
320 *
321 * ==== Estimate optimal workspace. ====
322 *
323  jw = min( nw, kbot-ktop+1 )
324  IF( jw.LE.2 ) THEN
325  lwkopt = 1
326  ELSE
327 *
328 * ==== Workspace query call to ZGEHRD ====
329 *
330  CALL zgehrd( jw, 1, jw-1, t, ldt, work, work, -1, info )
331  lwk1 = int( work( 1 ) )
332 *
333 * ==== Workspace query call to ZUNMHR ====
334 *
335  CALL zunmhr( 'R', 'N', jw, jw, 1, jw-1, t, ldt, work, v, ldv,
336  $ work, -1, info )
337  lwk2 = int( work( 1 ) )
338 *
339 * ==== Workspace query call to ZLAQR4 ====
340 *
341  CALL zlaqr4( .true., .true., jw, 1, jw, t, ldt, sh, 1, jw, v,
342  $ ldv, work, -1, infqr )
343  lwk3 = int( work( 1 ) )
344 *
345 * ==== Optimal workspace ====
346 *
347  lwkopt = max( jw+max( lwk1, lwk2 ), lwk3 )
348  END IF
349 *
350 * ==== Quick return in case of workspace query. ====
351 *
352  IF( lwork.EQ.-1 ) THEN
353  work( 1 ) = dcmplx( lwkopt, 0 )
354  RETURN
355  END IF
356 *
357 * ==== Nothing to do ...
358 * ... for an empty active block ... ====
359  ns = 0
360  nd = 0
361  work( 1 ) = one
362  IF( ktop.GT.kbot )
363  $ RETURN
364 * ... nor for an empty deflation window. ====
365  IF( nw.LT.1 )
366  $ RETURN
367 *
368 * ==== Machine constants ====
369 *
370  safmin = dlamch( 'SAFE MINIMUM' )
371  safmax = rone / safmin
372  CALL dlabad( safmin, safmax )
373  ulp = dlamch( 'PRECISION' )
374  smlnum = safmin*( dble( n ) / ulp )
375 *
376 * ==== Setup deflation window ====
377 *
378  jw = min( nw, kbot-ktop+1 )
379  kwtop = kbot - jw + 1
380  IF( kwtop.EQ.ktop ) THEN
381  s = zero
382  ELSE
383  s = h( kwtop, kwtop-1 )
384  END IF
385 *
386  IF( kbot.EQ.kwtop ) THEN
387 *
388 * ==== 1-by-1 deflation window: not much to do ====
389 *
390  sh( kwtop ) = h( kwtop, kwtop )
391  ns = 1
392  nd = 0
393  IF( cabs1( s ).LE.max( smlnum, ulp*cabs1( h( kwtop,
394  $ kwtop ) ) ) ) THEN
395  ns = 0
396  nd = 1
397  IF( kwtop.GT.ktop )
398  $ h( kwtop, kwtop-1 ) = zero
399  END IF
400  work( 1 ) = one
401  RETURN
402  END IF
403 *
404 * ==== Convert to spike-triangular form. (In case of a
405 * . rare QR failure, this routine continues to do
406 * . aggressive early deflation using that part of
407 * . the deflation window that converged using INFQR
408 * . here and there to keep track.) ====
409 *
410  CALL zlacpy( 'U', jw, jw, h( kwtop, kwtop ), ldh, t, ldt )
411  CALL zcopy( jw-1, h( kwtop+1, kwtop ), ldh+1, t( 2, 1 ), ldt+1 )
412 *
413  CALL zlaset( 'A', jw, jw, zero, one, v, ldv )
414  nmin = ilaenv( 12, 'ZLAQR3', 'SV', jw, 1, jw, lwork )
415  IF( jw.GT.nmin ) THEN
416  CALL zlaqr4( .true., .true., jw, 1, jw, t, ldt, sh( kwtop ), 1,
417  $ jw, v, ldv, work, lwork, infqr )
418  ELSE
419  CALL zlahqr( .true., .true., jw, 1, jw, t, ldt, sh( kwtop ), 1,
420  $ jw, v, ldv, infqr )
421  END IF
422 *
423 * ==== Deflation detection loop ====
424 *
425  ns = jw
426  ilst = infqr + 1
427  DO 10 knt = infqr + 1, jw
428 *
429 * ==== Small spike tip deflation test ====
430 *
431  foo = cabs1( t( ns, ns ) )
432  IF( foo.EQ.rzero )
433  $ foo = cabs1( s )
434  IF( cabs1( s )*cabs1( v( 1, ns ) ).LE.max( smlnum, ulp*foo ) )
435  $ THEN
436 *
437 * ==== One more converged eigenvalue ====
438 *
439  ns = ns - 1
440  ELSE
441 *
442 * ==== One undeflatable eigenvalue. Move it up out of the
443 * . way. (ZTREXC can not fail in this case.) ====
444 *
445  ifst = ns
446  CALL ztrexc( 'V', jw, t, ldt, v, ldv, ifst, ilst, info )
447  ilst = ilst + 1
448  END IF
449  10 CONTINUE
450 *
451 * ==== Return to Hessenberg form ====
452 *
453  IF( ns.EQ.0 )
454  $ s = zero
455 *
456  IF( ns.LT.jw ) THEN
457 *
458 * ==== sorting the diagonal of T improves accuracy for
459 * . graded matrices. ====
460 *
461  DO 30 i = infqr + 1, ns
462  ifst = i
463  DO 20 j = i + 1, ns
464  IF( cabs1( t( j, j ) ).GT.cabs1( t( ifst, ifst ) ) )
465  $ ifst = j
466  20 CONTINUE
467  ilst = i
468  IF( ifst.NE.ilst )
469  $ CALL ztrexc( 'V', jw, t, ldt, v, ldv, ifst, ilst, info )
470  30 CONTINUE
471  END IF
472 *
473 * ==== Restore shift/eigenvalue array from T ====
474 *
475  DO 40 i = infqr + 1, jw
476  sh( kwtop+i-1 ) = t( i, i )
477  40 CONTINUE
478 *
479 *
480  IF( ns.LT.jw .OR. s.EQ.zero ) THEN
481  IF( ns.GT.1 .AND. s.NE.zero ) THEN
482 *
483 * ==== Reflect spike back into lower triangle ====
484 *
485  CALL zcopy( ns, v, ldv, work, 1 )
486  DO 50 i = 1, ns
487  work( i ) = dconjg( work( i ) )
488  50 CONTINUE
489  beta = work( 1 )
490  CALL zlarfg( ns, beta, work( 2 ), 1, tau )
491  work( 1 ) = one
492 *
493  CALL zlaset( 'L', jw-2, jw-2, zero, zero, t( 3, 1 ), ldt )
494 *
495  CALL zlarf( 'L', ns, jw, work, 1, dconjg( tau ), t, ldt,
496  $ work( jw+1 ) )
497  CALL zlarf( 'R', ns, ns, work, 1, tau, t, ldt,
498  $ work( jw+1 ) )
499  CALL zlarf( 'R', jw, ns, work, 1, tau, v, ldv,
500  $ work( jw+1 ) )
501 *
502  CALL zgehrd( jw, 1, ns, t, ldt, work, work( jw+1 ),
503  $ lwork-jw, info )
504  END IF
505 *
506 * ==== Copy updated reduced window into place ====
507 *
508  IF( kwtop.GT.1 )
509  $ h( kwtop, kwtop-1 ) = s*dconjg( v( 1, 1 ) )
510  CALL zlacpy( 'U', jw, jw, t, ldt, h( kwtop, kwtop ), ldh )
511  CALL zcopy( jw-1, t( 2, 1 ), ldt+1, h( kwtop+1, kwtop ),
512  $ ldh+1 )
513 *
514 * ==== Accumulate orthogonal matrix in order update
515 * . H and Z, if requested. ====
516 *
517  IF( ns.GT.1 .AND. s.NE.zero )
518  $ CALL zunmhr( 'R', 'N', jw, ns, 1, ns, t, ldt, work, v, ldv,
519  $ work( jw+1 ), lwork-jw, info )
520 *
521 * ==== Update vertical slab in H ====
522 *
523  IF( wantt ) THEN
524  ltop = 1
525  ELSE
526  ltop = ktop
527  END IF
528  DO 60 krow = ltop, kwtop - 1, nv
529  kln = min( nv, kwtop-krow )
530  CALL zgemm( 'N', 'N', kln, jw, jw, one, h( krow, kwtop ),
531  $ ldh, v, ldv, zero, wv, ldwv )
532  CALL zlacpy( 'A', kln, jw, wv, ldwv, h( krow, kwtop ), ldh )
533  60 CONTINUE
534 *
535 * ==== Update horizontal slab in H ====
536 *
537  IF( wantt ) THEN
538  DO 70 kcol = kbot + 1, n, nh
539  kln = min( nh, n-kcol+1 )
540  CALL zgemm( 'C', 'N', jw, kln, jw, one, v, ldv,
541  $ h( kwtop, kcol ), ldh, zero, t, ldt )
542  CALL zlacpy( 'A', jw, kln, t, ldt, h( kwtop, kcol ),
543  $ ldh )
544  70 CONTINUE
545  END IF
546 *
547 * ==== Update vertical slab in Z ====
548 *
549  IF( wantz ) THEN
550  DO 80 krow = iloz, ihiz, nv
551  kln = min( nv, ihiz-krow+1 )
552  CALL zgemm( 'N', 'N', kln, jw, jw, one, z( krow, kwtop ),
553  $ ldz, v, ldv, zero, wv, ldwv )
554  CALL zlacpy( 'A', kln, jw, wv, ldwv, z( krow, kwtop ),
555  $ ldz )
556  80 CONTINUE
557  END IF
558  END IF
559 *
560 * ==== Return the number of deflations ... ====
561 *
562  nd = jw - ns
563 *
564 * ==== ... and the number of shifts. (Subtracting
565 * . INFQR from the spike length takes care
566 * . of the case of a rare QR failure while
567 * . calculating eigenvalues of the deflation
568 * . window.) ====
569 *
570  ns = ns - infqr
571 *
572 * ==== Return optimal workspace. ====
573 *
574  work( 1 ) = dcmplx( lwkopt, 0 )
575 *
576 * ==== End of ZLAQR3 ====
577 *
578  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 zcopy(N, ZX, INCX, ZY, INCY)
ZCOPY
Definition: zcopy.f:52
subroutine zlarfg(N, ALPHA, X, INCX, TAU)
ZLARFG generates an elementary reflector (Householder matrix).
Definition: zlarfg.f:108
subroutine zgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
ZGEHRD
Definition: zgehrd.f:169
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:189
subroutine zunmhr(SIDE, TRANS, M, N, ILO, IHI, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
ZUNMHR
Definition: zunmhr.f:180
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
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:76
subroutine zlaqr4(WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ, IHIZ, Z, LDZ, WORK, LWORK, INFO)
ZLAQR4 computes the eigenvalues of a Hessenberg matrix, and optionally the matrices from the Schur de...
Definition: zlaqr4.f:249
subroutine zlahqr(WANTT, WANTZ, N, ILO, IHI, H, LDH, W, ILOZ, IHIZ, Z, LDZ, INFO)
ZLAHQR computes the eigenvalues and Schur factorization of an upper Hessenberg matrix, using the double-shift/single-shift QR algorithm.
Definition: zlahqr.f:197
subroutine ztrexc(COMPQ, N, T, LDT, Q, LDQ, IFST, ILST, INFO)
ZTREXC
Definition: ztrexc.f:126
subroutine zlarf(SIDE, M, N, V, INCV, TAU, C, LDC, WORK)
ZLARF applies an elementary reflector to a general rectangular matrix.
Definition: zlarf.f:130
subroutine zlaqr3(WANTT, WANTZ, N, KTOP, KBOT, NW, H, LDH, ILOZ, IHIZ, Z, LDZ, NS, ND, SH, V, LDV, NH, T, LDT, NV, WV, LDWV, WORK, LWORK)
ZLAQR3 performs the unitary similarity transformation of a Hessenberg matrix to detect and deflate fu...
Definition: zlaqr3.f:269