![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLAHQR 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DLAHQR + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlahqr.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlahqr.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlahqr.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DLAHQR( WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI, 00022 * ILOZ, IHIZ, Z, LDZ, INFO ) 00023 * 00024 * .. Scalar Arguments .. 00025 * INTEGER IHI, IHIZ, ILO, ILOZ, INFO, LDH, LDZ, N 00026 * LOGICAL WANTT, WANTZ 00027 * .. 00028 * .. Array Arguments .. 00029 * DOUBLE PRECISION H( LDH, * ), WI( * ), WR( * ), Z( LDZ, * ) 00030 * .. 00031 * 00032 * 00033 *> \par Purpose: 00034 * ============= 00035 *> 00036 *> \verbatim 00037 *> 00038 *> DLAHQR is an auxiliary routine called by DHSEQR to update the 00039 *> eigenvalues and Schur decomposition already computed by DHSEQR, by 00040 *> dealing with the Hessenberg submatrix in rows and columns ILO to 00041 *> IHI. 00042 *> \endverbatim 00043 * 00044 * Arguments: 00045 * ========== 00046 * 00047 *> \param[in] WANTT 00048 *> \verbatim 00049 *> WANTT is LOGICAL 00050 *> = .TRUE. : the full Schur form T is required; 00051 *> = .FALSE.: only eigenvalues are required. 00052 *> \endverbatim 00053 *> 00054 *> \param[in] WANTZ 00055 *> \verbatim 00056 *> WANTZ is LOGICAL 00057 *> = .TRUE. : the matrix of Schur vectors Z is required; 00058 *> = .FALSE.: Schur vectors are not required. 00059 *> \endverbatim 00060 *> 00061 *> \param[in] N 00062 *> \verbatim 00063 *> N is INTEGER 00064 *> The order of the matrix H. N >= 0. 00065 *> \endverbatim 00066 *> 00067 *> \param[in] ILO 00068 *> \verbatim 00069 *> ILO is INTEGER 00070 *> \endverbatim 00071 *> 00072 *> \param[in] IHI 00073 *> \verbatim 00074 *> IHI is INTEGER 00075 *> It is assumed that H is already upper quasi-triangular in 00076 *> rows and columns IHI+1:N, and that H(ILO,ILO-1) = 0 (unless 00077 *> ILO = 1). DLAHQR works primarily with the Hessenberg 00078 *> submatrix in rows and columns ILO to IHI, but applies 00079 *> transformations to all of H if WANTT is .TRUE.. 00080 *> 1 <= ILO <= max(1,IHI); IHI <= N. 00081 *> \endverbatim 00082 *> 00083 *> \param[in,out] H 00084 *> \verbatim 00085 *> H is DOUBLE PRECISION array, dimension (LDH,N) 00086 *> On entry, the upper Hessenberg matrix H. 00087 *> On exit, if INFO is zero and if WANTT is .TRUE., H is upper 00088 *> quasi-triangular in rows and columns ILO:IHI, with any 00089 *> 2-by-2 diagonal blocks in standard form. If INFO is zero 00090 *> and WANTT is .FALSE., the contents of H are unspecified on 00091 *> exit. The output state of H if INFO is nonzero is given 00092 *> below under the description of INFO. 00093 *> \endverbatim 00094 *> 00095 *> \param[in] LDH 00096 *> \verbatim 00097 *> LDH is INTEGER 00098 *> The leading dimension of the array H. LDH >= max(1,N). 00099 *> \endverbatim 00100 *> 00101 *> \param[out] WR 00102 *> \verbatim 00103 *> WR is DOUBLE PRECISION array, dimension (N) 00104 *> \endverbatim 00105 *> 00106 *> \param[out] WI 00107 *> \verbatim 00108 *> WI is DOUBLE PRECISION array, dimension (N) 00109 *> The real and imaginary parts, respectively, of the computed 00110 *> eigenvalues ILO to IHI are stored in the corresponding 00111 *> elements of WR and WI. If two eigenvalues are computed as a 00112 *> complex conjugate pair, they are stored in consecutive 00113 *> elements of WR and WI, say the i-th and (i+1)th, with 00114 *> WI(i) > 0 and WI(i+1) < 0. If WANTT is .TRUE., the 00115 *> eigenvalues are stored in the same order as on the diagonal 00116 *> of the Schur form returned in H, with WR(i) = H(i,i), and, if 00117 *> H(i:i+1,i:i+1) is a 2-by-2 diagonal block, 00118 *> WI(i) = sqrt(H(i+1,i)*H(i,i+1)) and WI(i+1) = -WI(i). 00119 *> \endverbatim 00120 *> 00121 *> \param[in] ILOZ 00122 *> \verbatim 00123 *> ILOZ is INTEGER 00124 *> \endverbatim 00125 *> 00126 *> \param[in] IHIZ 00127 *> \verbatim 00128 *> IHIZ is INTEGER 00129 *> Specify the rows of Z to which transformations must be 00130 *> applied if WANTZ is .TRUE.. 00131 *> 1 <= ILOZ <= ILO; IHI <= IHIZ <= N. 00132 *> \endverbatim 00133 *> 00134 *> \param[in,out] Z 00135 *> \verbatim 00136 *> Z is DOUBLE PRECISION array, dimension (LDZ,N) 00137 *> If WANTZ is .TRUE., on entry Z must contain the current 00138 *> matrix Z of transformations accumulated by DHSEQR, and on 00139 *> exit Z has been updated; transformations are applied only to 00140 *> the submatrix Z(ILOZ:IHIZ,ILO:IHI). 00141 *> If WANTZ is .FALSE., Z is not referenced. 00142 *> \endverbatim 00143 *> 00144 *> \param[in] LDZ 00145 *> \verbatim 00146 *> LDZ is INTEGER 00147 *> The leading dimension of the array Z. LDZ >= max(1,N). 00148 *> \endverbatim 00149 *> 00150 *> \param[out] INFO 00151 *> \verbatim 00152 *> INFO is INTEGER 00153 *> = 0: successful exit 00154 *> .GT. 0: If INFO = i, DLAHQR failed to compute all the 00155 *> eigenvalues ILO to IHI in a total of 30 iterations 00156 *> per eigenvalue; elements i+1:ihi of WR and WI 00157 *> contain those eigenvalues which have been 00158 *> successfully computed. 00159 *> 00160 *> If INFO .GT. 0 and WANTT is .FALSE., then on exit, 00161 *> the remaining unconverged eigenvalues are the 00162 *> eigenvalues of the upper Hessenberg matrix rows 00163 *> and columns ILO thorugh INFO of the final, output 00164 *> value of H. 00165 *> 00166 *> If INFO .GT. 0 and WANTT is .TRUE., then on exit 00167 *> (*) (initial value of H)*U = U*(final value of H) 00168 *> where U is an orthognal matrix. The final 00169 *> value of H is upper Hessenberg and triangular in 00170 *> rows and columns INFO+1 through IHI. 00171 *> 00172 *> If INFO .GT. 0 and WANTZ is .TRUE., then on exit 00173 *> (final value of Z) = (initial value of Z)*U 00174 *> where U is the orthogonal matrix in (*) 00175 *> (regardless of the value of WANTT.) 00176 *> \endverbatim 00177 * 00178 * Authors: 00179 * ======== 00180 * 00181 *> \author Univ. of Tennessee 00182 *> \author Univ. of California Berkeley 00183 *> \author Univ. of Colorado Denver 00184 *> \author NAG Ltd. 00185 * 00186 *> \date November 2011 00187 * 00188 *> \ingroup doubleOTHERauxiliary 00189 * 00190 *> \par Further Details: 00191 * ===================== 00192 *> 00193 *> \verbatim 00194 *> 00195 *> 02-96 Based on modifications by 00196 *> David Day, Sandia National Laboratory, USA 00197 *> 00198 *> 12-04 Further modifications by 00199 *> Ralph Byers, University of Kansas, USA 00200 *> This is a modified version of DLAHQR from LAPACK version 3.0. 00201 *> It is (1) more robust against overflow and underflow and 00202 *> (2) adopts the more conservative Ahues & Tisseur stopping 00203 *> criterion (LAWN 122, 1997). 00204 *> \endverbatim 00205 *> 00206 * ===================================================================== 00207 SUBROUTINE DLAHQR( WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI, 00208 $ ILOZ, IHIZ, Z, LDZ, INFO ) 00209 * 00210 * -- LAPACK auxiliary routine (version 3.4.0) -- 00211 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00212 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00213 * November 2011 00214 * 00215 * .. Scalar Arguments .. 00216 INTEGER IHI, IHIZ, ILO, ILOZ, INFO, LDH, LDZ, N 00217 LOGICAL WANTT, WANTZ 00218 * .. 00219 * .. Array Arguments .. 00220 DOUBLE PRECISION H( LDH, * ), WI( * ), WR( * ), Z( LDZ, * ) 00221 * .. 00222 * 00223 * ========================================================= 00224 * 00225 * .. Parameters .. 00226 INTEGER ITMAX 00227 PARAMETER ( ITMAX = 30 ) 00228 DOUBLE PRECISION ZERO, ONE, TWO 00229 PARAMETER ( ZERO = 0.0d0, ONE = 1.0d0, TWO = 2.0d0 ) 00230 DOUBLE PRECISION DAT1, DAT2 00231 PARAMETER ( DAT1 = 3.0d0 / 4.0d0, DAT2 = -0.4375d0 ) 00232 * .. 00233 * .. Local Scalars .. 00234 DOUBLE PRECISION AA, AB, BA, BB, CS, DET, H11, H12, H21, H21S, 00235 $ H22, RT1I, RT1R, RT2I, RT2R, RTDISC, S, SAFMAX, 00236 $ SAFMIN, SMLNUM, SN, SUM, T1, T2, T3, TR, TST, 00237 $ ULP, V2, V3 00238 INTEGER I, I1, I2, ITS, J, K, L, M, NH, NR, NZ 00239 * .. 00240 * .. Local Arrays .. 00241 DOUBLE PRECISION V( 3 ) 00242 * .. 00243 * .. External Functions .. 00244 DOUBLE PRECISION DLAMCH 00245 EXTERNAL DLAMCH 00246 * .. 00247 * .. External Subroutines .. 00248 EXTERNAL DCOPY, DLABAD, DLANV2, DLARFG, DROT 00249 * .. 00250 * .. Intrinsic Functions .. 00251 INTRINSIC ABS, DBLE, MAX, MIN, SQRT 00252 * .. 00253 * .. Executable Statements .. 00254 * 00255 INFO = 0 00256 * 00257 * Quick return if possible 00258 * 00259 IF( N.EQ.0 ) 00260 $ RETURN 00261 IF( ILO.EQ.IHI ) THEN 00262 WR( ILO ) = H( ILO, ILO ) 00263 WI( ILO ) = ZERO 00264 RETURN 00265 END IF 00266 * 00267 * ==== clear out the trash ==== 00268 DO 10 J = ILO, IHI - 3 00269 H( J+2, J ) = ZERO 00270 H( J+3, J ) = ZERO 00271 10 CONTINUE 00272 IF( ILO.LE.IHI-2 ) 00273 $ H( IHI, IHI-2 ) = ZERO 00274 * 00275 NH = IHI - ILO + 1 00276 NZ = IHIZ - ILOZ + 1 00277 * 00278 * Set machine-dependent constants for the stopping criterion. 00279 * 00280 SAFMIN = DLAMCH( 'SAFE MINIMUM' ) 00281 SAFMAX = ONE / SAFMIN 00282 CALL DLABAD( SAFMIN, SAFMAX ) 00283 ULP = DLAMCH( 'PRECISION' ) 00284 SMLNUM = SAFMIN*( DBLE( NH ) / ULP ) 00285 * 00286 * I1 and I2 are the indices of the first row and last column of H 00287 * to which transformations must be applied. If eigenvalues only are 00288 * being computed, I1 and I2 are set inside the main loop. 00289 * 00290 IF( WANTT ) THEN 00291 I1 = 1 00292 I2 = N 00293 END IF 00294 * 00295 * The main loop begins here. I is the loop index and decreases from 00296 * IHI to ILO in steps of 1 or 2. Each iteration of the loop works 00297 * with the active submatrix in rows and columns L to I. 00298 * Eigenvalues I+1 to IHI have already converged. Either L = ILO or 00299 * H(L,L-1) is negligible so that the matrix splits. 00300 * 00301 I = IHI 00302 20 CONTINUE 00303 L = ILO 00304 IF( I.LT.ILO ) 00305 $ GO TO 160 00306 * 00307 * Perform QR iterations on rows and columns ILO to I until a 00308 * submatrix of order 1 or 2 splits off at the bottom because a 00309 * subdiagonal element has become negligible. 00310 * 00311 DO 140 ITS = 0, ITMAX 00312 * 00313 * Look for a single small subdiagonal element. 00314 * 00315 DO 30 K = I, L + 1, -1 00316 IF( ABS( H( K, K-1 ) ).LE.SMLNUM ) 00317 $ GO TO 40 00318 TST = ABS( H( K-1, K-1 ) ) + ABS( H( K, K ) ) 00319 IF( TST.EQ.ZERO ) THEN 00320 IF( K-2.GE.ILO ) 00321 $ TST = TST + ABS( H( K-1, K-2 ) ) 00322 IF( K+1.LE.IHI ) 00323 $ TST = TST + ABS( H( K+1, K ) ) 00324 END IF 00325 * ==== The following is a conservative small subdiagonal 00326 * . deflation criterion due to Ahues & Tisseur (LAWN 122, 00327 * . 1997). It has better mathematical foundation and 00328 * . improves accuracy in some cases. ==== 00329 IF( ABS( H( K, K-1 ) ).LE.ULP*TST ) THEN 00330 AB = MAX( ABS( H( K, K-1 ) ), ABS( H( K-1, K ) ) ) 00331 BA = MIN( ABS( H( K, K-1 ) ), ABS( H( K-1, K ) ) ) 00332 AA = MAX( ABS( H( K, K ) ), 00333 $ ABS( H( K-1, K-1 )-H( K, K ) ) ) 00334 BB = MIN( ABS( H( K, K ) ), 00335 $ ABS( H( K-1, K-1 )-H( K, K ) ) ) 00336 S = AA + AB 00337 IF( BA*( AB / S ).LE.MAX( SMLNUM, 00338 $ ULP*( BB*( AA / S ) ) ) )GO TO 40 00339 END IF 00340 30 CONTINUE 00341 40 CONTINUE 00342 L = K 00343 IF( L.GT.ILO ) THEN 00344 * 00345 * H(L,L-1) is negligible 00346 * 00347 H( L, L-1 ) = ZERO 00348 END IF 00349 * 00350 * Exit from loop if a submatrix of order 1 or 2 has split off. 00351 * 00352 IF( L.GE.I-1 ) 00353 $ GO TO 150 00354 * 00355 * Now the active submatrix is in rows and columns L to I. If 00356 * eigenvalues only are being computed, only the active submatrix 00357 * need be transformed. 00358 * 00359 IF( .NOT.WANTT ) THEN 00360 I1 = L 00361 I2 = I 00362 END IF 00363 * 00364 IF( ITS.EQ.10 ) THEN 00365 * 00366 * Exceptional shift. 00367 * 00368 S = ABS( H( L+1, L ) ) + ABS( H( L+2, L+1 ) ) 00369 H11 = DAT1*S + H( L, L ) 00370 H12 = DAT2*S 00371 H21 = S 00372 H22 = H11 00373 ELSE IF( ITS.EQ.20 ) THEN 00374 * 00375 * Exceptional shift. 00376 * 00377 S = ABS( H( I, I-1 ) ) + ABS( H( I-1, I-2 ) ) 00378 H11 = DAT1*S + H( I, I ) 00379 H12 = DAT2*S 00380 H21 = S 00381 H22 = H11 00382 ELSE 00383 * 00384 * Prepare to use Francis' double shift 00385 * (i.e. 2nd degree generalized Rayleigh quotient) 00386 * 00387 H11 = H( I-1, I-1 ) 00388 H21 = H( I, I-1 ) 00389 H12 = H( I-1, I ) 00390 H22 = H( I, I ) 00391 END IF 00392 S = ABS( H11 ) + ABS( H12 ) + ABS( H21 ) + ABS( H22 ) 00393 IF( S.EQ.ZERO ) THEN 00394 RT1R = ZERO 00395 RT1I = ZERO 00396 RT2R = ZERO 00397 RT2I = ZERO 00398 ELSE 00399 H11 = H11 / S 00400 H21 = H21 / S 00401 H12 = H12 / S 00402 H22 = H22 / S 00403 TR = ( H11+H22 ) / TWO 00404 DET = ( H11-TR )*( H22-TR ) - H12*H21 00405 RTDISC = SQRT( ABS( DET ) ) 00406 IF( DET.GE.ZERO ) THEN 00407 * 00408 * ==== complex conjugate shifts ==== 00409 * 00410 RT1R = TR*S 00411 RT2R = RT1R 00412 RT1I = RTDISC*S 00413 RT2I = -RT1I 00414 ELSE 00415 * 00416 * ==== real shifts (use only one of them) ==== 00417 * 00418 RT1R = TR + RTDISC 00419 RT2R = TR - RTDISC 00420 IF( ABS( RT1R-H22 ).LE.ABS( RT2R-H22 ) ) THEN 00421 RT1R = RT1R*S 00422 RT2R = RT1R 00423 ELSE 00424 RT2R = RT2R*S 00425 RT1R = RT2R 00426 END IF 00427 RT1I = ZERO 00428 RT2I = ZERO 00429 END IF 00430 END IF 00431 * 00432 * Look for two consecutive small subdiagonal elements. 00433 * 00434 DO 50 M = I - 2, L, -1 00435 * Determine the effect of starting the double-shift QR 00436 * iteration at row M, and see if this would make H(M,M-1) 00437 * negligible. (The following uses scaling to avoid 00438 * overflows and most underflows.) 00439 * 00440 H21S = H( M+1, M ) 00441 S = ABS( H( M, M )-RT2R ) + ABS( RT2I ) + ABS( H21S ) 00442 H21S = H( M+1, M ) / S 00443 V( 1 ) = H21S*H( M, M+1 ) + ( H( M, M )-RT1R )* 00444 $ ( ( H( M, M )-RT2R ) / S ) - RT1I*( RT2I / S ) 00445 V( 2 ) = H21S*( H( M, M )+H( M+1, M+1 )-RT1R-RT2R ) 00446 V( 3 ) = H21S*H( M+2, M+1 ) 00447 S = ABS( V( 1 ) ) + ABS( V( 2 ) ) + ABS( V( 3 ) ) 00448 V( 1 ) = V( 1 ) / S 00449 V( 2 ) = V( 2 ) / S 00450 V( 3 ) = V( 3 ) / S 00451 IF( M.EQ.L ) 00452 $ GO TO 60 00453 IF( ABS( H( M, M-1 ) )*( ABS( V( 2 ) )+ABS( V( 3 ) ) ).LE. 00454 $ ULP*ABS( V( 1 ) )*( ABS( H( M-1, M-1 ) )+ABS( H( M, 00455 $ M ) )+ABS( H( M+1, M+1 ) ) ) )GO TO 60 00456 50 CONTINUE 00457 60 CONTINUE 00458 * 00459 * Double-shift QR step 00460 * 00461 DO 130 K = M, I - 1 00462 * 00463 * The first iteration of this loop determines a reflection G 00464 * from the vector V and applies it from left and right to H, 00465 * thus creating a nonzero bulge below the subdiagonal. 00466 * 00467 * Each subsequent iteration determines a reflection G to 00468 * restore the Hessenberg form in the (K-1)th column, and thus 00469 * chases the bulge one step toward the bottom of the active 00470 * submatrix. NR is the order of G. 00471 * 00472 NR = MIN( 3, I-K+1 ) 00473 IF( K.GT.M ) 00474 $ CALL DCOPY( NR, H( K, K-1 ), 1, V, 1 ) 00475 CALL DLARFG( NR, V( 1 ), V( 2 ), 1, T1 ) 00476 IF( K.GT.M ) THEN 00477 H( K, K-1 ) = V( 1 ) 00478 H( K+1, K-1 ) = ZERO 00479 IF( K.LT.I-1 ) 00480 $ H( K+2, K-1 ) = ZERO 00481 ELSE IF( M.GT.L ) THEN 00482 * ==== Use the following instead of 00483 * . H( K, K-1 ) = -H( K, K-1 ) to 00484 * . avoid a bug when v(2) and v(3) 00485 * . underflow. ==== 00486 H( K, K-1 ) = H( K, K-1 )*( ONE-T1 ) 00487 END IF 00488 V2 = V( 2 ) 00489 T2 = T1*V2 00490 IF( NR.EQ.3 ) THEN 00491 V3 = V( 3 ) 00492 T3 = T1*V3 00493 * 00494 * Apply G from the left to transform the rows of the matrix 00495 * in columns K to I2. 00496 * 00497 DO 70 J = K, I2 00498 SUM = H( K, J ) + V2*H( K+1, J ) + V3*H( K+2, J ) 00499 H( K, J ) = H( K, J ) - SUM*T1 00500 H( K+1, J ) = H( K+1, J ) - SUM*T2 00501 H( K+2, J ) = H( K+2, J ) - SUM*T3 00502 70 CONTINUE 00503 * 00504 * Apply G from the right to transform the columns of the 00505 * matrix in rows I1 to min(K+3,I). 00506 * 00507 DO 80 J = I1, MIN( K+3, I ) 00508 SUM = H( J, K ) + V2*H( J, K+1 ) + V3*H( J, K+2 ) 00509 H( J, K ) = H( J, K ) - SUM*T1 00510 H( J, K+1 ) = H( J, K+1 ) - SUM*T2 00511 H( J, K+2 ) = H( J, K+2 ) - SUM*T3 00512 80 CONTINUE 00513 * 00514 IF( WANTZ ) THEN 00515 * 00516 * Accumulate transformations in the matrix Z 00517 * 00518 DO 90 J = ILOZ, IHIZ 00519 SUM = Z( J, K ) + V2*Z( J, K+1 ) + V3*Z( J, K+2 ) 00520 Z( J, K ) = Z( J, K ) - SUM*T1 00521 Z( J, K+1 ) = Z( J, K+1 ) - SUM*T2 00522 Z( J, K+2 ) = Z( J, K+2 ) - SUM*T3 00523 90 CONTINUE 00524 END IF 00525 ELSE IF( NR.EQ.2 ) THEN 00526 * 00527 * Apply G from the left to transform the rows of the matrix 00528 * in columns K to I2. 00529 * 00530 DO 100 J = K, I2 00531 SUM = H( K, J ) + V2*H( K+1, J ) 00532 H( K, J ) = H( K, J ) - SUM*T1 00533 H( K+1, J ) = H( K+1, J ) - SUM*T2 00534 100 CONTINUE 00535 * 00536 * Apply G from the right to transform the columns of the 00537 * matrix in rows I1 to min(K+3,I). 00538 * 00539 DO 110 J = I1, I 00540 SUM = H( J, K ) + V2*H( J, K+1 ) 00541 H( J, K ) = H( J, K ) - SUM*T1 00542 H( J, K+1 ) = H( J, K+1 ) - SUM*T2 00543 110 CONTINUE 00544 * 00545 IF( WANTZ ) THEN 00546 * 00547 * Accumulate transformations in the matrix Z 00548 * 00549 DO 120 J = ILOZ, IHIZ 00550 SUM = Z( J, K ) + V2*Z( J, K+1 ) 00551 Z( J, K ) = Z( J, K ) - SUM*T1 00552 Z( J, K+1 ) = Z( J, K+1 ) - SUM*T2 00553 120 CONTINUE 00554 END IF 00555 END IF 00556 130 CONTINUE 00557 * 00558 140 CONTINUE 00559 * 00560 * Failure to converge in remaining number of iterations 00561 * 00562 INFO = I 00563 RETURN 00564 * 00565 150 CONTINUE 00566 * 00567 IF( L.EQ.I ) THEN 00568 * 00569 * H(I,I-1) is negligible: one eigenvalue has converged. 00570 * 00571 WR( I ) = H( I, I ) 00572 WI( I ) = ZERO 00573 ELSE IF( L.EQ.I-1 ) THEN 00574 * 00575 * H(I-1,I-2) is negligible: a pair of eigenvalues have converged. 00576 * 00577 * Transform the 2-by-2 submatrix to standard Schur form, 00578 * and compute and store the eigenvalues. 00579 * 00580 CALL DLANV2( H( I-1, I-1 ), H( I-1, I ), H( I, I-1 ), 00581 $ H( I, I ), WR( I-1 ), WI( I-1 ), WR( I ), WI( I ), 00582 $ CS, SN ) 00583 * 00584 IF( WANTT ) THEN 00585 * 00586 * Apply the transformation to the rest of H. 00587 * 00588 IF( I2.GT.I ) 00589 $ CALL DROT( I2-I, H( I-1, I+1 ), LDH, H( I, I+1 ), LDH, 00590 $ CS, SN ) 00591 CALL DROT( I-I1-1, H( I1, I-1 ), 1, H( I1, I ), 1, CS, SN ) 00592 END IF 00593 IF( WANTZ ) THEN 00594 * 00595 * Apply the transformation to Z. 00596 * 00597 CALL DROT( NZ, Z( ILOZ, I-1 ), 1, Z( ILOZ, I ), 1, CS, SN ) 00598 END IF 00599 END IF 00600 * 00601 * return to start of the main loop with new value of I. 00602 * 00603 I = L - 1 00604 GO TO 20 00605 * 00606 160 CONTINUE 00607 RETURN 00608 * 00609 * End of DLAHQR 00610 * 00611 END