![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief <b> ZGGES computes the eigenvalues, the Schur form, and, optionally, the matrix of Schur vectors for GE matrices</b> 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download ZGGES + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zgges.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zgges.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zgges.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE ZGGES( JOBVSL, JOBVSR, SORT, SELCTG, N, A, LDA, B, LDB, 00022 * SDIM, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, WORK, 00023 * LWORK, RWORK, BWORK, INFO ) 00024 * 00025 * .. Scalar Arguments .. 00026 * CHARACTER JOBVSL, JOBVSR, SORT 00027 * INTEGER INFO, LDA, LDB, LDVSL, LDVSR, LWORK, N, SDIM 00028 * .. 00029 * .. Array Arguments .. 00030 * LOGICAL BWORK( * ) 00031 * DOUBLE PRECISION RWORK( * ) 00032 * COMPLEX*16 A( LDA, * ), ALPHA( * ), B( LDB, * ), 00033 * $ BETA( * ), VSL( LDVSL, * ), VSR( LDVSR, * ), 00034 * $ WORK( * ) 00035 * .. 00036 * .. Function Arguments .. 00037 * LOGICAL SELCTG 00038 * EXTERNAL SELCTG 00039 * .. 00040 * 00041 * 00042 *> \par Purpose: 00043 * ============= 00044 *> 00045 *> \verbatim 00046 *> 00047 *> ZGGES computes for a pair of N-by-N complex nonsymmetric matrices 00048 *> (A,B), the generalized eigenvalues, the generalized complex Schur 00049 *> form (S, T), and optionally left and/or right Schur vectors (VSL 00050 *> and VSR). This gives the generalized Schur factorization 00051 *> 00052 *> (A,B) = ( (VSL)*S*(VSR)**H, (VSL)*T*(VSR)**H ) 00053 *> 00054 *> where (VSR)**H is the conjugate-transpose of VSR. 00055 *> 00056 *> Optionally, it also orders the eigenvalues so that a selected cluster 00057 *> of eigenvalues appears in the leading diagonal blocks of the upper 00058 *> triangular matrix S and the upper triangular matrix T. The leading 00059 *> columns of VSL and VSR then form an unitary basis for the 00060 *> corresponding left and right eigenspaces (deflating subspaces). 00061 *> 00062 *> (If only the generalized eigenvalues are needed, use the driver 00063 *> ZGGEV instead, which is faster.) 00064 *> 00065 *> A generalized eigenvalue for a pair of matrices (A,B) is a scalar w 00066 *> or a ratio alpha/beta = w, such that A - w*B is singular. It is 00067 *> usually represented as the pair (alpha,beta), as there is a 00068 *> reasonable interpretation for beta=0, and even for both being zero. 00069 *> 00070 *> A pair of matrices (S,T) is in generalized complex Schur form if S 00071 *> and T are upper triangular and, in addition, the diagonal elements 00072 *> of T are non-negative real numbers. 00073 *> \endverbatim 00074 * 00075 * Arguments: 00076 * ========== 00077 * 00078 *> \param[in] JOBVSL 00079 *> \verbatim 00080 *> JOBVSL is CHARACTER*1 00081 *> = 'N': do not compute the left Schur vectors; 00082 *> = 'V': compute the left Schur vectors. 00083 *> \endverbatim 00084 *> 00085 *> \param[in] JOBVSR 00086 *> \verbatim 00087 *> JOBVSR is CHARACTER*1 00088 *> = 'N': do not compute the right Schur vectors; 00089 *> = 'V': compute the right Schur vectors. 00090 *> \endverbatim 00091 *> 00092 *> \param[in] SORT 00093 *> \verbatim 00094 *> SORT is CHARACTER*1 00095 *> Specifies whether or not to order the eigenvalues on the 00096 *> diagonal of the generalized Schur form. 00097 *> = 'N': Eigenvalues are not ordered; 00098 *> = 'S': Eigenvalues are ordered (see SELCTG). 00099 *> \endverbatim 00100 *> 00101 *> \param[in] SELCTG 00102 *> \verbatim 00103 *> SELCTG is a LOGICAL FUNCTION of two COMPLEX*16 arguments 00104 *> SELCTG must be declared EXTERNAL in the calling subroutine. 00105 *> If SORT = 'N', SELCTG is not referenced. 00106 *> If SORT = 'S', SELCTG is used to select eigenvalues to sort 00107 *> to the top left of the Schur form. 00108 *> An eigenvalue ALPHA(j)/BETA(j) is selected if 00109 *> SELCTG(ALPHA(j),BETA(j)) is true. 00110 *> 00111 *> Note that a selected complex eigenvalue may no longer satisfy 00112 *> SELCTG(ALPHA(j),BETA(j)) = .TRUE. after ordering, since 00113 *> ordering may change the value of complex eigenvalues 00114 *> (especially if the eigenvalue is ill-conditioned), in this 00115 *> case INFO is set to N+2 (See INFO below). 00116 *> \endverbatim 00117 *> 00118 *> \param[in] N 00119 *> \verbatim 00120 *> N is INTEGER 00121 *> The order of the matrices A, B, VSL, and VSR. N >= 0. 00122 *> \endverbatim 00123 *> 00124 *> \param[in,out] A 00125 *> \verbatim 00126 *> A is COMPLEX*16 array, dimension (LDA, N) 00127 *> On entry, the first of the pair of matrices. 00128 *> On exit, A has been overwritten by its generalized Schur 00129 *> form S. 00130 *> \endverbatim 00131 *> 00132 *> \param[in] LDA 00133 *> \verbatim 00134 *> LDA is INTEGER 00135 *> The leading dimension of A. LDA >= max(1,N). 00136 *> \endverbatim 00137 *> 00138 *> \param[in,out] B 00139 *> \verbatim 00140 *> B is COMPLEX*16 array, dimension (LDB, N) 00141 *> On entry, the second of the pair of matrices. 00142 *> On exit, B has been overwritten by its generalized Schur 00143 *> form T. 00144 *> \endverbatim 00145 *> 00146 *> \param[in] LDB 00147 *> \verbatim 00148 *> LDB is INTEGER 00149 *> The leading dimension of B. LDB >= max(1,N). 00150 *> \endverbatim 00151 *> 00152 *> \param[out] SDIM 00153 *> \verbatim 00154 *> SDIM is INTEGER 00155 *> If SORT = 'N', SDIM = 0. 00156 *> If SORT = 'S', SDIM = number of eigenvalues (after sorting) 00157 *> for which SELCTG is true. 00158 *> \endverbatim 00159 *> 00160 *> \param[out] ALPHA 00161 *> \verbatim 00162 *> ALPHA is COMPLEX*16 array, dimension (N) 00163 *> \endverbatim 00164 *> 00165 *> \param[out] BETA 00166 *> \verbatim 00167 *> BETA is COMPLEX*16 array, dimension (N) 00168 *> On exit, ALPHA(j)/BETA(j), j=1,...,N, will be the 00169 *> generalized eigenvalues. ALPHA(j), j=1,...,N and BETA(j), 00170 *> j=1,...,N are the diagonals of the complex Schur form (A,B) 00171 *> output by ZGGES. The BETA(j) will be non-negative real. 00172 *> 00173 *> Note: the quotients ALPHA(j)/BETA(j) may easily over- or 00174 *> underflow, and BETA(j) may even be zero. Thus, the user 00175 *> should avoid naively computing the ratio alpha/beta. 00176 *> However, ALPHA will be always less than and usually 00177 *> comparable with norm(A) in magnitude, and BETA always less 00178 *> than and usually comparable with norm(B). 00179 *> \endverbatim 00180 *> 00181 *> \param[out] VSL 00182 *> \verbatim 00183 *> VSL is COMPLEX*16 array, dimension (LDVSL,N) 00184 *> If JOBVSL = 'V', VSL will contain the left Schur vectors. 00185 *> Not referenced if JOBVSL = 'N'. 00186 *> \endverbatim 00187 *> 00188 *> \param[in] LDVSL 00189 *> \verbatim 00190 *> LDVSL is INTEGER 00191 *> The leading dimension of the matrix VSL. LDVSL >= 1, and 00192 *> if JOBVSL = 'V', LDVSL >= N. 00193 *> \endverbatim 00194 *> 00195 *> \param[out] VSR 00196 *> \verbatim 00197 *> VSR is COMPLEX*16 array, dimension (LDVSR,N) 00198 *> If JOBVSR = 'V', VSR will contain the right Schur vectors. 00199 *> Not referenced if JOBVSR = 'N'. 00200 *> \endverbatim 00201 *> 00202 *> \param[in] LDVSR 00203 *> \verbatim 00204 *> LDVSR is INTEGER 00205 *> The leading dimension of the matrix VSR. LDVSR >= 1, and 00206 *> if JOBVSR = 'V', LDVSR >= N. 00207 *> \endverbatim 00208 *> 00209 *> \param[out] WORK 00210 *> \verbatim 00211 *> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK)) 00212 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK. 00213 *> \endverbatim 00214 *> 00215 *> \param[in] LWORK 00216 *> \verbatim 00217 *> LWORK is INTEGER 00218 *> The dimension of the array WORK. LWORK >= max(1,2*N). 00219 *> For good performance, LWORK must generally be larger. 00220 *> 00221 *> If LWORK = -1, then a workspace query is assumed; the routine 00222 *> only calculates the optimal size of the WORK array, returns 00223 *> this value as the first entry of the WORK array, and no error 00224 *> message related to LWORK is issued by XERBLA. 00225 *> \endverbatim 00226 *> 00227 *> \param[out] RWORK 00228 *> \verbatim 00229 *> RWORK is DOUBLE PRECISION array, dimension (8*N) 00230 *> \endverbatim 00231 *> 00232 *> \param[out] BWORK 00233 *> \verbatim 00234 *> BWORK is LOGICAL array, dimension (N) 00235 *> Not referenced if SORT = 'N'. 00236 *> \endverbatim 00237 *> 00238 *> \param[out] INFO 00239 *> \verbatim 00240 *> INFO is INTEGER 00241 *> = 0: successful exit 00242 *> < 0: if INFO = -i, the i-th argument had an illegal value. 00243 *> =1,...,N: 00244 *> The QZ iteration failed. (A,B) are not in Schur 00245 *> form, but ALPHA(j) and BETA(j) should be correct for 00246 *> j=INFO+1,...,N. 00247 *> > N: =N+1: other than QZ iteration failed in ZHGEQZ 00248 *> =N+2: after reordering, roundoff changed values of 00249 *> some complex eigenvalues so that leading 00250 *> eigenvalues in the Generalized Schur form no 00251 *> longer satisfy SELCTG=.TRUE. This could also 00252 *> be caused due to scaling. 00253 *> =N+3: reordering falied in ZTGSEN. 00254 *> \endverbatim 00255 * 00256 * Authors: 00257 * ======== 00258 * 00259 *> \author Univ. of Tennessee 00260 *> \author Univ. of California Berkeley 00261 *> \author Univ. of Colorado Denver 00262 *> \author NAG Ltd. 00263 * 00264 *> \date November 2011 00265 * 00266 *> \ingroup complex16GEeigen 00267 * 00268 * ===================================================================== 00269 SUBROUTINE ZGGES( JOBVSL, JOBVSR, SORT, SELCTG, N, A, LDA, B, LDB, 00270 $ SDIM, ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, WORK, 00271 $ LWORK, RWORK, BWORK, INFO ) 00272 * 00273 * -- LAPACK driver routine (version 3.4.0) -- 00274 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00275 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00276 * November 2011 00277 * 00278 * .. Scalar Arguments .. 00279 CHARACTER JOBVSL, JOBVSR, SORT 00280 INTEGER INFO, LDA, LDB, LDVSL, LDVSR, LWORK, N, SDIM 00281 * .. 00282 * .. Array Arguments .. 00283 LOGICAL BWORK( * ) 00284 DOUBLE PRECISION RWORK( * ) 00285 COMPLEX*16 A( LDA, * ), ALPHA( * ), B( LDB, * ), 00286 $ BETA( * ), VSL( LDVSL, * ), VSR( LDVSR, * ), 00287 $ WORK( * ) 00288 * .. 00289 * .. Function Arguments .. 00290 LOGICAL SELCTG 00291 EXTERNAL SELCTG 00292 * .. 00293 * 00294 * ===================================================================== 00295 * 00296 * .. Parameters .. 00297 DOUBLE PRECISION ZERO, ONE 00298 PARAMETER ( ZERO = 0.0D0, ONE = 1.0D0 ) 00299 COMPLEX*16 CZERO, CONE 00300 PARAMETER ( CZERO = ( 0.0D0, 0.0D0 ), 00301 $ CONE = ( 1.0D0, 0.0D0 ) ) 00302 * .. 00303 * .. Local Scalars .. 00304 LOGICAL CURSL, ILASCL, ILBSCL, ILVSL, ILVSR, LASTSL, 00305 $ LQUERY, WANTST 00306 INTEGER I, ICOLS, IERR, IHI, IJOBVL, IJOBVR, ILEFT, 00307 $ ILO, IRIGHT, IROWS, IRWRK, ITAU, IWRK, LWKMIN, 00308 $ LWKOPT 00309 DOUBLE PRECISION ANRM, ANRMTO, BIGNUM, BNRM, BNRMTO, EPS, PVSL, 00310 $ PVSR, SMLNUM 00311 * .. 00312 * .. Local Arrays .. 00313 INTEGER IDUM( 1 ) 00314 DOUBLE PRECISION DIF( 2 ) 00315 * .. 00316 * .. External Subroutines .. 00317 EXTERNAL DLABAD, XERBLA, ZGEQRF, ZGGBAK, ZGGBAL, ZGGHRD, 00318 $ ZHGEQZ, ZLACPY, ZLASCL, ZLASET, ZTGSEN, ZUNGQR, 00319 $ ZUNMQR 00320 * .. 00321 * .. External Functions .. 00322 LOGICAL LSAME 00323 INTEGER ILAENV 00324 DOUBLE PRECISION DLAMCH, ZLANGE 00325 EXTERNAL LSAME, ILAENV, DLAMCH, ZLANGE 00326 * .. 00327 * .. Intrinsic Functions .. 00328 INTRINSIC MAX, SQRT 00329 * .. 00330 * .. Executable Statements .. 00331 * 00332 * Decode the input arguments 00333 * 00334 IF( LSAME( JOBVSL, 'N' ) ) THEN 00335 IJOBVL = 1 00336 ILVSL = .FALSE. 00337 ELSE IF( LSAME( JOBVSL, 'V' ) ) THEN 00338 IJOBVL = 2 00339 ILVSL = .TRUE. 00340 ELSE 00341 IJOBVL = -1 00342 ILVSL = .FALSE. 00343 END IF 00344 * 00345 IF( LSAME( JOBVSR, 'N' ) ) THEN 00346 IJOBVR = 1 00347 ILVSR = .FALSE. 00348 ELSE IF( LSAME( JOBVSR, 'V' ) ) THEN 00349 IJOBVR = 2 00350 ILVSR = .TRUE. 00351 ELSE 00352 IJOBVR = -1 00353 ILVSR = .FALSE. 00354 END IF 00355 * 00356 WANTST = LSAME( SORT, 'S' ) 00357 * 00358 * Test the input arguments 00359 * 00360 INFO = 0 00361 LQUERY = ( LWORK.EQ.-1 ) 00362 IF( IJOBVL.LE.0 ) THEN 00363 INFO = -1 00364 ELSE IF( IJOBVR.LE.0 ) THEN 00365 INFO = -2 00366 ELSE IF( ( .NOT.WANTST ) .AND. ( .NOT.LSAME( SORT, 'N' ) ) ) THEN 00367 INFO = -3 00368 ELSE IF( N.LT.0 ) THEN 00369 INFO = -5 00370 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 00371 INFO = -7 00372 ELSE IF( LDB.LT.MAX( 1, N ) ) THEN 00373 INFO = -9 00374 ELSE IF( LDVSL.LT.1 .OR. ( ILVSL .AND. LDVSL.LT.N ) ) THEN 00375 INFO = -14 00376 ELSE IF( LDVSR.LT.1 .OR. ( ILVSR .AND. LDVSR.LT.N ) ) THEN 00377 INFO = -16 00378 END IF 00379 * 00380 * Compute workspace 00381 * (Note: Comments in the code beginning "Workspace:" describe the 00382 * minimal amount of workspace needed at that point in the code, 00383 * as well as the preferred amount for good performance. 00384 * NB refers to the optimal block size for the immediately 00385 * following subroutine, as returned by ILAENV.) 00386 * 00387 IF( INFO.EQ.0 ) THEN 00388 LWKMIN = MAX( 1, 2*N ) 00389 LWKOPT = MAX( 1, N + N*ILAENV( 1, 'ZGEQRF', ' ', N, 1, N, 0 ) ) 00390 LWKOPT = MAX( LWKOPT, N + 00391 $ N*ILAENV( 1, 'ZUNMQR', ' ', N, 1, N, -1 ) ) 00392 IF( ILVSL ) THEN 00393 LWKOPT = MAX( LWKOPT, N + 00394 $ N*ILAENV( 1, 'ZUNGQR', ' ', N, 1, N, -1 ) ) 00395 END IF 00396 WORK( 1 ) = LWKOPT 00397 * 00398 IF( LWORK.LT.LWKMIN .AND. .NOT.LQUERY ) 00399 $ INFO = -18 00400 END IF 00401 * 00402 IF( INFO.NE.0 ) THEN 00403 CALL XERBLA( 'ZGGES ', -INFO ) 00404 RETURN 00405 ELSE IF( LQUERY ) THEN 00406 RETURN 00407 END IF 00408 * 00409 * Quick return if possible 00410 * 00411 IF( N.EQ.0 ) THEN 00412 SDIM = 0 00413 RETURN 00414 END IF 00415 * 00416 * Get machine constants 00417 * 00418 EPS = DLAMCH( 'P' ) 00419 SMLNUM = DLAMCH( 'S' ) 00420 BIGNUM = ONE / SMLNUM 00421 CALL DLABAD( SMLNUM, BIGNUM ) 00422 SMLNUM = SQRT( SMLNUM ) / EPS 00423 BIGNUM = ONE / SMLNUM 00424 * 00425 * Scale A if max element outside range [SMLNUM,BIGNUM] 00426 * 00427 ANRM = ZLANGE( 'M', N, N, A, LDA, RWORK ) 00428 ILASCL = .FALSE. 00429 IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN 00430 ANRMTO = SMLNUM 00431 ILASCL = .TRUE. 00432 ELSE IF( ANRM.GT.BIGNUM ) THEN 00433 ANRMTO = BIGNUM 00434 ILASCL = .TRUE. 00435 END IF 00436 * 00437 IF( ILASCL ) 00438 $ CALL ZLASCL( 'G', 0, 0, ANRM, ANRMTO, N, N, A, LDA, IERR ) 00439 * 00440 * Scale B if max element outside range [SMLNUM,BIGNUM] 00441 * 00442 BNRM = ZLANGE( 'M', N, N, B, LDB, RWORK ) 00443 ILBSCL = .FALSE. 00444 IF( BNRM.GT.ZERO .AND. BNRM.LT.SMLNUM ) THEN 00445 BNRMTO = SMLNUM 00446 ILBSCL = .TRUE. 00447 ELSE IF( BNRM.GT.BIGNUM ) THEN 00448 BNRMTO = BIGNUM 00449 ILBSCL = .TRUE. 00450 END IF 00451 * 00452 IF( ILBSCL ) 00453 $ CALL ZLASCL( 'G', 0, 0, BNRM, BNRMTO, N, N, B, LDB, IERR ) 00454 * 00455 * Permute the matrix to make it more nearly triangular 00456 * (Real Workspace: need 6*N) 00457 * 00458 ILEFT = 1 00459 IRIGHT = N + 1 00460 IRWRK = IRIGHT + N 00461 CALL ZGGBAL( 'P', N, A, LDA, B, LDB, ILO, IHI, RWORK( ILEFT ), 00462 $ RWORK( IRIGHT ), RWORK( IRWRK ), IERR ) 00463 * 00464 * Reduce B to triangular form (QR decomposition of B) 00465 * (Complex Workspace: need N, prefer N*NB) 00466 * 00467 IROWS = IHI + 1 - ILO 00468 ICOLS = N + 1 - ILO 00469 ITAU = 1 00470 IWRK = ITAU + IROWS 00471 CALL ZGEQRF( IROWS, ICOLS, B( ILO, ILO ), LDB, WORK( ITAU ), 00472 $ WORK( IWRK ), LWORK+1-IWRK, IERR ) 00473 * 00474 * Apply the orthogonal transformation to matrix A 00475 * (Complex Workspace: need N, prefer N*NB) 00476 * 00477 CALL ZUNMQR( 'L', 'C', IROWS, ICOLS, IROWS, B( ILO, ILO ), LDB, 00478 $ WORK( ITAU ), A( ILO, ILO ), LDA, WORK( IWRK ), 00479 $ LWORK+1-IWRK, IERR ) 00480 * 00481 * Initialize VSL 00482 * (Complex Workspace: need N, prefer N*NB) 00483 * 00484 IF( ILVSL ) THEN 00485 CALL ZLASET( 'Full', N, N, CZERO, CONE, VSL, LDVSL ) 00486 IF( IROWS.GT.1 ) THEN 00487 CALL ZLACPY( 'L', IROWS-1, IROWS-1, B( ILO+1, ILO ), LDB, 00488 $ VSL( ILO+1, ILO ), LDVSL ) 00489 END IF 00490 CALL ZUNGQR( IROWS, IROWS, IROWS, VSL( ILO, ILO ), LDVSL, 00491 $ WORK( ITAU ), WORK( IWRK ), LWORK+1-IWRK, IERR ) 00492 END IF 00493 * 00494 * Initialize VSR 00495 * 00496 IF( ILVSR ) 00497 $ CALL ZLASET( 'Full', N, N, CZERO, CONE, VSR, LDVSR ) 00498 * 00499 * Reduce to generalized Hessenberg form 00500 * (Workspace: none needed) 00501 * 00502 CALL ZGGHRD( JOBVSL, JOBVSR, N, ILO, IHI, A, LDA, B, LDB, VSL, 00503 $ LDVSL, VSR, LDVSR, IERR ) 00504 * 00505 SDIM = 0 00506 * 00507 * Perform QZ algorithm, computing Schur vectors if desired 00508 * (Complex Workspace: need N) 00509 * (Real Workspace: need N) 00510 * 00511 IWRK = ITAU 00512 CALL ZHGEQZ( 'S', JOBVSL, JOBVSR, N, ILO, IHI, A, LDA, B, LDB, 00513 $ ALPHA, BETA, VSL, LDVSL, VSR, LDVSR, WORK( IWRK ), 00514 $ LWORK+1-IWRK, RWORK( IRWRK ), IERR ) 00515 IF( IERR.NE.0 ) THEN 00516 IF( IERR.GT.0 .AND. IERR.LE.N ) THEN 00517 INFO = IERR 00518 ELSE IF( IERR.GT.N .AND. IERR.LE.2*N ) THEN 00519 INFO = IERR - N 00520 ELSE 00521 INFO = N + 1 00522 END IF 00523 GO TO 30 00524 END IF 00525 * 00526 * Sort eigenvalues ALPHA/BETA if desired 00527 * (Workspace: none needed) 00528 * 00529 IF( WANTST ) THEN 00530 * 00531 * Undo scaling on eigenvalues before selecting 00532 * 00533 IF( ILASCL ) 00534 $ CALL ZLASCL( 'G', 0, 0, ANRM, ANRMTO, N, 1, ALPHA, N, IERR ) 00535 IF( ILBSCL ) 00536 $ CALL ZLASCL( 'G', 0, 0, BNRM, BNRMTO, N, 1, BETA, N, IERR ) 00537 * 00538 * Select eigenvalues 00539 * 00540 DO 10 I = 1, N 00541 BWORK( I ) = SELCTG( ALPHA( I ), BETA( I ) ) 00542 10 CONTINUE 00543 * 00544 CALL ZTGSEN( 0, ILVSL, ILVSR, BWORK, N, A, LDA, B, LDB, ALPHA, 00545 $ BETA, VSL, LDVSL, VSR, LDVSR, SDIM, PVSL, PVSR, 00546 $ DIF, WORK( IWRK ), LWORK-IWRK+1, IDUM, 1, IERR ) 00547 IF( IERR.EQ.1 ) 00548 $ INFO = N + 3 00549 * 00550 END IF 00551 * 00552 * Apply back-permutation to VSL and VSR 00553 * (Workspace: none needed) 00554 * 00555 IF( ILVSL ) 00556 $ CALL ZGGBAK( 'P', 'L', N, ILO, IHI, RWORK( ILEFT ), 00557 $ RWORK( IRIGHT ), N, VSL, LDVSL, IERR ) 00558 IF( ILVSR ) 00559 $ CALL ZGGBAK( 'P', 'R', N, ILO, IHI, RWORK( ILEFT ), 00560 $ RWORK( IRIGHT ), N, VSR, LDVSR, IERR ) 00561 * 00562 * Undo scaling 00563 * 00564 IF( ILASCL ) THEN 00565 CALL ZLASCL( 'U', 0, 0, ANRMTO, ANRM, N, N, A, LDA, IERR ) 00566 CALL ZLASCL( 'G', 0, 0, ANRMTO, ANRM, N, 1, ALPHA, N, IERR ) 00567 END IF 00568 * 00569 IF( ILBSCL ) THEN 00570 CALL ZLASCL( 'U', 0, 0, BNRMTO, BNRM, N, N, B, LDB, IERR ) 00571 CALL ZLASCL( 'G', 0, 0, BNRMTO, BNRM, N, 1, BETA, N, IERR ) 00572 END IF 00573 * 00574 IF( WANTST ) THEN 00575 * 00576 * Check if reordering is correct 00577 * 00578 LASTSL = .TRUE. 00579 SDIM = 0 00580 DO 20 I = 1, N 00581 CURSL = SELCTG( ALPHA( I ), BETA( I ) ) 00582 IF( CURSL ) 00583 $ SDIM = SDIM + 1 00584 IF( CURSL .AND. .NOT.LASTSL ) 00585 $ INFO = N + 2 00586 LASTSL = CURSL 00587 20 CONTINUE 00588 * 00589 END IF 00590 * 00591 30 CONTINUE 00592 * 00593 WORK( 1 ) = LWKOPT 00594 * 00595 RETURN 00596 * 00597 * End of ZGGES 00598 * 00599 END