LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
zggev3.f
Go to the documentation of this file.
1 *> \brief <b> ZGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices (blocked algorithm)</b>
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download ZGGEV3 + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zggev3.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zggev3.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zggev3.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE ZGGEV3( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA,
22 * VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO )
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER JOBVL, JOBVR
26 * INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
27 * ..
28 * .. Array Arguments ..
29 * DOUBLE PRECISION RWORK( * )
30 * COMPLEX*16 A( LDA, * ), ALPHA( * ), B( LDB, * ),
31 * $ BETA( * ), VL( LDVL, * ), VR( LDVR, * ),
32 * $ WORK( * )
33 * ..
34 *
35 *
36 *> \par Purpose:
37 * =============
38 *>
39 *> \verbatim
40 *>
41 *> ZGGEV3 computes for a pair of N-by-N complex nonsymmetric matrices
42 *> (A,B), the generalized eigenvalues, and optionally, the left and/or
43 *> right generalized eigenvectors.
44 *>
45 *> A generalized eigenvalue for a pair of matrices (A,B) is a scalar
46 *> lambda or a ratio alpha/beta = lambda, such that A - lambda*B is
47 *> singular. It is usually represented as the pair (alpha,beta), as
48 *> there is a reasonable interpretation for beta=0, and even for both
49 *> being zero.
50 *>
51 *> The right generalized eigenvector v(j) corresponding to the
52 *> generalized eigenvalue lambda(j) of (A,B) satisfies
53 *>
54 *> A * v(j) = lambda(j) * B * v(j).
55 *>
56 *> The left generalized eigenvector u(j) corresponding to the
57 *> generalized eigenvalues lambda(j) of (A,B) satisfies
58 *>
59 *> u(j)**H * A = lambda(j) * u(j)**H * B
60 *>
61 *> where u(j)**H is the conjugate-transpose of u(j).
62 *> \endverbatim
63 *
64 * Arguments:
65 * ==========
66 *
67 *> \param[in] JOBVL
68 *> \verbatim
69 *> JOBVL is CHARACTER*1
70 *> = 'N': do not compute the left generalized eigenvectors;
71 *> = 'V': compute the left generalized eigenvectors.
72 *> \endverbatim
73 *>
74 *> \param[in] JOBVR
75 *> \verbatim
76 *> JOBVR is CHARACTER*1
77 *> = 'N': do not compute the right generalized eigenvectors;
78 *> = 'V': compute the right generalized eigenvectors.
79 *> \endverbatim
80 *>
81 *> \param[in] N
82 *> \verbatim
83 *> N is INTEGER
84 *> The order of the matrices A, B, VL, and VR. N >= 0.
85 *> \endverbatim
86 *>
87 *> \param[in,out] A
88 *> \verbatim
89 *> A is COMPLEX*16 array, dimension (LDA, N)
90 *> On entry, the matrix A in the pair (A,B).
91 *> On exit, A has been overwritten.
92 *> \endverbatim
93 *>
94 *> \param[in] LDA
95 *> \verbatim
96 *> LDA is INTEGER
97 *> The leading dimension of A. LDA >= max(1,N).
98 *> \endverbatim
99 *>
100 *> \param[in,out] B
101 *> \verbatim
102 *> B is COMPLEX*16 array, dimension (LDB, N)
103 *> On entry, the matrix B in the pair (A,B).
104 *> On exit, B has been overwritten.
105 *> \endverbatim
106 *>
107 *> \param[in] LDB
108 *> \verbatim
109 *> LDB is INTEGER
110 *> The leading dimension of B. LDB >= max(1,N).
111 *> \endverbatim
112 *>
113 *> \param[out] ALPHA
114 *> \verbatim
115 *> ALPHA is COMPLEX*16 array, dimension (N)
116 *> \endverbatim
117 *>
118 *> \param[out] BETA
119 *> \verbatim
120 *> BETA is COMPLEX*16 array, dimension (N)
121 *> On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the
122 *> generalized eigenvalues.
123 *>
124 *> Note: the quotients ALPHA(j)/BETA(j) may easily over- or
125 *> underflow, and BETA(j) may even be zero. Thus, the user
126 *> should avoid naively computing the ratio alpha/beta.
127 *> However, ALPHA will be always less than and usually
128 *> comparable with norm(A) in magnitude, and BETA always less
129 *> than and usually comparable with norm(B).
130 *> \endverbatim
131 *>
132 *> \param[out] VL
133 *> \verbatim
134 *> VL is COMPLEX*16 array, dimension (LDVL,N)
135 *> If JOBVL = 'V', the left generalized eigenvectors u(j) are
136 *> stored one after another in the columns of VL, in the same
137 *> order as their eigenvalues.
138 *> Each eigenvector is scaled so the largest component has
139 *> abs(real part) + abs(imag. part) = 1.
140 *> Not referenced if JOBVL = 'N'.
141 *> \endverbatim
142 *>
143 *> \param[in] LDVL
144 *> \verbatim
145 *> LDVL is INTEGER
146 *> The leading dimension of the matrix VL. LDVL >= 1, and
147 *> if JOBVL = 'V', LDVL >= N.
148 *> \endverbatim
149 *>
150 *> \param[out] VR
151 *> \verbatim
152 *> VR is COMPLEX*16 array, dimension (LDVR,N)
153 *> If JOBVR = 'V', the right generalized eigenvectors v(j) are
154 *> stored one after another in the columns of VR, in the same
155 *> order as their eigenvalues.
156 *> Each eigenvector is scaled so the largest component has
157 *> abs(real part) + abs(imag. part) = 1.
158 *> Not referenced if JOBVR = 'N'.
159 *> \endverbatim
160 *>
161 *> \param[in] LDVR
162 *> \verbatim
163 *> LDVR is INTEGER
164 *> The leading dimension of the matrix VR. LDVR >= 1, and
165 *> if JOBVR = 'V', LDVR >= N.
166 *> \endverbatim
167 *>
168 *> \param[out] WORK
169 *> \verbatim
170 *> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
171 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
172 *> \endverbatim
173 *>
174 *> \param[in] LWORK
175 *> \verbatim
176 *> LWORK is INTEGER
177 *> The dimension of the array WORK.
178 *>
179 *> If LWORK = -1, then a workspace query is assumed; the routine
180 *> only calculates the optimal size of the WORK array, returns
181 *> this value as the first entry of the WORK array, and no error
182 *> message related to LWORK is issued by XERBLA.
183 *> \endverbatim
184 *>
185 *> \param[out] RWORK
186 *> \verbatim
187 *> RWORK is DOUBLE PRECISION array, dimension (8*N)
188 *> \endverbatim
189 *>
190 *> \param[out] INFO
191 *> \verbatim
192 *> INFO is INTEGER
193 *> = 0: successful exit
194 *> < 0: if INFO = -i, the i-th argument had an illegal value.
195 *> =1,...,N:
196 *> The QZ iteration failed. No eigenvectors have been
197 *> calculated, but ALPHA(j) and BETA(j) should be
198 *> correct for j=INFO+1,...,N.
199 *> > N: =N+1: other then QZ iteration failed in DHGEQZ,
200 *> =N+2: error return from DTGEVC.
201 *> \endverbatim
202 *
203 * Authors:
204 * ========
205 *
206 *> \author Univ. of Tennessee
207 *> \author Univ. of California Berkeley
208 *> \author Univ. of Colorado Denver
209 *> \author NAG Ltd.
210 *
211 *> \date January 2015
212 *
213 *> \ingroup complex16GEeigen
214 *
215 * =====================================================================
216  SUBROUTINE zggev3( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA,
217  $ vl, ldvl, vr, ldvr, work, lwork, rwork, info )
218 *
219 * -- LAPACK driver routine (version 3.6.1) --
220 * -- LAPACK is a software package provided by Univ. of Tennessee, --
221 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
222 * January 2015
223 *
224 * .. Scalar Arguments ..
225  CHARACTER JOBVL, JOBVR
226  INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
227 * ..
228 * .. Array Arguments ..
229  DOUBLE PRECISION RWORK( * )
230  COMPLEX*16 A( lda, * ), ALPHA( * ), B( ldb, * ),
231  $ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
232  $ work( * )
233 * ..
234 *
235 * =====================================================================
236 *
237 * .. Parameters ..
238  DOUBLE PRECISION ZERO, ONE
239  parameter ( zero = 0.0d0, one = 1.0d0 )
240  COMPLEX*16 CZERO, CONE
241  parameter ( czero = ( 0.0d0, 0.0d0 ),
242  $ cone = ( 1.0d0, 0.0d0 ) )
243 * ..
244 * .. Local Scalars ..
245  LOGICAL ILASCL, ILBSCL, ILV, ILVL, ILVR, LQUERY
246  CHARACTER CHTEMP
247  INTEGER ICOLS, IERR, IHI, IJOBVL, IJOBVR, ILEFT, ILO,
248  $ in, iright, irows, irwrk, itau, iwrk, jc, jr,
249  $ lwkopt
250  DOUBLE PRECISION ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS,
251  $ smlnum, temp
252  COMPLEX*16 X
253 * ..
254 * .. Local Arrays ..
255  LOGICAL LDUMMA( 1 )
256 * ..
257 * .. External Subroutines ..
258  EXTERNAL dlabad, xerbla, zgeqrf, zggbak, zggbal, zgghd3,
260  $ zunmqr
261 * ..
262 * .. External Functions ..
263  LOGICAL LSAME
264  DOUBLE PRECISION DLAMCH, ZLANGE
265  EXTERNAL lsame, dlamch, zlange
266 * ..
267 * .. Intrinsic Functions ..
268  INTRINSIC abs, dble, dimag, max, sqrt
269 * ..
270 * .. Statement Functions ..
271  DOUBLE PRECISION ABS1
272 * ..
273 * .. Statement Function definitions ..
274  abs1( x ) = abs( dble( x ) ) + abs( dimag( x ) )
275 * ..
276 * .. Executable Statements ..
277 *
278 * Decode the input arguments
279 *
280  IF( lsame( jobvl, 'N' ) ) THEN
281  ijobvl = 1
282  ilvl = .false.
283  ELSE IF( lsame( jobvl, 'V' ) ) THEN
284  ijobvl = 2
285  ilvl = .true.
286  ELSE
287  ijobvl = -1
288  ilvl = .false.
289  END IF
290 *
291  IF( lsame( jobvr, 'N' ) ) THEN
292  ijobvr = 1
293  ilvr = .false.
294  ELSE IF( lsame( jobvr, 'V' ) ) THEN
295  ijobvr = 2
296  ilvr = .true.
297  ELSE
298  ijobvr = -1
299  ilvr = .false.
300  END IF
301  ilv = ilvl .OR. ilvr
302 *
303 * Test the input arguments
304 *
305  info = 0
306  lquery = ( lwork.EQ.-1 )
307  IF( ijobvl.LE.0 ) THEN
308  info = -1
309  ELSE IF( ijobvr.LE.0 ) THEN
310  info = -2
311  ELSE IF( n.LT.0 ) THEN
312  info = -3
313  ELSE IF( lda.LT.max( 1, n ) ) THEN
314  info = -5
315  ELSE IF( ldb.LT.max( 1, n ) ) THEN
316  info = -7
317  ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
318  info = -11
319  ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
320  info = -13
321  ELSE IF( lwork.LT.max( 1, 2*n ) .AND. .NOT.lquery ) THEN
322  info = -15
323  END IF
324 *
325 * Compute workspace
326 *
327  IF( info.EQ.0 ) THEN
328  CALL zgeqrf( n, n, b, ldb, work, work, -1, ierr )
329  lwkopt = max( 1, n+int( work( 1 ) ) )
330  CALL zunmqr( 'L', 'C', n, n, n, b, ldb, work, a, lda, work,
331  $ -1, ierr )
332  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
333  IF( ilvl ) THEN
334  CALL zungqr( n, n, n, vl, ldvl, work, work, -1, ierr )
335  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
336  END IF
337  IF( ilv ) THEN
338  CALL zgghd3( jobvl, jobvr, n, 1, n, a, lda, b, ldb, vl,
339  $ ldvl, vr, ldvr, work, -1, ierr )
340  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
341  CALL zhgeqz( 'S', jobvl, jobvr, n, 1, n, a, lda, b, ldb,
342  $ alpha, beta, vl, ldvl, vr, ldvr, work, -1,
343  $ rwork, ierr )
344  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
345  ELSE
346  CALL zgghd3( jobvl, jobvr, n, 1, n, a, lda, b, ldb, vl,
347  $ ldvl, vr, ldvr, work, -1, ierr )
348  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
349  CALL zhgeqz( 'E', jobvl, jobvr, n, 1, n, a, lda, b, ldb,
350  $ alpha, beta, vl, ldvl, vr, ldvr, work, -1,
351  $ rwork, ierr )
352  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
353  END IF
354  work( 1 ) = dcmplx( lwkopt )
355  END IF
356 *
357  IF( info.NE.0 ) THEN
358  CALL xerbla( 'ZGGEV3 ', -info )
359  RETURN
360  ELSE IF( lquery ) THEN
361  RETURN
362  END IF
363 *
364 * Quick return if possible
365 *
366  IF( n.EQ.0 )
367  $ RETURN
368 *
369 * Get machine constants
370 *
371  eps = dlamch( 'E' )*dlamch( 'B' )
372  smlnum = dlamch( 'S' )
373  bignum = one / smlnum
374  CALL dlabad( smlnum, bignum )
375  smlnum = sqrt( smlnum ) / eps
376  bignum = one / smlnum
377 *
378 * Scale A if max element outside range [SMLNUM,BIGNUM]
379 *
380  anrm = zlange( 'M', n, n, a, lda, rwork )
381  ilascl = .false.
382  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
383  anrmto = smlnum
384  ilascl = .true.
385  ELSE IF( anrm.GT.bignum ) THEN
386  anrmto = bignum
387  ilascl = .true.
388  END IF
389  IF( ilascl )
390  $ CALL zlascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
391 *
392 * Scale B if max element outside range [SMLNUM,BIGNUM]
393 *
394  bnrm = zlange( 'M', n, n, b, ldb, rwork )
395  ilbscl = .false.
396  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
397  bnrmto = smlnum
398  ilbscl = .true.
399  ELSE IF( bnrm.GT.bignum ) THEN
400  bnrmto = bignum
401  ilbscl = .true.
402  END IF
403  IF( ilbscl )
404  $ CALL zlascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
405 *
406 * Permute the matrices A, B to isolate eigenvalues if possible
407 *
408  ileft = 1
409  iright = n + 1
410  irwrk = iright + n
411  CALL zggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
412  $ rwork( iright ), rwork( irwrk ), ierr )
413 *
414 * Reduce B to triangular form (QR decomposition of B)
415 *
416  irows = ihi + 1 - ilo
417  IF( ilv ) THEN
418  icols = n + 1 - ilo
419  ELSE
420  icols = irows
421  END IF
422  itau = 1
423  iwrk = itau + irows
424  CALL zgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
425  $ work( iwrk ), lwork+1-iwrk, ierr )
426 *
427 * Apply the orthogonal transformation to matrix A
428 *
429  CALL zunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
430  $ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
431  $ lwork+1-iwrk, ierr )
432 *
433 * Initialize VL
434 *
435  IF( ilvl ) THEN
436  CALL zlaset( 'Full', n, n, czero, cone, vl, ldvl )
437  IF( irows.GT.1 ) THEN
438  CALL zlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
439  $ vl( ilo+1, ilo ), ldvl )
440  END IF
441  CALL zungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
442  $ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
443  END IF
444 *
445 * Initialize VR
446 *
447  IF( ilvr )
448  $ CALL zlaset( 'Full', n, n, czero, cone, vr, ldvr )
449 *
450 * Reduce to generalized Hessenberg form
451 *
452  IF( ilv ) THEN
453 *
454 * Eigenvectors requested -- work on whole matrix.
455 *
456  CALL zgghd3( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
457  $ ldvl, vr, ldvr, work( iwrk ), lwork+1-iwrk, ierr )
458  ELSE
459  CALL zgghd3( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
460  $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr,
461  $ work( iwrk ), lwork+1-iwrk, ierr )
462  END IF
463 *
464 * Perform QZ algorithm (Compute eigenvalues, and optionally, the
465 * Schur form and Schur vectors)
466 *
467  iwrk = itau
468  IF( ilv ) THEN
469  chtemp = 'S'
470  ELSE
471  chtemp = 'E'
472  END IF
473  CALL zhgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
474  $ alpha, beta, vl, ldvl, vr, ldvr, work( iwrk ),
475  $ lwork+1-iwrk, rwork( irwrk ), ierr )
476  IF( ierr.NE.0 ) THEN
477  IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
478  info = ierr
479  ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
480  info = ierr - n
481  ELSE
482  info = n + 1
483  END IF
484  GO TO 70
485  END IF
486 *
487 * Compute Eigenvectors
488 *
489  IF( ilv ) THEN
490  IF( ilvl ) THEN
491  IF( ilvr ) THEN
492  chtemp = 'B'
493  ELSE
494  chtemp = 'L'
495  END IF
496  ELSE
497  chtemp = 'R'
498  END IF
499 *
500  CALL ztgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
501  $ vr, ldvr, n, in, work( iwrk ), rwork( irwrk ),
502  $ ierr )
503  IF( ierr.NE.0 ) THEN
504  info = n + 2
505  GO TO 70
506  END IF
507 *
508 * Undo balancing on VL and VR and normalization
509 *
510  IF( ilvl ) THEN
511  CALL zggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
512  $ rwork( iright ), n, vl, ldvl, ierr )
513  DO 30 jc = 1, n
514  temp = zero
515  DO 10 jr = 1, n
516  temp = max( temp, abs1( vl( jr, jc ) ) )
517  10 CONTINUE
518  IF( temp.LT.smlnum )
519  $ GO TO 30
520  temp = one / temp
521  DO 20 jr = 1, n
522  vl( jr, jc ) = vl( jr, jc )*temp
523  20 CONTINUE
524  30 CONTINUE
525  END IF
526  IF( ilvr ) THEN
527  CALL zggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
528  $ rwork( iright ), n, vr, ldvr, ierr )
529  DO 60 jc = 1, n
530  temp = zero
531  DO 40 jr = 1, n
532  temp = max( temp, abs1( vr( jr, jc ) ) )
533  40 CONTINUE
534  IF( temp.LT.smlnum )
535  $ GO TO 60
536  temp = one / temp
537  DO 50 jr = 1, n
538  vr( jr, jc ) = vr( jr, jc )*temp
539  50 CONTINUE
540  60 CONTINUE
541  END IF
542  END IF
543 *
544 * Undo scaling if necessary
545 *
546  70 CONTINUE
547 *
548  IF( ilascl )
549  $ CALL zlascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
550 *
551  IF( ilbscl )
552  $ CALL zlascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
553 *
554  work( 1 ) = dcmplx( lwkopt )
555  RETURN
556 *
557 * End of ZGGEV3
558 *
559  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 zgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
ZGEQRF VARIANT: left-looking Level 3 BLAS of the algorithm.
Definition: zgeqrf.f:151
subroutine zggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
ZGGBAK
Definition: zggbak.f:150
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 zggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
ZGGBAL
Definition: zggbal.f:179
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:76
subroutine zunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
ZUNMQR
Definition: zunmqr.f:169
subroutine zgghd3(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
ZGGHD3
Definition: zgghd3.f:229
subroutine zggev3(JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO)
ZGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices ...
Definition: zggev3.f:218
subroutine zungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
ZUNGQR
Definition: zungqr.f:130
subroutine ztgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
ZTGEVC
Definition: ztgevc.f:221
subroutine zhgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
ZHGEQZ
Definition: zhgeqz.f:286
subroutine zlascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
ZLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: zlascl.f:145