LAPACK  3.6.0 LAPACK: Linear Algebra PACKage
complex
Collaboration diagram for complex:


This browser is not able to show SVG: try Firefox, Chrome, Safari, or Opera instead.

## Functions

subroutine cgees (JOBVS, SORT, SELECT, N, A, LDA, SDIM, W, VS, LDVS, WORK, LWORK, RWORK, BWORK, INFO)
CGEES computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices More...

subroutine cgeesx (JOBVS, SORT, SELECT, SENSE, N, A, LDA, SDIM, W, VS, LDVS, RCONDE, RCONDV, WORK, LWORK, RWORK, BWORK, INFO)
CGEESX computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices More...

subroutine cgeev (JOBVL, JOBVR, N, A, LDA, W, VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO)
CGEEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices More...

subroutine cgeevx (BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, W, VL, LDVL, VR, LDVR, ILO, IHI, SCALE, ABNRM, RCONDE, RCONDV, WORK, LWORK, RWORK, INFO)
CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices More...

subroutine cgges (JOBVSL, JOBVSR, SORT, SELCTG, N, A, LDA, B, LDB, SDIM, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK, RWORK, BWORK, INFO)
CGGES computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices More...

subroutine cgges3 (JOBVSL, JOBVSR, SORT, SELCTG, N, A, LDA, B, LDB, SDIM, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK, RWORK, BWORK, INFO)
CGGES3 computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices (blocked algorithm) More...

subroutine cggesx (JOBVSL, JOBVSR, SORT, SELCTG, SENSE, N, A, LDA, B, LDB, SDIM, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, RCONDE, RCONDV, WORK, LWORK, RWORK, IWORK, LIWORK, BWORK, INFO)
CGGESX computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices More...

subroutine cggev (JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO)
CGGEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices More...

subroutine cggev3 (JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO)
CGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices (blocked algorithm) More...

subroutine cggevx (BALANC, JOBVL, JOBVR, SENSE, N, A, LDA, B, LDB, ALPHA, BETA, VL, LDVL, VR, LDVR, ILO, IHI, LSCALE, RSCALE, ABNRM, BBNRM, RCONDE, RCONDV, WORK, LWORK, RWORK, IWORK, BWORK, INFO)
CGGEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices More...

subroutine cgegs (JOBVSL, JOBVSR, N, A, LDA, B, LDB, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, WORK, LWORK, RWORK, INFO)
CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices More...

subroutine cgegv (JOBVL, JOBVR, N, A, LDA, B, LDB, ALPHA, BETA, VL, LDVL, VR, LDVR, WORK, LWORK, RWORK, INFO)
CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices More...

## Detailed Description

This is the group of complex eigenvalue driver functions for GE matrices

## Function Documentation

 subroutine cgees ( character JOBVS, character SORT, external SELECT, integer N, complex, dimension( lda, * ) A, integer LDA, integer SDIM, complex, dimension( * ) W, complex, dimension( ldvs, * ) VS, integer LDVS, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, logical, dimension( * ) BWORK, integer INFO )

CGEES computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices

Purpose:
``` CGEES computes for an N-by-N complex nonsymmetric matrix A, the
eigenvalues, the Schur form T, and, optionally, the matrix of Schur
vectors Z.  This gives the Schur factorization A = Z*T*(Z**H).

Optionally, it also orders the eigenvalues on the diagonal of the
Schur form so that selected eigenvalues are at the top left.
The leading columns of Z then form an orthonormal basis for the
invariant subspace corresponding to the selected eigenvalues.

A complex matrix is in Schur form if it is upper triangular.```
Parameters
 [in] JOBVS ``` JOBVS is CHARACTER*1 = 'N': Schur vectors are not computed; = 'V': Schur vectors are computed.``` [in] SORT ``` SORT is CHARACTER*1 Specifies whether or not to order the eigenvalues on the diagonal of the Schur form. = 'N': Eigenvalues are not ordered: = 'S': Eigenvalues are ordered (see SELECT).``` [in] SELECT ``` SELECT is a LOGICAL FUNCTION of one COMPLEX argument SELECT must be declared EXTERNAL in the calling subroutine. If SORT = 'S', SELECT is used to select eigenvalues to order to the top left of the Schur form. IF SORT = 'N', SELECT is not referenced. The eigenvalue W(j) is selected if SELECT(W(j)) is true.``` [in] N ``` N is INTEGER The order of the matrix A. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA,N) On entry, the N-by-N matrix A. On exit, A has been overwritten by its Schur form T.``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= max(1,N).``` [out] SDIM ``` SDIM is INTEGER If SORT = 'N', SDIM = 0. If SORT = 'S', SDIM = number of eigenvalues for which SELECT is true.``` [out] W ``` W is COMPLEX array, dimension (N) W contains the computed eigenvalues, in the same order that they appear on the diagonal of the output Schur form T.``` [out] VS ``` VS is COMPLEX array, dimension (LDVS,N) If JOBVS = 'V', VS contains the unitary matrix Z of Schur vectors. If JOBVS = 'N', VS is not referenced.``` [in] LDVS ``` LDVS is INTEGER The leading dimension of the array VS. LDVS >= 1; if JOBVS = 'V', LDVS >= N.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,2*N). For good performance, LWORK must generally be larger. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (N)` [out] BWORK ``` BWORK is LOGICAL array, dimension (N) Not referenced if SORT = 'N'.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. > 0: if INFO = i, and i is <= N: the QR algorithm failed to compute all the eigenvalues; elements 1:ILO-1 and i+1:N of W contain those eigenvalues which have converged; if JOBVS = 'V', VS contains the matrix which reduces A to its partially converged Schur form. = N+1: the eigenvalues could not be reordered because some eigenvalues were too close to separate (the problem is very ill-conditioned); = N+2: after reordering, roundoff changed values of some complex eigenvalues so that leading eigenvalues in the Schur form no longer satisfy SELECT = .TRUE.. This could also be caused by underflow due to scaling.```
Date
November 2011

Definition at line 199 of file cgees.f.

199 *
200 * -- LAPACK driver routine (version 3.4.0) --
201 * -- LAPACK is a software package provided by Univ. of Tennessee, --
202 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
203 * November 2011
204 *
205 * .. Scalar Arguments ..
206  CHARACTER jobvs, sort
207  INTEGER info, lda, ldvs, lwork, n, sdim
208 * ..
209 * .. Array Arguments ..
210  LOGICAL bwork( * )
211  REAL rwork( * )
212  COMPLEX a( lda, * ), vs( ldvs, * ), w( * ), work( * )
213 * ..
214 * .. Function Arguments ..
215  LOGICAL select
216  EXTERNAL SELECT
217 * ..
218 *
219 * =====================================================================
220 *
221 * .. Parameters ..
222  REAL zero, one
223  parameter( zero = 0.0e0, one = 1.0e0 )
224 * ..
225 * .. Local Scalars ..
226  LOGICAL lquery, scalea, wantst, wantvs
227  INTEGER hswork, i, ibal, icond, ierr, ieval, ihi, ilo,
228  \$ itau, iwrk, maxwrk, minwrk
229  REAL anrm, bignum, cscale, eps, s, sep, smlnum
230 * ..
231 * .. Local Arrays ..
232  REAL dum( 1 )
233 * ..
234 * .. External Subroutines ..
235  EXTERNAL ccopy, cgebak, cgebal, cgehrd, chseqr, clacpy,
237 * ..
238 * .. External Functions ..
239  LOGICAL lsame
240  INTEGER ilaenv
241  REAL clange, slamch
242  EXTERNAL lsame, ilaenv, clange, slamch
243 * ..
244 * .. Intrinsic Functions ..
245  INTRINSIC max, sqrt
246 * ..
247 * .. Executable Statements ..
248 *
249 * Test the input arguments
250 *
251  info = 0
252  lquery = ( lwork.EQ.-1 )
253  wantvs = lsame( jobvs, 'V' )
254  wantst = lsame( sort, 'S' )
255  IF( ( .NOT.wantvs ) .AND. ( .NOT.lsame( jobvs, 'N' ) ) ) THEN
256  info = -1
257  ELSE IF( ( .NOT.wantst ) .AND. ( .NOT.lsame( sort, 'N' ) ) ) THEN
258  info = -2
259  ELSE IF( n.LT.0 ) THEN
260  info = -4
261  ELSE IF( lda.LT.max( 1, n ) ) THEN
262  info = -6
263  ELSE IF( ldvs.LT.1 .OR. ( wantvs .AND. ldvs.LT.n ) ) THEN
264  info = -10
265  END IF
266 *
267 * Compute workspace
268 * (Note: Comments in the code beginning "Workspace:" describe the
269 * minimal amount of workspace needed at that point in the code,
270 * as well as the preferred amount for good performance.
271 * CWorkspace refers to complex workspace, and RWorkspace to real
272 * workspace. NB refers to the optimal block size for the
273 * immediately following subroutine, as returned by ILAENV.
274 * HSWORK refers to the workspace preferred by CHSEQR, as
275 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
276 * the worst case.)
277 *
278  IF( info.EQ.0 ) THEN
279  IF( n.EQ.0 ) THEN
280  minwrk = 1
281  maxwrk = 1
282  ELSE
283  maxwrk = n + n*ilaenv( 1, 'CGEHRD', ' ', n, 1, n, 0 )
284  minwrk = 2*n
285 *
286  CALL chseqr( 'S', jobvs, n, 1, n, a, lda, w, vs, ldvs,
287  \$ work, -1, ieval )
288  hswork = work( 1 )
289 *
290  IF( .NOT.wantvs ) THEN
291  maxwrk = max( maxwrk, hswork )
292  ELSE
293  maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
294  \$ ' ', n, 1, n, -1 ) )
295  maxwrk = max( maxwrk, hswork )
296  END IF
297  END IF
298  work( 1 ) = maxwrk
299 *
300  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
301  info = -12
302  END IF
303  END IF
304 *
305  IF( info.NE.0 ) THEN
306  CALL xerbla( 'CGEES ', -info )
307  RETURN
308  ELSE IF( lquery ) THEN
309  RETURN
310  END IF
311 *
312 * Quick return if possible
313 *
314  IF( n.EQ.0 ) THEN
315  sdim = 0
316  RETURN
317  END IF
318 *
319 * Get machine constants
320 *
321  eps = slamch( 'P' )
322  smlnum = slamch( 'S' )
323  bignum = one / smlnum
324  CALL slabad( smlnum, bignum )
325  smlnum = sqrt( smlnum ) / eps
326  bignum = one / smlnum
327 *
328 * Scale A if max element outside range [SMLNUM,BIGNUM]
329 *
330  anrm = clange( 'M', n, n, a, lda, dum )
331  scalea = .false.
332  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
333  scalea = .true.
334  cscale = smlnum
335  ELSE IF( anrm.GT.bignum ) THEN
336  scalea = .true.
337  cscale = bignum
338  END IF
339  IF( scalea )
340  \$ CALL clascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
341 *
342 * Permute the matrix to make it more nearly triangular
343 * (CWorkspace: none)
344 * (RWorkspace: need N)
345 *
346  ibal = 1
347  CALL cgebal( 'P', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
348 *
349 * Reduce to upper Hessenberg form
350 * (CWorkspace: need 2*N, prefer N+N*NB)
351 * (RWorkspace: none)
352 *
353  itau = 1
354  iwrk = n + itau
355  CALL cgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
356  \$ lwork-iwrk+1, ierr )
357 *
358  IF( wantvs ) THEN
359 *
360 * Copy Householder vectors to VS
361 *
362  CALL clacpy( 'L', n, n, a, lda, vs, ldvs )
363 *
364 * Generate unitary matrix in VS
365 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
366 * (RWorkspace: none)
367 *
368  CALL cunghr( n, ilo, ihi, vs, ldvs, work( itau ), work( iwrk ),
369  \$ lwork-iwrk+1, ierr )
370  END IF
371 *
372  sdim = 0
373 *
374 * Perform QR iteration, accumulating Schur vectors in VS if desired
375 * (CWorkspace: need 1, prefer HSWORK (see comments) )
376 * (RWorkspace: none)
377 *
378  iwrk = itau
379  CALL chseqr( 'S', jobvs, n, ilo, ihi, a, lda, w, vs, ldvs,
380  \$ work( iwrk ), lwork-iwrk+1, ieval )
381  IF( ieval.GT.0 )
382  \$ info = ieval
383 *
384 * Sort eigenvalues if desired
385 *
386  IF( wantst .AND. info.EQ.0 ) THEN
387  IF( scalea )
388  \$ CALL clascl( 'G', 0, 0, cscale, anrm, n, 1, w, n, ierr )
389  DO 10 i = 1, n
390  bwork( i ) = SELECT( w( i ) )
391  10 CONTINUE
392 *
393 * Reorder eigenvalues and transform Schur vectors
394 * (CWorkspace: none)
395 * (RWorkspace: none)
396 *
397  CALL ctrsen( 'N', jobvs, bwork, n, a, lda, vs, ldvs, w, sdim,
398  \$ s, sep, work( iwrk ), lwork-iwrk+1, icond )
399  END IF
400 *
401  IF( wantvs ) THEN
402 *
403 * Undo balancing
404 * (CWorkspace: none)
405 * (RWorkspace: need N)
406 *
407  CALL cgebak( 'P', 'R', n, ilo, ihi, rwork( ibal ), n, vs, ldvs,
408  \$ ierr )
409  END IF
410 *
411  IF( scalea ) THEN
412 *
413 * Undo scaling for the Schur form of A
414 *
415  CALL clascl( 'U', 0, 0, cscale, anrm, n, n, a, lda, ierr )
416  CALL ccopy( n, a, lda+1, w, 1 )
417  END IF
418 *
419  work( 1 ) = maxwrk
420  RETURN
421 *
422 * End of CGEES
423 *
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine ctrsen(JOB, COMPQ, SELECT, N, T, LDT, Q, LDQ, W, M, S, SEP, WORK, LWORK, INFO)
CTRSEN
Definition: ctrsen.f:266
subroutine cgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CGEHRD
Definition: cgehrd.f:169
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cgebal(JOB, N, A, LDA, ILO, IHI, SCALE, INFO)
CGEBAL
Definition: cgebal.f:163
subroutine cgebak(JOB, SIDE, N, ILO, IHI, SCALE, M, V, LDV, INFO)
CGEBAK
Definition: cgebak.f:133
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine chseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, W, Z, LDZ, WORK, LWORK, INFO)
CHSEQR
Definition: chseqr.f:301
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:52
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cunghr(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CUNGHR
Definition: cunghr.f:128

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cgeesx ( character JOBVS, character SORT, external SELECT, character SENSE, integer N, complex, dimension( lda, * ) A, integer LDA, integer SDIM, complex, dimension( * ) W, complex, dimension( ldvs, * ) VS, integer LDVS, real RCONDE, real RCONDV, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, logical, dimension( * ) BWORK, integer INFO )

CGEESX computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices

Purpose:
``` CGEESX computes for an N-by-N complex nonsymmetric matrix A, the
eigenvalues, the Schur form T, and, optionally, the matrix of Schur
vectors Z.  This gives the Schur factorization A = Z*T*(Z**H).

Optionally, it also orders the eigenvalues on the diagonal of the
Schur form so that selected eigenvalues are at the top left;
computes a reciprocal condition number for the average of the
selected eigenvalues (RCONDE); and computes a reciprocal condition
number for the right invariant subspace corresponding to the
selected eigenvalues (RCONDV).  The leading columns of Z form an
orthonormal basis for this invariant subspace.

For further explanation of the reciprocal condition numbers RCONDE
and RCONDV, see Section 4.10 of the LAPACK Users' Guide (where
these quantities are called s and sep respectively).

A complex matrix is in Schur form if it is upper triangular.```
Parameters
 [in] JOBVS ``` JOBVS is CHARACTER*1 = 'N': Schur vectors are not computed; = 'V': Schur vectors are computed.``` [in] SORT ``` SORT is CHARACTER*1 Specifies whether or not to order the eigenvalues on the diagonal of the Schur form. = 'N': Eigenvalues are not ordered; = 'S': Eigenvalues are ordered (see SELECT).``` [in] SELECT ``` SELECT is procedure) LOGICAL FUNCTION of one COMPLEX argument SELECT must be declared EXTERNAL in the calling subroutine. If SORT = 'S', SELECT is used to select eigenvalues to order to the top left of the Schur form. If SORT = 'N', SELECT is not referenced. An eigenvalue W(j) is selected if SELECT(W(j)) is true.``` [in] SENSE ``` SENSE is CHARACTER*1 Determines which reciprocal condition numbers are computed. = 'N': None are computed; = 'E': Computed for average of selected eigenvalues only; = 'V': Computed for selected right invariant subspace only; = 'B': Computed for both. If SENSE = 'E', 'V' or 'B', SORT must equal 'S'.``` [in] N ``` N is INTEGER The order of the matrix A. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the N-by-N matrix A. On exit, A is overwritten by its Schur form T.``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= max(1,N).``` [out] SDIM ``` SDIM is INTEGER If SORT = 'N', SDIM = 0. If SORT = 'S', SDIM = number of eigenvalues for which SELECT is true.``` [out] W ``` W is COMPLEX array, dimension (N) W contains the computed eigenvalues, in the same order that they appear on the diagonal of the output Schur form T.``` [out] VS ``` VS is COMPLEX array, dimension (LDVS,N) If JOBVS = 'V', VS contains the unitary matrix Z of Schur vectors. If JOBVS = 'N', VS is not referenced.``` [in] LDVS ``` LDVS is INTEGER The leading dimension of the array VS. LDVS >= 1, and if JOBVS = 'V', LDVS >= N.``` [out] RCONDE ``` RCONDE is REAL If SENSE = 'E' or 'B', RCONDE contains the reciprocal condition number for the average of the selected eigenvalues. Not referenced if SENSE = 'N' or 'V'.``` [out] RCONDV ``` RCONDV is REAL If SENSE = 'V' or 'B', RCONDV contains the reciprocal condition number for the selected right invariant subspace. Not referenced if SENSE = 'N' or 'E'.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,2*N). Also, if SENSE = 'E' or 'V' or 'B', LWORK >= 2*SDIM*(N-SDIM), where SDIM is the number of selected eigenvalues computed by this routine. Note that 2*SDIM*(N-SDIM) <= N*N/2. Note also that an error is only returned if LWORK < max(1,2*N), but if SENSE = 'E' or 'V' or 'B' this may not be large enough. For good performance, LWORK must generally be larger. If LWORK = -1, then a workspace query is assumed; the routine only calculates upper bound on the optimal size of the array WORK, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (N)` [out] BWORK ``` BWORK is LOGICAL array, dimension (N) Not referenced if SORT = 'N'.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. > 0: if INFO = i, and i is <= N: the QR algorithm failed to compute all the eigenvalues; elements 1:ILO-1 and i+1:N of W contain those eigenvalues which have converged; if JOBVS = 'V', VS contains the transformation which reduces A to its partially converged Schur form. = N+1: the eigenvalues could not be reordered because some eigenvalues were too close to separate (the problem is very ill-conditioned); = N+2: after reordering, roundoff changed values of some complex eigenvalues so that leading eigenvalues in the Schur form no longer satisfy SELECT=.TRUE. This could also be caused by underflow due to scaling.```
Date
November 2011

Definition at line 241 of file cgeesx.f.

241 *
242 * -- LAPACK driver routine (version 3.4.0) --
243 * -- LAPACK is a software package provided by Univ. of Tennessee, --
244 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
245 * November 2011
246 *
247 * .. Scalar Arguments ..
248  CHARACTER jobvs, sense, sort
249  INTEGER info, lda, ldvs, lwork, n, sdim
250  REAL rconde, rcondv
251 * ..
252 * .. Array Arguments ..
253  LOGICAL bwork( * )
254  REAL rwork( * )
255  COMPLEX a( lda, * ), vs( ldvs, * ), w( * ), work( * )
256 * ..
257 * .. Function Arguments ..
258  LOGICAL select
259  EXTERNAL SELECT
260 * ..
261 *
262 * =====================================================================
263 *
264 * .. Parameters ..
265  REAL zero, one
266  parameter( zero = 0.0e0, one = 1.0e0 )
267 * ..
268 * .. Local Scalars ..
269  LOGICAL lquery, scalea, wantsb, wantse, wantsn, wantst,
270  \$ wantsv, wantvs
271  INTEGER hswork, i, ibal, icond, ierr, ieval, ihi, ilo,
272  \$ itau, iwrk, lwrk, maxwrk, minwrk
273  REAL anrm, bignum, cscale, eps, smlnum
274 * ..
275 * .. Local Arrays ..
276  REAL dum( 1 )
277 * ..
278 * .. External Subroutines ..
279  EXTERNAL ccopy, cgebak, cgebal, cgehrd, chseqr, clacpy,
281 * ..
282 * .. External Functions ..
283  LOGICAL lsame
284  INTEGER ilaenv
285  REAL clange, slamch
286  EXTERNAL lsame, ilaenv, clange, slamch
287 * ..
288 * .. Intrinsic Functions ..
289  INTRINSIC max, sqrt
290 * ..
291 * .. Executable Statements ..
292 *
293 * Test the input arguments
294 *
295  info = 0
296  wantvs = lsame( jobvs, 'V' )
297  wantst = lsame( sort, 'S' )
298  wantsn = lsame( sense, 'N' )
299  wantse = lsame( sense, 'E' )
300  wantsv = lsame( sense, 'V' )
301  wantsb = lsame( sense, 'B' )
302  lquery = ( lwork.EQ.-1 )
303 *
304  IF( ( .NOT.wantvs ) .AND. ( .NOT.lsame( jobvs, 'N' ) ) ) THEN
305  info = -1
306  ELSE IF( ( .NOT.wantst ) .AND. ( .NOT.lsame( sort, 'N' ) ) ) THEN
307  info = -2
308  ELSE IF( .NOT.( wantsn .OR. wantse .OR. wantsv .OR. wantsb ) .OR.
309  \$ ( .NOT.wantst .AND. .NOT.wantsn ) ) THEN
310  info = -4
311  ELSE IF( n.LT.0 ) THEN
312  info = -5
313  ELSE IF( lda.LT.max( 1, n ) ) THEN
314  info = -7
315  ELSE IF( ldvs.LT.1 .OR. ( wantvs .AND. ldvs.LT.n ) ) THEN
316  info = -11
317  END IF
318 *
319 * Compute workspace
320 * (Note: Comments in the code beginning "Workspace:" describe the
321 * minimal amount of real workspace needed at that point in the
322 * code, as well as the preferred amount for good performance.
323 * CWorkspace refers to complex workspace, and RWorkspace to real
324 * workspace. NB refers to the optimal block size for the
325 * immediately following subroutine, as returned by ILAENV.
326 * HSWORK refers to the workspace preferred by CHSEQR, as
327 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
328 * the worst case.
329 * If SENSE = 'E', 'V' or 'B', then the amount of workspace needed
330 * depends on SDIM, which is computed by the routine CTRSEN later
331 * in the code.)
332 *
333  IF( info.EQ.0 ) THEN
334  IF( n.EQ.0 ) THEN
335  minwrk = 1
336  lwrk = 1
337  ELSE
338  maxwrk = n + n*ilaenv( 1, 'CGEHRD', ' ', n, 1, n, 0 )
339  minwrk = 2*n
340 *
341  CALL chseqr( 'S', jobvs, n, 1, n, a, lda, w, vs, ldvs,
342  \$ work, -1, ieval )
343  hswork = work( 1 )
344 *
345  IF( .NOT.wantvs ) THEN
346  maxwrk = max( maxwrk, hswork )
347  ELSE
348  maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
349  \$ ' ', n, 1, n, -1 ) )
350  maxwrk = max( maxwrk, hswork )
351  END IF
352  lwrk = maxwrk
353  IF( .NOT.wantsn )
354  \$ lwrk = max( lwrk, ( n*n )/2 )
355  END IF
356  work( 1 ) = lwrk
357 *
358  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
359  info = -15
360  END IF
361  END IF
362 *
363  IF( info.NE.0 ) THEN
364  CALL xerbla( 'CGEESX', -info )
365  RETURN
366  ELSE IF( lquery ) THEN
367  RETURN
368  END IF
369 *
370 * Quick return if possible
371 *
372  IF( n.EQ.0 ) THEN
373  sdim = 0
374  RETURN
375  END IF
376 *
377 * Get machine constants
378 *
379  eps = slamch( 'P' )
380  smlnum = slamch( 'S' )
381  bignum = one / smlnum
382  CALL slabad( smlnum, bignum )
383  smlnum = sqrt( smlnum ) / eps
384  bignum = one / smlnum
385 *
386 * Scale A if max element outside range [SMLNUM,BIGNUM]
387 *
388  anrm = clange( 'M', n, n, a, lda, dum )
389  scalea = .false.
390  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
391  scalea = .true.
392  cscale = smlnum
393  ELSE IF( anrm.GT.bignum ) THEN
394  scalea = .true.
395  cscale = bignum
396  END IF
397  IF( scalea )
398  \$ CALL clascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
399 *
400 *
401 * Permute the matrix to make it more nearly triangular
402 * (CWorkspace: none)
403 * (RWorkspace: need N)
404 *
405  ibal = 1
406  CALL cgebal( 'P', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
407 *
408 * Reduce to upper Hessenberg form
409 * (CWorkspace: need 2*N, prefer N+N*NB)
410 * (RWorkspace: none)
411 *
412  itau = 1
413  iwrk = n + itau
414  CALL cgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
415  \$ lwork-iwrk+1, ierr )
416 *
417  IF( wantvs ) THEN
418 *
419 * Copy Householder vectors to VS
420 *
421  CALL clacpy( 'L', n, n, a, lda, vs, ldvs )
422 *
423 * Generate unitary matrix in VS
424 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
425 * (RWorkspace: none)
426 *
427  CALL cunghr( n, ilo, ihi, vs, ldvs, work( itau ), work( iwrk ),
428  \$ lwork-iwrk+1, ierr )
429  END IF
430 *
431  sdim = 0
432 *
433 * Perform QR iteration, accumulating Schur vectors in VS if desired
434 * (CWorkspace: need 1, prefer HSWORK (see comments) )
435 * (RWorkspace: none)
436 *
437  iwrk = itau
438  CALL chseqr( 'S', jobvs, n, ilo, ihi, a, lda, w, vs, ldvs,
439  \$ work( iwrk ), lwork-iwrk+1, ieval )
440  IF( ieval.GT.0 )
441  \$ info = ieval
442 *
443 * Sort eigenvalues if desired
444 *
445  IF( wantst .AND. info.EQ.0 ) THEN
446  IF( scalea )
447  \$ CALL clascl( 'G', 0, 0, cscale, anrm, n, 1, w, n, ierr )
448  DO 10 i = 1, n
449  bwork( i ) = SELECT( w( i ) )
450  10 CONTINUE
451 *
452 * Reorder eigenvalues, transform Schur vectors, and compute
453 * reciprocal condition numbers
454 * (CWorkspace: if SENSE is not 'N', need 2*SDIM*(N-SDIM)
455 * otherwise, need none )
456 * (RWorkspace: none)
457 *
458  CALL ctrsen( sense, jobvs, bwork, n, a, lda, vs, ldvs, w, sdim,
459  \$ rconde, rcondv, work( iwrk ), lwork-iwrk+1,
460  \$ icond )
461  IF( .NOT.wantsn )
462  \$ maxwrk = max( maxwrk, 2*sdim*( n-sdim ) )
463  IF( icond.EQ.-14 ) THEN
464 *
465 * Not enough complex workspace
466 *
467  info = -15
468  END IF
469  END IF
470 *
471  IF( wantvs ) THEN
472 *
473 * Undo balancing
474 * (CWorkspace: none)
475 * (RWorkspace: need N)
476 *
477  CALL cgebak( 'P', 'R', n, ilo, ihi, rwork( ibal ), n, vs, ldvs,
478  \$ ierr )
479  END IF
480 *
481  IF( scalea ) THEN
482 *
483 * Undo scaling for the Schur form of A
484 *
485  CALL clascl( 'U', 0, 0, cscale, anrm, n, n, a, lda, ierr )
486  CALL ccopy( n, a, lda+1, w, 1 )
487  IF( ( wantsv .OR. wantsb ) .AND. info.EQ.0 ) THEN
488  dum( 1 ) = rcondv
489  CALL slascl( 'G', 0, 0, cscale, anrm, 1, 1, dum, 1, ierr )
490  rcondv = dum( 1 )
491  END IF
492  END IF
493 *
494  work( 1 ) = maxwrk
495  RETURN
496 *
497 * End of CGEESX
498 *
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine ctrsen(JOB, COMPQ, SELECT, N, T, LDT, Q, LDQ, W, M, S, SEP, WORK, LWORK, INFO)
CTRSEN
Definition: ctrsen.f:266
subroutine cgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CGEHRD
Definition: cgehrd.f:169
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cgebal(JOB, N, A, LDA, ILO, IHI, SCALE, INFO)
CGEBAL
Definition: cgebal.f:163
subroutine cgebak(JOB, SIDE, N, ILO, IHI, SCALE, M, V, LDV, INFO)
CGEBAK
Definition: cgebak.f:133
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine chseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, W, Z, LDZ, WORK, LWORK, INFO)
CHSEQR
Definition: chseqr.f:301
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine ccopy(N, CX, INCX, CY, INCY)
CCOPY
Definition: ccopy.f:52
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: slascl.f:141
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cunghr(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CUNGHR
Definition: cunghr.f:128

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cgeev ( character JOBVL, character JOBVR, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( * ) W, complex, dimension( ldvl, * ) VL, integer LDVL, complex, dimension( ldvr, * ) VR, integer LDVR, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, integer INFO )

CGEEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

Purpose:
``` CGEEV computes for an N-by-N complex nonsymmetric matrix A, the
eigenvalues and, optionally, the left and/or right eigenvectors.

The right eigenvector v(j) of A satisfies
A * v(j) = lambda(j) * v(j)
where lambda(j) is its eigenvalue.
The left eigenvector u(j) of A satisfies
u(j)**H * A = lambda(j) * u(j)**H
where u(j)**H denotes the conjugate transpose of u(j).

The computed eigenvectors are normalized to have Euclidean norm
equal to 1 and largest component real.```
Parameters
 [in] JOBVL ``` JOBVL is CHARACTER*1 = 'N': left eigenvectors of A are not computed; = 'V': left eigenvectors of are computed.``` [in] JOBVR ``` JOBVR is CHARACTER*1 = 'N': right eigenvectors of A are not computed; = 'V': right eigenvectors of A are computed.``` [in] N ``` N is INTEGER The order of the matrix A. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA,N) On entry, the N-by-N matrix A. On exit, A has been overwritten.``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= max(1,N).``` [out] W ``` W is COMPLEX array, dimension (N) W contains the computed eigenvalues.``` [out] VL ``` VL is COMPLEX array, dimension (LDVL,N) If JOBVL = 'V', the left eigenvectors u(j) are stored one after another in the columns of VL, in the same order as their eigenvalues. If JOBVL = 'N', VL is not referenced. u(j) = VL(:,j), the j-th column of VL.``` [in] LDVL ``` LDVL is INTEGER The leading dimension of the array VL. LDVL >= 1; if JOBVL = 'V', LDVL >= N.``` [out] VR ``` VR is COMPLEX array, dimension (LDVR,N) If JOBVR = 'V', the right eigenvectors v(j) are stored one after another in the columns of VR, in the same order as their eigenvalues. If JOBVR = 'N', VR is not referenced. v(j) = VR(:,j), the j-th column of VR.``` [in] LDVR ``` LDVR is INTEGER The leading dimension of the array VR. LDVR >= 1; if JOBVR = 'V', LDVR >= N.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,2*N). For good performance, LWORK must generally be larger. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (2*N)` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. > 0: if INFO = i, the QR algorithm failed to compute all the eigenvalues, and no eigenvectors have been computed; elements and i+1:N of W contain eigenvalues which have converged.```
Date
November 2011

Definition at line 179 of file cgeev.f.

179 *
180 * -- LAPACK driver routine (version 3.4.0) --
181 * -- LAPACK is a software package provided by Univ. of Tennessee, --
182 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
183 * November 2011
184 *
185 * .. Scalar Arguments ..
186  CHARACTER jobvl, jobvr
187  INTEGER info, lda, ldvl, ldvr, lwork, n
188 * ..
189 * .. Array Arguments ..
190  REAL rwork( * )
191  COMPLEX a( lda, * ), vl( ldvl, * ), vr( ldvr, * ),
192  \$ w( * ), work( * )
193 * ..
194 *
195 * =====================================================================
196 *
197 * .. Parameters ..
198  REAL zero, one
199  parameter( zero = 0.0e0, one = 1.0e0 )
200 * ..
201 * .. Local Scalars ..
202  LOGICAL lquery, scalea, wantvl, wantvr
203  CHARACTER side
204  INTEGER hswork, i, ibal, ierr, ihi, ilo, irwork, itau,
205  \$ iwrk, k, maxwrk, minwrk, nout
206  REAL anrm, bignum, cscale, eps, scl, smlnum
207  COMPLEX tmp
208 * ..
209 * .. Local Arrays ..
210  LOGICAL select( 1 )
211  REAL dum( 1 )
212 * ..
213 * .. External Subroutines ..
214  EXTERNAL cgebak, cgebal, cgehrd, chseqr, clacpy, clascl,
216 * ..
217 * .. External Functions ..
218  LOGICAL lsame
219  INTEGER ilaenv, isamax
220  REAL clange, scnrm2, slamch
221  EXTERNAL lsame, ilaenv, isamax, clange, scnrm2, slamch
222 * ..
223 * .. Intrinsic Functions ..
224  INTRINSIC aimag, cmplx, conjg, max, REAL, sqrt
225 * ..
226 * .. Executable Statements ..
227 *
228 * Test the input arguments
229 *
230  info = 0
231  lquery = ( lwork.EQ.-1 )
232  wantvl = lsame( jobvl, 'V' )
233  wantvr = lsame( jobvr, 'V' )
234  IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl, 'N' ) ) ) THEN
235  info = -1
236  ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr, 'N' ) ) ) THEN
237  info = -2
238  ELSE IF( n.LT.0 ) THEN
239  info = -3
240  ELSE IF( lda.LT.max( 1, n ) ) THEN
241  info = -5
242  ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) ) THEN
243  info = -8
244  ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) ) THEN
245  info = -10
246  END IF
247
248 *
249 * Compute workspace
250 * (Note: Comments in the code beginning "Workspace:" describe the
251 * minimal amount of workspace needed at that point in the code,
252 * as well as the preferred amount for good performance.
253 * CWorkspace refers to complex workspace, and RWorkspace to real
254 * workspace. NB refers to the optimal block size for the
255 * immediately following subroutine, as returned by ILAENV.
256 * HSWORK refers to the workspace preferred by CHSEQR, as
257 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
258 * the worst case.)
259 *
260  IF( info.EQ.0 ) THEN
261  IF( n.EQ.0 ) THEN
262  minwrk = 1
263  maxwrk = 1
264  ELSE
265  maxwrk = n + n*ilaenv( 1, 'CGEHRD', ' ', n, 1, n, 0 )
266  minwrk = 2*n
267  IF( wantvl ) THEN
268  maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
269  \$ ' ', n, 1, n, -1 ) )
270  CALL chseqr( 'S', 'V', n, 1, n, a, lda, w, vl, ldvl,
271  \$ work, -1, info )
272  ELSE IF( wantvr ) THEN
273  maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
274  \$ ' ', n, 1, n, -1 ) )
275  CALL chseqr( 'S', 'V', n, 1, n, a, lda, w, vr, ldvr,
276  \$ work, -1, info )
277  ELSE
278  CALL chseqr( 'E', 'N', n, 1, n, a, lda, w, vr, ldvr,
279  \$ work, -1, info )
280  END IF
281  hswork = work( 1 )
282  maxwrk = max( maxwrk, hswork, minwrk )
283  END IF
284  work( 1 ) = maxwrk
285 *
286  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
287  info = -12
288  END IF
289  END IF
290 *
291  IF( info.NE.0 ) THEN
292  CALL xerbla( 'CGEEV ', -info )
293  RETURN
294  ELSE IF( lquery ) THEN
295  RETURN
296  END IF
297 *
298 * Quick return if possible
299 *
300  IF( n.EQ.0 )
301  \$ RETURN
302 *
303 * Get machine constants
304 *
305  eps = slamch( 'P' )
306  smlnum = slamch( 'S' )
307  bignum = one / smlnum
308  CALL slabad( smlnum, bignum )
309  smlnum = sqrt( smlnum ) / eps
310  bignum = one / smlnum
311 *
312 * Scale A if max element outside range [SMLNUM,BIGNUM]
313 *
314  anrm = clange( 'M', n, n, a, lda, dum )
315  scalea = .false.
316  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
317  scalea = .true.
318  cscale = smlnum
319  ELSE IF( anrm.GT.bignum ) THEN
320  scalea = .true.
321  cscale = bignum
322  END IF
323  IF( scalea )
324  \$ CALL clascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
325 *
326 * Balance the matrix
327 * (CWorkspace: none)
328 * (RWorkspace: need N)
329 *
330  ibal = 1
331  CALL cgebal( 'B', n, a, lda, ilo, ihi, rwork( ibal ), ierr )
332 *
333 * Reduce to upper Hessenberg form
334 * (CWorkspace: need 2*N, prefer N+N*NB)
335 * (RWorkspace: none)
336 *
337  itau = 1
338  iwrk = itau + n
339  CALL cgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
340  \$ lwork-iwrk+1, ierr )
341 *
342  IF( wantvl ) THEN
343 *
344 * Want left eigenvectors
345 * Copy Householder vectors to VL
346 *
347  side = 'L'
348  CALL clacpy( 'L', n, n, a, lda, vl, ldvl )
349 *
350 * Generate unitary matrix in VL
351 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
352 * (RWorkspace: none)
353 *
354  CALL cunghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
355  \$ lwork-iwrk+1, ierr )
356 *
357 * Perform QR iteration, accumulating Schur vectors in VL
358 * (CWorkspace: need 1, prefer HSWORK (see comments) )
359 * (RWorkspace: none)
360 *
361  iwrk = itau
362  CALL chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,
363  \$ work( iwrk ), lwork-iwrk+1, info )
364 *
365  IF( wantvr ) THEN
366 *
367 * Want left and right eigenvectors
368 * Copy Schur vectors to VR
369 *
370  side = 'B'
371  CALL clacpy( 'F', n, n, vl, ldvl, vr, ldvr )
372  END IF
373 *
374  ELSE IF( wantvr ) THEN
375 *
376 * Want right eigenvectors
377 * Copy Householder vectors to VR
378 *
379  side = 'R'
380  CALL clacpy( 'L', n, n, a, lda, vr, ldvr )
381 *
382 * Generate unitary matrix in VR
383 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
384 * (RWorkspace: none)
385 *
386  CALL cunghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
387  \$ lwork-iwrk+1, ierr )
388 *
389 * Perform QR iteration, accumulating Schur vectors in VR
390 * (CWorkspace: need 1, prefer HSWORK (see comments) )
391 * (RWorkspace: none)
392 *
393  iwrk = itau
394  CALL chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,
395  \$ work( iwrk ), lwork-iwrk+1, info )
396 *
397  ELSE
398 *
399 * Compute eigenvalues only
400 * (CWorkspace: need 1, prefer HSWORK (see comments) )
401 * (RWorkspace: none)
402 *
403  iwrk = itau
404  CALL chseqr( 'E', 'N', n, ilo, ihi, a, lda, w, vr, ldvr,
405  \$ work( iwrk ), lwork-iwrk+1, info )
406  END IF
407 *
408 * If INFO > 0 from CHSEQR, then quit
409 *
410  IF( info.GT.0 )
411  \$ GO TO 50
412 *
413  IF( wantvl .OR. wantvr ) THEN
414 *
415 * Compute left and/or right eigenvectors
416 * (CWorkspace: need 2*N)
417 * (RWorkspace: need 2*N)
418 *
419  irwork = ibal + n
420  CALL ctrevc( side, 'B', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
421  \$ n, nout, work( iwrk ), rwork( irwork ), ierr )
422  END IF
423 *
424  IF( wantvl ) THEN
425 *
426 * Undo balancing of left eigenvectors
427 * (CWorkspace: none)
428 * (RWorkspace: need N)
429 *
430  CALL cgebak( 'B', 'L', n, ilo, ihi, rwork( ibal ), n, vl, ldvl,
431  \$ ierr )
432 *
433 * Normalize left eigenvectors and make largest component real
434 *
435  DO 20 i = 1, n
436  scl = one / scnrm2( n, vl( 1, i ), 1 )
437  CALL csscal( n, scl, vl( 1, i ), 1 )
438  DO 10 k = 1, n
439  rwork( irwork+k-1 ) = REAL( VL( K, I ) )**2 +
440  \$ aimag( vl( k, i ) )**2
441  10 CONTINUE
442  k = isamax( n, rwork( irwork ), 1 )
443  tmp = conjg( vl( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
444  CALL cscal( n, tmp, vl( 1, i ), 1 )
445  vl( k, i ) = cmplx( REAL( VL( K, I ) ), zero )
446  20 CONTINUE
447  END IF
448 *
449  IF( wantvr ) THEN
450 *
451 * Undo balancing of right eigenvectors
452 * (CWorkspace: none)
453 * (RWorkspace: need N)
454 *
455  CALL cgebak( 'B', 'R', n, ilo, ihi, rwork( ibal ), n, vr, ldvr,
456  \$ ierr )
457 *
458 * Normalize right eigenvectors and make largest component real
459 *
460  DO 40 i = 1, n
461  scl = one / scnrm2( n, vr( 1, i ), 1 )
462  CALL csscal( n, scl, vr( 1, i ), 1 )
463  DO 30 k = 1, n
464  rwork( irwork+k-1 ) = REAL( VR( K, I ) )**2 +
465  \$ aimag( vr( k, i ) )**2
466  30 CONTINUE
467  k = isamax( n, rwork( irwork ), 1 )
468  tmp = conjg( vr( k, i ) ) / sqrt( rwork( irwork+k-1 ) )
469  CALL cscal( n, tmp, vr( 1, i ), 1 )
470  vr( k, i ) = cmplx( REAL( VR( K, I ) ), zero )
471  40 CONTINUE
472  END IF
473 *
474 * Undo scaling if necessary
475 *
476  50 CONTINUE
477  IF( scalea ) THEN
478  CALL clascl( 'G', 0, 0, cscale, anrm, n-info, 1, w( info+1 ),
479  \$ max( n-info, 1 ), ierr )
480  IF( info.GT.0 ) THEN
481  CALL clascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, w, n, ierr )
482  END IF
483  END IF
484 *
485  work( 1 ) = maxwrk
486  RETURN
487 *
488 * End of CGEEV
489 *
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
integer function isamax(N, SX, INCX)
ISAMAX
Definition: isamax.f:53
subroutine cgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CGEHRD
Definition: cgehrd.f:169
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cgebal(JOB, N, A, LDA, ILO, IHI, SCALE, INFO)
CGEBAL
Definition: cgebal.f:163
subroutine cgebak(JOB, SIDE, N, ILO, IHI, SCALE, M, V, LDV, INFO)
CGEBAK
Definition: cgebak.f:133
subroutine csscal(N, SA, CX, INCX)
CSSCAL
Definition: csscal.f:54
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
real function scnrm2(N, X, INCX)
SCNRM2
Definition: scnrm2.f:56
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine chseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, W, Z, LDZ, WORK, LWORK, INFO)
CHSEQR
Definition: chseqr.f:301
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine ctrevc(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTREVC
Definition: ctrevc.f:220
subroutine cscal(N, CA, CX, INCX)
CSCAL
Definition: cscal.f:54
subroutine cunghr(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CUNGHR
Definition: cunghr.f:128

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cgeevx ( character BALANC, character JOBVL, character JOBVR, character SENSE, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( * ) W, complex, dimension( ldvl, * ) VL, integer LDVL, complex, dimension( ldvr, * ) VR, integer LDVR, integer ILO, integer IHI, real, dimension( * ) SCALE, real ABNRM, real, dimension( * ) RCONDE, real, dimension( * ) RCONDV, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, integer INFO )

CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

Purpose:
``` CGEEVX computes for an N-by-N complex nonsymmetric matrix A, the
eigenvalues and, optionally, the left and/or right eigenvectors.

Optionally also, it computes a balancing transformation to improve
the conditioning of the eigenvalues and eigenvectors (ILO, IHI,
SCALE, and ABNRM), reciprocal condition numbers for the eigenvalues
(RCONDE), and reciprocal condition numbers for the right
eigenvectors (RCONDV).

The right eigenvector v(j) of A satisfies
A * v(j) = lambda(j) * v(j)
where lambda(j) is its eigenvalue.
The left eigenvector u(j) of A satisfies
u(j)**H * A = lambda(j) * u(j)**H
where u(j)**H denotes the conjugate transpose of u(j).

The computed eigenvectors are normalized to have Euclidean norm
equal to 1 and largest component real.

Balancing a matrix means permuting the rows and columns to make it
more nearly upper triangular, and applying a diagonal similarity
transformation D * A * D**(-1), where D is a diagonal matrix, to
make its rows and columns closer in norm and the condition numbers
of its eigenvalues and eigenvectors smaller.  The computed
reciprocal condition numbers correspond to the balanced matrix.
Permuting rows and columns will not change the condition numbers
(in exact arithmetic) but diagonal scaling will.  For further
explanation of balancing, see section 4.10.2 of the LAPACK
Users' Guide.```
Parameters
 [in] BALANC ``` BALANC is CHARACTER*1 Indicates how the input matrix should be diagonally scaled and/or permuted to improve the conditioning of its eigenvalues. = 'N': Do not diagonally scale or permute; = 'P': Perform permutations to make the matrix more nearly upper triangular. Do not diagonally scale; = 'S': Diagonally scale the matrix, ie. replace A by D*A*D**(-1), where D is a diagonal matrix chosen to make the rows and columns of A more equal in norm. Do not permute; = 'B': Both diagonally scale and permute A. Computed reciprocal condition numbers will be for the matrix after balancing and/or permuting. Permuting does not change condition numbers (in exact arithmetic), but balancing does.``` [in] JOBVL ``` JOBVL is CHARACTER*1 = 'N': left eigenvectors of A are not computed; = 'V': left eigenvectors of A are computed. If SENSE = 'E' or 'B', JOBVL must = 'V'.``` [in] JOBVR ``` JOBVR is CHARACTER*1 = 'N': right eigenvectors of A are not computed; = 'V': right eigenvectors of A are computed. If SENSE = 'E' or 'B', JOBVR must = 'V'.``` [in] SENSE ``` SENSE is CHARACTER*1 Determines which reciprocal condition numbers are computed. = 'N': None are computed; = 'E': Computed for eigenvalues only; = 'V': Computed for right eigenvectors only; = 'B': Computed for eigenvalues and right eigenvectors. If SENSE = 'E' or 'B', both left and right eigenvectors must also be computed (JOBVL = 'V' and JOBVR = 'V').``` [in] N ``` N is INTEGER The order of the matrix A. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA,N) On entry, the N-by-N matrix A. On exit, A has been overwritten. If JOBVL = 'V' or JOBVR = 'V', A contains the Schur form of the balanced version of the matrix A.``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= max(1,N).``` [out] W ``` W is COMPLEX array, dimension (N) W contains the computed eigenvalues.``` [out] VL ``` VL is COMPLEX array, dimension (LDVL,N) If JOBVL = 'V', the left eigenvectors u(j) are stored one after another in the columns of VL, in the same order as their eigenvalues. If JOBVL = 'N', VL is not referenced. u(j) = VL(:,j), the j-th column of VL.``` [in] LDVL ``` LDVL is INTEGER The leading dimension of the array VL. LDVL >= 1; if JOBVL = 'V', LDVL >= N.``` [out] VR ``` VR is COMPLEX array, dimension (LDVR,N) If JOBVR = 'V', the right eigenvectors v(j) are stored one after another in the columns of VR, in the same order as their eigenvalues. If JOBVR = 'N', VR is not referenced. v(j) = VR(:,j), the j-th column of VR.``` [in] LDVR ``` LDVR is INTEGER The leading dimension of the array VR. LDVR >= 1; if JOBVR = 'V', LDVR >= N.``` [out] ILO ` ILO is INTEGER` [out] IHI ``` IHI is INTEGER ILO and IHI are integer values determined when A was balanced. The balanced A(i,j) = 0 if I > J and J = 1,...,ILO-1 or I = IHI+1,...,N.``` [out] SCALE ``` SCALE is REAL array, dimension (N) Details of the permutations and scaling factors applied when balancing A. If P(j) is the index of the row and column interchanged with row and column j, and D(j) is the scaling factor applied to row and column j, then SCALE(J) = P(J), for J = 1,...,ILO-1 = D(J), for J = ILO,...,IHI = P(J) for J = IHI+1,...,N. The order in which the interchanges are made is N to IHI+1, then 1 to ILO-1.``` [out] ABNRM ``` ABNRM is REAL The one-norm of the balanced matrix (the maximum of the sum of absolute values of elements of any column).``` [out] RCONDE ``` RCONDE is REAL array, dimension (N) RCONDE(j) is the reciprocal condition number of the j-th eigenvalue.``` [out] RCONDV ``` RCONDV is REAL array, dimension (N) RCONDV(j) is the reciprocal condition number of the j-th right eigenvector.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. If SENSE = 'N' or 'E', LWORK >= max(1,2*N), and if SENSE = 'V' or 'B', LWORK >= N*N+2*N. For good performance, LWORK must generally be larger. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (2*N)` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. > 0: if INFO = i, the QR algorithm failed to compute all the eigenvalues, and no eigenvectors or condition numbers have been computed; elements 1:ILO-1 and i+1:N of W contain eigenvalues which have converged.```
Date
November 2011

Definition at line 287 of file cgeevx.f.

287 *
288 * -- LAPACK driver routine (version 3.4.0) --
289 * -- LAPACK is a software package provided by Univ. of Tennessee, --
290 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
291 * November 2011
292 *
293 * .. Scalar Arguments ..
294  CHARACTER balanc, jobvl, jobvr, sense
295  INTEGER ihi, ilo, info, lda, ldvl, ldvr, lwork, n
296  REAL abnrm
297 * ..
298 * .. Array Arguments ..
299  REAL rconde( * ), rcondv( * ), rwork( * ),
300  \$ scale( * )
301  COMPLEX a( lda, * ), vl( ldvl, * ), vr( ldvr, * ),
302  \$ w( * ), work( * )
303 * ..
304 *
305 * =====================================================================
306 *
307 * .. Parameters ..
308  REAL zero, one
309  parameter( zero = 0.0e0, one = 1.0e0 )
310 * ..
311 * .. Local Scalars ..
312  LOGICAL lquery, scalea, wantvl, wantvr, wntsnb, wntsne,
313  \$ wntsnn, wntsnv
314  CHARACTER job, side
315  INTEGER hswork, i, icond, ierr, itau, iwrk, k, maxwrk,
316  \$ minwrk, nout
317  REAL anrm, bignum, cscale, eps, scl, smlnum
318  COMPLEX tmp
319 * ..
320 * .. Local Arrays ..
321  LOGICAL select( 1 )
322  REAL dum( 1 )
323 * ..
324 * .. External Subroutines ..
325  EXTERNAL cgebak, cgebal, cgehrd, chseqr, clacpy, clascl,
327  \$ slascl, xerbla
328 * ..
329 * .. External Functions ..
330  LOGICAL lsame
331  INTEGER ilaenv, isamax
332  REAL clange, scnrm2, slamch
333  EXTERNAL lsame, ilaenv, isamax, clange, scnrm2, slamch
334 * ..
335 * .. Intrinsic Functions ..
336  INTRINSIC aimag, cmplx, conjg, max, REAL, sqrt
337 * ..
338 * .. Executable Statements ..
339 *
340 * Test the input arguments
341 *
342  info = 0
343  lquery = ( lwork.EQ.-1 )
344  wantvl = lsame( jobvl, 'V' )
345  wantvr = lsame( jobvr, 'V' )
346  wntsnn = lsame( sense, 'N' )
347  wntsne = lsame( sense, 'E' )
348  wntsnv = lsame( sense, 'V' )
349  wntsnb = lsame( sense, 'B' )
350  IF( .NOT.( lsame( balanc, 'N' ) .OR. lsame( balanc, 'S' ) .OR.
351  \$ lsame( balanc, 'P' ) .OR. lsame( balanc, 'B' ) ) ) THEN
352  info = -1
353  ELSE IF( ( .NOT.wantvl ) .AND. ( .NOT.lsame( jobvl, 'N' ) ) ) THEN
354  info = -2
355  ELSE IF( ( .NOT.wantvr ) .AND. ( .NOT.lsame( jobvr, 'N' ) ) ) THEN
356  info = -3
357  ELSE IF( .NOT.( wntsnn .OR. wntsne .OR. wntsnb .OR. wntsnv ) .OR.
358  \$ ( ( wntsne .OR. wntsnb ) .AND. .NOT.( wantvl .AND.
359  \$ wantvr ) ) ) THEN
360  info = -4
361  ELSE IF( n.LT.0 ) THEN
362  info = -5
363  ELSE IF( lda.LT.max( 1, n ) ) THEN
364  info = -7
365  ELSE IF( ldvl.LT.1 .OR. ( wantvl .AND. ldvl.LT.n ) ) THEN
366  info = -10
367  ELSE IF( ldvr.LT.1 .OR. ( wantvr .AND. ldvr.LT.n ) ) THEN
368  info = -12
369  END IF
370 *
371 * Compute workspace
372 * (Note: Comments in the code beginning "Workspace:" describe the
373 * minimal amount of workspace needed at that point in the code,
374 * as well as the preferred amount for good performance.
375 * CWorkspace refers to complex workspace, and RWorkspace to real
376 * workspace. NB refers to the optimal block size for the
377 * immediately following subroutine, as returned by ILAENV.
378 * HSWORK refers to the workspace preferred by CHSEQR, as
379 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N,
380 * the worst case.)
381 *
382  IF( info.EQ.0 ) THEN
383  IF( n.EQ.0 ) THEN
384  minwrk = 1
385  maxwrk = 1
386  ELSE
387  maxwrk = n + n*ilaenv( 1, 'CGEHRD', ' ', n, 1, n, 0 )
388 *
389  IF( wantvl ) THEN
390  CALL chseqr( 'S', 'V', n, 1, n, a, lda, w, vl, ldvl,
391  \$ work, -1, info )
392  ELSE IF( wantvr ) THEN
393  CALL chseqr( 'S', 'V', n, 1, n, a, lda, w, vr, ldvr,
394  \$ work, -1, info )
395  ELSE
396  IF( wntsnn ) THEN
397  CALL chseqr( 'E', 'N', n, 1, n, a, lda, w, vr, ldvr,
398  \$ work, -1, info )
399  ELSE
400  CALL chseqr( 'S', 'N', n, 1, n, a, lda, w, vr, ldvr,
401  \$ work, -1, info )
402  END IF
403  END IF
404  hswork = work( 1 )
405 *
406  IF( ( .NOT.wantvl ) .AND. ( .NOT.wantvr ) ) THEN
407  minwrk = 2*n
408  IF( .NOT.( wntsnn .OR. wntsne ) )
409  \$ minwrk = max( minwrk, n*n + 2*n )
410  maxwrk = max( maxwrk, hswork )
411  IF( .NOT.( wntsnn .OR. wntsne ) )
412  \$ maxwrk = max( maxwrk, n*n + 2*n )
413  ELSE
414  minwrk = 2*n
415  IF( .NOT.( wntsnn .OR. wntsne ) )
416  \$ minwrk = max( minwrk, n*n + 2*n )
417  maxwrk = max( maxwrk, hswork )
418  maxwrk = max( maxwrk, n + ( n - 1 )*ilaenv( 1, 'CUNGHR',
419  \$ ' ', n, 1, n, -1 ) )
420  IF( .NOT.( wntsnn .OR. wntsne ) )
421  \$ maxwrk = max( maxwrk, n*n + 2*n )
422  maxwrk = max( maxwrk, 2*n )
423  END IF
424  maxwrk = max( maxwrk, minwrk )
425  END IF
426  work( 1 ) = maxwrk
427 *
428  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
429  info = -20
430  END IF
431  END IF
432 *
433  IF( info.NE.0 ) THEN
434  CALL xerbla( 'CGEEVX', -info )
435  RETURN
436  ELSE IF( lquery ) THEN
437  RETURN
438  END IF
439 *
440 * Quick return if possible
441 *
442  IF( n.EQ.0 )
443  \$ RETURN
444 *
445 * Get machine constants
446 *
447  eps = slamch( 'P' )
448  smlnum = slamch( 'S' )
449  bignum = one / smlnum
450  CALL slabad( smlnum, bignum )
451  smlnum = sqrt( smlnum ) / eps
452  bignum = one / smlnum
453 *
454 * Scale A if max element outside range [SMLNUM,BIGNUM]
455 *
456  icond = 0
457  anrm = clange( 'M', n, n, a, lda, dum )
458  scalea = .false.
459  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
460  scalea = .true.
461  cscale = smlnum
462  ELSE IF( anrm.GT.bignum ) THEN
463  scalea = .true.
464  cscale = bignum
465  END IF
466  IF( scalea )
467  \$ CALL clascl( 'G', 0, 0, anrm, cscale, n, n, a, lda, ierr )
468 *
469 * Balance the matrix and compute ABNRM
470 *
471  CALL cgebal( balanc, n, a, lda, ilo, ihi, scale, ierr )
472  abnrm = clange( '1', n, n, a, lda, dum )
473  IF( scalea ) THEN
474  dum( 1 ) = abnrm
475  CALL slascl( 'G', 0, 0, cscale, anrm, 1, 1, dum, 1, ierr )
476  abnrm = dum( 1 )
477  END IF
478 *
479 * Reduce to upper Hessenberg form
480 * (CWorkspace: need 2*N, prefer N+N*NB)
481 * (RWorkspace: none)
482 *
483  itau = 1
484  iwrk = itau + n
485  CALL cgehrd( n, ilo, ihi, a, lda, work( itau ), work( iwrk ),
486  \$ lwork-iwrk+1, ierr )
487 *
488  IF( wantvl ) THEN
489 *
490 * Want left eigenvectors
491 * Copy Householder vectors to VL
492 *
493  side = 'L'
494  CALL clacpy( 'L', n, n, a, lda, vl, ldvl )
495 *
496 * Generate unitary matrix in VL
497 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
498 * (RWorkspace: none)
499 *
500  CALL cunghr( n, ilo, ihi, vl, ldvl, work( itau ), work( iwrk ),
501  \$ lwork-iwrk+1, ierr )
502 *
503 * Perform QR iteration, accumulating Schur vectors in VL
504 * (CWorkspace: need 1, prefer HSWORK (see comments) )
505 * (RWorkspace: none)
506 *
507  iwrk = itau
508  CALL chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vl, ldvl,
509  \$ work( iwrk ), lwork-iwrk+1, info )
510 *
511  IF( wantvr ) THEN
512 *
513 * Want left and right eigenvectors
514 * Copy Schur vectors to VR
515 *
516  side = 'B'
517  CALL clacpy( 'F', n, n, vl, ldvl, vr, ldvr )
518  END IF
519 *
520  ELSE IF( wantvr ) THEN
521 *
522 * Want right eigenvectors
523 * Copy Householder vectors to VR
524 *
525  side = 'R'
526  CALL clacpy( 'L', n, n, a, lda, vr, ldvr )
527 *
528 * Generate unitary matrix in VR
529 * (CWorkspace: need 2*N-1, prefer N+(N-1)*NB)
530 * (RWorkspace: none)
531 *
532  CALL cunghr( n, ilo, ihi, vr, ldvr, work( itau ), work( iwrk ),
533  \$ lwork-iwrk+1, ierr )
534 *
535 * Perform QR iteration, accumulating Schur vectors in VR
536 * (CWorkspace: need 1, prefer HSWORK (see comments) )
537 * (RWorkspace: none)
538 *
539  iwrk = itau
540  CALL chseqr( 'S', 'V', n, ilo, ihi, a, lda, w, vr, ldvr,
541  \$ work( iwrk ), lwork-iwrk+1, info )
542 *
543  ELSE
544 *
545 * Compute eigenvalues only
546 * If condition numbers desired, compute Schur form
547 *
548  IF( wntsnn ) THEN
549  job = 'E'
550  ELSE
551  job = 'S'
552  END IF
553 *
554 * (CWorkspace: need 1, prefer HSWORK (see comments) )
555 * (RWorkspace: none)
556 *
557  iwrk = itau
558  CALL chseqr( job, 'N', n, ilo, ihi, a, lda, w, vr, ldvr,
559  \$ work( iwrk ), lwork-iwrk+1, info )
560  END IF
561 *
562 * If INFO > 0 from CHSEQR, then quit
563 *
564  IF( info.GT.0 )
565  \$ GO TO 50
566 *
567  IF( wantvl .OR. wantvr ) THEN
568 *
569 * Compute left and/or right eigenvectors
570 * (CWorkspace: need 2*N)
571 * (RWorkspace: need N)
572 *
573  CALL ctrevc( side, 'B', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
574  \$ n, nout, work( iwrk ), rwork, ierr )
575  END IF
576 *
577 * Compute condition numbers if desired
578 * (CWorkspace: need N*N+2*N unless SENSE = 'E')
579 * (RWorkspace: need 2*N unless SENSE = 'E')
580 *
581  IF( .NOT.wntsnn ) THEN
582  CALL ctrsna( sense, 'A', SELECT, n, a, lda, vl, ldvl, vr, ldvr,
583  \$ rconde, rcondv, n, nout, work( iwrk ), n, rwork,
584  \$ icond )
585  END IF
586 *
587  IF( wantvl ) THEN
588 *
589 * Undo balancing of left eigenvectors
590 *
591  CALL cgebak( balanc, 'L', n, ilo, ihi, scale, n, vl, ldvl,
592  \$ ierr )
593 *
594 * Normalize left eigenvectors and make largest component real
595 *
596  DO 20 i = 1, n
597  scl = one / scnrm2( n, vl( 1, i ), 1 )
598  CALL csscal( n, scl, vl( 1, i ), 1 )
599  DO 10 k = 1, n
600  rwork( k ) = REAL( VL( K, I ) )**2 +
601  \$ aimag( vl( k, i ) )**2
602  10 CONTINUE
603  k = isamax( n, rwork, 1 )
604  tmp = conjg( vl( k, i ) ) / sqrt( rwork( k ) )
605  CALL cscal( n, tmp, vl( 1, i ), 1 )
606  vl( k, i ) = cmplx( REAL( VL( K, I ) ), zero )
607  20 CONTINUE
608  END IF
609 *
610  IF( wantvr ) THEN
611 *
612 * Undo balancing of right eigenvectors
613 *
614  CALL cgebak( balanc, 'R', n, ilo, ihi, scale, n, vr, ldvr,
615  \$ ierr )
616 *
617 * Normalize right eigenvectors and make largest component real
618 *
619  DO 40 i = 1, n
620  scl = one / scnrm2( n, vr( 1, i ), 1 )
621  CALL csscal( n, scl, vr( 1, i ), 1 )
622  DO 30 k = 1, n
623  rwork( k ) = REAL( VR( K, I ) )**2 +
624  \$ aimag( vr( k, i ) )**2
625  30 CONTINUE
626  k = isamax( n, rwork, 1 )
627  tmp = conjg( vr( k, i ) ) / sqrt( rwork( k ) )
628  CALL cscal( n, tmp, vr( 1, i ), 1 )
629  vr( k, i ) = cmplx( REAL( VR( K, I ) ), zero )
630  40 CONTINUE
631  END IF
632 *
633 * Undo scaling if necessary
634 *
635  50 CONTINUE
636  IF( scalea ) THEN
637  CALL clascl( 'G', 0, 0, cscale, anrm, n-info, 1, w( info+1 ),
638  \$ max( n-info, 1 ), ierr )
639  IF( info.EQ.0 ) THEN
640  IF( ( wntsnv .OR. wntsnb ) .AND. icond.EQ.0 )
641  \$ CALL slascl( 'G', 0, 0, cscale, anrm, n, 1, rcondv, n,
642  \$ ierr )
643  ELSE
644  CALL clascl( 'G', 0, 0, cscale, anrm, ilo-1, 1, w, n, ierr )
645  END IF
646  END IF
647 *
648  work( 1 ) = maxwrk
649  RETURN
650 *
651 * End of CGEEVX
652 *
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
integer function isamax(N, SX, INCX)
ISAMAX
Definition: isamax.f:53
subroutine cgehrd(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CGEHRD
Definition: cgehrd.f:169
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cgebal(JOB, N, A, LDA, ILO, IHI, SCALE, INFO)
CGEBAL
Definition: cgebal.f:163
subroutine cgebak(JOB, SIDE, N, ILO, IHI, SCALE, M, V, LDV, INFO)
CGEBAK
Definition: cgebak.f:133
subroutine csscal(N, SA, CX, INCX)
CSSCAL
Definition: csscal.f:54
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
real function scnrm2(N, X, INCX)
SCNRM2
Definition: scnrm2.f:56
subroutine ctrsna(JOB, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, S, SEP, MM, M, WORK, LDWORK, RWORK, INFO)
CTRSNA
Definition: ctrsna.f:251
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine chseqr(JOB, COMPZ, N, ILO, IHI, H, LDH, W, Z, LDZ, WORK, LWORK, INFO)
CHSEQR
Definition: chseqr.f:301
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: slascl.f:141
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine ctrevc(SIDE, HOWMNY, SELECT, N, T, LDT, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTREVC
Definition: ctrevc.f:220
subroutine cscal(N, CA, CX, INCX)
CSCAL
Definition: cscal.f:54
subroutine cunghr(N, ILO, IHI, A, LDA, TAU, WORK, LWORK, INFO)
CUNGHR
Definition: cunghr.f:128

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cgegs ( character JOBVSL, character JOBVSR, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( ldb, * ) B, integer LDB, complex, dimension( * ) ALPHA, complex, dimension( * ) BETA, complex, dimension( ldvsl, * ) VSL, integer LDVSL, complex, dimension( ldvsr, * ) VSR, integer LDVSR, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, integer INFO )

CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

Purpose:
``` This routine is deprecated and has been replaced by routine CGGES.

CGEGS computes the eigenvalues, Schur form, and, optionally, the
left and or/right Schur vectors of a complex matrix pair (A,B).
Given two square matrices A and B, the generalized Schur
factorization has the form

A = Q*S*Z**H,  B = Q*T*Z**H

where Q and Z are unitary matrices and S and T are upper triangular.
The columns of Q are the left Schur vectors
and the columns of Z are the right Schur vectors.

If only the eigenvalues of (A,B) are needed, the driver routine
CGEGV should be used instead.  See CGEGV for a description of the
eigenvalues of the generalized nonsymmetric eigenvalue problem
(GNEP).```
Parameters
 [in] JOBVSL ``` JOBVSL is CHARACTER*1 = 'N': do not compute the left Schur vectors; = 'V': compute the left Schur vectors (returned in VSL).``` [in] JOBVSR ``` JOBVSR is CHARACTER*1 = 'N': do not compute the right Schur vectors; = 'V': compute the right Schur vectors (returned in VSR).``` [in] N ``` N is INTEGER The order of the matrices A, B, VSL, and VSR. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the matrix A. On exit, the upper triangular matrix S from the generalized Schur factorization.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. LDA >= max(1,N).``` [in,out] B ``` B is COMPLEX array, dimension (LDB, N) On entry, the matrix B. On exit, the upper triangular matrix T from the generalized Schur factorization.``` [in] LDB ``` LDB is INTEGER The leading dimension of B. LDB >= max(1,N).``` [out] ALPHA ``` ALPHA is COMPLEX array, dimension (N) The complex scalars alpha that define the eigenvalues of GNEP. ALPHA(j) = S(j,j), the diagonal element of the Schur form of A.``` [out] BETA ``` BETA is COMPLEX array, dimension (N) The non-negative real scalars beta that define the eigenvalues of GNEP. BETA(j) = T(j,j), the diagonal element of the triangular factor T. Together, the quantities alpha = ALPHA(j) and beta = BETA(j) represent the j-th eigenvalue of the matrix pair (A,B), in one of the forms lambda = alpha/beta or mu = beta/alpha. Since either lambda or mu may overflow, they should not, in general, be computed.``` [out] VSL ``` VSL is COMPLEX array, dimension (LDVSL,N) If JOBVSL = 'V', the matrix of left Schur vectors Q. Not referenced if JOBVSL = 'N'.``` [in] LDVSL ``` LDVSL is INTEGER The leading dimension of the matrix VSL. LDVSL >= 1, and if JOBVSL = 'V', LDVSL >= N.``` [out] VSR ``` VSR is COMPLEX array, dimension (LDVSR,N) If JOBVSR = 'V', the matrix of right Schur vectors Z. Not referenced if JOBVSR = 'N'.``` [in] LDVSR ``` LDVSR is INTEGER The leading dimension of the matrix VSR. LDVSR >= 1, and if JOBVSR = 'V', LDVSR >= N.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,2*N). For good performance, LWORK must generally be larger. To compute the optimal value of LWORK, call ILAENV to get blocksizes (for CGEQRF, CUNMQR, and CUNGQR.) Then compute: NB -- MAX of the blocksizes for CGEQRF, CUNMQR, and CUNGQR; the optimal LWORK is N*(NB+1). If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (3*N)` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. =1,...,N: The QZ iteration failed. (A,B) are not in Schur form, but ALPHA(j) and BETA(j) should be correct for j=INFO+1,...,N. > N: errors that usually indicate LAPACK problems: =N+1: error return from CGGBAL =N+2: error return from CGEQRF =N+3: error return from CUNMQR =N+4: error return from CUNGQR =N+5: error return from CGGHRD =N+6: error return from CHGEQZ (other than failed iteration) =N+7: error return from CGGBAK (computing VSL) =N+8: error return from CGGBAK (computing VSR) =N+9: error return from CLASCL (various places)```
Date
November 2011

Definition at line 227 of file cgegs.f.

227 *
228 * -- LAPACK driver routine (version 3.4.0) --
229 * -- LAPACK is a software package provided by Univ. of Tennessee, --
230 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
231 * November 2011
232 *
233 * .. Scalar Arguments ..
234  CHARACTER jobvsl, jobvsr
235  INTEGER info, lda, ldb, ldvsl, ldvsr, lwork, n
236 * ..
237 * .. Array Arguments ..
238  REAL rwork( * )
239  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
240  \$ beta( * ), vsl( ldvsl, * ), vsr( ldvsr, * ),
241  \$ work( * )
242 * ..
243 *
244 * =====================================================================
245 *
246 * .. Parameters ..
247  REAL zero, one
248  parameter( zero = 0.0e0, one = 1.0e0 )
249  COMPLEX czero, cone
250  parameter( czero = ( 0.0e0, 0.0e0 ),
251  \$ cone = ( 1.0e0, 0.0e0 ) )
252 * ..
253 * .. Local Scalars ..
254  LOGICAL ilascl, ilbscl, ilvsl, ilvsr, lquery
255  INTEGER icols, ihi, iinfo, ijobvl, ijobvr, ileft,
256  \$ ilo, iright, irows, irwork, itau, iwork,
257  \$ lopt, lwkmin, lwkopt, nb, nb1, nb2, nb3
258  REAL anrm, anrmto, bignum, bnrm, bnrmto, eps,
259  \$ safmin, smlnum
260 * ..
261 * .. External Subroutines ..
262  EXTERNAL cgeqrf, cggbak, cggbal, cgghrd, chgeqz, clacpy,
264 * ..
265 * .. External Functions ..
266  LOGICAL lsame
267  INTEGER ilaenv
268  REAL clange, slamch
269  EXTERNAL ilaenv, lsame, clange, slamch
270 * ..
271 * .. Intrinsic Functions ..
272  INTRINSIC int, max
273 * ..
274 * .. Executable Statements ..
275 *
276 * Decode the input arguments
277 *
278  IF( lsame( jobvsl, 'N' ) ) THEN
279  ijobvl = 1
280  ilvsl = .false.
281  ELSE IF( lsame( jobvsl, 'V' ) ) THEN
282  ijobvl = 2
283  ilvsl = .true.
284  ELSE
285  ijobvl = -1
286  ilvsl = .false.
287  END IF
288 *
289  IF( lsame( jobvsr, 'N' ) ) THEN
290  ijobvr = 1
291  ilvsr = .false.
292  ELSE IF( lsame( jobvsr, 'V' ) ) THEN
293  ijobvr = 2
294  ilvsr = .true.
295  ELSE
296  ijobvr = -1
297  ilvsr = .false.
298  END IF
299 *
300 * Test the input arguments
301 *
302  lwkmin = max( 2*n, 1 )
303  lwkopt = lwkmin
304  work( 1 ) = lwkopt
305  lquery = ( lwork.EQ.-1 )
306  info = 0
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( ldvsl.LT.1 .OR. ( ilvsl .AND. ldvsl.LT.n ) ) THEN
318  info = -11
319  ELSE IF( ldvsr.LT.1 .OR. ( ilvsr .AND. ldvsr.LT.n ) ) THEN
320  info = -13
321  ELSE IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
322  info = -15
323  END IF
324 *
325  IF( info.EQ.0 ) THEN
326  nb1 = ilaenv( 1, 'CGEQRF', ' ', n, n, -1, -1 )
327  nb2 = ilaenv( 1, 'CUNMQR', ' ', n, n, n, -1 )
328  nb3 = ilaenv( 1, 'CUNGQR', ' ', n, n, n, -1 )
329  nb = max( nb1, nb2, nb3 )
330  lopt = n*(nb+1)
331  work( 1 ) = lopt
332  END IF
333 *
334  IF( info.NE.0 ) THEN
335  CALL xerbla( 'CGEGS ', -info )
336  RETURN
337  ELSE IF( lquery ) THEN
338  RETURN
339  END IF
340 *
341 * Quick return if possible
342 *
343  IF( n.EQ.0 )
344  \$ RETURN
345 *
346 * Get machine constants
347 *
348  eps = slamch( 'E' )*slamch( 'B' )
349  safmin = slamch( 'S' )
350  smlnum = n*safmin / eps
351  bignum = one / smlnum
352 *
353 * Scale A if max element outside range [SMLNUM,BIGNUM]
354 *
355  anrm = clange( 'M', n, n, a, lda, rwork )
356  ilascl = .false.
357  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
358  anrmto = smlnum
359  ilascl = .true.
360  ELSE IF( anrm.GT.bignum ) THEN
361  anrmto = bignum
362  ilascl = .true.
363  END IF
364 *
365  IF( ilascl ) THEN
366  CALL clascl( 'G', -1, -1, anrm, anrmto, n, n, a, lda, iinfo )
367  IF( iinfo.NE.0 ) THEN
368  info = n + 9
369  RETURN
370  END IF
371  END IF
372 *
373 * Scale B if max element outside range [SMLNUM,BIGNUM]
374 *
375  bnrm = clange( 'M', n, n, b, ldb, rwork )
376  ilbscl = .false.
377  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
378  bnrmto = smlnum
379  ilbscl = .true.
380  ELSE IF( bnrm.GT.bignum ) THEN
381  bnrmto = bignum
382  ilbscl = .true.
383  END IF
384 *
385  IF( ilbscl ) THEN
386  CALL clascl( 'G', -1, -1, bnrm, bnrmto, n, n, b, ldb, iinfo )
387  IF( iinfo.NE.0 ) THEN
388  info = n + 9
389  RETURN
390  END IF
391  END IF
392 *
393 * Permute the matrix to make it more nearly triangular
394 *
395  ileft = 1
396  iright = n + 1
397  irwork = iright + n
398  iwork = 1
399  CALL cggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
400  \$ rwork( iright ), rwork( irwork ), iinfo )
401  IF( iinfo.NE.0 ) THEN
402  info = n + 1
403  GO TO 10
404  END IF
405 *
406 * Reduce B to triangular form, and initialize VSL and/or VSR
407 *
408  irows = ihi + 1 - ilo
409  icols = n + 1 - ilo
410  itau = iwork
411  iwork = itau + irows
412  CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
413  \$ work( iwork ), lwork+1-iwork, iinfo )
414  IF( iinfo.GE.0 )
415  \$ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
416  IF( iinfo.NE.0 ) THEN
417  info = n + 2
418  GO TO 10
419  END IF
420 *
421  CALL cunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
422  \$ work( itau ), a( ilo, ilo ), lda, work( iwork ),
423  \$ lwork+1-iwork, iinfo )
424  IF( iinfo.GE.0 )
425  \$ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
426  IF( iinfo.NE.0 ) THEN
427  info = n + 3
428  GO TO 10
429  END IF
430 *
431  IF( ilvsl ) THEN
432  CALL claset( 'Full', n, n, czero, cone, vsl, ldvsl )
433  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
434  \$ vsl( ilo+1, ilo ), ldvsl )
435  CALL cungqr( irows, irows, irows, vsl( ilo, ilo ), ldvsl,
436  \$ work( itau ), work( iwork ), lwork+1-iwork,
437  \$ iinfo )
438  IF( iinfo.GE.0 )
439  \$ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
440  IF( iinfo.NE.0 ) THEN
441  info = n + 4
442  GO TO 10
443  END IF
444  END IF
445 *
446  IF( ilvsr )
447  \$ CALL claset( 'Full', n, n, czero, cone, vsr, ldvsr )
448 *
449 * Reduce to generalized Hessenberg form
450 *
451  CALL cgghrd( jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb, vsl,
452  \$ ldvsl, vsr, ldvsr, iinfo )
453  IF( iinfo.NE.0 ) THEN
454  info = n + 5
455  GO TO 10
456  END IF
457 *
458 * Perform QZ algorithm, computing Schur vectors if desired
459 *
460  iwork = itau
461  CALL chgeqz( 'S', jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb,
462  \$ alpha, beta, vsl, ldvsl, vsr, ldvsr, work( iwork ),
463  \$ lwork+1-iwork, rwork( irwork ), iinfo )
464  IF( iinfo.GE.0 )
465  \$ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
466  IF( iinfo.NE.0 ) THEN
467  IF( iinfo.GT.0 .AND. iinfo.LE.n ) THEN
468  info = iinfo
469  ELSE IF( iinfo.GT.n .AND. iinfo.LE.2*n ) THEN
470  info = iinfo - n
471  ELSE
472  info = n + 6
473  END IF
474  GO TO 10
475  END IF
476 *
477 * Apply permutation to VSL and VSR
478 *
479  IF( ilvsl ) THEN
480  CALL cggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
481  \$ rwork( iright ), n, vsl, ldvsl, iinfo )
482  IF( iinfo.NE.0 ) THEN
483  info = n + 7
484  GO TO 10
485  END IF
486  END IF
487  IF( ilvsr ) THEN
488  CALL cggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
489  \$ rwork( iright ), n, vsr, ldvsr, iinfo )
490  IF( iinfo.NE.0 ) THEN
491  info = n + 8
492  GO TO 10
493  END IF
494  END IF
495 *
496 * Undo scaling
497 *
498  IF( ilascl ) THEN
499  CALL clascl( 'U', -1, -1, anrmto, anrm, n, n, a, lda, iinfo )
500  IF( iinfo.NE.0 ) THEN
501  info = n + 9
502  RETURN
503  END IF
504  CALL clascl( 'G', -1, -1, anrmto, anrm, n, 1, alpha, n, iinfo )
505  IF( iinfo.NE.0 ) THEN
506  info = n + 9
507  RETURN
508  END IF
509  END IF
510 *
511  IF( ilbscl ) THEN
512  CALL clascl( 'U', -1, -1, bnrmto, bnrm, n, n, b, ldb, iinfo )
513  IF( iinfo.NE.0 ) THEN
514  info = n + 9
515  RETURN
516  END IF
517  CALL clascl( 'G', -1, -1, bnrmto, bnrm, n, 1, beta, n, iinfo )
518  IF( iinfo.NE.0 ) THEN
519  info = n + 9
520  RETURN
521  END IF
522  END IF
523 *
524  10 CONTINUE
525  work( 1 ) = lwkopt
526 *
527  RETURN
528 *
529 * End of CGEGS
530 *
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:170
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:179
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:150
subroutine cgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
CGGHRD
Definition: cgghrd.f:206
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:130

Here is the call graph for this function:

 subroutine cgegv ( character JOBVL, character JOBVR, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( ldb, * ) B, integer LDB, complex, dimension( * ) ALPHA, complex, dimension( * ) BETA, complex, dimension( ldvl, * ) VL, integer LDVL, complex, dimension( ldvr, * ) VR, integer LDVR, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, integer INFO )

CGEEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

Purpose:
``` This routine is deprecated and has been replaced by routine CGGEV.

CGEGV computes the eigenvalues and, optionally, the left and/or right
eigenvectors of a complex matrix pair (A,B).
Given two square matrices A and B,
the generalized nonsymmetric eigenvalue problem (GNEP) is to find the
eigenvalues lambda and corresponding (non-zero) eigenvectors x such
that
A*x = lambda*B*x.

An alternate form is to find the eigenvalues mu and corresponding
eigenvectors y such that
mu*A*y = B*y.

These two forms are equivalent with mu = 1/lambda and x = y if
neither lambda nor mu is zero.  In order to deal with the case that
lambda or mu is zero or small, two values alpha and beta are returned
for each eigenvalue, such that lambda = alpha/beta and
mu = beta/alpha.

The vectors x and y in the above equations are right eigenvectors of
the matrix pair (A,B).  Vectors u and v satisfying
u**H*A = lambda*u**H*B  or  mu*v**H*A = v**H*B
are left eigenvectors of (A,B).

Note: this routine performs "full balancing" on A and B```
Parameters
 [in] JOBVL ``` JOBVL is CHARACTER*1 = 'N': do not compute the left generalized eigenvectors; = 'V': compute the left generalized eigenvectors (returned in VL).``` [in] JOBVR ``` JOBVR is CHARACTER*1 = 'N': do not compute the right generalized eigenvectors; = 'V': compute the right generalized eigenvectors (returned in VR).``` [in] N ``` N is INTEGER The order of the matrices A, B, VL, and VR. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the matrix A. If JOBVL = 'V' or JOBVR = 'V', then on exit A contains the Schur form of A from the generalized Schur factorization of the pair (A,B) after balancing. If no eigenvectors were computed, then only the diagonal elements of the Schur form will be correct. See CGGHRD and CHGEQZ for details.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. LDA >= max(1,N).``` [in,out] B ``` B is COMPLEX array, dimension (LDB, N) On entry, the matrix B. If JOBVL = 'V' or JOBVR = 'V', then on exit B contains the upper triangular matrix obtained from B in the generalized Schur factorization of the pair (A,B) after balancing. If no eigenvectors were computed, then only the diagonal elements of B will be correct. See CGGHRD and CHGEQZ for details.``` [in] LDB ``` LDB is INTEGER The leading dimension of B. LDB >= max(1,N).``` [out] ALPHA ``` ALPHA is COMPLEX array, dimension (N) The complex scalars alpha that define the eigenvalues of GNEP.``` [out] BETA ``` BETA is COMPLEX array, dimension (N) The complex scalars beta that define the eigenvalues of GNEP. Together, the quantities alpha = ALPHA(j) and beta = BETA(j) represent the j-th eigenvalue of the matrix pair (A,B), in one of the forms lambda = alpha/beta or mu = beta/alpha. Since either lambda or mu may overflow, they should not, in general, be computed.``` [out] VL ``` VL is COMPLEX array, dimension (LDVL,N) If JOBVL = 'V', the left eigenvectors u(j) are stored in the columns of VL, in the same order as their eigenvalues. Each eigenvector is scaled so that its largest component has abs(real part) + abs(imag. part) = 1, except for eigenvectors corresponding to an eigenvalue with alpha = beta = 0, which are set to zero. Not referenced if JOBVL = 'N'.``` [in] LDVL ``` LDVL is INTEGER The leading dimension of the matrix VL. LDVL >= 1, and if JOBVL = 'V', LDVL >= N.``` [out] VR ``` VR is COMPLEX array, dimension (LDVR,N) If JOBVR = 'V', the right eigenvectors x(j) are stored in the columns of VR, in the same order as their eigenvalues. Each eigenvector is scaled so that its largest component has abs(real part) + abs(imag. part) = 1, except for eigenvectors corresponding to an eigenvalue with alpha = beta = 0, which are set to zero. Not referenced if JOBVR = 'N'.``` [in] LDVR ``` LDVR is INTEGER The leading dimension of the matrix VR. LDVR >= 1, and if JOBVR = 'V', LDVR >= N.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,2*N). For good performance, LWORK must generally be larger. To compute the optimal value of LWORK, call ILAENV to get blocksizes (for CGEQRF, CUNMQR, and CUNGQR.) Then compute: NB -- MAX of the blocksizes for CGEQRF, CUNMQR, and CUNGQR; The optimal LWORK is MAX( 2*N, N*(NB+1) ). If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (8*N)` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. =1,...,N: The QZ iteration failed. No eigenvectors have been calculated, but ALPHA(j) and BETA(j) should be correct for j=INFO+1,...,N. > N: errors that usually indicate LAPACK problems: =N+1: error return from CGGBAL =N+2: error return from CGEQRF =N+3: error return from CUNMQR =N+4: error return from CUNGQR =N+5: error return from CGGHRD =N+6: error return from CHGEQZ (other than failed iteration) =N+7: error return from CTGEVC =N+8: error return from CGGBAK (computing VL) =N+9: error return from CGGBAK (computing VR) =N+10: error return from CLASCL (various calls)```
Date
November 2011
Further Details:
```  Balancing
---------

This driver calls CGGBAL to both permute and scale rows and columns
of A and B.  The permutations PL and PR are chosen so that PL*A*PR
and PL*B*R will be upper triangular except for the diagonal blocks
A(i:j,i:j) and B(i:j,i:j), with i and j as close together as
possible.  The diagonal scaling matrices DL and DR are chosen so
that the pair  DL*PL*A*PR*DR, DL*PL*B*PR*DR have elements close to
one (except for the elements that start out zero.)

After the eigenvalues and eigenvectors of the balanced matrices
have been computed, CGGBAK transforms the eigenvectors back to what
they would have been (in perfect arithmetic) if they had not been
balanced.

Contents of A and B on Exit
-------- -- - --- - -- ----

If any eigenvectors are computed (either JOBVL='V' or JOBVR='V' or
both), then on exit the arrays A and B will contain the complex Schur
form[*] of the "balanced" versions of A and B.  If no eigenvectors
are computed, then only the diagonal blocks will be correct.

[*] In other words, upper triangular form.```

Definition at line 284 of file cgegv.f.

284 *
285 * -- LAPACK driver routine (version 3.4.0) --
286 * -- LAPACK is a software package provided by Univ. of Tennessee, --
287 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
288 * November 2011
289 *
290 * .. Scalar Arguments ..
291  CHARACTER jobvl, jobvr
292  INTEGER info, lda, ldb, ldvl, ldvr, lwork, n
293 * ..
294 * .. Array Arguments ..
295  REAL rwork( * )
296  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
297  \$ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
298  \$ work( * )
299 * ..
300 *
301 * =====================================================================
302 *
303 * .. Parameters ..
304  REAL zero, one
305  parameter( zero = 0.0e0, one = 1.0e0 )
306  COMPLEX czero, cone
307  parameter( czero = ( 0.0e0, 0.0e0 ),
308  \$ cone = ( 1.0e0, 0.0e0 ) )
309 * ..
310 * .. Local Scalars ..
311  LOGICAL ilimit, ilv, ilvl, ilvr, lquery
312  CHARACTER chtemp
313  INTEGER icols, ihi, iinfo, ijobvl, ijobvr, ileft, ilo,
314  \$ in, iright, irows, irwork, itau, iwork, jc, jr,
315  \$ lopt, lwkmin, lwkopt, nb, nb1, nb2, nb3
316  REAL absai, absar, absb, anrm, anrm1, anrm2, bnrm,
317  \$ bnrm1, bnrm2, eps, safmax, safmin, salfai,
318  \$ salfar, sbeta, scale, temp
319  COMPLEX x
320 * ..
321 * .. Local Arrays ..
322  LOGICAL ldumma( 1 )
323 * ..
324 * .. External Subroutines ..
325  EXTERNAL cgeqrf, cggbak, cggbal, cgghrd, chgeqz, clacpy,
327 * ..
328 * .. External Functions ..
329  LOGICAL lsame
330  INTEGER ilaenv
331  REAL clange, slamch
332  EXTERNAL ilaenv, lsame, clange, slamch
333 * ..
334 * .. Intrinsic Functions ..
335  INTRINSIC abs, aimag, cmplx, int, max, real
336 * ..
337 * .. Statement Functions ..
338  REAL abs1
339 * ..
340 * .. Statement Function definitions ..
341  abs1( x ) = abs( REAL( X ) ) + abs( aimag( x ) )
342 * ..
343 * .. Executable Statements ..
344 *
345 * Decode the input arguments
346 *
347  IF( lsame( jobvl, 'N' ) ) THEN
348  ijobvl = 1
349  ilvl = .false.
350  ELSE IF( lsame( jobvl, 'V' ) ) THEN
351  ijobvl = 2
352  ilvl = .true.
353  ELSE
354  ijobvl = -1
355  ilvl = .false.
356  END IF
357 *
358  IF( lsame( jobvr, 'N' ) ) THEN
359  ijobvr = 1
360  ilvr = .false.
361  ELSE IF( lsame( jobvr, 'V' ) ) THEN
362  ijobvr = 2
363  ilvr = .true.
364  ELSE
365  ijobvr = -1
366  ilvr = .false.
367  END IF
368  ilv = ilvl .OR. ilvr
369 *
370 * Test the input arguments
371 *
372  lwkmin = max( 2*n, 1 )
373  lwkopt = lwkmin
374  work( 1 ) = lwkopt
375  lquery = ( lwork.EQ.-1 )
376  info = 0
377  IF( ijobvl.LE.0 ) THEN
378  info = -1
379  ELSE IF( ijobvr.LE.0 ) THEN
380  info = -2
381  ELSE IF( n.LT.0 ) THEN
382  info = -3
383  ELSE IF( lda.LT.max( 1, n ) ) THEN
384  info = -5
385  ELSE IF( ldb.LT.max( 1, n ) ) THEN
386  info = -7
387  ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
388  info = -11
389  ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
390  info = -13
391  ELSE IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
392  info = -15
393  END IF
394 *
395  IF( info.EQ.0 ) THEN
396  nb1 = ilaenv( 1, 'CGEQRF', ' ', n, n, -1, -1 )
397  nb2 = ilaenv( 1, 'CUNMQR', ' ', n, n, n, -1 )
398  nb3 = ilaenv( 1, 'CUNGQR', ' ', n, n, n, -1 )
399  nb = max( nb1, nb2, nb3 )
400  lopt = max( 2*n, n*(nb+1) )
401  work( 1 ) = lopt
402  END IF
403 *
404  IF( info.NE.0 ) THEN
405  CALL xerbla( 'CGEGV ', -info )
406  RETURN
407  ELSE IF( lquery ) THEN
408  RETURN
409  END IF
410 *
411 * Quick return if possible
412 *
413  IF( n.EQ.0 )
414  \$ RETURN
415 *
416 * Get machine constants
417 *
418  eps = slamch( 'E' )*slamch( 'B' )
419  safmin = slamch( 'S' )
420  safmin = safmin + safmin
421  safmax = one / safmin
422 *
423 * Scale A
424 *
425  anrm = clange( 'M', n, n, a, lda, rwork )
426  anrm1 = anrm
427  anrm2 = one
428  IF( anrm.LT.one ) THEN
429  IF( safmax*anrm.LT.one ) THEN
430  anrm1 = safmin
431  anrm2 = safmax*anrm
432  END IF
433  END IF
434 *
435  IF( anrm.GT.zero ) THEN
436  CALL clascl( 'G', -1, -1, anrm, one, n, n, a, lda, iinfo )
437  IF( iinfo.NE.0 ) THEN
438  info = n + 10
439  RETURN
440  END IF
441  END IF
442 *
443 * Scale B
444 *
445  bnrm = clange( 'M', n, n, b, ldb, rwork )
446  bnrm1 = bnrm
447  bnrm2 = one
448  IF( bnrm.LT.one ) THEN
449  IF( safmax*bnrm.LT.one ) THEN
450  bnrm1 = safmin
451  bnrm2 = safmax*bnrm
452  END IF
453  END IF
454 *
455  IF( bnrm.GT.zero ) THEN
456  CALL clascl( 'G', -1, -1, bnrm, one, n, n, b, ldb, iinfo )
457  IF( iinfo.NE.0 ) THEN
458  info = n + 10
459  RETURN
460  END IF
461  END IF
462 *
463 * Permute the matrix to make it more nearly triangular
464 * Also "balance" the matrix.
465 *
466  ileft = 1
467  iright = n + 1
468  irwork = iright + n
469  CALL cggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
470  \$ rwork( iright ), rwork( irwork ), iinfo )
471  IF( iinfo.NE.0 ) THEN
472  info = n + 1
473  GO TO 80
474  END IF
475 *
476 * Reduce B to triangular form, and initialize VL and/or VR
477 *
478  irows = ihi + 1 - ilo
479  IF( ilv ) THEN
480  icols = n + 1 - ilo
481  ELSE
482  icols = irows
483  END IF
484  itau = 1
485  iwork = itau + irows
486  CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
487  \$ work( iwork ), lwork+1-iwork, iinfo )
488  IF( iinfo.GE.0 )
489  \$ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
490  IF( iinfo.NE.0 ) THEN
491  info = n + 2
492  GO TO 80
493  END IF
494 *
495  CALL cunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
496  \$ work( itau ), a( ilo, ilo ), lda, work( iwork ),
497  \$ lwork+1-iwork, iinfo )
498  IF( iinfo.GE.0 )
499  \$ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
500  IF( iinfo.NE.0 ) THEN
501  info = n + 3
502  GO TO 80
503  END IF
504 *
505  IF( ilvl ) THEN
506  CALL claset( 'Full', n, n, czero, cone, vl, ldvl )
507  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
508  \$ vl( ilo+1, ilo ), ldvl )
509  CALL cungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
510  \$ work( itau ), work( iwork ), lwork+1-iwork,
511  \$ iinfo )
512  IF( iinfo.GE.0 )
513  \$ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
514  IF( iinfo.NE.0 ) THEN
515  info = n + 4
516  GO TO 80
517  END IF
518  END IF
519 *
520  IF( ilvr )
521  \$ CALL claset( 'Full', n, n, czero, cone, vr, ldvr )
522 *
523 * Reduce to generalized Hessenberg form
524 *
525  IF( ilv ) THEN
526 *
527 * Eigenvectors requested -- work on whole matrix.
528 *
529  CALL cgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
530  \$ ldvl, vr, ldvr, iinfo )
531  ELSE
532  CALL cgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
533  \$ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, iinfo )
534  END IF
535  IF( iinfo.NE.0 ) THEN
536  info = n + 5
537  GO TO 80
538  END IF
539 *
540 * Perform QZ algorithm
541 *
542  iwork = itau
543  IF( ilv ) THEN
544  chtemp = 'S'
545  ELSE
546  chtemp = 'E'
547  END IF
548  CALL chgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
549  \$ alpha, beta, vl, ldvl, vr, ldvr, work( iwork ),
550  \$ lwork+1-iwork, rwork( irwork ), iinfo )
551  IF( iinfo.GE.0 )
552  \$ lwkopt = max( lwkopt, int( work( iwork ) )+iwork-1 )
553  IF( iinfo.NE.0 ) THEN
554  IF( iinfo.GT.0 .AND. iinfo.LE.n ) THEN
555  info = iinfo
556  ELSE IF( iinfo.GT.n .AND. iinfo.LE.2*n ) THEN
557  info = iinfo - n
558  ELSE
559  info = n + 6
560  END IF
561  GO TO 80
562  END IF
563 *
564  IF( ilv ) THEN
565 *
566 * Compute Eigenvectors
567 *
568  IF( ilvl ) THEN
569  IF( ilvr ) THEN
570  chtemp = 'B'
571  ELSE
572  chtemp = 'L'
573  END IF
574  ELSE
575  chtemp = 'R'
576  END IF
577 *
578  CALL ctgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
579  \$ vr, ldvr, n, in, work( iwork ), rwork( irwork ),
580  \$ iinfo )
581  IF( iinfo.NE.0 ) THEN
582  info = n + 7
583  GO TO 80
584  END IF
585 *
586 * Undo balancing on VL and VR, rescale
587 *
588  IF( ilvl ) THEN
589  CALL cggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
590  \$ rwork( iright ), n, vl, ldvl, iinfo )
591  IF( iinfo.NE.0 ) THEN
592  info = n + 8
593  GO TO 80
594  END IF
595  DO 30 jc = 1, n
596  temp = zero
597  DO 10 jr = 1, n
598  temp = max( temp, abs1( vl( jr, jc ) ) )
599  10 CONTINUE
600  IF( temp.LT.safmin )
601  \$ GO TO 30
602  temp = one / temp
603  DO 20 jr = 1, n
604  vl( jr, jc ) = vl( jr, jc )*temp
605  20 CONTINUE
606  30 CONTINUE
607  END IF
608  IF( ilvr ) THEN
609  CALL cggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
610  \$ rwork( iright ), n, vr, ldvr, iinfo )
611  IF( iinfo.NE.0 ) THEN
612  info = n + 9
613  GO TO 80
614  END IF
615  DO 60 jc = 1, n
616  temp = zero
617  DO 40 jr = 1, n
618  temp = max( temp, abs1( vr( jr, jc ) ) )
619  40 CONTINUE
620  IF( temp.LT.safmin )
621  \$ GO TO 60
622  temp = one / temp
623  DO 50 jr = 1, n
624  vr( jr, jc ) = vr( jr, jc )*temp
625  50 CONTINUE
626  60 CONTINUE
627  END IF
628 *
629 * End of eigenvector calculation
630 *
631  END IF
632 *
633 * Undo scaling in alpha, beta
634 *
635 * Note: this does not give the alpha and beta for the unscaled
636 * problem.
637 *
638 * Un-scaling is limited to avoid underflow in alpha and beta
639 * if they are significant.
640 *
641  DO 70 jc = 1, n
642  absar = abs( REAL( ALPHA( JC ) ) )
643  absai = abs( aimag( alpha( jc ) ) )
644  absb = abs( REAL( BETA( JC ) ) )
645  salfar = anrm*REAL( ALPHA( JC ) )
646  salfai = anrm*aimag( alpha( jc ) )
647  sbeta = bnrm*REAL( BETA( JC ) )
648  ilimit = .false.
649  scale = one
650 *
651 * Check for significant underflow in imaginary part of ALPHA
652 *
653  IF( abs( salfai ).LT.safmin .AND. absai.GE.
654  \$ max( safmin, eps*absar, eps*absb ) ) THEN
655  ilimit = .true.
656  scale = ( safmin / anrm1 ) / max( safmin, anrm2*absai )
657  END IF
658 *
659 * Check for significant underflow in real part of ALPHA
660 *
661  IF( abs( salfar ).LT.safmin .AND. absar.GE.
662  \$ max( safmin, eps*absai, eps*absb ) ) THEN
663  ilimit = .true.
664  scale = max( scale, ( safmin / anrm1 ) /
665  \$ max( safmin, anrm2*absar ) )
666  END IF
667 *
668 * Check for significant underflow in BETA
669 *
670  IF( abs( sbeta ).LT.safmin .AND. absb.GE.
671  \$ max( safmin, eps*absar, eps*absai ) ) THEN
672  ilimit = .true.
673  scale = max( scale, ( safmin / bnrm1 ) /
674  \$ max( safmin, bnrm2*absb ) )
675  END IF
676 *
677 * Check for possible overflow when limiting scaling
678 *
679  IF( ilimit ) THEN
680  temp = ( scale*safmin )*max( abs( salfar ), abs( salfai ),
681  \$ abs( sbeta ) )
682  IF( temp.GT.one )
683  \$ scale = scale / temp
684  IF( scale.LT.one )
685  \$ ilimit = .false.
686  END IF
687 *
688 * Recompute un-scaled ALPHA, BETA if necessary.
689 *
690  IF( ilimit ) THEN
691  salfar = ( scale*REAL( ALPHA( JC ) ) )*anrm
692  salfai = ( scale*aimag( alpha( jc ) ) )*anrm
693  sbeta = ( scale*beta( jc ) )*bnrm
694  END IF
695  alpha( jc ) = cmplx( salfar, salfai )
696  beta( jc ) = sbeta
697  70 CONTINUE
698 *
699  80 CONTINUE
700  work( 1 ) = lwkopt
701 *
702  RETURN
703 *
704 * End of CGEGV
705 *
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:170
subroutine ctgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTGEVC
Definition: ctgevc.f:221
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:179
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:150
subroutine cgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
CGGHRD
Definition: cgghrd.f:206
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:130

Here is the call graph for this function:

 subroutine cgges ( character JOBVSL, character JOBVSR, character SORT, external SELCTG, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( ldb, * ) B, integer LDB, integer SDIM, complex, dimension( * ) ALPHA, complex, dimension( * ) BETA, complex, dimension( ldvsl, * ) VSL, integer LDVSL, complex, dimension( ldvsr, * ) VSR, integer LDVSR, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, logical, dimension( * ) BWORK, integer INFO )

CGGES computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices

Purpose:
``` CGGES computes for a pair of N-by-N complex nonsymmetric matrices
(A,B), the generalized eigenvalues, the generalized complex Schur
form (S, T), and optionally left and/or right Schur vectors (VSL
and VSR). This gives the generalized Schur factorization

(A,B) = ( (VSL)*S*(VSR)**H, (VSL)*T*(VSR)**H )

where (VSR)**H is the conjugate-transpose of VSR.

Optionally, it also orders the eigenvalues so that a selected cluster
of eigenvalues appears in the leading diagonal blocks of the upper
triangular matrix S and the upper triangular matrix T. The leading
columns of VSL and VSR then form an unitary basis for the
corresponding left and right eigenspaces (deflating subspaces).

(If only the generalized eigenvalues are needed, use the driver

A generalized eigenvalue for a pair of matrices (A,B) is a scalar w
or a ratio alpha/beta = w, such that  A - w*B is singular.  It is
usually represented as the pair (alpha,beta), as there is a
reasonable interpretation for beta=0, and even for both being zero.

A pair of matrices (S,T) is in generalized complex Schur form if S
and T are upper triangular and, in addition, the diagonal elements
of T are non-negative real numbers.```
Parameters
 [in] JOBVSL ``` JOBVSL is CHARACTER*1 = 'N': do not compute the left Schur vectors; = 'V': compute the left Schur vectors.``` [in] JOBVSR ``` JOBVSR is CHARACTER*1 = 'N': do not compute the right Schur vectors; = 'V': compute the right Schur vectors.``` [in] SORT ``` SORT is CHARACTER*1 Specifies whether or not to order the eigenvalues on the diagonal of the generalized Schur form. = 'N': Eigenvalues are not ordered; = 'S': Eigenvalues are ordered (see SELCTG).``` [in] SELCTG ``` SELCTG is a LOGICAL FUNCTION of two COMPLEX arguments SELCTG must be declared EXTERNAL in the calling subroutine. If SORT = 'N', SELCTG is not referenced. If SORT = 'S', SELCTG is used to select eigenvalues to sort to the top left of the Schur form. An eigenvalue ALPHA(j)/BETA(j) is selected if SELCTG(ALPHA(j),BETA(j)) is true. Note that a selected complex eigenvalue may no longer satisfy SELCTG(ALPHA(j),BETA(j)) = .TRUE. after ordering, since ordering may change the value of complex eigenvalues (especially if the eigenvalue is ill-conditioned), in this case INFO is set to N+2 (See INFO below).``` [in] N ``` N is INTEGER The order of the matrices A, B, VSL, and VSR. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the first of the pair of matrices. On exit, A has been overwritten by its generalized Schur form S.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. LDA >= max(1,N).``` [in,out] B ``` B is COMPLEX array, dimension (LDB, N) On entry, the second of the pair of matrices. On exit, B has been overwritten by its generalized Schur form T.``` [in] LDB ``` LDB is INTEGER The leading dimension of B. LDB >= max(1,N).``` [out] SDIM ``` SDIM is INTEGER If SORT = 'N', SDIM = 0. If SORT = 'S', SDIM = number of eigenvalues (after sorting) for which SELCTG is true.``` [out] ALPHA ` ALPHA is COMPLEX array, dimension (N)` [out] BETA ``` BETA is COMPLEX array, dimension (N) On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the generalized eigenvalues. ALPHA(j), j=1,...,N and BETA(j), j=1,...,N are the diagonals of the complex Schur form (A,B) output by CGGES. The BETA(j) will be non-negative real. Note: the quotients ALPHA(j)/BETA(j) may easily over- or underflow, and BETA(j) may even be zero. Thus, the user should avoid naively computing the ratio alpha/beta. However, ALPHA will be always less than and usually comparable with norm(A) in magnitude, and BETA always less than and usually comparable with norm(B).``` [out] VSL ``` VSL is COMPLEX array, dimension (LDVSL,N) If JOBVSL = 'V', VSL will contain the left Schur vectors. Not referenced if JOBVSL = 'N'.``` [in] LDVSL ``` LDVSL is INTEGER The leading dimension of the matrix VSL. LDVSL >= 1, and if JOBVSL = 'V', LDVSL >= N.``` [out] VSR ``` VSR is COMPLEX array, dimension (LDVSR,N) If JOBVSR = 'V', VSR will contain the right Schur vectors. Not referenced if JOBVSR = 'N'.``` [in] LDVSR ``` LDVSR is INTEGER The leading dimension of the matrix VSR. LDVSR >= 1, and if JOBVSR = 'V', LDVSR >= N.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,2*N). For good performance, LWORK must generally be larger. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (8*N)` [out] BWORK ``` BWORK is LOGICAL array, dimension (N) Not referenced if SORT = 'N'.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. =1,...,N: The QZ iteration failed. (A,B) are not in Schur form, but ALPHA(j) and BETA(j) should be correct for j=INFO+1,...,N. > N: =N+1: other than QZ iteration failed in CHGEQZ =N+2: after reordering, roundoff changed values of some complex eigenvalues so that leading eigenvalues in the Generalized Schur form no longer satisfy SELCTG=.TRUE. This could also be caused due to scaling. =N+3: reordering failed in CTGSEN.```
Date
November 2015

Definition at line 272 of file cgges.f.

272 *
273 * -- LAPACK driver routine (version 3.6.0) --
274 * -- LAPACK is a software package provided by Univ. of Tennessee, --
275 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
276 * November 2015
277 *
278 * .. Scalar Arguments ..
279  CHARACTER jobvsl, jobvsr, sort
280  INTEGER info, lda, ldb, ldvsl, ldvsr, lwork, n, sdim
281 * ..
282 * .. Array Arguments ..
283  LOGICAL bwork( * )
284  REAL rwork( * )
285  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
286  \$ beta( * ), vsl( ldvsl, * ), vsr( ldvsr, * ),
287  \$ work( * )
288 * ..
289 * .. Function Arguments ..
290  LOGICAL selctg
291  EXTERNAL selctg
292 * ..
293 *
294 * =====================================================================
295 *
296 * .. Parameters ..
297  REAL zero, one
298  parameter( zero = 0.0e0, one = 1.0e0 )
299  COMPLEX czero, cone
300  parameter( czero = ( 0.0e0, 0.0e0 ),
301  \$ cone = ( 1.0e0, 0.0e0 ) )
302 * ..
303 * .. Local Scalars ..
304  LOGICAL cursl, ilascl, ilbscl, ilvsl, ilvsr, lastsl,
305  \$ lquery, wantst
306  INTEGER i, icols, ierr, ihi, ijobvl, ijobvr, ileft,
307  \$ ilo, iright, irows, irwrk, itau, iwrk, lwkmin,
308  \$ lwkopt
309  REAL anrm, anrmto, bignum, bnrm, bnrmto, eps, pvsl,
310  \$ pvsr, smlnum
311 * ..
312 * .. Local Arrays ..
313  INTEGER idum( 1 )
314  REAL dif( 2 )
315 * ..
316 * .. External Subroutines ..
317  EXTERNAL cgeqrf, cggbak, cggbal, cgghrd, chgeqz, clacpy,
319  \$ xerbla
320 * ..
321 * .. External Functions ..
322  LOGICAL lsame
323  INTEGER ilaenv
324  REAL clange, slamch
325  EXTERNAL lsame, ilaenv, clange, slamch
326 * ..
327 * .. Intrinsic Functions ..
328  INTRINSIC max, sqrt
329 * ..
330 * .. Executable Statements ..
331 *
332 * Decode the input arguments
333 *
334  IF( lsame( jobvsl, 'N' ) ) THEN
335  ijobvl = 1
336  ilvsl = .false.
337  ELSE IF( lsame( jobvsl, 'V' ) ) THEN
338  ijobvl = 2
339  ilvsl = .true.
340  ELSE
341  ijobvl = -1
342  ilvsl = .false.
343  END IF
344 *
345  IF( lsame( jobvsr, 'N' ) ) THEN
346  ijobvr = 1
347  ilvsr = .false.
348  ELSE IF( lsame( jobvsr, 'V' ) ) THEN
349  ijobvr = 2
350  ilvsr = .true.
351  ELSE
352  ijobvr = -1
353  ilvsr = .false.
354  END IF
355 *
356  wantst = lsame( sort, 'S' )
357 *
358 * Test the input arguments
359 *
360  info = 0
361  lquery = ( lwork.EQ.-1 )
362  IF( ijobvl.LE.0 ) THEN
363  info = -1
364  ELSE IF( ijobvr.LE.0 ) THEN
365  info = -2
366  ELSE IF( ( .NOT.wantst ) .AND. ( .NOT.lsame( sort, 'N' ) ) ) THEN
367  info = -3
368  ELSE IF( n.LT.0 ) THEN
369  info = -5
370  ELSE IF( lda.LT.max( 1, n ) ) THEN
371  info = -7
372  ELSE IF( ldb.LT.max( 1, n ) ) THEN
373  info = -9
374  ELSE IF( ldvsl.LT.1 .OR. ( ilvsl .AND. ldvsl.LT.n ) ) THEN
375  info = -14
376  ELSE IF( ldvsr.LT.1 .OR. ( ilvsr .AND. ldvsr.LT.n ) ) THEN
377  info = -16
378  END IF
379 *
380 * Compute workspace
381 * (Note: Comments in the code beginning "Workspace:" describe the
382 * minimal amount of workspace needed at that point in the code,
383 * as well as the preferred amount for good performance.
384 * NB refers to the optimal block size for the immediately
385 * following subroutine, as returned by ILAENV.)
386 *
387  IF( info.EQ.0 ) THEN
388  lwkmin = max( 1, 2*n )
389  lwkopt = max( 1, n + n*ilaenv( 1, 'CGEQRF', ' ', n, 1, n, 0 ) )
390  lwkopt = max( lwkopt, n +
391  \$ n*ilaenv( 1, 'CUNMQR', ' ', n, 1, n, -1 ) )
392  IF( ilvsl ) THEN
393  lwkopt = max( lwkopt, n +
394  \$ n*ilaenv( 1, 'CUNGQR', ' ', n, 1, n, -1 ) )
395  END IF
396  work( 1 ) = lwkopt
397 *
398  IF( lwork.LT.lwkmin .AND. .NOT.lquery )
399  \$ info = -18
400  END IF
401 *
402  IF( info.NE.0 ) THEN
403  CALL xerbla( 'CGGES ', -info )
404  RETURN
405  ELSE IF( lquery ) THEN
406  RETURN
407  END IF
408 *
409 * Quick return if possible
410 *
411  IF( n.EQ.0 ) THEN
412  sdim = 0
413  RETURN
414  END IF
415 *
416 * Get machine constants
417 *
418  eps = slamch( 'P' )
419  smlnum = slamch( 'S' )
420  bignum = one / smlnum
421  CALL slabad( smlnum, bignum )
422  smlnum = sqrt( smlnum ) / eps
423  bignum = one / smlnum
424 *
425 * Scale A if max element outside range [SMLNUM,BIGNUM]
426 *
427  anrm = clange( 'M', n, n, a, lda, rwork )
428  ilascl = .false.
429  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
430  anrmto = smlnum
431  ilascl = .true.
432  ELSE IF( anrm.GT.bignum ) THEN
433  anrmto = bignum
434  ilascl = .true.
435  END IF
436 *
437  IF( ilascl )
438  \$ CALL clascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
439 *
440 * Scale B if max element outside range [SMLNUM,BIGNUM]
441 *
442  bnrm = clange( 'M', n, n, b, ldb, rwork )
443  ilbscl = .false.
444  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
445  bnrmto = smlnum
446  ilbscl = .true.
447  ELSE IF( bnrm.GT.bignum ) THEN
448  bnrmto = bignum
449  ilbscl = .true.
450  END IF
451 *
452  IF( ilbscl )
453  \$ CALL clascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
454 *
455 * Permute the matrix to make it more nearly triangular
456 * (Real Workspace: need 6*N)
457 *
458  ileft = 1
459  iright = n + 1
460  irwrk = iright + n
461  CALL cggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
462  \$ rwork( iright ), rwork( irwrk ), ierr )
463 *
464 * Reduce B to triangular form (QR decomposition of B)
465 * (Complex Workspace: need N, prefer N*NB)
466 *
467  irows = ihi + 1 - ilo
468  icols = n + 1 - ilo
469  itau = 1
470  iwrk = itau + irows
471  CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
472  \$ work( iwrk ), lwork+1-iwrk, ierr )
473 *
474 * Apply the orthogonal transformation to matrix A
475 * (Complex Workspace: need N, prefer N*NB)
476 *
477  CALL cunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
478  \$ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
479  \$ lwork+1-iwrk, ierr )
480 *
481 * Initialize VSL
482 * (Complex Workspace: need N, prefer N*NB)
483 *
484  IF( ilvsl ) THEN
485  CALL claset( 'Full', n, n, czero, cone, vsl, ldvsl )
486  IF( irows.GT.1 ) THEN
487  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
488  \$ vsl( ilo+1, ilo ), ldvsl )
489  END IF
490  CALL cungqr( irows, irows, irows, vsl( ilo, ilo ), ldvsl,
491  \$ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
492  END IF
493 *
494 * Initialize VSR
495 *
496  IF( ilvsr )
497  \$ CALL claset( 'Full', n, n, czero, cone, vsr, ldvsr )
498 *
499 * Reduce to generalized Hessenberg form
500 * (Workspace: none needed)
501 *
502  CALL cgghrd( jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb, vsl,
503  \$ ldvsl, vsr, ldvsr, ierr )
504 *
505  sdim = 0
506 *
507 * Perform QZ algorithm, computing Schur vectors if desired
508 * (Complex Workspace: need N)
509 * (Real Workspace: need N)
510 *
511  iwrk = itau
512  CALL chgeqz( 'S', jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb,
513  \$ alpha, beta, vsl, ldvsl, vsr, ldvsr, work( iwrk ),
514  \$ lwork+1-iwrk, rwork( irwrk ), ierr )
515  IF( ierr.NE.0 ) THEN
516  IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
517  info = ierr
518  ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
519  info = ierr - n
520  ELSE
521  info = n + 1
522  END IF
523  GO TO 30
524  END IF
525 *
526 * Sort eigenvalues ALPHA/BETA if desired
527 * (Workspace: none needed)
528 *
529  IF( wantst ) THEN
530 *
531 * Undo scaling on eigenvalues before selecting
532 *
533  IF( ilascl )
534  \$ CALL clascl( 'G', 0, 0, anrm, anrmto, n, 1, alpha, n, ierr )
535  IF( ilbscl )
536  \$ CALL clascl( 'G', 0, 0, bnrm, bnrmto, n, 1, beta, n, ierr )
537 *
538 * Select eigenvalues
539 *
540  DO 10 i = 1, n
541  bwork( i ) = selctg( alpha( i ), beta( i ) )
542  10 CONTINUE
543 *
544  CALL ctgsen( 0, ilvsl, ilvsr, bwork, n, a, lda, b, ldb, alpha,
545  \$ beta, vsl, ldvsl, vsr, ldvsr, sdim, pvsl, pvsr,
546  \$ dif, work( iwrk ), lwork-iwrk+1, idum, 1, ierr )
547  IF( ierr.EQ.1 )
548  \$ info = n + 3
549 *
550  END IF
551 *
552 * Apply back-permutation to VSL and VSR
553 * (Workspace: none needed)
554 *
555  IF( ilvsl )
556  \$ CALL cggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
557  \$ rwork( iright ), n, vsl, ldvsl, ierr )
558  IF( ilvsr )
559  \$ CALL cggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
560  \$ rwork( iright ), n, vsr, ldvsr, ierr )
561 *
562 * Undo scaling
563 *
564  IF( ilascl ) THEN
565  CALL clascl( 'U', 0, 0, anrmto, anrm, n, n, a, lda, ierr )
566  CALL clascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
567  END IF
568 *
569  IF( ilbscl ) THEN
570  CALL clascl( 'U', 0, 0, bnrmto, bnrm, n, n, b, ldb, ierr )
571  CALL clascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
572  END IF
573 *
574  IF( wantst ) THEN
575 *
576 * Check if reordering is correct
577 *
578  lastsl = .true.
579  sdim = 0
580  DO 20 i = 1, n
581  cursl = selctg( alpha( i ), beta( i ) )
582  IF( cursl )
583  \$ sdim = sdim + 1
584  IF( cursl .AND. .NOT.lastsl )
585  \$ info = n + 2
586  lastsl = cursl
587  20 CONTINUE
588 *
589  END IF
590 *
591  30 CONTINUE
592 *
593  work( 1 ) = lwkopt
594 *
595  RETURN
596 *
597 * End of CGGES
598 *
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:170
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:179
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:150
subroutine cgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
CGGHRD
Definition: cgghrd.f:206
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine ctgsen(IJOB, WANTQ, WANTZ, SELECT, N, A, LDA, B, LDB, ALPHA, BETA, Q, LDQ, Z, LDZ, M, PL, PR, DIF, WORK, LWORK, IWORK, LIWORK, INFO)
CTGSEN
Definition: ctgsen.f:435
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:130

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cgges3 ( character JOBVSL, character JOBVSR, character SORT, external SELCTG, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( ldb, * ) B, integer LDB, integer SDIM, complex, dimension( * ) ALPHA, complex, dimension( * ) BETA, complex, dimension( ldvsl, * ) VSL, integer LDVSL, complex, dimension( ldvsr, * ) VSR, integer LDVSR, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, logical, dimension( * ) BWORK, integer INFO )

CGGES3 computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices (blocked algorithm)

Purpose:
``` CGGES3 computes for a pair of N-by-N complex nonsymmetric matrices
(A,B), the generalized eigenvalues, the generalized complex Schur
form (S, T), and optionally left and/or right Schur vectors (VSL
and VSR). This gives the generalized Schur factorization

(A,B) = ( (VSL)*S*(VSR)**H, (VSL)*T*(VSR)**H )

where (VSR)**H is the conjugate-transpose of VSR.

Optionally, it also orders the eigenvalues so that a selected cluster
of eigenvalues appears in the leading diagonal blocks of the upper
triangular matrix S and the upper triangular matrix T. The leading
columns of VSL and VSR then form an unitary basis for the
corresponding left and right eigenspaces (deflating subspaces).

(If only the generalized eigenvalues are needed, use the driver

A generalized eigenvalue for a pair of matrices (A,B) is a scalar w
or a ratio alpha/beta = w, such that  A - w*B is singular.  It is
usually represented as the pair (alpha,beta), as there is a
reasonable interpretation for beta=0, and even for both being zero.

A pair of matrices (S,T) is in generalized complex Schur form if S
and T are upper triangular and, in addition, the diagonal elements
of T are non-negative real numbers.```
Parameters
 [in] JOBVSL ``` JOBVSL is CHARACTER*1 = 'N': do not compute the left Schur vectors; = 'V': compute the left Schur vectors.``` [in] JOBVSR ``` JOBVSR is CHARACTER*1 = 'N': do not compute the right Schur vectors; = 'V': compute the right Schur vectors.``` [in] SORT ``` SORT is CHARACTER*1 Specifies whether or not to order the eigenvalues on the diagonal of the generalized Schur form. = 'N': Eigenvalues are not ordered; = 'S': Eigenvalues are ordered (see SELCTG).``` [in] SELCTG ``` SELCTG is a LOGICAL FUNCTION of two COMPLEX arguments SELCTG must be declared EXTERNAL in the calling subroutine. If SORT = 'N', SELCTG is not referenced. If SORT = 'S', SELCTG is used to select eigenvalues to sort to the top left of the Schur form. An eigenvalue ALPHA(j)/BETA(j) is selected if SELCTG(ALPHA(j),BETA(j)) is true. Note that a selected complex eigenvalue may no longer satisfy SELCTG(ALPHA(j),BETA(j)) = .TRUE. after ordering, since ordering may change the value of complex eigenvalues (especially if the eigenvalue is ill-conditioned), in this case INFO is set to N+2 (See INFO below).``` [in] N ``` N is INTEGER The order of the matrices A, B, VSL, and VSR. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the first of the pair of matrices. On exit, A has been overwritten by its generalized Schur form S.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. LDA >= max(1,N).``` [in,out] B ``` B is COMPLEX array, dimension (LDB, N) On entry, the second of the pair of matrices. On exit, B has been overwritten by its generalized Schur form T.``` [in] LDB ``` LDB is INTEGER The leading dimension of B. LDB >= max(1,N).``` [out] SDIM ``` SDIM is INTEGER If SORT = 'N', SDIM = 0. If SORT = 'S', SDIM = number of eigenvalues (after sorting) for which SELCTG is true.``` [out] ALPHA ` ALPHA is COMPLEX array, dimension (N)` [out] BETA ``` BETA is COMPLEX array, dimension (N) On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the generalized eigenvalues. ALPHA(j), j=1,...,N and BETA(j), j=1,...,N are the diagonals of the complex Schur form (A,B) output by CGGES3. The BETA(j) will be non-negative real. Note: the quotients ALPHA(j)/BETA(j) may easily over- or underflow, and BETA(j) may even be zero. Thus, the user should avoid naively computing the ratio alpha/beta. However, ALPHA will be always less than and usually comparable with norm(A) in magnitude, and BETA always less than and usually comparable with norm(B).``` [out] VSL ``` VSL is COMPLEX array, dimension (LDVSL,N) If JOBVSL = 'V', VSL will contain the left Schur vectors. Not referenced if JOBVSL = 'N'.``` [in] LDVSL ``` LDVSL is INTEGER The leading dimension of the matrix VSL. LDVSL >= 1, and if JOBVSL = 'V', LDVSL >= N.``` [out] VSR ``` VSR is COMPLEX array, dimension (LDVSR,N) If JOBVSR = 'V', VSR will contain the right Schur vectors. Not referenced if JOBVSR = 'N'.``` [in] LDVSR ``` LDVSR is INTEGER The leading dimension of the matrix VSR. LDVSR >= 1, and if JOBVSR = 'V', LDVSR >= N.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (8*N)` [out] BWORK ``` BWORK is LOGICAL array, dimension (N) Not referenced if SORT = 'N'.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. =1,...,N: The QZ iteration failed. (A,B) are not in Schur form, but ALPHA(j) and BETA(j) should be correct for j=INFO+1,...,N. > N: =N+1: other than QZ iteration failed in CHGEQZ =N+2: after reordering, roundoff changed values of some complex eigenvalues so that leading eigenvalues in the Generalized Schur form no longer satisfy SELCTG=.TRUE. This could also be caused due to scaling. =N+3: reordering failed in CTGSEN.```
Date
January 2015

Definition at line 271 of file cgges3.f.

271 *
272 * -- LAPACK driver routine (version 3.6.0) --
273 * -- LAPACK is a software package provided by Univ. of Tennessee, --
274 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
275 * January 2015
276 *
277 * .. Scalar Arguments ..
278  CHARACTER jobvsl, jobvsr, sort
279  INTEGER info, lda, ldb, ldvsl, ldvsr, lwork, n, sdim
280 * ..
281 * .. Array Arguments ..
282  LOGICAL bwork( * )
283  REAL rwork( * )
284  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
285  \$ beta( * ), vsl( ldvsl, * ), vsr( ldvsr, * ),
286  \$ work( * )
287 * ..
288 * .. Function Arguments ..
289  LOGICAL selctg
290  EXTERNAL selctg
291 * ..
292 *
293 * =====================================================================
294 *
295 * .. Parameters ..
296  REAL zero, one
297  parameter( zero = 0.0e0, one = 1.0e0 )
298  COMPLEX czero, cone
299  parameter( czero = ( 0.0e0, 0.0e0 ),
300  \$ cone = ( 1.0e0, 0.0e0 ) )
301 * ..
302 * .. Local Scalars ..
303  LOGICAL cursl, ilascl, ilbscl, ilvsl, ilvsr, lastsl,
304  \$ lquery, wantst
305  INTEGER i, icols, ierr, ihi, ijobvl, ijobvr, ileft,
306  \$ ilo, iright, irows, irwrk, itau, iwrk, lwkopt
307  REAL anrm, anrmto, bignum, bnrm, bnrmto, eps, pvsl,
308  \$ pvsr, smlnum
309 * ..
310 * .. Local Arrays ..
311  INTEGER idum( 1 )
312  REAL dif( 2 )
313 * ..
314 * .. External Subroutines ..
315  EXTERNAL cgeqrf, cggbak, cggbal, cgghd3, chgeqz, clacpy,
317  \$ xerbla
318 * ..
319 * .. External Functions ..
320  LOGICAL lsame
321  REAL clange, slamch
322  EXTERNAL lsame, clange, slamch
323 * ..
324 * .. Intrinsic Functions ..
325  INTRINSIC max, sqrt
326 * ..
327 * .. Executable Statements ..
328 *
329 * Decode the input arguments
330 *
331  IF( lsame( jobvsl, 'N' ) ) THEN
332  ijobvl = 1
333  ilvsl = .false.
334  ELSE IF( lsame( jobvsl, 'V' ) ) THEN
335  ijobvl = 2
336  ilvsl = .true.
337  ELSE
338  ijobvl = -1
339  ilvsl = .false.
340  END IF
341 *
342  IF( lsame( jobvsr, 'N' ) ) THEN
343  ijobvr = 1
344  ilvsr = .false.
345  ELSE IF( lsame( jobvsr, 'V' ) ) THEN
346  ijobvr = 2
347  ilvsr = .true.
348  ELSE
349  ijobvr = -1
350  ilvsr = .false.
351  END IF
352 *
353  wantst = lsame( sort, 'S' )
354 *
355 * Test the input arguments
356 *
357  info = 0
358  lquery = ( lwork.EQ.-1 )
359  IF( ijobvl.LE.0 ) THEN
360  info = -1
361  ELSE IF( ijobvr.LE.0 ) THEN
362  info = -2
363  ELSE IF( ( .NOT.wantst ) .AND. ( .NOT.lsame( sort, 'N' ) ) ) THEN
364  info = -3
365  ELSE IF( n.LT.0 ) THEN
366  info = -5
367  ELSE IF( lda.LT.max( 1, n ) ) THEN
368  info = -7
369  ELSE IF( ldb.LT.max( 1, n ) ) THEN
370  info = -9
371  ELSE IF( ldvsl.LT.1 .OR. ( ilvsl .AND. ldvsl.LT.n ) ) THEN
372  info = -14
373  ELSE IF( ldvsr.LT.1 .OR. ( ilvsr .AND. ldvsr.LT.n ) ) THEN
374  info = -16
375  ELSE IF( lwork.LT.max( 1, 2*n ) .AND. .NOT.lquery ) THEN
376  info = -18
377  END IF
378 *
379 * Compute workspace
380 *
381  IF( info.EQ.0 ) THEN
382  CALL cgeqrf( n, n, b, ldb, work, work, -1, ierr )
383  lwkopt = max( 1, n + int( work( 1 ) ) )
384  CALL cunmqr( 'L', 'C', n, n, n, b, ldb, work, a, lda, work,
385  \$ -1, ierr )
386  lwkopt = max( lwkopt, n + int( work( 1 ) ) )
387  IF( ilvsl ) THEN
388  CALL cungqr( n, n, n, vsl, ldvsl, work, work, -1,
389  \$ ierr )
390  lwkopt = max( lwkopt, n + int( work( 1 ) ) )
391  END IF
392  CALL cgghd3( jobvsl, jobvsr, n, 1, n, a, lda, b, ldb, vsl,
393  \$ ldvsl, vsr, ldvsr, work, -1, ierr )
394  lwkopt = max( lwkopt, n + int( work( 1 ) ) )
395  CALL chgeqz( 'S', jobvsl, jobvsr, n, 1, n, a, lda, b, ldb,
396  \$ alpha, beta, vsl, ldvsl, vsr, ldvsr, work, -1,
397  \$ work, ierr )
398  lwkopt = max( lwkopt, int( work( 1 ) ) )
399  IF( wantst ) THEN
400  CALL ctgsen( 0, ilvsl, ilvsr, bwork, n, a, lda, b, ldb,
401  \$ alpha, beta, vsl, ldvsl, vsr, ldvsr, sdim,
402  \$ pvsl, pvsr, dif, work, -1, idum, 1, ierr )
403  lwkopt = max( lwkopt, int( work( 1 ) ) )
404  END IF
405  work( 1 ) = cmplx( lwkopt )
406  END IF
407
408 *
409  IF( info.NE.0 ) THEN
410  CALL xerbla( 'CGGES3 ', -info )
411  RETURN
412  ELSE IF( lquery ) THEN
413  RETURN
414  END IF
415 *
416 * Quick return if possible
417 *
418  IF( n.EQ.0 ) THEN
419  sdim = 0
420  RETURN
421  END IF
422 *
423 * Get machine constants
424 *
425  eps = slamch( 'P' )
426  smlnum = slamch( 'S' )
427  bignum = one / smlnum
428  CALL slabad( smlnum, bignum )
429  smlnum = sqrt( smlnum ) / eps
430  bignum = one / smlnum
431 *
432 * Scale A if max element outside range [SMLNUM,BIGNUM]
433 *
434  anrm = clange( 'M', n, n, a, lda, rwork )
435  ilascl = .false.
436  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
437  anrmto = smlnum
438  ilascl = .true.
439  ELSE IF( anrm.GT.bignum ) THEN
440  anrmto = bignum
441  ilascl = .true.
442  END IF
443 *
444  IF( ilascl )
445  \$ CALL clascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
446 *
447 * Scale B if max element outside range [SMLNUM,BIGNUM]
448 *
449  bnrm = clange( 'M', n, n, b, ldb, rwork )
450  ilbscl = .false.
451  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
452  bnrmto = smlnum
453  ilbscl = .true.
454  ELSE IF( bnrm.GT.bignum ) THEN
455  bnrmto = bignum
456  ilbscl = .true.
457  END IF
458 *
459  IF( ilbscl )
460  \$ CALL clascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
461 *
462 * Permute the matrix to make it more nearly triangular
463 *
464  ileft = 1
465  iright = n + 1
466  irwrk = iright + n
467  CALL cggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
468  \$ rwork( iright ), rwork( irwrk ), ierr )
469 *
470 * Reduce B to triangular form (QR decomposition of B)
471 *
472  irows = ihi + 1 - ilo
473  icols = n + 1 - ilo
474  itau = 1
475  iwrk = itau + irows
476  CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
477  \$ work( iwrk ), lwork+1-iwrk, ierr )
478 *
479 * Apply the orthogonal transformation to matrix A
480 *
481  CALL cunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
482  \$ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
483  \$ lwork+1-iwrk, ierr )
484 *
485 * Initialize VSL
486 *
487  IF( ilvsl ) THEN
488  CALL claset( 'Full', n, n, czero, cone, vsl, ldvsl )
489  IF( irows.GT.1 ) THEN
490  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
491  \$ vsl( ilo+1, ilo ), ldvsl )
492  END IF
493  CALL cungqr( irows, irows, irows, vsl( ilo, ilo ), ldvsl,
494  \$ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
495  END IF
496 *
497 * Initialize VSR
498 *
499  IF( ilvsr )
500  \$ CALL claset( 'Full', n, n, czero, cone, vsr, ldvsr )
501 *
502 * Reduce to generalized Hessenberg form
503 *
504  CALL cgghd3( jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb, vsl,
505  \$ ldvsl, vsr, ldvsr, work( iwrk ), lwork+1-iwrk, ierr )
506 *
507  sdim = 0
508 *
509 * Perform QZ algorithm, computing Schur vectors if desired
510 *
511  iwrk = itau
512  CALL chgeqz( 'S', jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb,
513  \$ alpha, beta, vsl, ldvsl, vsr, ldvsr, work( iwrk ),
514  \$ lwork+1-iwrk, rwork( irwrk ), ierr )
515  IF( ierr.NE.0 ) THEN
516  IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
517  info = ierr
518  ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
519  info = ierr - n
520  ELSE
521  info = n + 1
522  END IF
523  GO TO 30
524  END IF
525 *
526 * Sort eigenvalues ALPHA/BETA if desired
527 *
528  IF( wantst ) THEN
529 *
530 * Undo scaling on eigenvalues before selecting
531 *
532  IF( ilascl )
533  \$ CALL clascl( 'G', 0, 0, anrm, anrmto, n, 1, alpha, n, ierr )
534  IF( ilbscl )
535  \$ CALL clascl( 'G', 0, 0, bnrm, bnrmto, n, 1, beta, n, ierr )
536 *
537 * Select eigenvalues
538 *
539  DO 10 i = 1, n
540  bwork( i ) = selctg( alpha( i ), beta( i ) )
541  10 CONTINUE
542 *
543  CALL ctgsen( 0, ilvsl, ilvsr, bwork, n, a, lda, b, ldb, alpha,
544  \$ beta, vsl, ldvsl, vsr, ldvsr, sdim, pvsl, pvsr,
545  \$ dif, work( iwrk ), lwork-iwrk+1, idum, 1, ierr )
546  IF( ierr.EQ.1 )
547  \$ info = n + 3
548 *
549  END IF
550 *
551 * Apply back-permutation to VSL and VSR
552 *
553  IF( ilvsl )
554  \$ CALL cggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
555  \$ rwork( iright ), n, vsl, ldvsl, ierr )
556  IF( ilvsr )
557  \$ CALL cggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
558  \$ rwork( iright ), n, vsr, ldvsr, ierr )
559 *
560 * Undo scaling
561 *
562  IF( ilascl ) THEN
563  CALL clascl( 'U', 0, 0, anrmto, anrm, n, n, a, lda, ierr )
564  CALL clascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
565  END IF
566 *
567  IF( ilbscl ) THEN
568  CALL clascl( 'U', 0, 0, bnrmto, bnrm, n, n, b, ldb, ierr )
569  CALL clascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
570  END IF
571 *
572  IF( wantst ) THEN
573 *
574 * Check if reordering is correct
575 *
576  lastsl = .true.
577  sdim = 0
578  DO 20 i = 1, n
579  cursl = selctg( alpha( i ), beta( i ) )
580  IF( cursl )
581  \$ sdim = sdim + 1
582  IF( cursl .AND. .NOT.lastsl )
583  \$ info = n + 2
584  lastsl = cursl
585  20 CONTINUE
586 *
587  END IF
588 *
589  30 CONTINUE
590 *
591  work( 1 ) = cmplx( lwkopt )
592 *
593  RETURN
594 *
595 * End of CGGES3
596 *
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:170
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cgghd3(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
CGGHD3
Definition: cgghd3.f:233
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:179
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:150
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine ctgsen(IJOB, WANTQ, WANTZ, SELECT, N, A, LDA, B, LDB, ALPHA, BETA, Q, LDQ, Z, LDZ, M, PL, PR, DIF, WORK, LWORK, IWORK, LIWORK, INFO)
CTGSEN
Definition: ctgsen.f:435
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:130

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cggesx ( character JOBVSL, character JOBVSR, character SORT, external SELCTG, character SENSE, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( ldb, * ) B, integer LDB, integer SDIM, complex, dimension( * ) ALPHA, complex, dimension( * ) BETA, complex, dimension( ldvsl, * ) VSL, integer LDVSL, complex, dimension( ldvsr, * ) VSR, integer LDVSR, real, dimension( 2 ) RCONDE, real, dimension( 2 ) RCONDV, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, integer, dimension( * ) IWORK, integer LIWORK, logical, dimension( * ) BWORK, integer INFO )

CGGESX computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices

Purpose:
``` CGGESX computes for a pair of N-by-N complex nonsymmetric matrices
(A,B), the generalized eigenvalues, the complex Schur form (S,T),
and, optionally, the left and/or right matrices of Schur vectors (VSL
and VSR).  This gives the generalized Schur factorization

(A,B) = ( (VSL) S (VSR)**H, (VSL) T (VSR)**H )

where (VSR)**H is the conjugate-transpose of VSR.

Optionally, it also orders the eigenvalues so that a selected cluster
of eigenvalues appears in the leading diagonal blocks of the upper
triangular matrix S and the upper triangular matrix T; computes
a reciprocal condition number for the average of the selected
eigenvalues (RCONDE); and computes a reciprocal condition number for
the right and left deflating subspaces corresponding to the selected
eigenvalues (RCONDV). The leading columns of VSL and VSR then form
an orthonormal basis for the corresponding left and right eigenspaces
(deflating subspaces).

A generalized eigenvalue for a pair of matrices (A,B) is a scalar w
or a ratio alpha/beta = w, such that  A - w*B is singular.  It is
usually represented as the pair (alpha,beta), as there is a
reasonable interpretation for beta=0 or for both being zero.

A pair of matrices (S,T) is in generalized complex Schur form if T is
upper triangular with non-negative diagonal and S is upper
triangular.```
Parameters
 [in] JOBVSL ``` JOBVSL is CHARACTER*1 = 'N': do not compute the left Schur vectors; = 'V': compute the left Schur vectors.``` [in] JOBVSR ``` JOBVSR is CHARACTER*1 = 'N': do not compute the right Schur vectors; = 'V': compute the right Schur vectors.``` [in] SORT ``` SORT is CHARACTER*1 Specifies whether or not to order the eigenvalues on the diagonal of the generalized Schur form. = 'N': Eigenvalues are not ordered; = 'S': Eigenvalues are ordered (see SELCTG).``` [in] SELCTG ``` SELCTG is procedure) LOGICAL FUNCTION of two COMPLEX arguments SELCTG must be declared EXTERNAL in the calling subroutine. If SORT = 'N', SELCTG is not referenced. If SORT = 'S', SELCTG is used to select eigenvalues to sort to the top left of the Schur form. Note that a selected complex eigenvalue may no longer satisfy SELCTG(ALPHA(j),BETA(j)) = .TRUE. after ordering, since ordering may change the value of complex eigenvalues (especially if the eigenvalue is ill-conditioned), in this case INFO is set to N+3 see INFO below).``` [in] SENSE ``` SENSE is CHARACTER*1 Determines which reciprocal condition numbers are computed. = 'N' : None are computed; = 'E' : Computed for average of selected eigenvalues only; = 'V' : Computed for selected deflating subspaces only; = 'B' : Computed for both. If SENSE = 'E', 'V', or 'B', SORT must equal 'S'.``` [in] N ``` N is INTEGER The order of the matrices A, B, VSL, and VSR. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the first of the pair of matrices. On exit, A has been overwritten by its generalized Schur form S.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. LDA >= max(1,N).``` [in,out] B ``` B is COMPLEX array, dimension (LDB, N) On entry, the second of the pair of matrices. On exit, B has been overwritten by its generalized Schur form T.``` [in] LDB ``` LDB is INTEGER The leading dimension of B. LDB >= max(1,N).``` [out] SDIM ``` SDIM is INTEGER If SORT = 'N', SDIM = 0. If SORT = 'S', SDIM = number of eigenvalues (after sorting) for which SELCTG is true.``` [out] ALPHA ` ALPHA is COMPLEX array, dimension (N)` [out] BETA ``` BETA is COMPLEX array, dimension (N) On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the generalized eigenvalues. ALPHA(j) and BETA(j),j=1,...,N are the diagonals of the complex Schur form (S,T). BETA(j) will be non-negative real. Note: the quotients ALPHA(j)/BETA(j) may easily over- or underflow, and BETA(j) may even be zero. Thus, the user should avoid naively computing the ratio alpha/beta. However, ALPHA will be always less than and usually comparable with norm(A) in magnitude, and BETA always less than and usually comparable with norm(B).``` [out] VSL ``` VSL is COMPLEX array, dimension (LDVSL,N) If JOBVSL = 'V', VSL will contain the left Schur vectors. Not referenced if JOBVSL = 'N'.``` [in] LDVSL ``` LDVSL is INTEGER The leading dimension of the matrix VSL. LDVSL >=1, and if JOBVSL = 'V', LDVSL >= N.``` [out] VSR ``` VSR is COMPLEX array, dimension (LDVSR,N) If JOBVSR = 'V', VSR will contain the right Schur vectors. Not referenced if JOBVSR = 'N'.``` [in] LDVSR ``` LDVSR is INTEGER The leading dimension of the matrix VSR. LDVSR >= 1, and if JOBVSR = 'V', LDVSR >= N.``` [out] RCONDE ``` RCONDE is REAL array, dimension ( 2 ) If SENSE = 'E' or 'B', RCONDE(1) and RCONDE(2) contain the reciprocal condition numbers for the average of the selected eigenvalues. Not referenced if SENSE = 'N' or 'V'.``` [out] RCONDV ``` RCONDV is REAL array, dimension ( 2 ) If SENSE = 'V' or 'B', RCONDV(1) and RCONDV(2) contain the reciprocal condition number for the selected deflating subspaces. Not referenced if SENSE = 'N' or 'E'.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. If N = 0, LWORK >= 1, else if SENSE = 'E', 'V', or 'B', LWORK >= MAX(1,2*N,2*SDIM*(N-SDIM)), else LWORK >= MAX(1,2*N). Note that 2*SDIM*(N-SDIM) <= N*N/2. Note also that an error is only returned if LWORK < MAX(1,2*N), but if SENSE = 'E' or 'V' or 'B' this may not be large enough. If LWORK = -1, then a workspace query is assumed; the routine only calculates the bound on the optimal size of the WORK array and the minimum size of the IWORK array, returns these values as the first entries of the WORK and IWORK arrays, and no error message related to LWORK or LIWORK is issued by XERBLA.``` [out] RWORK ``` RWORK is REAL array, dimension ( 8*N ) Real workspace.``` [out] IWORK ``` IWORK is INTEGER array, dimension (MAX(1,LIWORK)) On exit, if INFO = 0, IWORK(1) returns the minimum LIWORK.``` [in] LIWORK ``` LIWORK is INTEGER The dimension of the array WORK. If SENSE = 'N' or N = 0, LIWORK >= 1, otherwise LIWORK >= N+2. If LIWORK = -1, then a workspace query is assumed; the routine only calculates the bound on the optimal size of the WORK array and the minimum size of the IWORK array, returns these values as the first entries of the WORK and IWORK arrays, and no error message related to LWORK or LIWORK is issued by XERBLA.``` [out] BWORK ``` BWORK is LOGICAL array, dimension (N) Not referenced if SORT = 'N'.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. = 1,...,N: The QZ iteration failed. (A,B) are not in Schur form, but ALPHA(j) and BETA(j) should be correct for j=INFO+1,...,N. > N: =N+1: other than QZ iteration failed in CHGEQZ =N+2: after reordering, roundoff changed values of some complex eigenvalues so that leading eigenvalues in the Generalized Schur form no longer satisfy SELCTG=.TRUE. This could also be caused due to scaling. =N+3: reordering failed in CTGSEN.```
Date
November 2011

Definition at line 332 of file cggesx.f.

332 *
333 * -- LAPACK driver routine (version 3.4.0) --
334 * -- LAPACK is a software package provided by Univ. of Tennessee, --
335 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
336 * November 2011
337 *
338 * .. Scalar Arguments ..
339  CHARACTER jobvsl, jobvsr, sense, sort
340  INTEGER info, lda, ldb, ldvsl, ldvsr, liwork, lwork, n,
341  \$ sdim
342 * ..
343 * .. Array Arguments ..
344  LOGICAL bwork( * )
345  INTEGER iwork( * )
346  REAL rconde( 2 ), rcondv( 2 ), rwork( * )
347  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
348  \$ beta( * ), vsl( ldvsl, * ), vsr( ldvsr, * ),
349  \$ work( * )
350 * ..
351 * .. Function Arguments ..
352  LOGICAL selctg
353  EXTERNAL selctg
354 * ..
355 *
356 * =====================================================================
357 *
358 * .. Parameters ..
359  REAL zero, one
360  parameter( zero = 0.0e+0, one = 1.0e+0 )
361  COMPLEX czero, cone
362  parameter( czero = ( 0.0e+0, 0.0e+0 ),
363  \$ cone = ( 1.0e+0, 0.0e+0 ) )
364 * ..
365 * .. Local Scalars ..
366  LOGICAL cursl, ilascl, ilbscl, ilvsl, ilvsr, lastsl,
367  \$ lquery, wantsb, wantse, wantsn, wantst, wantsv
368  INTEGER i, icols, ierr, ihi, ijob, ijobvl, ijobvr,
369  \$ ileft, ilo, iright, irows, irwrk, itau, iwrk,
370  \$ liwmin, lwrk, maxwrk, minwrk
371  REAL anrm, anrmto, bignum, bnrm, bnrmto, eps, pl,
372  \$ pr, smlnum
373 * ..
374 * .. Local Arrays ..
375  REAL dif( 2 )
376 * ..
377 * .. External Subroutines ..
378  EXTERNAL cgeqrf, cggbak, cggbal, cgghrd, chgeqz, clacpy,
380  \$ xerbla
381 * ..
382 * .. External Functions ..
383  LOGICAL lsame
384  INTEGER ilaenv
385  REAL clange, slamch
386  EXTERNAL lsame, ilaenv, clange, slamch
387 * ..
388 * .. Intrinsic Functions ..
389  INTRINSIC max, sqrt
390 * ..
391 * .. Executable Statements ..
392 *
393 * Decode the input arguments
394 *
395  IF( lsame( jobvsl, 'N' ) ) THEN
396  ijobvl = 1
397  ilvsl = .false.
398  ELSE IF( lsame( jobvsl, 'V' ) ) THEN
399  ijobvl = 2
400  ilvsl = .true.
401  ELSE
402  ijobvl = -1
403  ilvsl = .false.
404  END IF
405 *
406  IF( lsame( jobvsr, 'N' ) ) THEN
407  ijobvr = 1
408  ilvsr = .false.
409  ELSE IF( lsame( jobvsr, 'V' ) ) THEN
410  ijobvr = 2
411  ilvsr = .true.
412  ELSE
413  ijobvr = -1
414  ilvsr = .false.
415  END IF
416 *
417  wantst = lsame( sort, 'S' )
418  wantsn = lsame( sense, 'N' )
419  wantse = lsame( sense, 'E' )
420  wantsv = lsame( sense, 'V' )
421  wantsb = lsame( sense, 'B' )
422  lquery = ( lwork.EQ.-1 .OR. liwork.EQ.-1 )
423  IF( wantsn ) THEN
424  ijob = 0
425  ELSE IF( wantse ) THEN
426  ijob = 1
427  ELSE IF( wantsv ) THEN
428  ijob = 2
429  ELSE IF( wantsb ) THEN
430  ijob = 4
431  END IF
432 *
433 * Test the input arguments
434 *
435  info = 0
436  IF( ijobvl.LE.0 ) THEN
437  info = -1
438  ELSE IF( ijobvr.LE.0 ) THEN
439  info = -2
440  ELSE IF( ( .NOT.wantst ) .AND. ( .NOT.lsame( sort, 'N' ) ) ) THEN
441  info = -3
442  ELSE IF( .NOT.( wantsn .OR. wantse .OR. wantsv .OR. wantsb ) .OR.
443  \$ ( .NOT.wantst .AND. .NOT.wantsn ) ) THEN
444  info = -5
445  ELSE IF( n.LT.0 ) THEN
446  info = -6
447  ELSE IF( lda.LT.max( 1, n ) ) THEN
448  info = -8
449  ELSE IF( ldb.LT.max( 1, n ) ) THEN
450  info = -10
451  ELSE IF( ldvsl.LT.1 .OR. ( ilvsl .AND. ldvsl.LT.n ) ) THEN
452  info = -15
453  ELSE IF( ldvsr.LT.1 .OR. ( ilvsr .AND. ldvsr.LT.n ) ) THEN
454  info = -17
455  END IF
456 *
457 * Compute workspace
458 * (Note: Comments in the code beginning "Workspace:" describe the
459 * minimal amount of workspace needed at that point in the code,
460 * as well as the preferred amount for good performance.
461 * NB refers to the optimal block size for the immediately
462 * following subroutine, as returned by ILAENV.)
463 *
464  IF( info.EQ.0 ) THEN
465  IF( n.GT.0) THEN
466  minwrk = 2*n
467  maxwrk = n*(1 + ilaenv( 1, 'CGEQRF', ' ', n, 1, n, 0 ) )
468  maxwrk = max( maxwrk, n*( 1 +
469  \$ ilaenv( 1, 'CUNMQR', ' ', n, 1, n, -1 ) ) )
470  IF( ilvsl ) THEN
471  maxwrk = max( maxwrk, n*( 1 +
472  \$ ilaenv( 1, 'CUNGQR', ' ', n, 1, n, -1 ) ) )
473  END IF
474  lwrk = maxwrk
475  IF( ijob.GE.1 )
476  \$ lwrk = max( lwrk, n*n/2 )
477  ELSE
478  minwrk = 1
479  maxwrk = 1
480  lwrk = 1
481  END IF
482  work( 1 ) = lwrk
483  IF( wantsn .OR. n.EQ.0 ) THEN
484  liwmin = 1
485  ELSE
486  liwmin = n + 2
487  END IF
488  iwork( 1 ) = liwmin
489 *
490  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
491  info = -21
492  ELSE IF( liwork.LT.liwmin .AND. .NOT.lquery) THEN
493  info = -24
494  END IF
495  END IF
496 *
497  IF( info.NE.0 ) THEN
498  CALL xerbla( 'CGGESX', -info )
499  RETURN
500  ELSE IF (lquery) THEN
501  RETURN
502  END IF
503 *
504 * Quick return if possible
505 *
506  IF( n.EQ.0 ) THEN
507  sdim = 0
508  RETURN
509  END IF
510 *
511 * Get machine constants
512 *
513  eps = slamch( 'P' )
514  smlnum = slamch( 'S' )
515  bignum = one / smlnum
516  CALL slabad( smlnum, bignum )
517  smlnum = sqrt( smlnum ) / eps
518  bignum = one / smlnum
519 *
520 * Scale A if max element outside range [SMLNUM,BIGNUM]
521 *
522  anrm = clange( 'M', n, n, a, lda, rwork )
523  ilascl = .false.
524  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
525  anrmto = smlnum
526  ilascl = .true.
527  ELSE IF( anrm.GT.bignum ) THEN
528  anrmto = bignum
529  ilascl = .true.
530  END IF
531  IF( ilascl )
532  \$ CALL clascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
533 *
534 * Scale B if max element outside range [SMLNUM,BIGNUM]
535 *
536  bnrm = clange( 'M', n, n, b, ldb, rwork )
537  ilbscl = .false.
538  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
539  bnrmto = smlnum
540  ilbscl = .true.
541  ELSE IF( bnrm.GT.bignum ) THEN
542  bnrmto = bignum
543  ilbscl = .true.
544  END IF
545  IF( ilbscl )
546  \$ CALL clascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
547 *
548 * Permute the matrix to make it more nearly triangular
549 * (Real Workspace: need 6*N)
550 *
551  ileft = 1
552  iright = n + 1
553  irwrk = iright + n
554  CALL cggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
555  \$ rwork( iright ), rwork( irwrk ), ierr )
556 *
557 * Reduce B to triangular form (QR decomposition of B)
558 * (Complex Workspace: need N, prefer N*NB)
559 *
560  irows = ihi + 1 - ilo
561  icols = n + 1 - ilo
562  itau = 1
563  iwrk = itau + irows
564  CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
565  \$ work( iwrk ), lwork+1-iwrk, ierr )
566 *
567 * Apply the unitary transformation to matrix A
568 * (Complex Workspace: need N, prefer N*NB)
569 *
570  CALL cunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
571  \$ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
572  \$ lwork+1-iwrk, ierr )
573 *
574 * Initialize VSL
575 * (Complex Workspace: need N, prefer N*NB)
576 *
577  IF( ilvsl ) THEN
578  CALL claset( 'Full', n, n, czero, cone, vsl, ldvsl )
579  IF( irows.GT.1 ) THEN
580  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
581  \$ vsl( ilo+1, ilo ), ldvsl )
582  END IF
583  CALL cungqr( irows, irows, irows, vsl( ilo, ilo ), ldvsl,
584  \$ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
585  END IF
586 *
587 * Initialize VSR
588 *
589  IF( ilvsr )
590  \$ CALL claset( 'Full', n, n, czero, cone, vsr, ldvsr )
591 *
592 * Reduce to generalized Hessenberg form
593 * (Workspace: none needed)
594 *
595  CALL cgghrd( jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb, vsl,
596  \$ ldvsl, vsr, ldvsr, ierr )
597 *
598  sdim = 0
599 *
600 * Perform QZ algorithm, computing Schur vectors if desired
601 * (Complex Workspace: need N)
602 * (Real Workspace: need N)
603 *
604  iwrk = itau
605  CALL chgeqz( 'S', jobvsl, jobvsr, n, ilo, ihi, a, lda, b, ldb,
606  \$ alpha, beta, vsl, ldvsl, vsr, ldvsr, work( iwrk ),
607  \$ lwork+1-iwrk, rwork( irwrk ), ierr )
608  IF( ierr.NE.0 ) THEN
609  IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
610  info = ierr
611  ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
612  info = ierr - n
613  ELSE
614  info = n + 1
615  END IF
616  GO TO 40
617  END IF
618 *
619 * Sort eigenvalues ALPHA/BETA and compute the reciprocal of
620 * condition number(s)
621 *
622  IF( wantst ) THEN
623 *
624 * Undo scaling on eigenvalues before SELCTGing
625 *
626  IF( ilascl )
627  \$ CALL clascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
628  IF( ilbscl )
629  \$ CALL clascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
630 *
631 * Select eigenvalues
632 *
633  DO 10 i = 1, n
634  bwork( i ) = selctg( alpha( i ), beta( i ) )
635  10 CONTINUE
636 *
637 * Reorder eigenvalues, transform Generalized Schur vectors, and
638 * compute reciprocal condition numbers
639 * (Complex Workspace: If IJOB >= 1, need MAX(1, 2*SDIM*(N-SDIM))
640 * otherwise, need 1 )
641 *
642  CALL ctgsen( ijob, ilvsl, ilvsr, bwork, n, a, lda, b, ldb,
643  \$ alpha, beta, vsl, ldvsl, vsr, ldvsr, sdim, pl, pr,
644  \$ dif, work( iwrk ), lwork-iwrk+1, iwork, liwork,
645  \$ ierr )
646 *
647  IF( ijob.GE.1 )
648  \$ maxwrk = max( maxwrk, 2*sdim*( n-sdim ) )
649  IF( ierr.EQ.-21 ) THEN
650 *
651 * not enough complex workspace
652 *
653  info = -21
654  ELSE
655  IF( ijob.EQ.1 .OR. ijob.EQ.4 ) THEN
656  rconde( 1 ) = pl
657  rconde( 2 ) = pr
658  END IF
659  IF( ijob.EQ.2 .OR. ijob.EQ.4 ) THEN
660  rcondv( 1 ) = dif( 1 )
661  rcondv( 2 ) = dif( 2 )
662  END IF
663  IF( ierr.EQ.1 )
664  \$ info = n + 3
665  END IF
666 *
667  END IF
668 *
669 * Apply permutation to VSL and VSR
670 * (Workspace: none needed)
671 *
672  IF( ilvsl )
673  \$ CALL cggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
674  \$ rwork( iright ), n, vsl, ldvsl, ierr )
675 *
676  IF( ilvsr )
677  \$ CALL cggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
678  \$ rwork( iright ), n, vsr, ldvsr, ierr )
679 *
680 * Undo scaling
681 *
682  IF( ilascl ) THEN
683  CALL clascl( 'U', 0, 0, anrmto, anrm, n, n, a, lda, ierr )
684  CALL clascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
685  END IF
686 *
687  IF( ilbscl ) THEN
688  CALL clascl( 'U', 0, 0, bnrmto, bnrm, n, n, b, ldb, ierr )
689  CALL clascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
690  END IF
691 *
692  IF( wantst ) THEN
693 *
694 * Check if reordering is correct
695 *
696  lastsl = .true.
697  sdim = 0
698  DO 30 i = 1, n
699  cursl = selctg( alpha( i ), beta( i ) )
700  IF( cursl )
701  \$ sdim = sdim + 1
702  IF( cursl .AND. .NOT.lastsl )
703  \$ info = n + 2
704  lastsl = cursl
705  30 CONTINUE
706 *
707  END IF
708 *
709  40 CONTINUE
710 *
711  work( 1 ) = maxwrk
712  iwork( 1 ) = liwmin
713 *
714  RETURN
715 *
716 * End of CGGESX
717 *
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:170
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:179
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:150
subroutine cgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
CGGHRD
Definition: cgghrd.f:206
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine ctgsen(IJOB, WANTQ, WANTZ, SELECT, N, A, LDA, B, LDB, ALPHA, BETA, Q, LDQ, Z, LDZ, M, PL, PR, DIF, WORK, LWORK, IWORK, LIWORK, INFO)
CTGSEN
Definition: ctgsen.f:435
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:130

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cggev ( character JOBVL, character JOBVR, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( ldb, * ) B, integer LDB, complex, dimension( * ) ALPHA, complex, dimension( * ) BETA, complex, dimension( ldvl, * ) VL, integer LDVL, complex, dimension( ldvr, * ) VR, integer LDVR, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, integer INFO )

CGGEV computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

Purpose:
``` CGGEV computes for a pair of N-by-N complex nonsymmetric matrices
(A,B), the generalized eigenvalues, and optionally, the left and/or
right generalized eigenvectors.

A generalized eigenvalue for a pair of matrices (A,B) is a scalar
lambda or a ratio alpha/beta = lambda, such that A - lambda*B is
singular. It is usually represented as the pair (alpha,beta), as
there is a reasonable interpretation for beta=0, and even for both
being zero.

The right generalized eigenvector v(j) corresponding to the
generalized eigenvalue lambda(j) of (A,B) satisfies

A * v(j) = lambda(j) * B * v(j).

The left generalized eigenvector u(j) corresponding to the
generalized eigenvalues lambda(j) of (A,B) satisfies

u(j)**H * A = lambda(j) * u(j)**H * B

where u(j)**H is the conjugate-transpose of u(j).```
Parameters
 [in] JOBVL ``` JOBVL is CHARACTER*1 = 'N': do not compute the left generalized eigenvectors; = 'V': compute the left generalized eigenvectors.``` [in] JOBVR ``` JOBVR is CHARACTER*1 = 'N': do not compute the right generalized eigenvectors; = 'V': compute the right generalized eigenvectors.``` [in] N ``` N is INTEGER The order of the matrices A, B, VL, and VR. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the matrix A in the pair (A,B). On exit, A has been overwritten.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. LDA >= max(1,N).``` [in,out] B ``` B is COMPLEX array, dimension (LDB, N) On entry, the matrix B in the pair (A,B). On exit, B has been overwritten.``` [in] LDB ``` LDB is INTEGER The leading dimension of B. LDB >= max(1,N).``` [out] ALPHA ` ALPHA is COMPLEX array, dimension (N)` [out] BETA ``` BETA is COMPLEX array, dimension (N) On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the generalized eigenvalues. Note: the quotients ALPHA(j)/BETA(j) may easily over- or underflow, and BETA(j) may even be zero. Thus, the user should avoid naively computing the ratio alpha/beta. However, ALPHA will be always less than and usually comparable with norm(A) in magnitude, and BETA always less than and usually comparable with norm(B).``` [out] VL ``` VL is COMPLEX array, dimension (LDVL,N) If JOBVL = 'V', the left generalized eigenvectors u(j) are stored one after another in the columns of VL, in the same order as their eigenvalues. Each eigenvector is scaled so the largest component has abs(real part) + abs(imag. part) = 1. Not referenced if JOBVL = 'N'.``` [in] LDVL ``` LDVL is INTEGER The leading dimension of the matrix VL. LDVL >= 1, and if JOBVL = 'V', LDVL >= N.``` [out] VR ``` VR is COMPLEX array, dimension (LDVR,N) If JOBVR = 'V', the right generalized eigenvectors v(j) are stored one after another in the columns of VR, in the same order as their eigenvalues. Each eigenvector is scaled so the largest component has abs(real part) + abs(imag. part) = 1. Not referenced if JOBVR = 'N'.``` [in] LDVR ``` LDVR is INTEGER The leading dimension of the matrix VR. LDVR >= 1, and if JOBVR = 'V', LDVR >= N.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,2*N). For good performance, LWORK must generally be larger. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (8*N)` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. =1,...,N: The QZ iteration failed. No eigenvectors have been calculated, but ALPHA(j) and BETA(j) should be correct for j=INFO+1,...,N. > N: =N+1: other then QZ iteration failed in SHGEQZ, =N+2: error return from STGEVC.```
Date
April 2012

Definition at line 219 of file cggev.f.

219 *
220 * -- LAPACK driver routine (version 3.4.1) --
221 * -- LAPACK is a software package provided by Univ. of Tennessee, --
222 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
223 * April 2012
224 *
225 * .. Scalar Arguments ..
226  CHARACTER jobvl, jobvr
227  INTEGER info, lda, ldb, ldvl, ldvr, lwork, n
228 * ..
229 * .. Array Arguments ..
230  REAL rwork( * )
231  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
232  \$ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
233  \$ work( * )
234 * ..
235 *
236 * =====================================================================
237 *
238 * .. Parameters ..
239  REAL zero, one
240  parameter( zero = 0.0e0, one = 1.0e0 )
241  COMPLEX czero, cone
242  parameter( czero = ( 0.0e0, 0.0e0 ),
243  \$ cone = ( 1.0e0, 0.0e0 ) )
244 * ..
245 * .. Local Scalars ..
246  LOGICAL ilascl, ilbscl, ilv, ilvl, ilvr, lquery
247  CHARACTER chtemp
248  INTEGER icols, ierr, ihi, ijobvl, ijobvr, ileft, ilo,
249  \$ in, iright, irows, irwrk, itau, iwrk, jc, jr,
250  \$ lwkmin, lwkopt
251  REAL anrm, anrmto, bignum, bnrm, bnrmto, eps,
252  \$ smlnum, temp
253  COMPLEX x
254 * ..
255 * .. Local Arrays ..
256  LOGICAL ldumma( 1 )
257 * ..
258 * .. External Subroutines ..
259  EXTERNAL cgeqrf, cggbak, cggbal, cgghrd, chgeqz, clacpy,
261  \$ xerbla
262 * ..
263 * .. External Functions ..
264  LOGICAL lsame
265  INTEGER ilaenv
266  REAL clange, slamch
267  EXTERNAL lsame, ilaenv, clange, slamch
268 * ..
269 * .. Intrinsic Functions ..
270  INTRINSIC abs, aimag, max, REAL, sqrt
271 * ..
272 * .. Statement Functions ..
273  REAL abs1
274 * ..
275 * .. Statement Function definitions ..
276  abs1( x ) = abs( REAL( X ) ) + abs( aimag( x ) )
277 * ..
278 * .. Executable Statements ..
279 *
280 * Decode the input arguments
281 *
282  IF( lsame( jobvl, 'N' ) ) THEN
283  ijobvl = 1
284  ilvl = .false.
285  ELSE IF( lsame( jobvl, 'V' ) ) THEN
286  ijobvl = 2
287  ilvl = .true.
288  ELSE
289  ijobvl = -1
290  ilvl = .false.
291  END IF
292 *
293  IF( lsame( jobvr, 'N' ) ) THEN
294  ijobvr = 1
295  ilvr = .false.
296  ELSE IF( lsame( jobvr, 'V' ) ) THEN
297  ijobvr = 2
298  ilvr = .true.
299  ELSE
300  ijobvr = -1
301  ilvr = .false.
302  END IF
303  ilv = ilvl .OR. ilvr
304 *
305 * Test the input arguments
306 *
307  info = 0
308  lquery = ( lwork.EQ.-1 )
309  IF( ijobvl.LE.0 ) THEN
310  info = -1
311  ELSE IF( ijobvr.LE.0 ) THEN
312  info = -2
313  ELSE IF( n.LT.0 ) THEN
314  info = -3
315  ELSE IF( lda.LT.max( 1, n ) ) THEN
316  info = -5
317  ELSE IF( ldb.LT.max( 1, n ) ) THEN
318  info = -7
319  ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
320  info = -11
321  ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
322  info = -13
323  END IF
324 *
325 * Compute workspace
326 * (Note: Comments in the code beginning "Workspace:" describe the
327 * minimal amount of workspace needed at that point in the code,
328 * as well as the preferred amount for good performance.
329 * NB refers to the optimal block size for the immediately
330 * following subroutine, as returned by ILAENV. The workspace is
331 * computed assuming ILO = 1 and IHI = N, the worst case.)
332 *
333  IF( info.EQ.0 ) THEN
334  lwkmin = max( 1, 2*n )
335  lwkopt = max( 1, n + n*ilaenv( 1, 'CGEQRF', ' ', n, 1, n, 0 ) )
336  lwkopt = max( lwkopt, n +
337  \$ n*ilaenv( 1, 'CUNMQR', ' ', n, 1, n, 0 ) )
338  IF( ilvl ) THEN
339  lwkopt = max( lwkopt, n +
340  \$ n*ilaenv( 1, 'CUNGQR', ' ', n, 1, n, -1 ) )
341  END IF
342  work( 1 ) = lwkopt
343 *
344  IF( lwork.LT.lwkmin .AND. .NOT.lquery )
345  \$ info = -15
346  END IF
347 *
348  IF( info.NE.0 ) THEN
349  CALL xerbla( 'CGGEV ', -info )
350  RETURN
351  ELSE IF( lquery ) THEN
352  RETURN
353  END IF
354 *
355 * Quick return if possible
356 *
357  IF( n.EQ.0 )
358  \$ RETURN
359 *
360 * Get machine constants
361 *
362  eps = slamch( 'E' )*slamch( 'B' )
363  smlnum = slamch( 'S' )
364  bignum = one / smlnum
365  CALL slabad( smlnum, bignum )
366  smlnum = sqrt( smlnum ) / eps
367  bignum = one / smlnum
368 *
369 * Scale A if max element outside range [SMLNUM,BIGNUM]
370 *
371  anrm = clange( 'M', n, n, a, lda, rwork )
372  ilascl = .false.
373  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
374  anrmto = smlnum
375  ilascl = .true.
376  ELSE IF( anrm.GT.bignum ) THEN
377  anrmto = bignum
378  ilascl = .true.
379  END IF
380  IF( ilascl )
381  \$ CALL clascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
382 *
383 * Scale B if max element outside range [SMLNUM,BIGNUM]
384 *
385  bnrm = clange( 'M', n, n, b, ldb, rwork )
386  ilbscl = .false.
387  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
388  bnrmto = smlnum
389  ilbscl = .true.
390  ELSE IF( bnrm.GT.bignum ) THEN
391  bnrmto = bignum
392  ilbscl = .true.
393  END IF
394  IF( ilbscl )
395  \$ CALL clascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
396 *
397 * Permute the matrices A, B to isolate eigenvalues if possible
398 * (Real Workspace: need 6*N)
399 *
400  ileft = 1
401  iright = n + 1
402  irwrk = iright + n
403  CALL cggbal( 'P', n, a, lda, b, ldb, ilo, ihi, rwork( ileft ),
404  \$ rwork( iright ), rwork( irwrk ), ierr )
405 *
406 * Reduce B to triangular form (QR decomposition of B)
407 * (Complex Workspace: need N, prefer N*NB)
408 *
409  irows = ihi + 1 - ilo
410  IF( ilv ) THEN
411  icols = n + 1 - ilo
412  ELSE
413  icols = irows
414  END IF
415  itau = 1
416  iwrk = itau + irows
417  CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
418  \$ work( iwrk ), lwork+1-iwrk, ierr )
419 *
420 * Apply the orthogonal transformation to matrix A
421 * (Complex Workspace: need N, prefer N*NB)
422 *
423  CALL cunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
424  \$ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
425  \$ lwork+1-iwrk, ierr )
426 *
427 * Initialize VL
428 * (Complex Workspace: need N, prefer N*NB)
429 *
430  IF( ilvl ) THEN
431  CALL claset( 'Full', n, n, czero, cone, vl, ldvl )
432  IF( irows.GT.1 ) THEN
433  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
434  \$ vl( ilo+1, ilo ), ldvl )
435  END IF
436  CALL cungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
437  \$ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
438  END IF
439 *
440 * Initialize VR
441 *
442  IF( ilvr )
443  \$ CALL claset( 'Full', n, n, czero, cone, vr, ldvr )
444 *
445 * Reduce to generalized Hessenberg form
446 *
447  IF( ilv ) THEN
448 *
449 * Eigenvectors requested -- work on whole matrix.
450 *
451  CALL cgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
452  \$ ldvl, vr, ldvr, ierr )
453  ELSE
454  CALL cgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
455  \$ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, ierr )
456  END IF
457 *
458 * Perform QZ algorithm (Compute eigenvalues, and optionally, the
459 * Schur form and Schur vectors)
460 * (Complex Workspace: need N)
461 * (Real Workspace: need N)
462 *
463  iwrk = itau
464  IF( ilv ) THEN
465  chtemp = 'S'
466  ELSE
467  chtemp = 'E'
468  END IF
469  CALL chgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
470  \$ alpha, beta, vl, ldvl, vr, ldvr, work( iwrk ),
471  \$ lwork+1-iwrk, rwork( irwrk ), ierr )
472  IF( ierr.NE.0 ) THEN
473  IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
474  info = ierr
475  ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
476  info = ierr - n
477  ELSE
478  info = n + 1
479  END IF
480  GO TO 70
481  END IF
482 *
483 * Compute Eigenvectors
484 * (Real Workspace: need 2*N)
485 * (Complex Workspace: need 2*N)
486 *
487  IF( ilv ) THEN
488  IF( ilvl ) THEN
489  IF( ilvr ) THEN
490  chtemp = 'B'
491  ELSE
492  chtemp = 'L'
493  END IF
494  ELSE
495  chtemp = 'R'
496  END IF
497 *
498  CALL ctgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
499  \$ vr, ldvr, n, in, work( iwrk ), rwork( irwrk ),
500  \$ ierr )
501  IF( ierr.NE.0 ) THEN
502  info = n + 2
503  GO TO 70
504  END IF
505 *
506 * Undo balancing on VL and VR and normalization
507 * (Workspace: none needed)
508 *
509  IF( ilvl ) THEN
510  CALL cggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
511  \$ rwork( iright ), n, vl, ldvl, ierr )
512  DO 30 jc = 1, n
513  temp = zero
514  DO 10 jr = 1, n
515  temp = max( temp, abs1( vl( jr, jc ) ) )
516  10 CONTINUE
517  IF( temp.LT.smlnum )
518  \$ GO TO 30
519  temp = one / temp
520  DO 20 jr = 1, n
521  vl( jr, jc ) = vl( jr, jc )*temp
522  20 CONTINUE
523  30 CONTINUE
524  END IF
525  IF( ilvr ) THEN
526  CALL cggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
527  \$ rwork( iright ), n, vr, ldvr, ierr )
528  DO 60 jc = 1, n
529  temp = zero
530  DO 40 jr = 1, n
531  temp = max( temp, abs1( vr( jr, jc ) ) )
532  40 CONTINUE
533  IF( temp.LT.smlnum )
534  \$ GO TO 60
535  temp = one / temp
536  DO 50 jr = 1, n
537  vr( jr, jc ) = vr( jr, jc )*temp
538  50 CONTINUE
539  60 CONTINUE
540  END IF
541  END IF
542 *
543 * Undo scaling if necessary
544 *
545  70 CONTINUE
546 *
547  IF( ilascl )
548  \$ CALL clascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
549 *
550  IF( ilbscl )
551  \$ CALL clascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
552 *
553  work( 1 ) = lwkopt
554  RETURN
555 *
556 * End of CGGEV
557 *
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:170
subroutine ctgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTGEVC
Definition: ctgevc.f:221
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:179
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:150
subroutine cgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
CGGHRD
Definition: cgghrd.f:206
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:130

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cggev3 ( character JOBVL, character JOBVR, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( ldb, * ) B, integer LDB, complex, dimension( * ) ALPHA, complex, dimension( * ) BETA, complex, dimension( ldvl, * ) VL, integer LDVL, complex, dimension( ldvr, * ) VR, integer LDVR, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, integer INFO )

CGGEV3 computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices (blocked algorithm)

Purpose:
``` CGGEV3 computes for a pair of N-by-N complex nonsymmetric matrices
(A,B), the generalized eigenvalues, and optionally, the left and/or
right generalized eigenvectors.

A generalized eigenvalue for a pair of matrices (A,B) is a scalar
lambda or a ratio alpha/beta = lambda, such that A - lambda*B is
singular. It is usually represented as the pair (alpha,beta), as
there is a reasonable interpretation for beta=0, and even for both
being zero.

The right generalized eigenvector v(j) corresponding to the
generalized eigenvalue lambda(j) of (A,B) satisfies

A * v(j) = lambda(j) * B * v(j).

The left generalized eigenvector u(j) corresponding to the
generalized eigenvalues lambda(j) of (A,B) satisfies

u(j)**H * A = lambda(j) * u(j)**H * B

where u(j)**H is the conjugate-transpose of u(j).```
Parameters
 [in] JOBVL ``` JOBVL is CHARACTER*1 = 'N': do not compute the left generalized eigenvectors; = 'V': compute the left generalized eigenvectors.``` [in] JOBVR ``` JOBVR is CHARACTER*1 = 'N': do not compute the right generalized eigenvectors; = 'V': compute the right generalized eigenvectors.``` [in] N ``` N is INTEGER The order of the matrices A, B, VL, and VR. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the matrix A in the pair (A,B). On exit, A has been overwritten.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. LDA >= max(1,N).``` [in,out] B ``` B is COMPLEX array, dimension (LDB, N) On entry, the matrix B in the pair (A,B). On exit, B has been overwritten.``` [in] LDB ``` LDB is INTEGER The leading dimension of B. LDB >= max(1,N).``` [out] ALPHA ` ALPHA is COMPLEX array, dimension (N)` [out] BETA ``` BETA is COMPLEX array, dimension (N) On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the generalized eigenvalues. Note: the quotients ALPHA(j)/BETA(j) may easily over- or underflow, and BETA(j) may even be zero. Thus, the user should avoid naively computing the ratio alpha/beta. However, ALPHA will be always less than and usually comparable with norm(A) in magnitude, and BETA always less than and usually comparable with norm(B).``` [out] VL ``` VL is COMPLEX array, dimension (LDVL,N) If JOBVL = 'V', the left generalized eigenvectors u(j) are stored one after another in the columns of VL, in the same order as their eigenvalues. Each eigenvector is scaled so the largest component has abs(real part) + abs(imag. part) = 1. Not referenced if JOBVL = 'N'.``` [in] LDVL ``` LDVL is INTEGER The leading dimension of the matrix VL. LDVL >= 1, and if JOBVL = 'V', LDVL >= N.``` [out] VR ``` VR is COMPLEX array, dimension (LDVR,N) If JOBVR = 'V', the right generalized eigenvectors v(j) are stored one after another in the columns of VR, in the same order as their eigenvalues. Each eigenvector is scaled so the largest component has abs(real part) + abs(imag. part) = 1. Not referenced if JOBVR = 'N'.``` [in] LDVR ``` LDVR is INTEGER The leading dimension of the matrix VR. LDVR >= 1, and if JOBVR = 'V', LDVR >= N.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ` RWORK is REAL array, dimension (8*N)` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. =1,...,N: The QZ iteration failed. No eigenvectors have been calculated, but ALPHA(j) and BETA(j) should be correct for j=INFO+1,...,N. > N: =N+1: other then QZ iteration failed in SHGEQZ, =N+2: error return from STGEVC.```
Date
January 2015

Definition at line 218 of file cggev3.f.

218 *
219 * -- LAPACK driver routine (version 3.6.0) --
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  REAL rwork( * )
230  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
231  \$ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
232  \$ work( * )
233 * ..
234 *
235 * =====================================================================
236 *
237 * .. Parameters ..
238  REAL zero, one
239  parameter( zero = 0.0e0, one = 1.0e0 )
240  COMPLEX czero, cone
241  parameter( czero = ( 0.0e0, 0.0e0 ),
242  \$ cone = ( 1.0e0, 0.0e0 ) )
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  REAL anrm, anrmto, bignum, bnrm, bnrmto, eps,
251  \$ smlnum, temp
252  COMPLEX x
253 * ..
254 * .. Local Arrays ..
255  LOGICAL ldumma( 1 )
256 * ..
257 * .. External Subroutines ..
258  EXTERNAL cgeqrf, cggbak, cggbal, cgghd3, chgeqz, clacpy,
260  \$ xerbla
261 * ..
262 * .. External Functions ..
263  LOGICAL lsame
264  REAL clange, slamch
265  EXTERNAL lsame, clange, slamch
266 * ..
267 * .. Intrinsic Functions ..
268  INTRINSIC abs, aimag, max, REAL, sqrt
269 * ..
270 * .. Statement Functions ..
271  REAL abs1
272 * ..
273 * .. Statement Function definitions ..
274  abs1( x ) = abs( REAL( X ) ) + abs( aimag( 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 cgeqrf( n, n, b, ldb, work, work, -1, ierr )
329  lwkopt = max( n, n+int( work( 1 ) ) )
330  CALL cunmqr( '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 cungqr( 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 cgghd3( 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 chgeqz( 'S', jobvl, jobvr, n, 1, n, a, lda, b, ldb,
342  \$ alpha, beta, vl, ldvl, vr, ldvr, work,
343  \$ -1, work, ierr )
344  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
345  ELSE
346  CALL cgghd3( 'N', 'N', n, 1, n, a, lda, b, ldb, vl, ldvl,
347  \$ vr, ldvr, work, -1, ierr )
348  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
349  CALL chgeqz( 'E', jobvl, jobvr, n, 1, n, a, lda, b, ldb,
350  \$ alpha, beta, vl, ldvl, vr, ldvr, work,
351  \$ -1, work, ierr )
352  lwkopt = max( lwkopt, n+int( work( 1 ) ) )
353  END IF
354  work( 1 ) = cmplx( lwkopt )
355  END IF
356 *
357  IF( info.NE.0 ) THEN
358  CALL xerbla( 'CGGEV3 ', -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 = slamch( 'E' )*slamch( 'B' )
372  smlnum = slamch( 'S' )
373  bignum = one / smlnum
374  CALL slabad( 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 = clange( '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 clascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
391 *
392 * Scale B if max element outside range [SMLNUM,BIGNUM]
393 *
394  bnrm = clange( '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 clascl( '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 cggbal( '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 cgeqrf( 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 cunmqr( '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 claset( 'Full', n, n, czero, cone, vl, ldvl )
437  IF( irows.GT.1 ) THEN
438  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
439  \$ vl( ilo+1, ilo ), ldvl )
440  END IF
441  CALL cungqr( 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 claset( '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 cgghd3( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
457  \$ ldvl, vr, ldvr, work( iwrk ), lwork+1-iwrk,
458  \$ ierr )
459  ELSE
460  CALL cgghd3( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
461  \$ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr,
462  \$ work( iwrk ), lwork+1-iwrk, ierr )
463  END IF
464 *
465 * Perform QZ algorithm (Compute eigenvalues, and optionally, the
466 * Schur form and Schur vectors)
467 *
468  iwrk = itau
469  IF( ilv ) THEN
470  chtemp = 'S'
471  ELSE
472  chtemp = 'E'
473  END IF
474  CALL chgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
475  \$ alpha, beta, vl, ldvl, vr, ldvr, work( iwrk ),
476  \$ lwork+1-iwrk, rwork( irwrk ), ierr )
477  IF( ierr.NE.0 ) THEN
478  IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
479  info = ierr
480  ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
481  info = ierr - n
482  ELSE
483  info = n + 1
484  END IF
485  GO TO 70
486  END IF
487 *
488 * Compute Eigenvectors
489 *
490  IF( ilv ) THEN
491  IF( ilvl ) THEN
492  IF( ilvr ) THEN
493  chtemp = 'B'
494  ELSE
495  chtemp = 'L'
496  END IF
497  ELSE
498  chtemp = 'R'
499  END IF
500 *
501  CALL ctgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl, ldvl,
502  \$ vr, ldvr, n, in, work( iwrk ), rwork( irwrk ),
503  \$ ierr )
504  IF( ierr.NE.0 ) THEN
505  info = n + 2
506  GO TO 70
507  END IF
508 *
509 * Undo balancing on VL and VR and normalization
510 *
511  IF( ilvl ) THEN
512  CALL cggbak( 'P', 'L', n, ilo, ihi, rwork( ileft ),
513  \$ rwork( iright ), n, vl, ldvl, ierr )
514  DO 30 jc = 1, n
515  temp = zero
516  DO 10 jr = 1, n
517  temp = max( temp, abs1( vl( jr, jc ) ) )
518  10 CONTINUE
519  IF( temp.LT.smlnum )
520  \$ GO TO 30
521  temp = one / temp
522  DO 20 jr = 1, n
523  vl( jr, jc ) = vl( jr, jc )*temp
524  20 CONTINUE
525  30 CONTINUE
526  END IF
527  IF( ilvr ) THEN
528  CALL cggbak( 'P', 'R', n, ilo, ihi, rwork( ileft ),
529  \$ rwork( iright ), n, vr, ldvr, ierr )
530  DO 60 jc = 1, n
531  temp = zero
532  DO 40 jr = 1, n
533  temp = max( temp, abs1( vr( jr, jc ) ) )
534  40 CONTINUE
535  IF( temp.LT.smlnum )
536  \$ GO TO 60
537  temp = one / temp
538  DO 50 jr = 1, n
539  vr( jr, jc ) = vr( jr, jc )*temp
540  50 CONTINUE
541  60 CONTINUE
542  END IF
543  END IF
544 *
545 * Undo scaling if necessary
546 *
547  70 CONTINUE
548 *
549  IF( ilascl )
550  \$ CALL clascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
551 *
552  IF( ilbscl )
553  \$ CALL clascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
554 *
555  work( 1 ) = cmplx( lwkopt )
556  RETURN
557 *
558 * End of CGGEV3
559 *
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:170
subroutine ctgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTGEVC
Definition: ctgevc.f:221
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cgghd3(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, WORK, LWORK, INFO)
CGGHD3
Definition: cgghd3.f:233
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:179
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:150
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:130

Here is the call graph for this function:

Here is the caller graph for this function:

 subroutine cggevx ( character BALANC, character JOBVL, character JOBVR, character SENSE, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( ldb, * ) B, integer LDB, complex, dimension( * ) ALPHA, complex, dimension( * ) BETA, complex, dimension( ldvl, * ) VL, integer LDVL, complex, dimension( ldvr, * ) VR, integer LDVR, integer ILO, integer IHI, real, dimension( * ) LSCALE, real, dimension( * ) RSCALE, real ABNRM, real BBNRM, real, dimension( * ) RCONDE, real, dimension( * ) RCONDV, complex, dimension( * ) WORK, integer LWORK, real, dimension( * ) RWORK, integer, dimension( * ) IWORK, logical, dimension( * ) BWORK, integer INFO )

CGGEVX computes the eigenvalues and, optionally, the left and/or right eigenvectors for GE matrices

Purpose:
``` CGGEVX computes for a pair of N-by-N complex nonsymmetric matrices
(A,B) the generalized eigenvalues, and optionally, the left and/or
right generalized eigenvectors.

Optionally, it also computes a balancing transformation to improve
the conditioning of the eigenvalues and eigenvectors (ILO, IHI,
LSCALE, RSCALE, ABNRM, and BBNRM), reciprocal condition numbers for
the eigenvalues (RCONDE), and reciprocal condition numbers for the
right eigenvectors (RCONDV).

A generalized eigenvalue for a pair of matrices (A,B) is a scalar
lambda or a ratio alpha/beta = lambda, such that A - lambda*B is
singular. It is usually represented as the pair (alpha,beta), as
there is a reasonable interpretation for beta=0, and even for both
being zero.

The right eigenvector v(j) corresponding to the eigenvalue lambda(j)
of (A,B) satisfies
A * v(j) = lambda(j) * B * v(j) .
The left eigenvector u(j) corresponding to the eigenvalue lambda(j)
of (A,B) satisfies
u(j)**H * A  = lambda(j) * u(j)**H * B.
where u(j)**H is the conjugate-transpose of u(j).```
Parameters
 [in] BALANC ``` BALANC is CHARACTER*1 Specifies the balance option to be performed: = 'N': do not diagonally scale or permute; = 'P': permute only; = 'S': scale only; = 'B': both permute and scale. Computed reciprocal condition numbers will be for the matrices after permuting and/or balancing. Permuting does not change condition numbers (in exact arithmetic), but balancing does.``` [in] JOBVL ``` JOBVL is CHARACTER*1 = 'N': do not compute the left generalized eigenvectors; = 'V': compute the left generalized eigenvectors.``` [in] JOBVR ``` JOBVR is CHARACTER*1 = 'N': do not compute the right generalized eigenvectors; = 'V': compute the right generalized eigenvectors.``` [in] SENSE ``` SENSE is CHARACTER*1 Determines which reciprocal condition numbers are computed. = 'N': none are computed; = 'E': computed for eigenvalues only; = 'V': computed for eigenvectors only; = 'B': computed for eigenvalues and eigenvectors.``` [in] N ``` N is INTEGER The order of the matrices A, B, VL, and VR. N >= 0.``` [in,out] A ``` A is COMPLEX array, dimension (LDA, N) On entry, the matrix A in the pair (A,B). On exit, A has been overwritten. If JOBVL='V' or JOBVR='V' or both, then A contains the first part of the complex Schur form of the "balanced" versions of the input A and B.``` [in] LDA ``` LDA is INTEGER The leading dimension of A. LDA >= max(1,N).``` [in,out] B ``` B is COMPLEX array, dimension (LDB, N) On entry, the matrix B in the pair (A,B). On exit, B has been overwritten. If JOBVL='V' or JOBVR='V' or both, then B contains the second part of the complex Schur form of the "balanced" versions of the input A and B.``` [in] LDB ``` LDB is INTEGER The leading dimension of B. LDB >= max(1,N).``` [out] ALPHA ` ALPHA is COMPLEX array, dimension (N)` [out] BETA ``` BETA is COMPLEX array, dimension (N) On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the generalized eigenvalues. Note: the quotient ALPHA(j)/BETA(j) ) may easily over- or underflow, and BETA(j) may even be zero. Thus, the user should avoid naively computing the ratio ALPHA/BETA. However, ALPHA will be always less than and usually comparable with norm(A) in magnitude, and BETA always less than and usually comparable with norm(B).``` [out] VL ``` VL is COMPLEX array, dimension (LDVL,N) If JOBVL = 'V', the left generalized eigenvectors u(j) are stored one after another in the columns of VL, in the same order as their eigenvalues. Each eigenvector will be scaled so the largest component will have abs(real part) + abs(imag. part) = 1. Not referenced if JOBVL = 'N'.``` [in] LDVL ``` LDVL is INTEGER The leading dimension of the matrix VL. LDVL >= 1, and if JOBVL = 'V', LDVL >= N.``` [out] VR ``` VR is COMPLEX array, dimension (LDVR,N) If JOBVR = 'V', the right generalized eigenvectors v(j) are stored one after another in the columns of VR, in the same order as their eigenvalues. Each eigenvector will be scaled so the largest component will have abs(real part) + abs(imag. part) = 1. Not referenced if JOBVR = 'N'.``` [in] LDVR ``` LDVR is INTEGER The leading dimension of the matrix VR. LDVR >= 1, and if JOBVR = 'V', LDVR >= N.``` [out] ILO ` ILO is INTEGER` [out] IHI ``` IHI is INTEGER ILO and IHI are integer values such that on exit A(i,j) = 0 and B(i,j) = 0 if i > j and j = 1,...,ILO-1 or i = IHI+1,...,N. If BALANC = 'N' or 'S', ILO = 1 and IHI = N.``` [out] LSCALE ``` LSCALE is REAL array, dimension (N) Details of the permutations and scaling factors applied to the left side of A and B. If PL(j) is the index of the row interchanged with row j, and DL(j) is the scaling factor applied to row j, then LSCALE(j) = PL(j) for j = 1,...,ILO-1 = DL(j) for j = ILO,...,IHI = PL(j) for j = IHI+1,...,N. The order in which the interchanges are made is N to IHI+1, then 1 to ILO-1.``` [out] RSCALE ``` RSCALE is REAL array, dimension (N) Details of the permutations and scaling factors applied to the right side of A and B. If PR(j) is the index of the column interchanged with column j, and DR(j) is the scaling factor applied to column j, then RSCALE(j) = PR(j) for j = 1,...,ILO-1 = DR(j) for j = ILO,...,IHI = PR(j) for j = IHI+1,...,N The order in which the interchanges are made is N to IHI+1, then 1 to ILO-1.``` [out] ABNRM ``` ABNRM is REAL The one-norm of the balanced matrix A.``` [out] BBNRM ``` BBNRM is REAL The one-norm of the balanced matrix B.``` [out] RCONDE ``` RCONDE is REAL array, dimension (N) If SENSE = 'E' or 'B', the reciprocal condition numbers of the eigenvalues, stored in consecutive elements of the array. If SENSE = 'N' or 'V', RCONDE is not referenced.``` [out] RCONDV ``` RCONDV is REAL array, dimension (N) If SENSE = 'V' or 'B', the estimated reciprocal condition numbers of the eigenvectors, stored in consecutive elements of the array. If the eigenvalues cannot be reordered to compute RCONDV(j), RCONDV(j) is set to 0; this can only occur when the true value would be very small anyway. If SENSE = 'N' or 'E', RCONDV is not referenced.``` [out] WORK ``` WORK is COMPLEX array, dimension (MAX(1,LWORK)) On exit, if INFO = 0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= max(1,2*N). If SENSE = 'E', LWORK >= max(1,4*N). If SENSE = 'V' or 'B', LWORK >= max(1,2*N*N+2*N). If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] RWORK ``` RWORK is REAL array, dimension (lrwork) lrwork must be at least max(1,6*N) if BALANC = 'S' or 'B', and at least max(1,2*N) otherwise. Real workspace.``` [out] IWORK ``` IWORK is INTEGER array, dimension (N+2) If SENSE = 'E', IWORK is not referenced.``` [out] BWORK ``` BWORK is LOGICAL array, dimension (N) If SENSE = 'N', BWORK is not referenced.``` [out] INFO ``` INFO is INTEGER = 0: successful exit < 0: if INFO = -i, the i-th argument had an illegal value. = 1,...,N: The QZ iteration failed. No eigenvectors have been calculated, but ALPHA(j) and BETA(j) should be correct for j=INFO+1,...,N. > N: =N+1: other than QZ iteration failed in CHGEQZ. =N+2: error return from CTGEVC.```
Date
April 2012
Further Details:
```  Balancing a matrix pair (A,B) includes, first, permuting rows and
columns to isolate eigenvalues, second, applying diagonal similarity
transformation to the rows and columns to make the rows and columns
as close in norm as possible. The computed reciprocal condition
numbers correspond to the balanced matrix. Permuting rows and columns
will not change the condition numbers (in exact arithmetic) but
diagonal scaling will.  For further explanation of balancing, see
section 4.11.1.2 of LAPACK Users' Guide.

An approximate error bound on the chordal distance between the i-th
computed generalized eigenvalue w and the corresponding exact
eigenvalue lambda is

chord(w, lambda) <= EPS * norm(ABNRM, BBNRM) / RCONDE(I)

An approximate error bound for the angle between the i-th computed
eigenvector VL(i) or VR(i) is given by

EPS * norm(ABNRM, BBNRM) / DIF(i).

For further explanation of the reciprocal condition numbers RCONDE
and RCONDV, see section 4.11 of LAPACK User's Guide.```

Definition at line 376 of file cggevx.f.

376 *
377 * -- LAPACK driver routine (version 3.4.1) --
378 * -- LAPACK is a software package provided by Univ. of Tennessee, --
379 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
380 * April 2012
381 *
382 * .. Scalar Arguments ..
383  CHARACTER balanc, jobvl, jobvr, sense
384  INTEGER ihi, ilo, info, lda, ldb, ldvl, ldvr, lwork, n
385  REAL abnrm, bbnrm
386 * ..
387 * .. Array Arguments ..
388  LOGICAL bwork( * )
389  INTEGER iwork( * )
390  REAL lscale( * ), rconde( * ), rcondv( * ),
391  \$ rscale( * ), rwork( * )
392  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
393  \$ beta( * ), vl( ldvl, * ), vr( ldvr, * ),
394  \$ work( * )
395 * ..
396 *
397 * =====================================================================
398 *
399 * .. Parameters ..
400  REAL zero, one
401  parameter( zero = 0.0e+0, one = 1.0e+0 )
402  COMPLEX czero, cone
403  parameter( czero = ( 0.0e+0, 0.0e+0 ),
404  \$ cone = ( 1.0e+0, 0.0e+0 ) )
405 * ..
406 * .. Local Scalars ..
407  LOGICAL ilascl, ilbscl, ilv, ilvl, ilvr, lquery, noscl,
408  \$ wantsb, wantse, wantsn, wantsv
409  CHARACTER chtemp
410  INTEGER i, icols, ierr, ijobvl, ijobvr, in, irows,
411  \$ itau, iwrk, iwrk1, j, jc, jr, m, maxwrk, minwrk
412  REAL anrm, anrmto, bignum, bnrm, bnrmto, eps,
413  \$ smlnum, temp
414  COMPLEX x
415 * ..
416 * .. Local Arrays ..
417  LOGICAL ldumma( 1 )
418 * ..
419 * .. External Subroutines ..
420  EXTERNAL cgeqrf, cggbak, cggbal, cgghrd, chgeqz, clacpy,
423 * ..
424 * .. External Functions ..
425  LOGICAL lsame
426  INTEGER ilaenv
427  REAL clange, slamch
428  EXTERNAL lsame, ilaenv, clange, slamch
429 * ..
430 * .. Intrinsic Functions ..
431  INTRINSIC abs, aimag, max, REAL, sqrt
432 * ..
433 * .. Statement Functions ..
434  REAL abs1
435 * ..
436 * .. Statement Function definitions ..
437  abs1( x ) = abs( REAL( X ) ) + abs( aimag( x ) )
438 * ..
439 * .. Executable Statements ..
440 *
441 * Decode the input arguments
442 *
443  IF( lsame( jobvl, 'N' ) ) THEN
444  ijobvl = 1
445  ilvl = .false.
446  ELSE IF( lsame( jobvl, 'V' ) ) THEN
447  ijobvl = 2
448  ilvl = .true.
449  ELSE
450  ijobvl = -1
451  ilvl = .false.
452  END IF
453 *
454  IF( lsame( jobvr, 'N' ) ) THEN
455  ijobvr = 1
456  ilvr = .false.
457  ELSE IF( lsame( jobvr, 'V' ) ) THEN
458  ijobvr = 2
459  ilvr = .true.
460  ELSE
461  ijobvr = -1
462  ilvr = .false.
463  END IF
464  ilv = ilvl .OR. ilvr
465 *
466  noscl = lsame( balanc, 'N' ) .OR. lsame( balanc, 'P' )
467  wantsn = lsame( sense, 'N' )
468  wantse = lsame( sense, 'E' )
469  wantsv = lsame( sense, 'V' )
470  wantsb = lsame( sense, 'B' )
471 *
472 * Test the input arguments
473 *
474  info = 0
475  lquery = ( lwork.EQ.-1 )
476  IF( .NOT.( noscl .OR. lsame( balanc,'S' ) .OR.
477  \$ lsame( balanc, 'B' ) ) ) THEN
478  info = -1
479  ELSE IF( ijobvl.LE.0 ) THEN
480  info = -2
481  ELSE IF( ijobvr.LE.0 ) THEN
482  info = -3
483  ELSE IF( .NOT.( wantsn .OR. wantse .OR. wantsb .OR. wantsv ) )
484  \$ THEN
485  info = -4
486  ELSE IF( n.LT.0 ) THEN
487  info = -5
488  ELSE IF( lda.LT.max( 1, n ) ) THEN
489  info = -7
490  ELSE IF( ldb.LT.max( 1, n ) ) THEN
491  info = -9
492  ELSE IF( ldvl.LT.1 .OR. ( ilvl .AND. ldvl.LT.n ) ) THEN
493  info = -13
494  ELSE IF( ldvr.LT.1 .OR. ( ilvr .AND. ldvr.LT.n ) ) THEN
495  info = -15
496  END IF
497 *
498 * Compute workspace
499 * (Note: Comments in the code beginning "Workspace:" describe the
500 * minimal amount of workspace needed at that point in the code,
501 * as well as the preferred amount for good performance.
502 * NB refers to the optimal block size for the immediately
503 * following subroutine, as returned by ILAENV. The workspace is
504 * computed assuming ILO = 1 and IHI = N, the worst case.)
505 *
506  IF( info.EQ.0 ) THEN
507  IF( n.EQ.0 ) THEN
508  minwrk = 1
509  maxwrk = 1
510  ELSE
511  minwrk = 2*n
512  IF( wantse ) THEN
513  minwrk = 4*n
514  ELSE IF( wantsv .OR. wantsb ) THEN
515  minwrk = 2*n*( n + 1)
516  END IF
517  maxwrk = minwrk
518  maxwrk = max( maxwrk,
519  \$ n + n*ilaenv( 1, 'CGEQRF', ' ', n, 1, n, 0 ) )
520  maxwrk = max( maxwrk,
521  \$ n + n*ilaenv( 1, 'CUNMQR', ' ', n, 1, n, 0 ) )
522  IF( ilvl ) THEN
523  maxwrk = max( maxwrk, n +
524  \$ n*ilaenv( 1, 'CUNGQR', ' ', n, 1, n, 0 ) )
525  END IF
526  END IF
527  work( 1 ) = maxwrk
528 *
529  IF( lwork.LT.minwrk .AND. .NOT.lquery ) THEN
530  info = -25
531  END IF
532  END IF
533 *
534  IF( info.NE.0 ) THEN
535  CALL xerbla( 'CGGEVX', -info )
536  RETURN
537  ELSE IF( lquery ) THEN
538  RETURN
539  END IF
540 *
541 * Quick return if possible
542 *
543  IF( n.EQ.0 )
544  \$ RETURN
545 *
546 * Get machine constants
547 *
548  eps = slamch( 'P' )
549  smlnum = slamch( 'S' )
550  bignum = one / smlnum
551  CALL slabad( smlnum, bignum )
552  smlnum = sqrt( smlnum ) / eps
553  bignum = one / smlnum
554 *
555 * Scale A if max element outside range [SMLNUM,BIGNUM]
556 *
557  anrm = clange( 'M', n, n, a, lda, rwork )
558  ilascl = .false.
559  IF( anrm.GT.zero .AND. anrm.LT.smlnum ) THEN
560  anrmto = smlnum
561  ilascl = .true.
562  ELSE IF( anrm.GT.bignum ) THEN
563  anrmto = bignum
564  ilascl = .true.
565  END IF
566  IF( ilascl )
567  \$ CALL clascl( 'G', 0, 0, anrm, anrmto, n, n, a, lda, ierr )
568 *
569 * Scale B if max element outside range [SMLNUM,BIGNUM]
570 *
571  bnrm = clange( 'M', n, n, b, ldb, rwork )
572  ilbscl = .false.
573  IF( bnrm.GT.zero .AND. bnrm.LT.smlnum ) THEN
574  bnrmto = smlnum
575  ilbscl = .true.
576  ELSE IF( bnrm.GT.bignum ) THEN
577  bnrmto = bignum
578  ilbscl = .true.
579  END IF
580  IF( ilbscl )
581  \$ CALL clascl( 'G', 0, 0, bnrm, bnrmto, n, n, b, ldb, ierr )
582 *
583 * Permute and/or balance the matrix pair (A,B)
584 * (Real Workspace: need 6*N if BALANC = 'S' or 'B', 1 otherwise)
585 *
586  CALL cggbal( balanc, n, a, lda, b, ldb, ilo, ihi, lscale, rscale,
587  \$ rwork, ierr )
588 *
589 * Compute ABNRM and BBNRM
590 *
591  abnrm = clange( '1', n, n, a, lda, rwork( 1 ) )
592  IF( ilascl ) THEN
593  rwork( 1 ) = abnrm
594  CALL slascl( 'G', 0, 0, anrmto, anrm, 1, 1, rwork( 1 ), 1,
595  \$ ierr )
596  abnrm = rwork( 1 )
597  END IF
598 *
599  bbnrm = clange( '1', n, n, b, ldb, rwork( 1 ) )
600  IF( ilbscl ) THEN
601  rwork( 1 ) = bbnrm
602  CALL slascl( 'G', 0, 0, bnrmto, bnrm, 1, 1, rwork( 1 ), 1,
603  \$ ierr )
604  bbnrm = rwork( 1 )
605  END IF
606 *
607 * Reduce B to triangular form (QR decomposition of B)
608 * (Complex Workspace: need N, prefer N*NB )
609 *
610  irows = ihi + 1 - ilo
611  IF( ilv .OR. .NOT.wantsn ) THEN
612  icols = n + 1 - ilo
613  ELSE
614  icols = irows
615  END IF
616  itau = 1
617  iwrk = itau + irows
618  CALL cgeqrf( irows, icols, b( ilo, ilo ), ldb, work( itau ),
619  \$ work( iwrk ), lwork+1-iwrk, ierr )
620 *
621 * Apply the unitary transformation to A
622 * (Complex Workspace: need N, prefer N*NB)
623 *
624  CALL cunmqr( 'L', 'C', irows, icols, irows, b( ilo, ilo ), ldb,
625  \$ work( itau ), a( ilo, ilo ), lda, work( iwrk ),
626  \$ lwork+1-iwrk, ierr )
627 *
628 * Initialize VL and/or VR
629 * (Workspace: need N, prefer N*NB)
630 *
631  IF( ilvl ) THEN
632  CALL claset( 'Full', n, n, czero, cone, vl, ldvl )
633  IF( irows.GT.1 ) THEN
634  CALL clacpy( 'L', irows-1, irows-1, b( ilo+1, ilo ), ldb,
635  \$ vl( ilo+1, ilo ), ldvl )
636  END IF
637  CALL cungqr( irows, irows, irows, vl( ilo, ilo ), ldvl,
638  \$ work( itau ), work( iwrk ), lwork+1-iwrk, ierr )
639  END IF
640 *
641  IF( ilvr )
642  \$ CALL claset( 'Full', n, n, czero, cone, vr, ldvr )
643 *
644 * Reduce to generalized Hessenberg form
645 * (Workspace: none needed)
646 *
647  IF( ilv .OR. .NOT.wantsn ) THEN
648 *
649 * Eigenvectors requested -- work on whole matrix.
650 *
651  CALL cgghrd( jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb, vl,
652  \$ ldvl, vr, ldvr, ierr )
653  ELSE
654  CALL cgghrd( 'N', 'N', irows, 1, irows, a( ilo, ilo ), lda,
655  \$ b( ilo, ilo ), ldb, vl, ldvl, vr, ldvr, ierr )
656  END IF
657 *
658 * Perform QZ algorithm (Compute eigenvalues, and optionally, the
659 * Schur forms and Schur vectors)
660 * (Complex Workspace: need N)
661 * (Real Workspace: need N)
662 *
663  iwrk = itau
664  IF( ilv .OR. .NOT.wantsn ) THEN
665  chtemp = 'S'
666  ELSE
667  chtemp = 'E'
668  END IF
669 *
670  CALL chgeqz( chtemp, jobvl, jobvr, n, ilo, ihi, a, lda, b, ldb,
671  \$ alpha, beta, vl, ldvl, vr, ldvr, work( iwrk ),
672  \$ lwork+1-iwrk, rwork, ierr )
673  IF( ierr.NE.0 ) THEN
674  IF( ierr.GT.0 .AND. ierr.LE.n ) THEN
675  info = ierr
676  ELSE IF( ierr.GT.n .AND. ierr.LE.2*n ) THEN
677  info = ierr - n
678  ELSE
679  info = n + 1
680  END IF
681  GO TO 90
682  END IF
683 *
684 * Compute Eigenvectors and estimate condition numbers if desired
685 * CTGEVC: (Complex Workspace: need 2*N )
686 * (Real Workspace: need 2*N )
687 * CTGSNA: (Complex Workspace: need 2*N*N if SENSE='V' or 'B')
688 * (Integer Workspace: need N+2 )
689 *
690  IF( ilv .OR. .NOT.wantsn ) THEN
691  IF( ilv ) THEN
692  IF( ilvl ) THEN
693  IF( ilvr ) THEN
694  chtemp = 'B'
695  ELSE
696  chtemp = 'L'
697  END IF
698  ELSE
699  chtemp = 'R'
700  END IF
701 *
702  CALL ctgevc( chtemp, 'B', ldumma, n, a, lda, b, ldb, vl,
703  \$ ldvl, vr, ldvr, n, in, work( iwrk ), rwork,
704  \$ ierr )
705  IF( ierr.NE.0 ) THEN
706  info = n + 2
707  GO TO 90
708  END IF
709  END IF
710 *
711  IF( .NOT.wantsn ) THEN
712 *
713 * compute eigenvectors (STGEVC) and estimate condition
714 * numbers (STGSNA). Note that the definition of the condition
715 * number is not invariant under transformation (u,v) to
716 * (Q*u, Z*v), where (u,v) are eigenvectors of the generalized
717 * Schur form (S,T), Q and Z are orthogonal matrices. In order
718 * to avoid using extra 2*N*N workspace, we have to
719 * re-calculate eigenvectors and estimate the condition numbers
720 * one at a time.
721 *
722  DO 20 i = 1, n
723 *
724  DO 10 j = 1, n
725  bwork( j ) = .false.
726  10 CONTINUE
727  bwork( i ) = .true.
728 *
729  iwrk = n + 1
730  iwrk1 = iwrk + n
731 *
732  IF( wantse .OR. wantsb ) THEN
733  CALL ctgevc( 'B', 'S', bwork, n, a, lda, b, ldb,
734  \$ work( 1 ), n, work( iwrk ), n, 1, m,
735  \$ work( iwrk1 ), rwork, ierr )
736  IF( ierr.NE.0 ) THEN
737  info = n + 2
738  GO TO 90
739  END IF
740  END IF
741 *
742  CALL ctgsna( sense, 'S', bwork, n, a, lda, b, ldb,
743  \$ work( 1 ), n, work( iwrk ), n, rconde( i ),
744  \$ rcondv( i ), 1, m, work( iwrk1 ),
745  \$ lwork-iwrk1+1, iwork, ierr )
746 *
747  20 CONTINUE
748  END IF
749  END IF
750 *
751 * Undo balancing on VL and VR and normalization
752 * (Workspace: none needed)
753 *
754  IF( ilvl ) THEN
755  CALL cggbak( balanc, 'L', n, ilo, ihi, lscale, rscale, n, vl,
756  \$ ldvl, ierr )
757 *
758  DO 50 jc = 1, n
759  temp = zero
760  DO 30 jr = 1, n
761  temp = max( temp, abs1( vl( jr, jc ) ) )
762  30 CONTINUE
763  IF( temp.LT.smlnum )
764  \$ GO TO 50
765  temp = one / temp
766  DO 40 jr = 1, n
767  vl( jr, jc ) = vl( jr, jc )*temp
768  40 CONTINUE
769  50 CONTINUE
770  END IF
771 *
772  IF( ilvr ) THEN
773  CALL cggbak( balanc, 'R', n, ilo, ihi, lscale, rscale, n, vr,
774  \$ ldvr, ierr )
775  DO 80 jc = 1, n
776  temp = zero
777  DO 60 jr = 1, n
778  temp = max( temp, abs1( vr( jr, jc ) ) )
779  60 CONTINUE
780  IF( temp.LT.smlnum )
781  \$ GO TO 80
782  temp = one / temp
783  DO 70 jr = 1, n
784  vr( jr, jc ) = vr( jr, jc )*temp
785  70 CONTINUE
786  80 CONTINUE
787  END IF
788 *
789 * Undo scaling if necessary
790 *
791  90 CONTINUE
792 *
793  IF( ilascl )
794  \$ CALL clascl( 'G', 0, 0, anrmto, anrm, n, 1, alpha, n, ierr )
795 *
796  IF( ilbscl )
797  \$ CALL clascl( 'G', 0, 0, bnrmto, bnrm, n, 1, beta, n, ierr )
798 *
799  work( 1 ) = maxwrk
800  RETURN
801 *
802 * End of CGGEVX
803 *
subroutine claset(UPLO, M, N, ALPHA, BETA, A, LDA)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values...
Definition: claset.f:108
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:69
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine ctgsna(JOB, HOWMNY, SELECT, N, A, LDA, B, LDB, VL, LDVL, VR, LDVR, S, DIF, MM, M, WORK, LWORK, IWORK, INFO)
CTGSNA
Definition: ctgsna.f:313
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:170
subroutine ctgevc(SIDE, HOWMNY, SELECT, N, S, LDS, P, LDP, VL, LDVL, VR, LDVR, MM, M, WORK, RWORK, INFO)
CTGEVC
Definition: ctgevc.f:221
logical function lsame(CA, CB)
LSAME
Definition: lsame.f:55
subroutine cggbal(JOB, N, A, LDA, B, LDB, ILO, IHI, LSCALE, RSCALE, WORK, INFO)
CGGBAL
Definition: cggbal.f:179
subroutine clacpy(UPLO, M, N, A, LDA, B, LDB)
CLACPY copies all or part of one two-dimensional array to another.
Definition: clacpy.f:105
subroutine cggbak(JOB, SIDE, N, ILO, IHI, LSCALE, RSCALE, M, V, LDV, INFO)
CGGBAK
Definition: cggbak.f:150
subroutine cgghrd(COMPQ, COMPZ, N, ILO, IHI, A, LDA, B, LDB, Q, LDQ, Z, LDZ, INFO)
CGGHRD
Definition: cgghrd.f:206
subroutine clascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
CLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: clascl.f:141
subroutine chgeqz(JOB, COMPQ, COMPZ, N, ILO, IHI, H, LDH, T, LDT, ALPHA, BETA, Q, LDQ, Z, LDZ, WORK, LWORK, RWORK, INFO)
CHGEQZ
Definition: chgeqz.f:286
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine slascl(TYPE, KL, KU, CFROM, CTO, M, N, A, LDA, INFO)
SLASCL multiplies a general rectangular matrix by a real scalar defined as cto/cfrom.
Definition: slascl.f:141
real function clange(NORM, M, N, A, LDA, WORK)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition: clange.f:117
subroutine cungqr(M, N, K, A, LDA, TAU, WORK, LWORK, INFO)
CUNGQR
Definition: cungqr.f:130

Here is the call graph for this function:

Here is the caller graph for this function: