![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief <b> SGEESX 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 SGEESX + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/sgeesx.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/sgeesx.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/sgeesx.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE SGEESX( JOBVS, SORT, SELECT, SENSE, N, A, LDA, SDIM, 00022 * WR, WI, VS, LDVS, RCONDE, RCONDV, WORK, LWORK, 00023 * IWORK, LIWORK, BWORK, INFO ) 00024 * 00025 * .. Scalar Arguments .. 00026 * CHARACTER JOBVS, SENSE, SORT 00027 * INTEGER INFO, LDA, LDVS, LIWORK, LWORK, N, SDIM 00028 * REAL RCONDE, RCONDV 00029 * .. 00030 * .. Array Arguments .. 00031 * LOGICAL BWORK( * ) 00032 * INTEGER IWORK( * ) 00033 * REAL A( LDA, * ), VS( LDVS, * ), WI( * ), WORK( * ), 00034 * $ WR( * ) 00035 * .. 00036 * .. Function Arguments .. 00037 * LOGICAL SELECT 00038 * EXTERNAL SELECT 00039 * .. 00040 * 00041 * 00042 *> \par Purpose: 00043 * ============= 00044 *> 00045 *> \verbatim 00046 *> 00047 *> SGEESX computes for an N-by-N real nonsymmetric matrix A, the 00048 *> eigenvalues, the real Schur form T, and, optionally, the matrix of 00049 *> Schur vectors Z. This gives the Schur factorization A = Z*T*(Z**T). 00050 *> 00051 *> Optionally, it also orders the eigenvalues on the diagonal of the 00052 *> real Schur form so that selected eigenvalues are at the top left; 00053 *> computes a reciprocal condition number for the average of the 00054 *> selected eigenvalues (RCONDE); and computes a reciprocal condition 00055 *> number for the right invariant subspace corresponding to the 00056 *> selected eigenvalues (RCONDV). The leading columns of Z form an 00057 *> orthonormal basis for this invariant subspace. 00058 *> 00059 *> For further explanation of the reciprocal condition numbers RCONDE 00060 *> and RCONDV, see Section 4.10 of the LAPACK Users' Guide (where 00061 *> these quantities are called s and sep respectively). 00062 *> 00063 *> A real matrix is in real Schur form if it is upper quasi-triangular 00064 *> with 1-by-1 and 2-by-2 blocks. 2-by-2 blocks will be standardized in 00065 *> the form 00066 *> [ a b ] 00067 *> [ c a ] 00068 *> 00069 *> where b*c < 0. The eigenvalues of such a block are a +- sqrt(bc). 00070 *> \endverbatim 00071 * 00072 * Arguments: 00073 * ========== 00074 * 00075 *> \param[in] JOBVS 00076 *> \verbatim 00077 *> JOBVS is CHARACTER*1 00078 *> = 'N': Schur vectors are not computed; 00079 *> = 'V': Schur vectors are computed. 00080 *> \endverbatim 00081 *> 00082 *> \param[in] SORT 00083 *> \verbatim 00084 *> SORT is CHARACTER*1 00085 *> Specifies whether or not to order the eigenvalues on the 00086 *> diagonal of the Schur form. 00087 *> = 'N': Eigenvalues are not ordered; 00088 *> = 'S': Eigenvalues are ordered (see SELECT). 00089 *> \endverbatim 00090 *> 00091 *> \param[in] SELECT 00092 *> \verbatim 00093 *> SELECT is procedure) LOGICAL FUNCTION of two REAL arguments 00094 *> SELECT must be declared EXTERNAL in the calling subroutine. 00095 *> If SORT = 'S', SELECT is used to select eigenvalues to sort 00096 *> to the top left of the Schur form. 00097 *> If SORT = 'N', SELECT is not referenced. 00098 *> An eigenvalue WR(j)+sqrt(-1)*WI(j) is selected if 00099 *> SELECT(WR(j),WI(j)) is true; i.e., if either one of a 00100 *> complex conjugate pair of eigenvalues is selected, then both 00101 *> are. Note that a selected complex eigenvalue may no longer 00102 *> satisfy SELECT(WR(j),WI(j)) = .TRUE. after ordering, since 00103 *> ordering may change the value of complex eigenvalues 00104 *> (especially if the eigenvalue is ill-conditioned); in this 00105 *> case INFO may be set to N+3 (see INFO below). 00106 *> \endverbatim 00107 *> 00108 *> \param[in] SENSE 00109 *> \verbatim 00110 *> SENSE is CHARACTER*1 00111 *> Determines which reciprocal condition numbers are computed. 00112 *> = 'N': None are computed; 00113 *> = 'E': Computed for average of selected eigenvalues only; 00114 *> = 'V': Computed for selected right invariant subspace only; 00115 *> = 'B': Computed for both. 00116 *> If SENSE = 'E', 'V' or 'B', SORT must equal 'S'. 00117 *> \endverbatim 00118 *> 00119 *> \param[in] N 00120 *> \verbatim 00121 *> N is INTEGER 00122 *> The order of the matrix A. N >= 0. 00123 *> \endverbatim 00124 *> 00125 *> \param[in,out] A 00126 *> \verbatim 00127 *> A is REAL array, dimension (LDA, N) 00128 *> On entry, the N-by-N matrix A. 00129 *> On exit, A is overwritten by its real Schur form T. 00130 *> \endverbatim 00131 *> 00132 *> \param[in] LDA 00133 *> \verbatim 00134 *> LDA is INTEGER 00135 *> The leading dimension of the array A. LDA >= max(1,N). 00136 *> \endverbatim 00137 *> 00138 *> \param[out] SDIM 00139 *> \verbatim 00140 *> SDIM is INTEGER 00141 *> If SORT = 'N', SDIM = 0. 00142 *> If SORT = 'S', SDIM = number of eigenvalues (after sorting) 00143 *> for which SELECT is true. (Complex conjugate 00144 *> pairs for which SELECT is true for either 00145 *> eigenvalue count as 2.) 00146 *> \endverbatim 00147 *> 00148 *> \param[out] WR 00149 *> \verbatim 00150 *> WR is REAL array, dimension (N) 00151 *> \endverbatim 00152 *> 00153 *> \param[out] WI 00154 *> \verbatim 00155 *> WI is REAL array, dimension (N) 00156 *> WR and WI contain the real and imaginary parts, respectively, 00157 *> of the computed eigenvalues, in the same order that they 00158 *> appear on the diagonal of the output Schur form T. Complex 00159 *> conjugate pairs of eigenvalues appear consecutively with the 00160 *> eigenvalue having the positive imaginary part first. 00161 *> \endverbatim 00162 *> 00163 *> \param[out] VS 00164 *> \verbatim 00165 *> VS is REAL array, dimension (LDVS,N) 00166 *> If JOBVS = 'V', VS contains the orthogonal matrix Z of Schur 00167 *> vectors. 00168 *> If JOBVS = 'N', VS is not referenced. 00169 *> \endverbatim 00170 *> 00171 *> \param[in] LDVS 00172 *> \verbatim 00173 *> LDVS is INTEGER 00174 *> The leading dimension of the array VS. LDVS >= 1, and if 00175 *> JOBVS = 'V', LDVS >= N. 00176 *> \endverbatim 00177 *> 00178 *> \param[out] RCONDE 00179 *> \verbatim 00180 *> RCONDE is REAL 00181 *> If SENSE = 'E' or 'B', RCONDE contains the reciprocal 00182 *> condition number for the average of the selected eigenvalues. 00183 *> Not referenced if SENSE = 'N' or 'V'. 00184 *> \endverbatim 00185 *> 00186 *> \param[out] RCONDV 00187 *> \verbatim 00188 *> RCONDV is REAL 00189 *> If SENSE = 'V' or 'B', RCONDV contains the reciprocal 00190 *> condition number for the selected right invariant subspace. 00191 *> Not referenced if SENSE = 'N' or 'E'. 00192 *> \endverbatim 00193 *> 00194 *> \param[out] WORK 00195 *> \verbatim 00196 *> WORK is REAL array, dimension (MAX(1,LWORK)) 00197 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK. 00198 *> \endverbatim 00199 *> 00200 *> \param[in] LWORK 00201 *> \verbatim 00202 *> LWORK is INTEGER 00203 *> The dimension of the array WORK. LWORK >= max(1,3*N). 00204 *> Also, if SENSE = 'E' or 'V' or 'B', 00205 *> LWORK >= N+2*SDIM*(N-SDIM), where SDIM is the number of 00206 *> selected eigenvalues computed by this routine. Note that 00207 *> N+2*SDIM*(N-SDIM) <= N+N*N/2. Note also that an error is only 00208 *> returned if LWORK < max(1,3*N), but if SENSE = 'E' or 'V' or 00209 *> 'B' this may not be large enough. 00210 *> For good performance, LWORK must generally be larger. 00211 *> 00212 *> If LWORK = -1, then a workspace query is assumed; the routine 00213 *> only calculates upper bounds on the optimal sizes of the 00214 *> arrays WORK and IWORK, returns these values as the first 00215 *> entries of the WORK and IWORK arrays, and no error messages 00216 *> related to LWORK or LIWORK are issued by XERBLA. 00217 *> \endverbatim 00218 *> 00219 *> \param[out] IWORK 00220 *> \verbatim 00221 *> IWORK is INTEGER array, dimension (MAX(1,LIWORK)) 00222 *> On exit, if INFO = 0, IWORK(1) returns the optimal LIWORK. 00223 *> \endverbatim 00224 *> 00225 *> \param[in] LIWORK 00226 *> \verbatim 00227 *> LIWORK is INTEGER 00228 *> The dimension of the array IWORK. 00229 *> LIWORK >= 1; if SENSE = 'V' or 'B', LIWORK >= SDIM*(N-SDIM). 00230 *> Note that SDIM*(N-SDIM) <= N*N/4. Note also that an error is 00231 *> only returned if LIWORK < 1, but if SENSE = 'V' or 'B' this 00232 *> may not be large enough. 00233 *> 00234 *> If LIWORK = -1, then a workspace query is assumed; the 00235 *> routine only calculates upper bounds on the optimal sizes of 00236 *> the arrays WORK and IWORK, returns these values as the first 00237 *> entries of the WORK and IWORK arrays, and no error messages 00238 *> related to LWORK or LIWORK are issued by XERBLA. 00239 *> \endverbatim 00240 *> 00241 *> \param[out] BWORK 00242 *> \verbatim 00243 *> BWORK is LOGICAL array, dimension (N) 00244 *> Not referenced if SORT = 'N'. 00245 *> \endverbatim 00246 *> 00247 *> \param[out] INFO 00248 *> \verbatim 00249 *> INFO is INTEGER 00250 *> = 0: successful exit 00251 *> < 0: if INFO = -i, the i-th argument had an illegal value. 00252 *> > 0: if INFO = i, and i is 00253 *> <= N: the QR algorithm failed to compute all the 00254 *> eigenvalues; elements 1:ILO-1 and i+1:N of WR and WI 00255 *> contain those eigenvalues which have converged; if 00256 *> JOBVS = 'V', VS contains the transformation which 00257 *> reduces A to its partially converged Schur form. 00258 *> = N+1: the eigenvalues could not be reordered because some 00259 *> eigenvalues were too close to separate (the problem 00260 *> is very ill-conditioned); 00261 *> = N+2: after reordering, roundoff changed values of some 00262 *> complex eigenvalues so that leading eigenvalues in 00263 *> the Schur form no longer satisfy SELECT=.TRUE. This 00264 *> could also be caused by underflow due to scaling. 00265 *> \endverbatim 00266 * 00267 * Authors: 00268 * ======== 00269 * 00270 *> \author Univ. of Tennessee 00271 *> \author Univ. of California Berkeley 00272 *> \author Univ. of Colorado Denver 00273 *> \author NAG Ltd. 00274 * 00275 *> \date November 2011 00276 * 00277 *> \ingroup realGEeigen 00278 * 00279 * ===================================================================== 00280 SUBROUTINE SGEESX( JOBVS, SORT, SELECT, SENSE, N, A, LDA, SDIM, 00281 $ WR, WI, VS, LDVS, RCONDE, RCONDV, WORK, LWORK, 00282 $ IWORK, LIWORK, BWORK, INFO ) 00283 * 00284 * -- LAPACK driver routine (version 3.4.0) -- 00285 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00286 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00287 * November 2011 00288 * 00289 * .. Scalar Arguments .. 00290 CHARACTER JOBVS, SENSE, SORT 00291 INTEGER INFO, LDA, LDVS, LIWORK, LWORK, N, SDIM 00292 REAL RCONDE, RCONDV 00293 * .. 00294 * .. Array Arguments .. 00295 LOGICAL BWORK( * ) 00296 INTEGER IWORK( * ) 00297 REAL A( LDA, * ), VS( LDVS, * ), WI( * ), WORK( * ), 00298 $ WR( * ) 00299 * .. 00300 * .. Function Arguments .. 00301 LOGICAL SELECT 00302 EXTERNAL SELECT 00303 * .. 00304 * 00305 * ===================================================================== 00306 * 00307 * .. Parameters .. 00308 REAL ZERO, ONE 00309 PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0 ) 00310 * .. 00311 * .. Local Scalars .. 00312 LOGICAL CURSL, LASTSL, LQUERY, LST2SL, SCALEA, WANTSB, 00313 $ WANTSE, WANTSN, WANTST, WANTSV, WANTVS 00314 INTEGER HSWORK, I, I1, I2, IBAL, ICOND, IERR, IEVAL, 00315 $ IHI, ILO, INXT, IP, ITAU, IWRK, LWRK, LIWRK, 00316 $ MAXWRK, MINWRK 00317 REAL ANRM, BIGNUM, CSCALE, EPS, SMLNUM 00318 * .. 00319 * .. Local Arrays .. 00320 REAL DUM( 1 ) 00321 * .. 00322 * .. External Subroutines .. 00323 EXTERNAL SCOPY, SGEBAK, SGEBAL, SGEHRD, SHSEQR, SLABAD, 00324 $ SLACPY, SLASCL, SORGHR, SSWAP, STRSEN, XERBLA 00325 * .. 00326 * .. External Functions .. 00327 LOGICAL LSAME 00328 INTEGER ILAENV 00329 REAL SLAMCH, SLANGE 00330 EXTERNAL LSAME, ILAENV, SLAMCH, SLANGE 00331 * .. 00332 * .. Intrinsic Functions .. 00333 INTRINSIC MAX, SQRT 00334 * .. 00335 * .. Executable Statements .. 00336 * 00337 * Test the input arguments 00338 * 00339 INFO = 0 00340 WANTVS = LSAME( JOBVS, 'V' ) 00341 WANTST = LSAME( SORT, 'S' ) 00342 WANTSN = LSAME( SENSE, 'N' ) 00343 WANTSE = LSAME( SENSE, 'E' ) 00344 WANTSV = LSAME( SENSE, 'V' ) 00345 WANTSB = LSAME( SENSE, 'B' ) 00346 LQUERY = ( LWORK.EQ.-1 .OR. LIWORK.EQ.-1 ) 00347 * 00348 IF( ( .NOT.WANTVS ) .AND. ( .NOT.LSAME( JOBVS, 'N' ) ) ) THEN 00349 INFO = -1 00350 ELSE IF( ( .NOT.WANTST ) .AND. ( .NOT.LSAME( SORT, 'N' ) ) ) THEN 00351 INFO = -2 00352 ELSE IF( .NOT.( WANTSN .OR. WANTSE .OR. WANTSV .OR. WANTSB ) .OR. 00353 $ ( .NOT.WANTST .AND. .NOT.WANTSN ) ) THEN 00354 INFO = -4 00355 ELSE IF( N.LT.0 ) THEN 00356 INFO = -5 00357 ELSE IF( LDA.LT.MAX( 1, N ) ) THEN 00358 INFO = -7 00359 ELSE IF( LDVS.LT.1 .OR. ( WANTVS .AND. LDVS.LT.N ) ) THEN 00360 INFO = -12 00361 END IF 00362 * 00363 * Compute workspace 00364 * (Note: Comments in the code beginning "RWorkspace:" describe the 00365 * minimal amount of real workspace needed at that point in the 00366 * code, as well as the preferred amount for good performance. 00367 * IWorkspace refers to integer workspace. 00368 * NB refers to the optimal block size for the immediately 00369 * following subroutine, as returned by ILAENV. 00370 * HSWORK refers to the workspace preferred by SHSEQR, as 00371 * calculated below. HSWORK is computed assuming ILO=1 and IHI=N, 00372 * the worst case. 00373 * If SENSE = 'E', 'V' or 'B', then the amount of workspace needed 00374 * depends on SDIM, which is computed by the routine STRSEN later 00375 * in the code.) 00376 * 00377 IF( INFO.EQ.0 ) THEN 00378 LIWRK = 1 00379 IF( N.EQ.0 ) THEN 00380 MINWRK = 1 00381 LWRK = 1 00382 ELSE 00383 MAXWRK = 2*N + N*ILAENV( 1, 'SGEHRD', ' ', N, 1, N, 0 ) 00384 MINWRK = 3*N 00385 * 00386 CALL SHSEQR( 'S', JOBVS, N, 1, N, A, LDA, WR, WI, VS, LDVS, 00387 $ WORK, -1, IEVAL ) 00388 HSWORK = WORK( 1 ) 00389 * 00390 IF( .NOT.WANTVS ) THEN 00391 MAXWRK = MAX( MAXWRK, N + HSWORK ) 00392 ELSE 00393 MAXWRK = MAX( MAXWRK, 2*N + ( N - 1 )*ILAENV( 1, 00394 $ 'SORGHR', ' ', N, 1, N, -1 ) ) 00395 MAXWRK = MAX( MAXWRK, N + HSWORK ) 00396 END IF 00397 LWRK = MAXWRK 00398 IF( .NOT.WANTSN ) 00399 $ LWRK = MAX( LWRK, N + ( N*N )/2 ) 00400 IF( WANTSV .OR. WANTSB ) 00401 $ LIWRK = ( N*N )/4 00402 END IF 00403 IWORK( 1 ) = LIWRK 00404 WORK( 1 ) = LWRK 00405 * 00406 IF( LWORK.LT.MINWRK .AND. .NOT.LQUERY ) THEN 00407 INFO = -16 00408 ELSE IF( LIWORK.LT.1 .AND. .NOT.LQUERY ) THEN 00409 INFO = -18 00410 END IF 00411 END IF 00412 * 00413 IF( INFO.NE.0 ) THEN 00414 CALL XERBLA( 'SGEESX', -INFO ) 00415 RETURN 00416 ELSE IF( LQUERY ) THEN 00417 RETURN 00418 END IF 00419 * 00420 * Quick return if possible 00421 * 00422 IF( N.EQ.0 ) THEN 00423 SDIM = 0 00424 RETURN 00425 END IF 00426 * 00427 * Get machine constants 00428 * 00429 EPS = SLAMCH( 'P' ) 00430 SMLNUM = SLAMCH( 'S' ) 00431 BIGNUM = ONE / SMLNUM 00432 CALL SLABAD( SMLNUM, BIGNUM ) 00433 SMLNUM = SQRT( SMLNUM ) / EPS 00434 BIGNUM = ONE / SMLNUM 00435 * 00436 * Scale A if max element outside range [SMLNUM,BIGNUM] 00437 * 00438 ANRM = SLANGE( 'M', N, N, A, LDA, DUM ) 00439 SCALEA = .FALSE. 00440 IF( ANRM.GT.ZERO .AND. ANRM.LT.SMLNUM ) THEN 00441 SCALEA = .TRUE. 00442 CSCALE = SMLNUM 00443 ELSE IF( ANRM.GT.BIGNUM ) THEN 00444 SCALEA = .TRUE. 00445 CSCALE = BIGNUM 00446 END IF 00447 IF( SCALEA ) 00448 $ CALL SLASCL( 'G', 0, 0, ANRM, CSCALE, N, N, A, LDA, IERR ) 00449 * 00450 * Permute the matrix to make it more nearly triangular 00451 * (RWorkspace: need N) 00452 * 00453 IBAL = 1 00454 CALL SGEBAL( 'P', N, A, LDA, ILO, IHI, WORK( IBAL ), IERR ) 00455 * 00456 * Reduce to upper Hessenberg form 00457 * (RWorkspace: need 3*N, prefer 2*N+N*NB) 00458 * 00459 ITAU = N + IBAL 00460 IWRK = N + ITAU 00461 CALL SGEHRD( N, ILO, IHI, A, LDA, WORK( ITAU ), WORK( IWRK ), 00462 $ LWORK-IWRK+1, IERR ) 00463 * 00464 IF( WANTVS ) THEN 00465 * 00466 * Copy Householder vectors to VS 00467 * 00468 CALL SLACPY( 'L', N, N, A, LDA, VS, LDVS ) 00469 * 00470 * Generate orthogonal matrix in VS 00471 * (RWorkspace: need 3*N-1, prefer 2*N+(N-1)*NB) 00472 * 00473 CALL SORGHR( N, ILO, IHI, VS, LDVS, WORK( ITAU ), WORK( IWRK ), 00474 $ LWORK-IWRK+1, IERR ) 00475 END IF 00476 * 00477 SDIM = 0 00478 * 00479 * Perform QR iteration, accumulating Schur vectors in VS if desired 00480 * (RWorkspace: need N+1, prefer N+HSWORK (see comments) ) 00481 * 00482 IWRK = ITAU 00483 CALL SHSEQR( 'S', JOBVS, N, ILO, IHI, A, LDA, WR, WI, VS, LDVS, 00484 $ WORK( IWRK ), LWORK-IWRK+1, IEVAL ) 00485 IF( IEVAL.GT.0 ) 00486 $ INFO = IEVAL 00487 * 00488 * Sort eigenvalues if desired 00489 * 00490 IF( WANTST .AND. INFO.EQ.0 ) THEN 00491 IF( SCALEA ) THEN 00492 CALL SLASCL( 'G', 0, 0, CSCALE, ANRM, N, 1, WR, N, IERR ) 00493 CALL SLASCL( 'G', 0, 0, CSCALE, ANRM, N, 1, WI, N, IERR ) 00494 END IF 00495 DO 10 I = 1, N 00496 BWORK( I ) = SELECT( WR( I ), WI( I ) ) 00497 10 CONTINUE 00498 * 00499 * Reorder eigenvalues, transform Schur vectors, and compute 00500 * reciprocal condition numbers 00501 * (RWorkspace: if SENSE is not 'N', need N+2*SDIM*(N-SDIM) 00502 * otherwise, need N ) 00503 * (IWorkspace: if SENSE is 'V' or 'B', need SDIM*(N-SDIM) 00504 * otherwise, need 0 ) 00505 * 00506 CALL STRSEN( SENSE, JOBVS, BWORK, N, A, LDA, VS, LDVS, WR, WI, 00507 $ SDIM, RCONDE, RCONDV, WORK( IWRK ), LWORK-IWRK+1, 00508 $ IWORK, LIWORK, ICOND ) 00509 IF( .NOT.WANTSN ) 00510 $ MAXWRK = MAX( MAXWRK, N+2*SDIM*( N-SDIM ) ) 00511 IF( ICOND.EQ.-15 ) THEN 00512 * 00513 * Not enough real workspace 00514 * 00515 INFO = -16 00516 ELSE IF( ICOND.EQ.-17 ) THEN 00517 * 00518 * Not enough integer workspace 00519 * 00520 INFO = -18 00521 ELSE IF( ICOND.GT.0 ) THEN 00522 * 00523 * STRSEN failed to reorder or to restore standard Schur form 00524 * 00525 INFO = ICOND + N 00526 END IF 00527 END IF 00528 * 00529 IF( WANTVS ) THEN 00530 * 00531 * Undo balancing 00532 * (RWorkspace: need N) 00533 * 00534 CALL SGEBAK( 'P', 'R', N, ILO, IHI, WORK( IBAL ), N, VS, LDVS, 00535 $ IERR ) 00536 END IF 00537 * 00538 IF( SCALEA ) THEN 00539 * 00540 * Undo scaling for the Schur form of A 00541 * 00542 CALL SLASCL( 'H', 0, 0, CSCALE, ANRM, N, N, A, LDA, IERR ) 00543 CALL SCOPY( N, A, LDA+1, WR, 1 ) 00544 IF( ( WANTSV .OR. WANTSB ) .AND. INFO.EQ.0 ) THEN 00545 DUM( 1 ) = RCONDV 00546 CALL SLASCL( 'G', 0, 0, CSCALE, ANRM, 1, 1, DUM, 1, IERR ) 00547 RCONDV = DUM( 1 ) 00548 END IF 00549 IF( CSCALE.EQ.SMLNUM ) THEN 00550 * 00551 * If scaling back towards underflow, adjust WI if an 00552 * offdiagonal element of a 2-by-2 block in the Schur form 00553 * underflows. 00554 * 00555 IF( IEVAL.GT.0 ) THEN 00556 I1 = IEVAL + 1 00557 I2 = IHI - 1 00558 CALL SLASCL( 'G', 0, 0, CSCALE, ANRM, ILO-1, 1, WI, N, 00559 $ IERR ) 00560 ELSE IF( WANTST ) THEN 00561 I1 = 1 00562 I2 = N - 1 00563 ELSE 00564 I1 = ILO 00565 I2 = IHI - 1 00566 END IF 00567 INXT = I1 - 1 00568 DO 20 I = I1, I2 00569 IF( I.LT.INXT ) 00570 $ GO TO 20 00571 IF( WI( I ).EQ.ZERO ) THEN 00572 INXT = I + 1 00573 ELSE 00574 IF( A( I+1, I ).EQ.ZERO ) THEN 00575 WI( I ) = ZERO 00576 WI( I+1 ) = ZERO 00577 ELSE IF( A( I+1, I ).NE.ZERO .AND. A( I, I+1 ).EQ. 00578 $ ZERO ) THEN 00579 WI( I ) = ZERO 00580 WI( I+1 ) = ZERO 00581 IF( I.GT.1 ) 00582 $ CALL SSWAP( I-1, A( 1, I ), 1, A( 1, I+1 ), 1 ) 00583 IF( N.GT.I+1 ) 00584 $ CALL SSWAP( N-I-1, A( I, I+2 ), LDA, 00585 $ A( I+1, I+2 ), LDA ) 00586 CALL SSWAP( N, VS( 1, I ), 1, VS( 1, I+1 ), 1 ) 00587 A( I, I+1 ) = A( I+1, I ) 00588 A( I+1, I ) = ZERO 00589 END IF 00590 INXT = I + 2 00591 END IF 00592 20 CONTINUE 00593 END IF 00594 CALL SLASCL( 'G', 0, 0, CSCALE, ANRM, N-IEVAL, 1, 00595 $ WI( IEVAL+1 ), MAX( N-IEVAL, 1 ), IERR ) 00596 END IF 00597 * 00598 IF( WANTST .AND. INFO.EQ.0 ) THEN 00599 * 00600 * Check if reordering successful 00601 * 00602 LASTSL = .TRUE. 00603 LST2SL = .TRUE. 00604 SDIM = 0 00605 IP = 0 00606 DO 30 I = 1, N 00607 CURSL = SELECT( WR( I ), WI( I ) ) 00608 IF( WI( I ).EQ.ZERO ) THEN 00609 IF( CURSL ) 00610 $ SDIM = SDIM + 1 00611 IP = 0 00612 IF( CURSL .AND. .NOT.LASTSL ) 00613 $ INFO = N + 2 00614 ELSE 00615 IF( IP.EQ.1 ) THEN 00616 * 00617 * Last eigenvalue of conjugate pair 00618 * 00619 CURSL = CURSL .OR. LASTSL 00620 LASTSL = CURSL 00621 IF( CURSL ) 00622 $ SDIM = SDIM + 2 00623 IP = -1 00624 IF( CURSL .AND. .NOT.LST2SL ) 00625 $ INFO = N + 2 00626 ELSE 00627 * 00628 * First eigenvalue of conjugate pair 00629 * 00630 IP = 1 00631 END IF 00632 END IF 00633 LST2SL = LASTSL 00634 LASTSL = CURSL 00635 30 CONTINUE 00636 END IF 00637 * 00638 WORK( 1 ) = MAXWRK 00639 IF( WANTSV .OR. WANTSB ) THEN 00640 IWORK( 1 ) = SDIM*(N-SDIM) 00641 ELSE 00642 IWORK( 1 ) = 1 00643 END IF 00644 * 00645 RETURN 00646 * 00647 * End of SGEESX 00648 * 00649 END