![]() |
LAPACK
3.4.1
LAPACK: Linear Algebra PACKage
|
00001 *> \brief \b DLAEIN 00002 * 00003 * =========== DOCUMENTATION =========== 00004 * 00005 * Online html documentation available at 00006 * http://www.netlib.org/lapack/explore-html/ 00007 * 00008 *> \htmlonly 00009 *> Download DLAEIN + dependencies 00010 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlaein.f"> 00011 *> [TGZ]</a> 00012 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlaein.f"> 00013 *> [ZIP]</a> 00014 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlaein.f"> 00015 *> [TXT]</a> 00016 *> \endhtmlonly 00017 * 00018 * Definition: 00019 * =========== 00020 * 00021 * SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, B, 00022 * LDB, WORK, EPS3, SMLNUM, BIGNUM, INFO ) 00023 * 00024 * .. Scalar Arguments .. 00025 * LOGICAL NOINIT, RIGHTV 00026 * INTEGER INFO, LDB, LDH, N 00027 * DOUBLE PRECISION BIGNUM, EPS3, SMLNUM, WI, WR 00028 * .. 00029 * .. Array Arguments .. 00030 * DOUBLE PRECISION B( LDB, * ), H( LDH, * ), VI( * ), VR( * ), 00031 * $ WORK( * ) 00032 * .. 00033 * 00034 * 00035 *> \par Purpose: 00036 * ============= 00037 *> 00038 *> \verbatim 00039 *> 00040 *> DLAEIN uses inverse iteration to find a right or left eigenvector 00041 *> corresponding to the eigenvalue (WR,WI) of a real upper Hessenberg 00042 *> matrix H. 00043 *> \endverbatim 00044 * 00045 * Arguments: 00046 * ========== 00047 * 00048 *> \param[in] RIGHTV 00049 *> \verbatim 00050 *> RIGHTV is LOGICAL 00051 *> = .TRUE. : compute right eigenvector; 00052 *> = .FALSE.: compute left eigenvector. 00053 *> \endverbatim 00054 *> 00055 *> \param[in] NOINIT 00056 *> \verbatim 00057 *> NOINIT is LOGICAL 00058 *> = .TRUE. : no initial vector supplied in (VR,VI). 00059 *> = .FALSE.: initial vector supplied in (VR,VI). 00060 *> \endverbatim 00061 *> 00062 *> \param[in] N 00063 *> \verbatim 00064 *> N is INTEGER 00065 *> The order of the matrix H. N >= 0. 00066 *> \endverbatim 00067 *> 00068 *> \param[in] H 00069 *> \verbatim 00070 *> H is DOUBLE PRECISION array, dimension (LDH,N) 00071 *> The upper Hessenberg matrix H. 00072 *> \endverbatim 00073 *> 00074 *> \param[in] LDH 00075 *> \verbatim 00076 *> LDH is INTEGER 00077 *> The leading dimension of the array H. LDH >= max(1,N). 00078 *> \endverbatim 00079 *> 00080 *> \param[in] WR 00081 *> \verbatim 00082 *> WR is DOUBLE PRECISION 00083 *> \endverbatim 00084 *> 00085 *> \param[in] WI 00086 *> \verbatim 00087 *> WI is DOUBLE PRECISION 00088 *> The real and imaginary parts of the eigenvalue of H whose 00089 *> corresponding right or left eigenvector is to be computed. 00090 *> \endverbatim 00091 *> 00092 *> \param[in,out] VR 00093 *> \verbatim 00094 *> VR is DOUBLE PRECISION array, dimension (N) 00095 *> \endverbatim 00096 *> 00097 *> \param[in,out] VI 00098 *> \verbatim 00099 *> VI is DOUBLE PRECISION array, dimension (N) 00100 *> On entry, if NOINIT = .FALSE. and WI = 0.0, VR must contain 00101 *> a real starting vector for inverse iteration using the real 00102 *> eigenvalue WR; if NOINIT = .FALSE. and WI.ne.0.0, VR and VI 00103 *> must contain the real and imaginary parts of a complex 00104 *> starting vector for inverse iteration using the complex 00105 *> eigenvalue (WR,WI); otherwise VR and VI need not be set. 00106 *> On exit, if WI = 0.0 (real eigenvalue), VR contains the 00107 *> computed real eigenvector; if WI.ne.0.0 (complex eigenvalue), 00108 *> VR and VI contain the real and imaginary parts of the 00109 *> computed complex eigenvector. The eigenvector is normalized 00110 *> so that the component of largest magnitude has magnitude 1; 00111 *> here the magnitude of a complex number (x,y) is taken to be 00112 *> |x| + |y|. 00113 *> VI is not referenced if WI = 0.0. 00114 *> \endverbatim 00115 *> 00116 *> \param[out] B 00117 *> \verbatim 00118 *> B is DOUBLE PRECISION array, dimension (LDB,N) 00119 *> \endverbatim 00120 *> 00121 *> \param[in] LDB 00122 *> \verbatim 00123 *> LDB is INTEGER 00124 *> The leading dimension of the array B. LDB >= N+1. 00125 *> \endverbatim 00126 *> 00127 *> \param[out] WORK 00128 *> \verbatim 00129 *> WORK is DOUBLE PRECISION array, dimension (N) 00130 *> \endverbatim 00131 *> 00132 *> \param[in] EPS3 00133 *> \verbatim 00134 *> EPS3 is DOUBLE PRECISION 00135 *> A small machine-dependent value which is used to perturb 00136 *> close eigenvalues, and to replace zero pivots. 00137 *> \endverbatim 00138 *> 00139 *> \param[in] SMLNUM 00140 *> \verbatim 00141 *> SMLNUM is DOUBLE PRECISION 00142 *> A machine-dependent value close to the underflow threshold. 00143 *> \endverbatim 00144 *> 00145 *> \param[in] BIGNUM 00146 *> \verbatim 00147 *> BIGNUM is DOUBLE PRECISION 00148 *> A machine-dependent value close to the overflow threshold. 00149 *> \endverbatim 00150 *> 00151 *> \param[out] INFO 00152 *> \verbatim 00153 *> INFO is INTEGER 00154 *> = 0: successful exit 00155 *> = 1: inverse iteration did not converge; VR is set to the 00156 *> last iterate, and so is VI if WI.ne.0.0. 00157 *> \endverbatim 00158 * 00159 * Authors: 00160 * ======== 00161 * 00162 *> \author Univ. of Tennessee 00163 *> \author Univ. of California Berkeley 00164 *> \author Univ. of Colorado Denver 00165 *> \author NAG Ltd. 00166 * 00167 *> \date November 2011 00168 * 00169 *> \ingroup doubleOTHERauxiliary 00170 * 00171 * ===================================================================== 00172 SUBROUTINE DLAEIN( RIGHTV, NOINIT, N, H, LDH, WR, WI, VR, VI, B, 00173 $ LDB, WORK, EPS3, SMLNUM, BIGNUM, INFO ) 00174 * 00175 * -- LAPACK auxiliary routine (version 3.4.0) -- 00176 * -- LAPACK is a software package provided by Univ. of Tennessee, -- 00177 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..-- 00178 * November 2011 00179 * 00180 * .. Scalar Arguments .. 00181 LOGICAL NOINIT, RIGHTV 00182 INTEGER INFO, LDB, LDH, N 00183 DOUBLE PRECISION BIGNUM, EPS3, SMLNUM, WI, WR 00184 * .. 00185 * .. Array Arguments .. 00186 DOUBLE PRECISION B( LDB, * ), H( LDH, * ), VI( * ), VR( * ), 00187 $ WORK( * ) 00188 * .. 00189 * 00190 * ===================================================================== 00191 * 00192 * .. Parameters .. 00193 DOUBLE PRECISION ZERO, ONE, TENTH 00194 PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0, TENTH = 1.0D-1 ) 00195 * .. 00196 * .. Local Scalars .. 00197 CHARACTER NORMIN, TRANS 00198 INTEGER I, I1, I2, I3, IERR, ITS, J 00199 DOUBLE PRECISION ABSBII, ABSBJJ, EI, EJ, GROWTO, NORM, NRMSML, 00200 $ REC, ROOTN, SCALE, TEMP, VCRIT, VMAX, VNORM, W, 00201 $ W1, X, XI, XR, Y 00202 * .. 00203 * .. External Functions .. 00204 INTEGER IDAMAX 00205 DOUBLE PRECISION DASUM, DLAPY2, DNRM2 00206 EXTERNAL IDAMAX, DASUM, DLAPY2, DNRM2 00207 * .. 00208 * .. External Subroutines .. 00209 EXTERNAL DLADIV, DLATRS, DSCAL 00210 * .. 00211 * .. Intrinsic Functions .. 00212 INTRINSIC ABS, DBLE, MAX, SQRT 00213 * .. 00214 * .. Executable Statements .. 00215 * 00216 INFO = 0 00217 * 00218 * GROWTO is the threshold used in the acceptance test for an 00219 * eigenvector. 00220 * 00221 ROOTN = SQRT( DBLE( N ) ) 00222 GROWTO = TENTH / ROOTN 00223 NRMSML = MAX( ONE, EPS3*ROOTN )*SMLNUM 00224 * 00225 * Form B = H - (WR,WI)*I (except that the subdiagonal elements and 00226 * the imaginary parts of the diagonal elements are not stored). 00227 * 00228 DO 20 J = 1, N 00229 DO 10 I = 1, J - 1 00230 B( I, J ) = H( I, J ) 00231 10 CONTINUE 00232 B( J, J ) = H( J, J ) - WR 00233 20 CONTINUE 00234 * 00235 IF( WI.EQ.ZERO ) THEN 00236 * 00237 * Real eigenvalue. 00238 * 00239 IF( NOINIT ) THEN 00240 * 00241 * Set initial vector. 00242 * 00243 DO 30 I = 1, N 00244 VR( I ) = EPS3 00245 30 CONTINUE 00246 ELSE 00247 * 00248 * Scale supplied initial vector. 00249 * 00250 VNORM = DNRM2( N, VR, 1 ) 00251 CALL DSCAL( N, ( EPS3*ROOTN ) / MAX( VNORM, NRMSML ), VR, 00252 $ 1 ) 00253 END IF 00254 * 00255 IF( RIGHTV ) THEN 00256 * 00257 * LU decomposition with partial pivoting of B, replacing zero 00258 * pivots by EPS3. 00259 * 00260 DO 60 I = 1, N - 1 00261 EI = H( I+1, I ) 00262 IF( ABS( B( I, I ) ).LT.ABS( EI ) ) THEN 00263 * 00264 * Interchange rows and eliminate. 00265 * 00266 X = B( I, I ) / EI 00267 B( I, I ) = EI 00268 DO 40 J = I + 1, N 00269 TEMP = B( I+1, J ) 00270 B( I+1, J ) = B( I, J ) - X*TEMP 00271 B( I, J ) = TEMP 00272 40 CONTINUE 00273 ELSE 00274 * 00275 * Eliminate without interchange. 00276 * 00277 IF( B( I, I ).EQ.ZERO ) 00278 $ B( I, I ) = EPS3 00279 X = EI / B( I, I ) 00280 IF( X.NE.ZERO ) THEN 00281 DO 50 J = I + 1, N 00282 B( I+1, J ) = B( I+1, J ) - X*B( I, J ) 00283 50 CONTINUE 00284 END IF 00285 END IF 00286 60 CONTINUE 00287 IF( B( N, N ).EQ.ZERO ) 00288 $ B( N, N ) = EPS3 00289 * 00290 TRANS = 'N' 00291 * 00292 ELSE 00293 * 00294 * UL decomposition with partial pivoting of B, replacing zero 00295 * pivots by EPS3. 00296 * 00297 DO 90 J = N, 2, -1 00298 EJ = H( J, J-1 ) 00299 IF( ABS( B( J, J ) ).LT.ABS( EJ ) ) THEN 00300 * 00301 * Interchange columns and eliminate. 00302 * 00303 X = B( J, J ) / EJ 00304 B( J, J ) = EJ 00305 DO 70 I = 1, J - 1 00306 TEMP = B( I, J-1 ) 00307 B( I, J-1 ) = B( I, J ) - X*TEMP 00308 B( I, J ) = TEMP 00309 70 CONTINUE 00310 ELSE 00311 * 00312 * Eliminate without interchange. 00313 * 00314 IF( B( J, J ).EQ.ZERO ) 00315 $ B( J, J ) = EPS3 00316 X = EJ / B( J, J ) 00317 IF( X.NE.ZERO ) THEN 00318 DO 80 I = 1, J - 1 00319 B( I, J-1 ) = B( I, J-1 ) - X*B( I, J ) 00320 80 CONTINUE 00321 END IF 00322 END IF 00323 90 CONTINUE 00324 IF( B( 1, 1 ).EQ.ZERO ) 00325 $ B( 1, 1 ) = EPS3 00326 * 00327 TRANS = 'T' 00328 * 00329 END IF 00330 * 00331 NORMIN = 'N' 00332 DO 110 ITS = 1, N 00333 * 00334 * Solve U*x = scale*v for a right eigenvector 00335 * or U**T*x = scale*v for a left eigenvector, 00336 * overwriting x on v. 00337 * 00338 CALL DLATRS( 'Upper', TRANS, 'Nonunit', NORMIN, N, B, LDB, 00339 $ VR, SCALE, WORK, IERR ) 00340 NORMIN = 'Y' 00341 * 00342 * Test for sufficient growth in the norm of v. 00343 * 00344 VNORM = DASUM( N, VR, 1 ) 00345 IF( VNORM.GE.GROWTO*SCALE ) 00346 $ GO TO 120 00347 * 00348 * Choose new orthogonal starting vector and try again. 00349 * 00350 TEMP = EPS3 / ( ROOTN+ONE ) 00351 VR( 1 ) = EPS3 00352 DO 100 I = 2, N 00353 VR( I ) = TEMP 00354 100 CONTINUE 00355 VR( N-ITS+1 ) = VR( N-ITS+1 ) - EPS3*ROOTN 00356 110 CONTINUE 00357 * 00358 * Failure to find eigenvector in N iterations. 00359 * 00360 INFO = 1 00361 * 00362 120 CONTINUE 00363 * 00364 * Normalize eigenvector. 00365 * 00366 I = IDAMAX( N, VR, 1 ) 00367 CALL DSCAL( N, ONE / ABS( VR( I ) ), VR, 1 ) 00368 ELSE 00369 * 00370 * Complex eigenvalue. 00371 * 00372 IF( NOINIT ) THEN 00373 * 00374 * Set initial vector. 00375 * 00376 DO 130 I = 1, N 00377 VR( I ) = EPS3 00378 VI( I ) = ZERO 00379 130 CONTINUE 00380 ELSE 00381 * 00382 * Scale supplied initial vector. 00383 * 00384 NORM = DLAPY2( DNRM2( N, VR, 1 ), DNRM2( N, VI, 1 ) ) 00385 REC = ( EPS3*ROOTN ) / MAX( NORM, NRMSML ) 00386 CALL DSCAL( N, REC, VR, 1 ) 00387 CALL DSCAL( N, REC, VI, 1 ) 00388 END IF 00389 * 00390 IF( RIGHTV ) THEN 00391 * 00392 * LU decomposition with partial pivoting of B, replacing zero 00393 * pivots by EPS3. 00394 * 00395 * The imaginary part of the (i,j)-th element of U is stored in 00396 * B(j+1,i). 00397 * 00398 B( 2, 1 ) = -WI 00399 DO 140 I = 2, N 00400 B( I+1, 1 ) = ZERO 00401 140 CONTINUE 00402 * 00403 DO 170 I = 1, N - 1 00404 ABSBII = DLAPY2( B( I, I ), B( I+1, I ) ) 00405 EI = H( I+1, I ) 00406 IF( ABSBII.LT.ABS( EI ) ) THEN 00407 * 00408 * Interchange rows and eliminate. 00409 * 00410 XR = B( I, I ) / EI 00411 XI = B( I+1, I ) / EI 00412 B( I, I ) = EI 00413 B( I+1, I ) = ZERO 00414 DO 150 J = I + 1, N 00415 TEMP = B( I+1, J ) 00416 B( I+1, J ) = B( I, J ) - XR*TEMP 00417 B( J+1, I+1 ) = B( J+1, I ) - XI*TEMP 00418 B( I, J ) = TEMP 00419 B( J+1, I ) = ZERO 00420 150 CONTINUE 00421 B( I+2, I ) = -WI 00422 B( I+1, I+1 ) = B( I+1, I+1 ) - XI*WI 00423 B( I+2, I+1 ) = B( I+2, I+1 ) + XR*WI 00424 ELSE 00425 * 00426 * Eliminate without interchanging rows. 00427 * 00428 IF( ABSBII.EQ.ZERO ) THEN 00429 B( I, I ) = EPS3 00430 B( I+1, I ) = ZERO 00431 ABSBII = EPS3 00432 END IF 00433 EI = ( EI / ABSBII ) / ABSBII 00434 XR = B( I, I )*EI 00435 XI = -B( I+1, I )*EI 00436 DO 160 J = I + 1, N 00437 B( I+1, J ) = B( I+1, J ) - XR*B( I, J ) + 00438 $ XI*B( J+1, I ) 00439 B( J+1, I+1 ) = -XR*B( J+1, I ) - XI*B( I, J ) 00440 160 CONTINUE 00441 B( I+2, I+1 ) = B( I+2, I+1 ) - WI 00442 END IF 00443 * 00444 * Compute 1-norm of offdiagonal elements of i-th row. 00445 * 00446 WORK( I ) = DASUM( N-I, B( I, I+1 ), LDB ) + 00447 $ DASUM( N-I, B( I+2, I ), 1 ) 00448 170 CONTINUE 00449 IF( B( N, N ).EQ.ZERO .AND. B( N+1, N ).EQ.ZERO ) 00450 $ B( N, N ) = EPS3 00451 WORK( N ) = ZERO 00452 * 00453 I1 = N 00454 I2 = 1 00455 I3 = -1 00456 ELSE 00457 * 00458 * UL decomposition with partial pivoting of conjg(B), 00459 * replacing zero pivots by EPS3. 00460 * 00461 * The imaginary part of the (i,j)-th element of U is stored in 00462 * B(j+1,i). 00463 * 00464 B( N+1, N ) = WI 00465 DO 180 J = 1, N - 1 00466 B( N+1, J ) = ZERO 00467 180 CONTINUE 00468 * 00469 DO 210 J = N, 2, -1 00470 EJ = H( J, J-1 ) 00471 ABSBJJ = DLAPY2( B( J, J ), B( J+1, J ) ) 00472 IF( ABSBJJ.LT.ABS( EJ ) ) THEN 00473 * 00474 * Interchange columns and eliminate 00475 * 00476 XR = B( J, J ) / EJ 00477 XI = B( J+1, J ) / EJ 00478 B( J, J ) = EJ 00479 B( J+1, J ) = ZERO 00480 DO 190 I = 1, J - 1 00481 TEMP = B( I, J-1 ) 00482 B( I, J-1 ) = B( I, J ) - XR*TEMP 00483 B( J, I ) = B( J+1, I ) - XI*TEMP 00484 B( I, J ) = TEMP 00485 B( J+1, I ) = ZERO 00486 190 CONTINUE 00487 B( J+1, J-1 ) = WI 00488 B( J-1, J-1 ) = B( J-1, J-1 ) + XI*WI 00489 B( J, J-1 ) = B( J, J-1 ) - XR*WI 00490 ELSE 00491 * 00492 * Eliminate without interchange. 00493 * 00494 IF( ABSBJJ.EQ.ZERO ) THEN 00495 B( J, J ) = EPS3 00496 B( J+1, J ) = ZERO 00497 ABSBJJ = EPS3 00498 END IF 00499 EJ = ( EJ / ABSBJJ ) / ABSBJJ 00500 XR = B( J, J )*EJ 00501 XI = -B( J+1, J )*EJ 00502 DO 200 I = 1, J - 1 00503 B( I, J-1 ) = B( I, J-1 ) - XR*B( I, J ) + 00504 $ XI*B( J+1, I ) 00505 B( J, I ) = -XR*B( J+1, I ) - XI*B( I, J ) 00506 200 CONTINUE 00507 B( J, J-1 ) = B( J, J-1 ) + WI 00508 END IF 00509 * 00510 * Compute 1-norm of offdiagonal elements of j-th column. 00511 * 00512 WORK( J ) = DASUM( J-1, B( 1, J ), 1 ) + 00513 $ DASUM( J-1, B( J+1, 1 ), LDB ) 00514 210 CONTINUE 00515 IF( B( 1, 1 ).EQ.ZERO .AND. B( 2, 1 ).EQ.ZERO ) 00516 $ B( 1, 1 ) = EPS3 00517 WORK( 1 ) = ZERO 00518 * 00519 I1 = 1 00520 I2 = N 00521 I3 = 1 00522 END IF 00523 * 00524 DO 270 ITS = 1, N 00525 SCALE = ONE 00526 VMAX = ONE 00527 VCRIT = BIGNUM 00528 * 00529 * Solve U*(xr,xi) = scale*(vr,vi) for a right eigenvector, 00530 * or U**T*(xr,xi) = scale*(vr,vi) for a left eigenvector, 00531 * overwriting (xr,xi) on (vr,vi). 00532 * 00533 DO 250 I = I1, I2, I3 00534 * 00535 IF( WORK( I ).GT.VCRIT ) THEN 00536 REC = ONE / VMAX 00537 CALL DSCAL( N, REC, VR, 1 ) 00538 CALL DSCAL( N, REC, VI, 1 ) 00539 SCALE = SCALE*REC 00540 VMAX = ONE 00541 VCRIT = BIGNUM 00542 END IF 00543 * 00544 XR = VR( I ) 00545 XI = VI( I ) 00546 IF( RIGHTV ) THEN 00547 DO 220 J = I + 1, N 00548 XR = XR - B( I, J )*VR( J ) + B( J+1, I )*VI( J ) 00549 XI = XI - B( I, J )*VI( J ) - B( J+1, I )*VR( J ) 00550 220 CONTINUE 00551 ELSE 00552 DO 230 J = 1, I - 1 00553 XR = XR - B( J, I )*VR( J ) + B( I+1, J )*VI( J ) 00554 XI = XI - B( J, I )*VI( J ) - B( I+1, J )*VR( J ) 00555 230 CONTINUE 00556 END IF 00557 * 00558 W = ABS( B( I, I ) ) + ABS( B( I+1, I ) ) 00559 IF( W.GT.SMLNUM ) THEN 00560 IF( W.LT.ONE ) THEN 00561 W1 = ABS( XR ) + ABS( XI ) 00562 IF( W1.GT.W*BIGNUM ) THEN 00563 REC = ONE / W1 00564 CALL DSCAL( N, REC, VR, 1 ) 00565 CALL DSCAL( N, REC, VI, 1 ) 00566 XR = VR( I ) 00567 XI = VI( I ) 00568 SCALE = SCALE*REC 00569 VMAX = VMAX*REC 00570 END IF 00571 END IF 00572 * 00573 * Divide by diagonal element of B. 00574 * 00575 CALL DLADIV( XR, XI, B( I, I ), B( I+1, I ), VR( I ), 00576 $ VI( I ) ) 00577 VMAX = MAX( ABS( VR( I ) )+ABS( VI( I ) ), VMAX ) 00578 VCRIT = BIGNUM / VMAX 00579 ELSE 00580 DO 240 J = 1, N 00581 VR( J ) = ZERO 00582 VI( J ) = ZERO 00583 240 CONTINUE 00584 VR( I ) = ONE 00585 VI( I ) = ONE 00586 SCALE = ZERO 00587 VMAX = ONE 00588 VCRIT = BIGNUM 00589 END IF 00590 250 CONTINUE 00591 * 00592 * Test for sufficient growth in the norm of (VR,VI). 00593 * 00594 VNORM = DASUM( N, VR, 1 ) + DASUM( N, VI, 1 ) 00595 IF( VNORM.GE.GROWTO*SCALE ) 00596 $ GO TO 280 00597 * 00598 * Choose a new orthogonal starting vector and try again. 00599 * 00600 Y = EPS3 / ( ROOTN+ONE ) 00601 VR( 1 ) = EPS3 00602 VI( 1 ) = ZERO 00603 * 00604 DO 260 I = 2, N 00605 VR( I ) = Y 00606 VI( I ) = ZERO 00607 260 CONTINUE 00608 VR( N-ITS+1 ) = VR( N-ITS+1 ) - EPS3*ROOTN 00609 270 CONTINUE 00610 * 00611 * Failure to find eigenvector in N iterations 00612 * 00613 INFO = 1 00614 * 00615 280 CONTINUE 00616 * 00617 * Normalize eigenvector. 00618 * 00619 VNORM = ZERO 00620 DO 290 I = 1, N 00621 VNORM = MAX( VNORM, ABS( VR( I ) )+ABS( VI( I ) ) ) 00622 290 CONTINUE 00623 CALL DSCAL( N, ONE / VNORM, VR, 1 ) 00624 CALL DSCAL( N, ONE / VNORM, VI, 1 ) 00625 * 00626 END IF 00627 * 00628 RETURN 00629 * 00630 * End of DLAEIN 00631 * 00632 END