LAPACK  3.10.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 ZHGEQZ,
200 *> =N+2: error return from ZTGEVC.
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 *> \ingroup complex16GEeigen
212 *
213 * =====================================================================
214  SUBROUTINE zggev3( JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA,
215  $ VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO )
216 *
217 * -- LAPACK driver routine --
218 * -- LAPACK is a software package provided by Univ. of Tennessee, --
219 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
220 *
221 * .. Scalar Arguments ..
222  CHARACTER JOBVL, JOBVR
223  INTEGER INFO, LDA, LDB, LDVL, LDVR, LWORK, N
224 * ..
225 * .. Array Arguments ..
226  DOUBLE PRECISION RWORK( * )
227  COMPLEX*16 A( LDA, * ), ALPHA( * ), B( LDB, * ),
228  $ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
229  $ work( * )
230 * ..
231 *
232 * =====================================================================
233 *
234 * .. Parameters ..
235  DOUBLE PRECISION ZERO, ONE
236  parameter( zero = 0.0d0, one = 1.0d0 )
237  COMPLEX*16 CZERO, CONE
238  parameter( czero = ( 0.0d0, 0.0d0 ),
239  $ cone = ( 1.0d0, 0.0d0 ) )
240 * ..
241 * .. Local Scalars ..
242  LOGICAL ILASCL, ILBSCL, ILV, ILVL, ILVR, LQUERY
243  CHARACTER CHTEMP
244  INTEGER ICOLS, IERR, IHI, IJOBVL, IJOBVR, ILEFT, ILO,
245  $ in, iright, irows, irwrk, itau, iwrk, jc, jr,
246  $ lwkopt
247  DOUBLE PRECISION ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS,
248  $ smlnum, temp
249  COMPLEX*16 X
250 * ..
251 * .. Local Arrays ..
252  LOGICAL LDUMMA( 1 )
253 * ..
254 * .. External Subroutines ..
255  EXTERNAL dlabad, xerbla, zgeqrf, zggbak, zggbal, zgghd3,
257  $ zunmqr
258 * ..
259 * .. External Functions ..
260  LOGICAL LSAME
261  DOUBLE PRECISION DLAMCH, ZLANGE
262  EXTERNAL lsame, dlamch, zlange
263 * ..
264 * .. Intrinsic Functions ..
265  INTRINSIC abs, dble, dimag, max, sqrt
266 * ..
267 * .. Statement Functions ..
268  DOUBLE PRECISION ABS1
269 * ..
270 * .. Statement Function definitions ..
271  abs1( x ) = abs( dble( x ) ) + abs( dimag( x ) )
272 * ..
273 * .. Executable Statements ..
274 *
275 * Decode the input arguments
276 *
277  IF( lsame( jobvl, 'N' ) ) THEN
278  ijobvl = 1
279  ilvl = .false.
280  ELSE IF( lsame( jobvl, 'V' ) ) THEN
281  ijobvl = 2
282  ilvl = .true.
283  ELSE
284  ijobvl = -1
285  ilvl = .false.
286  END IF
287 *
288  IF( lsame( jobvr, 'N' ) ) THEN
289  ijobvr = 1
290  ilvr = .false.
291  ELSE IF( lsame( jobvr, 'V' ) ) THEN
292  ijobvr = 2
293  ilvr = .true.
294  ELSE
295  ijobvr = -1
296  ilvr = .false.
297  END IF
298  ilv = ilvl .OR. ilvr
299 *
300 * Test the input arguments
301 *
302  info = 0
303  lquery = ( lwork.EQ.-1 )
304  IF( ijobvl.LE.0 ) THEN
305  info = -1
306  ELSE IF( ijobvr.LE.0 ) THEN
307  info = -2
308  ELSE IF( n.LT.0 ) THEN
309  info = -3
310  ELSE IF( lda.LT.max( 1, n ) ) THEN
311  info = -5
312  ELSE IF( ldb.LT.max( 1, n ) ) THEN
313  info = -7
314  ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
315  info = -11
316  ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
317  info = -13
318  ELSE IF( lwork.LT.max( 1, 2*n ) .AND. .NOT.lquery ) THEN
319  info = -15
320  END IF
321 *
322 * Compute workspace
323 *
324  IF( info.EQ.0 ) THEN
325  CALL zgeqrf( n, n, b, ldb, work, work, -1, ierr )
326  lwkopt = max( 1, n+int( work( 1 ) ) )
327  CALL zunmqr( 'L', 'C', n, n, n, b, ldb, work, a, lda, work,
328  $ -1, ierr )
329  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
330  IF( ilvl ) THEN
331  CALL zungqr( n, n, n, vl, ldvl, work, work, -1, ierr )
332  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
333  END IF
334  IF( ilv ) THEN
335  CALL zgghd3( jobvl, jobvr, n, 1, n, a, lda, b, ldb, vl,
336  $ ldvl, vr, ldvr, work, -1, ierr )
337  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
338  CALL zlaqz0( 'S', jobvl, jobvr, n, 1, n, a, lda, b, ldb,
339  $ alpha, beta, vl, ldvl, vr, ldvr, work, -1,
340  $ rwork, 0, ierr )
341  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
342  ELSE
343  CALL zgghd3( jobvl, jobvr, n, 1, n, a, lda, b, ldb, vl,
344  $ ldvl, vr, ldvr, work, -1, ierr )
345  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
346  CALL zlaqz0( 'E', jobvl, jobvr, n, 1, n, a, lda, b, ldb,
347  $ alpha, beta, vl, ldvl, vr, ldvr, work, -1,
348  $ rwork, 0, ierr )
349  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
350  END IF
351  work( 1 ) = dcmplx( lwkopt )
352  END IF
353 *
354  IF( info.NE.0 ) THEN
355  CALL xerbla( 'ZGGEV3 ', -info )
356  RETURN
357  ELSE IF( lquery ) THEN
358  RETURN
359  END IF
360 *
361 * Quick return if possible
362 *
363  IF( n.EQ.0 )
364  $ RETURN
365 *
366 * Get machine constants
367 *
368  eps = dlamch( 'E' )*dlamch( 'B' )
369  smlnum = dlamch( 'S' )
370  bignum = one / smlnum
371  CALL dlabad( smlnum, bignum )
372  smlnum = sqrt( smlnum ) / eps
373  bignum = one / smlnum
374 *
375 * Scale A if max element outside range [SMLNUM,BIGNUM]
376 *
377  anrm = zlange( 'M', n, n, a, lda, rwork )
378  ilascl = .false.
379  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
380  anrmto = smlnum
381  ilascl = .true.
382  ELSE IF( anrm.GT.bignum ) THEN
383  anrmto = bignum
384  ilascl = .true.
385  END IF
386  IF( ilascl )
387  $ CALL zlascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
388 *
389 * Scale B if max element outside range [SMLNUM,BIGNUM]
390 *
391  bnrm = zlange( 'M', n, n, b, ldb, rwork )
392  ilbscl = .false.
393  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
394  bnrmto = smlnum
395  ilbscl = .true.
396  ELSE IF( bnrm.GT.bignum ) THEN
397  bnrmto = bignum
398  ilbscl = .true.
399  END IF
400  IF( ilbscl )
401  $ CALL zlascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
402 *
403 * Permute the matrices A, B to isolate eigenvalues if possible
404 *
405  ileft = 1
406  iright = n + 1
407  irwrk = iright + n
408  CALL zggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
409  $ rwork( iright ), rwork( irwrk ), ierr )
410 *
411 * Reduce B to triangular form (QR decomposition of B)
412 *
413  irows = ihi + 1 - ilo
414  IF( ilv ) THEN
415  icols = n + 1 - ilo
416  ELSE
417  icols = irows
418  END IF
419  itau = 1
420  iwrk = itau + irows
421  CALL zgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
422  $ work( iwrk ), lwork+1-iwrk, ierr )
423 *
424 * Apply the orthogonal transformation to matrix A
425 *
426  CALL zunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
427  $ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
428  $ lwork+1-iwrk, ierr )
429 *
430 * Initialize VL
431 *
432  IF( ilvl ) THEN
433  CALL zlaset( 'Full', n, n, czero, cone, vl, ldvl )
434  IF( irows.GT.1 ) THEN
435  CALL zlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
436  $ vl( ilo+1, ilo ), ldvl )
437  END IF
438  CALL zungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
439  $ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
440  END IF
441 *
442 * Initialize VR
443 *
444  IF( ilvr )
445  $ CALL zlaset( 'Full', n, n, czero, cone, vr, ldvr )
446 *
447 * Reduce to generalized Hessenberg form
448 *
449  IF( ilv ) THEN
450 *
451 * Eigenvectors requested -- work on whole matrix.
452 *
453  CALL zgghd3( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
454  $ ldvl, vr, ldvr, work( iwrk ), lwork+1-iwrk, ierr )
455  ELSE
456  CALL zgghd3( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
457  $ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr,
458  $ work( iwrk ), lwork+1-iwrk, ierr )
459  END IF
460 *
461 * Perform QZ algorithm (Compute eigenvalues, and optionally, the
462 * Schur form and Schur vectors)
463 *
464  iwrk = itau
465  IF( ilv ) THEN
466  chtemp = 'S'
467  ELSE
468  chtemp = 'E'
469  END IF
470  CALL zlaqz0( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
471  $ alpha, beta, vl, ldvl, vr, ldvr, work( iwrk ),
472  $ lwork+1-iwrk, rwork( irwrk ), 0, ierr )
473  IF( ierr.NE.0 ) THEN
474  IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
475  info = ierr
476  ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
477  info = ierr - n
478  ELSE
479  info = n + 1
480  END IF
481  GO TO 70
482  END IF
483 *
484 * Compute Eigenvectors
485 *
486  IF( ilv ) THEN
487  IF( ilvl ) THEN
488  IF( ilvr ) THEN
489  chtemp = 'B'
490  ELSE
491  chtemp = 'L'
492  END IF
493  ELSE
494  chtemp = 'R'
495  END IF
496 *
497  CALL ztgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
498  $ vr, ldvr, n, in, work( iwrk ), rwork( irwrk ),
499  $ ierr )
500  IF( ierr.NE.0 ) THEN
501  info = n + 2
502  GO TO 70
503  END IF
504 *
505 * Undo balancing on VL and VR and normalization
506 *
507  IF( ilvl ) THEN
508  CALL zggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
509  $ rwork( iright ), n, vl, ldvl, ierr )
510  DO 30 jc = 1, n
511  temp = zero
512  DO 10 jr = 1, n
513  temp = max( temp, abs1( vl( jr, jc ) ) )
514  10 CONTINUE
515  IF( temp.LT.smlnum )
516  $ GO TO 30
517  temp = one / temp
518  DO 20 jr = 1, n
519  vl( jr, jc ) = vl( jr, jc )*temp
520  20 CONTINUE
521  30 CONTINUE
522  END IF
523  IF( ilvr ) THEN
524  CALL zggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
525  $ rwork( iright ), n, vr, ldvr, ierr )
526  DO 60 jc = 1, n
527  temp = zero
528  DO 40 jr = 1, n
529  temp = max( temp, abs1( vr( jr, jc ) ) )
530  40 CONTINUE
531  IF( temp.LT.smlnum )
532  $ GO TO 60
533  temp = one / temp
534  DO 50 jr = 1, n
535  vr( jr, jc ) = vr( jr, jc )*temp
536  50 CONTINUE
537  60 CONTINUE
538  END IF
539  END IF
540 *
541 * Undo scaling if necessary
542 *
543  70 CONTINUE
544 *
545  IF( ilascl )
546  $ CALL zlascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
547 *
548  IF( ilbscl )
549  $ CALL zlascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
550 *
551  work( 1 ) = dcmplx( lwkopt )
552  RETURN
553 *
554 * End of ZGGEV3
555 *
556  END
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine zggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
ZGGBAL
Definition: zggbal.f:177
subroutine zggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
ZGGBAK
Definition: zggbak.f:148
recursive subroutine zlaqz0(WANTS, WANTQ, WANTZ, N, ILO, IHI, A, LDA, B, LDB, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, REC, INFO)
ZLAQZ0
Definition: zlaqz0.f:284
subroutine ztgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
ZTGEVC
Definition: ztgevc.f:219
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:216
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:143
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:103
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:106
subroutine zungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
ZUNGQR
Definition: zungqr.f:128
subroutine zgghd3(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
ZGGHD3
Definition: zgghd3.f:227
subroutine zunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
ZUNMQR
Definition: zunmqr.f:167
subroutine zgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
ZGEQRF VARIANT: left-looking Level 3 BLAS of the algorithm.
Definition: zgeqrf.f:151