LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zgegs.f
Go to the documentation of this file.
1*> \brief <b> ZGEGS computes the eigenvalues, Schur form, and, optionally, the left and or/right Schur vectors of a complex matrix pair (A,B)</b>
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download ZGEGS + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgegs.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgegs.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgegs.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE ZGEGS( JOBVSL, JOBVSR, N, A, LDA, B, LDB, ALPHA, BETA,
22* VSL, LDVSL, VSR, LDVSR, WORK, LWORK, RWORK,
23* INFO )
24*
25* .. Scalar Arguments ..
26* CHARACTER JOBVSL, JOBVSR
27* INTEGER INFO, LDA, LDB, LDVSL, LDVSR, LWORK, N
28* ..
29* .. Array Arguments ..
30* DOUBLE PRECISION RWORK( * )
31* COMPLEX*16 A( LDA, * ), ALPHA( * ), B( LDB, * ),
32* $ BETA( * ), VSL( LDVSL, * ), VSR( LDVSR, * ),
33* $ WORK( * )
34* ..
35*
36*
37*> \par Purpose:
38* =============
39*>
40*> \verbatim
41*>
42*> This routine is deprecated and has been replaced by routine ZGGES.
43*>
44*> ZGEGS computes the eigenvalues, Schur form, and, optionally, the
45*> left and or/right Schur vectors of a complex matrix pair (A,B).
46*> Given two square matrices A and B, the generalized Schur
47*> factorization has the form
48*>
49*> A = Q*S*Z**H, B = Q*T*Z**H
50*>
51*> where Q and Z are unitary matrices and S and T are upper triangular.
52*> The columns of Q are the left Schur vectors
53*> and the columns of Z are the right Schur vectors.
54*>
55*> If only the eigenvalues of (A,B) are needed, the driver routine
56*> ZGEGV should be used instead. See ZGEGV for a description of the
57*> eigenvalues of the generalized nonsymmetric eigenvalue problem
58*> (GNEP).
59*> \endverbatim
60*
61* Arguments:
62* ==========
63*
64*> \param[in] JOBVSL
65*> \verbatim
66*> JOBVSL is CHARACTER*1
67*> = 'N': do not compute the left Schur vectors;
68*> = 'V': compute the left Schur vectors (returned in VSL).
69*> \endverbatim
70*>
71*> \param[in] JOBVSR
72*> \verbatim
73*> JOBVSR is CHARACTER*1
74*> = 'N': do not compute the right Schur vectors;
75*> = 'V': compute the right Schur vectors (returned in VSR).
76*> \endverbatim
77*>
78*> \param[in] N
79*> \verbatim
80*> N is INTEGER
81*> The order of the matrices A, B, VSL, and VSR. N >= 0.
82*> \endverbatim
83*>
84*> \param[in,out] A
85*> \verbatim
86*> A is COMPLEX*16 array, dimension (LDA, N)
87*> On entry, the matrix A.
88*> On exit, the upper triangular matrix S from the generalized
89*> Schur factorization.
90*> \endverbatim
91*>
92*> \param[in] LDA
93*> \verbatim
94*> LDA is INTEGER
95*> The leading dimension of A. LDA >= max(1,N).
96*> \endverbatim
97*>
98*> \param[in,out] B
99*> \verbatim
100*> B is COMPLEX*16 array, dimension (LDB, N)
101*> On entry, the matrix B.
102*> On exit, the upper triangular matrix T from the generalized
103*> Schur factorization.
104*> \endverbatim
105*>
106*> \param[in] LDB
107*> \verbatim
108*> LDB is INTEGER
109*> The leading dimension of B. LDB >= max(1,N).
110*> \endverbatim
111*>
112*> \param[out] ALPHA
113*> \verbatim
114*> ALPHA is COMPLEX*16 array, dimension (N)
115*> The complex scalars alpha that define the eigenvalues of
116*> GNEP. ALPHA(j) = S(j,j), the diagonal element of the Schur
117*> form of A.
118*> \endverbatim
119*>
120*> \param[out] BETA
121*> \verbatim
122*> BETA is COMPLEX*16 array, dimension (N)
123*> The non-negative real scalars beta that define the
124*> eigenvalues of GNEP. BETA(j) = T(j,j), the diagonal element
125*> of the triangular factor T.
126*>
127*> Together, the quantities alpha = ALPHA(j) and beta = BETA(j)
128*> represent the j-th eigenvalue of the matrix pair (A,B), in
129*> one of the forms lambda = alpha/beta or mu = beta/alpha.
130*> Since either lambda or mu may overflow, they should not,
131*> in general, be computed.
132*> \endverbatim
133*>
134*> \param[out] VSL
135*> \verbatim
136*> VSL is COMPLEX*16 array, dimension (LDVSL,N)
137*> If JOBVSL = 'V', the matrix of left Schur vectors Q.
138*> Not referenced if JOBVSL = 'N'.
139*> \endverbatim
140*>
141*> \param[in] LDVSL
142*> \verbatim
143*> LDVSL is INTEGER
144*> The leading dimension of the matrix VSL. LDVSL >= 1, and
145*> if JOBVSL = 'V', LDVSL >= N.
146*> \endverbatim
147*>
148*> \param[out] VSR
149*> \verbatim
150*> VSR is COMPLEX*16 array, dimension (LDVSR,N)
151*> If JOBVSR = 'V', the matrix of right Schur vectors Z.
152*> Not referenced if JOBVSR = 'N'.
153*> \endverbatim
154*>
155*> \param[in] LDVSR
156*> \verbatim
157*> LDVSR is INTEGER
158*> The leading dimension of the matrix VSR. LDVSR >= 1, and
159*> if JOBVSR = 'V', LDVSR >= N.
160*> \endverbatim
161*>
162*> \param[out] WORK
163*> \verbatim
164*> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
165*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
166*> \endverbatim
167*>
168*> \param[in] LWORK
169*> \verbatim
170*> LWORK is INTEGER
171*> The dimension of the array WORK. LWORK >= max(1,2*N).
172*> For good performance, LWORK must generally be larger.
173*> To compute the optimal value of LWORK, call ILAENV to get
174*> blocksizes (for ZGEQRF, ZUNMQR, and CUNGQR.) Then compute:
175*> NB -- MAX of the blocksizes for ZGEQRF, ZUNMQR, and CUNGQR;
176*> the optimal LWORK is N*(NB+1).
177*>
178*> If LWORK = -1, then a workspace query is assumed; the routine
179*> only calculates the optimal size of the WORK array, returns
180*> this value as the first entry of the WORK array, and no error
181*> message related to LWORK is issued by XERBLA.
182*> \endverbatim
183*>
184*> \param[out] RWORK
185*> \verbatim
186*> RWORK is DOUBLE PRECISION array, dimension (3*N)
187*> \endverbatim
188*>
189*> \param[out] INFO
190*> \verbatim
191*> INFO is INTEGER
192*> = 0: successful exit
193*> < 0: if INFO = -i, the i-th argument had an illegal value.
194*> =1,...,N:
195*> The QZ iteration failed. (A,B) are not in Schur
196*> form, but ALPHA(j) and BETA(j) should be correct for
197*> j=INFO+1,...,N.
198*> > N: errors that usually indicate LAPACK problems:
199*> =N+1: error return from ZGGBAL
200*> =N+2: error return from ZGEQRF
201*> =N+3: error return from ZUNMQR
202*> =N+4: error return from ZUNGQR
203*> =N+5: error return from ZGGHRD
204*> =N+6: error return from ZHGEQZ (other than failed
205*> iteration)
206*> =N+7: error return from ZGGBAK (computing VSL)
207*> =N+8: error return from ZGGBAK (computing VSR)
208*> =N+9: error return from ZLASCL (various places)
209*> \endverbatim
210*
211* Authors:
212* ========
213*
214*> \author Univ. of Tennessee
215*> \author Univ. of California Berkeley
216*> \author Univ. of Colorado Denver
217*> \author NAG Ltd.
218*
219*> \ingroup complex16GEeigen
220*
221* =====================================================================
222 SUBROUTINE zgegs( JOBVSL, JOBVSR, N, A, LDA, B, LDB, ALPHA, BETA,
223 $ VSL, LDVSL, VSR, LDVSR, WORK, LWORK, RWORK,
224 $ INFO )
225*
226* -- LAPACK driver routine --
227* -- LAPACK is a software package provided by Univ. of Tennessee, --
228* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
229*
230* .. Scalar Arguments ..
231 CHARACTER JOBVSL, JOBVSR
232 INTEGER INFO, LDA, LDB, LDVSL, LDVSR, LWORK, N
233* ..
234* .. Array Arguments ..
235 DOUBLE PRECISION RWORK( * )
236 COMPLEX*16 A( LDA, * ), ALPHA( * ), B( LDB, * ),
237 $ beta( * ), vsl( ldvsl, * ), vsr( ldvsr, * ),
238 $ work( * )
239* ..
240*
241* =====================================================================
242*
243* .. Parameters ..
244 DOUBLE PRECISION ZERO, ONE
245 PARAMETER ( ZERO = 0.0d0, one = 1.0d0 )
246 COMPLEX*16 CZERO, CONE
247 parameter( czero = ( 0.0d0, 0.0d0 ),
248 $ cone = ( 1.0d0, 0.0d0 ) )
249* ..
250* .. Local Scalars ..
251 LOGICAL ILASCL, ILBSCL, ILVSL, ILVSR, LQUERY
252 INTEGER ICOLS, IHI, IINFO, IJOBVL, IJOBVR, ILEFT, ILO,
253 $ iright, irows, irwork, itau, iwork, lopt,
254 $ lwkmin, lwkopt, nb, nb1, nb2, nb3
255 DOUBLE PRECISION ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS,
256 $ SAFMIN, SMLNUM
257* ..
258* .. External Subroutines ..
259 EXTERNAL xerbla, zgeqrf, zggbak, zggbal, zgghrd, zhgeqz,
260 $ zlacpy, zlascl, zlaset, zungqr, zunmqr
261* ..
262* .. External Functions ..
263 LOGICAL LSAME
264 INTEGER ILAENV
265 DOUBLE PRECISION DLAMCH, ZLANGE
266 EXTERNAL lsame, ilaenv, dlamch, zlange
267* ..
268* .. Intrinsic Functions ..
269 INTRINSIC int, max
270* ..
271* .. Executable Statements ..
272*
273* Decode the input arguments
274*
275 IF( lsame( jobvsl, 'N' ) ) THEN
276 ijobvl = 1
277 ilvsl = .false.
278 ELSE IF( lsame( jobvsl, 'V' ) ) THEN
279 ijobvl = 2
280 ilvsl = .true.
281 ELSE
282 ijobvl = -1
283 ilvsl = .false.
284 END IF
285*
286 IF( lsame( jobvsr, 'N' ) ) THEN
287 ijobvr = 1
288 ilvsr = .false.
289 ELSE IF( lsame( jobvsr, 'V' ) ) THEN
290 ijobvr = 2
291 ilvsr = .true.
292 ELSE
293 ijobvr = -1
294 ilvsr = .false.
295 END IF
296*
297* Test the input arguments
298*
299 lwkmin = max( 2*n, 1 )
300 lwkopt = lwkmin
301 work( 1 ) = lwkopt
302 lquery = ( lwork.EQ.-1 )
303 info = 0
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( ldvsl.LT.1 .OR. ( ilvsl .AND. ldvsl.LT.n ) ) THEN
315 info = -11
316 ELSE IF( ldvsr.LT.1 .OR. ( ilvsr .AND. ldvsr.LT.n ) ) THEN
317 info = -13
318 ELSE IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
319 info = -15
320 END IF
321*
322 IF( info.EQ.0 ) THEN
323 nb1 = ilaenv( 1, 'ZGEQRF', ' ', n, n, -1, -1 )
324 nb2 = ilaenv( 1, 'ZUNMQR', ' ', n, n, n, -1 )
325 nb3 = ilaenv( 1, 'ZUNGQR', ' ', n, n, n, -1 )
326 nb = max( nb1, nb2, nb3 )
327 lopt = n*( nb+1 )
328 work( 1 ) = lopt
329 END IF
330*
331 IF( info.NE.0 ) THEN
332 CALL xerbla( 'ZGEGS ', -info )
333 RETURN
334 ELSE IF( lquery ) THEN
335 RETURN
336 END IF
337*
338* Quick return if possible
339*
340 IF( n.EQ.0 )
341 $ RETURN
342*
343* Get machine constants
344*
345 eps = dlamch( 'E' )*dlamch( 'B' )
346 safmin = dlamch( 'S' )
347 smlnum = n*safmin / eps
348 bignum = one / smlnum
349*
350* Scale A if max element outside range [SMLNUM,BIGNUM]
351*
352 anrm = zlange( 'M', n, n, a, lda, rwork )
353 ilascl = .false.
354 IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
355 anrmto = smlnum
356 ilascl = .true.
357 ELSE IF( anrm.GT.bignum ) THEN
358 anrmto = bignum
359 ilascl = .true.
360 END IF
361*
362 IF( ilascl ) THEN
363 CALL zlascl( 'G', -1, -1, anrm, anrmto, n, n, a, lda, iinfo )
364 IF( iinfo.NE.0 ) THEN
365 info = n + 9
366 RETURN
367 END IF
368 END IF
369*
370* Scale B if max element outside range [SMLNUM,BIGNUM]
371*
372 bnrm = zlange( 'M', n, n, b, ldb, rwork )
373 ilbscl = .false.
374 IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
375 bnrmto = smlnum
376 ilbscl = .true.
377 ELSE IF( bnrm.GT.bignum ) THEN
378 bnrmto = bignum
379 ilbscl = .true.
380 END IF
381*
382 IF( ilbscl ) THEN
383 CALL zlascl( 'G', -1, -1, bnrm, bnrmto, n, n, b, ldb, iinfo )
384 IF( iinfo.NE.0 ) THEN
385 info = n + 9
386 RETURN
387 END IF
388 END IF
389*
390* Permute the matrix to make it more nearly triangular
391*
392 ileft = 1
393 iright = n + 1
394 irwork = iright + n
395 iwork = 1
396 CALL zggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
397 $ rwork( iright ), rwork( irwork ), iinfo )
398 IF( iinfo.NE.0 ) THEN
399 info = n + 1
400 GO TO 10
401 END IF
402*
403* Reduce B to triangular form, and initialize VSL and/or VSR
404*
405 irows = ihi + 1 - ilo
406 icols = n + 1 - ilo
407 itau = iwork
408 iwork = itau + irows
409 CALL zgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
410 $ work( iwork ), lwork+1-iwork, iinfo )
411 IF( iinfo.GE.0 )
412 $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
413 IF( iinfo.NE.0 ) THEN
414 info = n + 2
415 GO TO 10
416 END IF
417*
418 CALL zunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
419 $ work( itau ), a( ilo, ilo ), lda, work( iwork ),
420 $ lwork+1-iwork, iinfo )
421 IF( iinfo.GE.0 )
422 $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
423 IF( iinfo.NE.0 ) THEN
424 info = n + 3
425 GO TO 10
426 END IF
427*
428 IF( ilvsl ) THEN
429 CALL zlaset( 'Full', n, n, czero, cone, vsl, ldvsl )
430 CALL zlacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
431 $ vsl( ilo+1, ilo ), ldvsl )
432 CALL zungqr( irows, irows, irows, vsl( ilo, ilo ), ldvsl,
433 $ work( itau ), work( iwork ), lwork+1-iwork,
434 $ iinfo )
435 IF( iinfo.GE.0 )
436 $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
437 IF( iinfo.NE.0 ) THEN
438 info = n + 4
439 GO TO 10
440 END IF
441 END IF
442*
443 IF( ilvsr )
444 $ CALL zlaset( 'Full', n, n, czero, cone, vsr, ldvsr )
445*
446* Reduce to generalized Hessenberg form
447*
448 CALL zgghrd( jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb, vsl,
449 $ ldvsl, vsr, ldvsr, iinfo )
450 IF( iinfo.NE.0 ) THEN
451 info = n + 5
452 GO TO 10
453 END IF
454*
455* Perform QZ algorithm, computing Schur vectors if desired
456*
457 iwork = itau
458 CALL zhgeqz( 'S', jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb,
459 $ alpha, beta, vsl, ldvsl, vsr, ldvsr, work( iwork ),
460 $ lwork+1-iwork, rwork( irwork ), iinfo )
461 IF( iinfo.GE.0 )
462 $ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
463 IF( iinfo.NE.0 ) THEN
464 IF( iinfo.GT.0 .AND. iinfo.LE.n ) THEN
465 info = iinfo
466 ELSE IF( iinfo.GT.n .AND. iinfo.LE.2*n ) THEN
467 info = iinfo - n
468 ELSE
469 info = n + 6
470 END IF
471 GO TO 10
472 END IF
473*
474* Apply permutation to VSL and VSR
475*
476 IF( ilvsl ) THEN
477 CALL zggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
478 $ rwork( iright ), n, vsl, ldvsl, iinfo )
479 IF( iinfo.NE.0 ) THEN
480 info = n + 7
481 GO TO 10
482 END IF
483 END IF
484 IF( ilvsr ) THEN
485 CALL zggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
486 $ rwork( iright ), n, vsr, ldvsr, iinfo )
487 IF( iinfo.NE.0 ) THEN
488 info = n + 8
489 GO TO 10
490 END IF
491 END IF
492*
493* Undo scaling
494*
495 IF( ilascl ) THEN
496 CALL zlascl( 'U', -1, -1, anrmto, anrm, n, n, a, lda, iinfo )
497 IF( iinfo.NE.0 ) THEN
498 info = n + 9
499 RETURN
500 END IF
501 CALL zlascl( 'G', -1, -1, anrmto, anrm, n, 1, alpha, n, iinfo )
502 IF( iinfo.NE.0 ) THEN
503 info = n + 9
504 RETURN
505 END IF
506 END IF
507*
508 IF( ilbscl ) THEN
509 CALL zlascl( 'U', -1, -1, bnrmto, bnrm, n, n, b, ldb, iinfo )
510 IF( iinfo.NE.0 ) THEN
511 info = n + 9
512 RETURN
513 END IF
514 CALL zlascl( 'G', -1, -1, bnrmto, bnrm, n, 1, beta, n, iinfo )
515 IF( iinfo.NE.0 ) THEN
516 info = n + 9
517 RETURN
518 END IF
519 END IF
520*
521 10 CONTINUE
522 work( 1 ) = lwkopt
523*
524 RETURN
525*
526* End of ZGEGS
527*
528 END
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
zgeqrf
subroutine zgeqrf(m, n, a, lda, tau, work, lwork, info)
ZGEQRF
Definition zgeqrf.f:146
zggbak
subroutine zggbak(job, side, n, ilo, ihi, lscale, rscale, m, v, ldv, info)
ZGGBAK
Definition zggbak.f:148
zggbal
subroutine zggbal(job, n, a, lda, b, ldb, ilo, ihi, lscale, rscale, work, info)
ZGGBAL
Definition zggbal.f:177
zgghrd
subroutine zgghrd(compq, compz, n, ilo, ihi, a, lda, b, ldb, q, ldq, z, ldz, info)
ZGGHRD
Definition zgghrd.f:204
zhgeqz
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:284
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:103
zlascl
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
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:106
zungqr
subroutine zungqr(m, n, k, a, lda, tau, work, lwork, info)
ZUNGQR
Definition zungqr.f:128
zunmqr
subroutine zunmqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
ZUNMQR
Definition zunmqr.f:167
zgegs
subroutine zgegs(jobvsl, jobvsr, n, a, lda, b, ldb, alpha, beta, vsl, ldvsl, vsr, ldvsr, work, lwork, rwork, info)
ZGEGS computes the eigenvalues, Schur form, and, optionally, the left and or/right Schur vectors of a...
Definition zgegs.f:225