Annotation of rpl/lapack/lapack/dlahqr.f, revision 1.8
1.8 ! bertrand 1: *> \brief \b DLAHQR
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: *> \htmlonly
! 9: *> Download DLAHQR + dependencies
! 10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlahqr.f">
! 11: *> [TGZ]</a>
! 12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlahqr.f">
! 13: *> [ZIP]</a>
! 14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlahqr.f">
! 15: *> [TXT]</a>
! 16: *> \endhtmlonly
! 17: *
! 18: * Definition:
! 19: * ===========
! 20: *
! 21: * SUBROUTINE DLAHQR( WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI,
! 22: * ILOZ, IHIZ, Z, LDZ, INFO )
! 23: *
! 24: * .. Scalar Arguments ..
! 25: * INTEGER IHI, IHIZ, ILO, ILOZ, INFO, LDH, LDZ, N
! 26: * LOGICAL WANTT, WANTZ
! 27: * ..
! 28: * .. Array Arguments ..
! 29: * DOUBLE PRECISION H( LDH, * ), WI( * ), WR( * ), Z( LDZ, * )
! 30: * ..
! 31: *
! 32: *
! 33: *> \par Purpose:
! 34: * =============
! 35: *>
! 36: *> \verbatim
! 37: *>
! 38: *> DLAHQR is an auxiliary routine called by DHSEQR to update the
! 39: *> eigenvalues and Schur decomposition already computed by DHSEQR, by
! 40: *> dealing with the Hessenberg submatrix in rows and columns ILO to
! 41: *> IHI.
! 42: *> \endverbatim
! 43: *
! 44: * Arguments:
! 45: * ==========
! 46: *
! 47: *> \param[in] WANTT
! 48: *> \verbatim
! 49: *> WANTT is LOGICAL
! 50: *> = .TRUE. : the full Schur form T is required;
! 51: *> = .FALSE.: only eigenvalues are required.
! 52: *> \endverbatim
! 53: *>
! 54: *> \param[in] WANTZ
! 55: *> \verbatim
! 56: *> WANTZ is LOGICAL
! 57: *> = .TRUE. : the matrix of Schur vectors Z is required;
! 58: *> = .FALSE.: Schur vectors are not required.
! 59: *> \endverbatim
! 60: *>
! 61: *> \param[in] N
! 62: *> \verbatim
! 63: *> N is INTEGER
! 64: *> The order of the matrix H. N >= 0.
! 65: *> \endverbatim
! 66: *>
! 67: *> \param[in] ILO
! 68: *> \verbatim
! 69: *> ILO is INTEGER
! 70: *> \endverbatim
! 71: *>
! 72: *> \param[in] IHI
! 73: *> \verbatim
! 74: *> IHI is INTEGER
! 75: *> It is assumed that H is already upper quasi-triangular in
! 76: *> rows and columns IHI+1:N, and that H(ILO,ILO-1) = 0 (unless
! 77: *> ILO = 1). DLAHQR works primarily with the Hessenberg
! 78: *> submatrix in rows and columns ILO to IHI, but applies
! 79: *> transformations to all of H if WANTT is .TRUE..
! 80: *> 1 <= ILO <= max(1,IHI); IHI <= N.
! 81: *> \endverbatim
! 82: *>
! 83: *> \param[in,out] H
! 84: *> \verbatim
! 85: *> H is DOUBLE PRECISION array, dimension (LDH,N)
! 86: *> On entry, the upper Hessenberg matrix H.
! 87: *> On exit, if INFO is zero and if WANTT is .TRUE., H is upper
! 88: *> quasi-triangular in rows and columns ILO:IHI, with any
! 89: *> 2-by-2 diagonal blocks in standard form. If INFO is zero
! 90: *> and WANTT is .FALSE., the contents of H are unspecified on
! 91: *> exit. The output state of H if INFO is nonzero is given
! 92: *> below under the description of INFO.
! 93: *> \endverbatim
! 94: *>
! 95: *> \param[in] LDH
! 96: *> \verbatim
! 97: *> LDH is INTEGER
! 98: *> The leading dimension of the array H. LDH >= max(1,N).
! 99: *> \endverbatim
! 100: *>
! 101: *> \param[out] WR
! 102: *> \verbatim
! 103: *> WR is DOUBLE PRECISION array, dimension (N)
! 104: *> \endverbatim
! 105: *>
! 106: *> \param[out] WI
! 107: *> \verbatim
! 108: *> WI is DOUBLE PRECISION array, dimension (N)
! 109: *> The real and imaginary parts, respectively, of the computed
! 110: *> eigenvalues ILO to IHI are stored in the corresponding
! 111: *> elements of WR and WI. If two eigenvalues are computed as a
! 112: *> complex conjugate pair, they are stored in consecutive
! 113: *> elements of WR and WI, say the i-th and (i+1)th, with
! 114: *> WI(i) > 0 and WI(i+1) < 0. If WANTT is .TRUE., the
! 115: *> eigenvalues are stored in the same order as on the diagonal
! 116: *> of the Schur form returned in H, with WR(i) = H(i,i), and, if
! 117: *> H(i:i+1,i:i+1) is a 2-by-2 diagonal block,
! 118: *> WI(i) = sqrt(H(i+1,i)*H(i,i+1)) and WI(i+1) = -WI(i).
! 119: *> \endverbatim
! 120: *>
! 121: *> \param[in] ILOZ
! 122: *> \verbatim
! 123: *> ILOZ is INTEGER
! 124: *> \endverbatim
! 125: *>
! 126: *> \param[in] IHIZ
! 127: *> \verbatim
! 128: *> IHIZ is INTEGER
! 129: *> Specify the rows of Z to which transformations must be
! 130: *> applied if WANTZ is .TRUE..
! 131: *> 1 <= ILOZ <= ILO; IHI <= IHIZ <= N.
! 132: *> \endverbatim
! 133: *>
! 134: *> \param[in,out] Z
! 135: *> \verbatim
! 136: *> Z is DOUBLE PRECISION array, dimension (LDZ,N)
! 137: *> If WANTZ is .TRUE., on entry Z must contain the current
! 138: *> matrix Z of transformations accumulated by DHSEQR, and on
! 139: *> exit Z has been updated; transformations are applied only to
! 140: *> the submatrix Z(ILOZ:IHIZ,ILO:IHI).
! 141: *> If WANTZ is .FALSE., Z is not referenced.
! 142: *> \endverbatim
! 143: *>
! 144: *> \param[in] LDZ
! 145: *> \verbatim
! 146: *> LDZ is INTEGER
! 147: *> The leading dimension of the array Z. LDZ >= max(1,N).
! 148: *> \endverbatim
! 149: *>
! 150: *> \param[out] INFO
! 151: *> \verbatim
! 152: *> INFO is INTEGER
! 153: *> = 0: successful exit
! 154: *> .GT. 0: If INFO = i, DLAHQR failed to compute all the
! 155: *> eigenvalues ILO to IHI in a total of 30 iterations
! 156: *> per eigenvalue; elements i+1:ihi of WR and WI
! 157: *> contain those eigenvalues which have been
! 158: *> successfully computed.
! 159: *>
! 160: *> If INFO .GT. 0 and WANTT is .FALSE., then on exit,
! 161: *> the remaining unconverged eigenvalues are the
! 162: *> eigenvalues of the upper Hessenberg matrix rows
! 163: *> and columns ILO thorugh INFO of the final, output
! 164: *> value of H.
! 165: *>
! 166: *> If INFO .GT. 0 and WANTT is .TRUE., then on exit
! 167: *> (*) (initial value of H)*U = U*(final value of H)
! 168: *> where U is an orthognal matrix. The final
! 169: *> value of H is upper Hessenberg and triangular in
! 170: *> rows and columns INFO+1 through IHI.
! 171: *>
! 172: *> If INFO .GT. 0 and WANTZ is .TRUE., then on exit
! 173: *> (final value of Z) = (initial value of Z)*U
! 174: *> where U is the orthogonal matrix in (*)
! 175: *> (regardless of the value of WANTT.)
! 176: *> \endverbatim
! 177: *
! 178: * Authors:
! 179: * ========
! 180: *
! 181: *> \author Univ. of Tennessee
! 182: *> \author Univ. of California Berkeley
! 183: *> \author Univ. of Colorado Denver
! 184: *> \author NAG Ltd.
! 185: *
! 186: *> \date November 2011
! 187: *
! 188: *> \ingroup doubleOTHERauxiliary
! 189: *
! 190: *> \par Further Details:
! 191: * =====================
! 192: *>
! 193: *> \verbatim
! 194: *>
! 195: *> 02-96 Based on modifications by
! 196: *> David Day, Sandia National Laboratory, USA
! 197: *>
! 198: *> 12-04 Further modifications by
! 199: *> Ralph Byers, University of Kansas, USA
! 200: *> This is a modified version of DLAHQR from LAPACK version 3.0.
! 201: *> It is (1) more robust against overflow and underflow and
! 202: *> (2) adopts the more conservative Ahues & Tisseur stopping
! 203: *> criterion (LAWN 122, 1997).
! 204: *> \endverbatim
! 205: *>
! 206: * =====================================================================
1.1 bertrand 207: SUBROUTINE DLAHQR( WANTT, WANTZ, N, ILO, IHI, H, LDH, WR, WI,
208: $ ILOZ, IHIZ, Z, LDZ, INFO )
209: *
1.8 ! bertrand 210: * -- LAPACK auxiliary routine (version 3.4.0) --
! 211: * -- LAPACK is a software package provided by Univ. of Tennessee, --
! 212: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
! 213: * November 2011
1.1 bertrand 214: *
215: * .. Scalar Arguments ..
216: INTEGER IHI, IHIZ, ILO, ILOZ, INFO, LDH, LDZ, N
217: LOGICAL WANTT, WANTZ
218: * ..
219: * .. Array Arguments ..
220: DOUBLE PRECISION H( LDH, * ), WI( * ), WR( * ), Z( LDZ, * )
221: * ..
222: *
1.8 ! bertrand 223: * =========================================================
1.1 bertrand 224: *
225: * .. Parameters ..
226: INTEGER ITMAX
227: PARAMETER ( ITMAX = 30 )
228: DOUBLE PRECISION ZERO, ONE, TWO
229: PARAMETER ( ZERO = 0.0d0, ONE = 1.0d0, TWO = 2.0d0 )
230: DOUBLE PRECISION DAT1, DAT2
231: PARAMETER ( DAT1 = 3.0d0 / 4.0d0, DAT2 = -0.4375d0 )
232: * ..
233: * .. Local Scalars ..
234: DOUBLE PRECISION AA, AB, BA, BB, CS, DET, H11, H12, H21, H21S,
235: $ H22, RT1I, RT1R, RT2I, RT2R, RTDISC, S, SAFMAX,
236: $ SAFMIN, SMLNUM, SN, SUM, T1, T2, T3, TR, TST,
237: $ ULP, V2, V3
238: INTEGER I, I1, I2, ITS, J, K, L, M, NH, NR, NZ
239: * ..
240: * .. Local Arrays ..
241: DOUBLE PRECISION V( 3 )
242: * ..
243: * .. External Functions ..
244: DOUBLE PRECISION DLAMCH
245: EXTERNAL DLAMCH
246: * ..
247: * .. External Subroutines ..
248: EXTERNAL DCOPY, DLABAD, DLANV2, DLARFG, DROT
249: * ..
250: * .. Intrinsic Functions ..
251: INTRINSIC ABS, DBLE, MAX, MIN, SQRT
252: * ..
253: * .. Executable Statements ..
254: *
255: INFO = 0
256: *
257: * Quick return if possible
258: *
259: IF( N.EQ.0 )
260: $ RETURN
261: IF( ILO.EQ.IHI ) THEN
262: WR( ILO ) = H( ILO, ILO )
263: WI( ILO ) = ZERO
264: RETURN
265: END IF
266: *
267: * ==== clear out the trash ====
268: DO 10 J = ILO, IHI - 3
269: H( J+2, J ) = ZERO
270: H( J+3, J ) = ZERO
271: 10 CONTINUE
272: IF( ILO.LE.IHI-2 )
273: $ H( IHI, IHI-2 ) = ZERO
274: *
275: NH = IHI - ILO + 1
276: NZ = IHIZ - ILOZ + 1
277: *
278: * Set machine-dependent constants for the stopping criterion.
279: *
280: SAFMIN = DLAMCH( 'SAFE MINIMUM' )
281: SAFMAX = ONE / SAFMIN
282: CALL DLABAD( SAFMIN, SAFMAX )
283: ULP = DLAMCH( 'PRECISION' )
284: SMLNUM = SAFMIN*( DBLE( NH ) / ULP )
285: *
286: * I1 and I2 are the indices of the first row and last column of H
287: * to which transformations must be applied. If eigenvalues only are
288: * being computed, I1 and I2 are set inside the main loop.
289: *
290: IF( WANTT ) THEN
291: I1 = 1
292: I2 = N
293: END IF
294: *
295: * The main loop begins here. I is the loop index and decreases from
296: * IHI to ILO in steps of 1 or 2. Each iteration of the loop works
297: * with the active submatrix in rows and columns L to I.
298: * Eigenvalues I+1 to IHI have already converged. Either L = ILO or
299: * H(L,L-1) is negligible so that the matrix splits.
300: *
301: I = IHI
302: 20 CONTINUE
303: L = ILO
304: IF( I.LT.ILO )
305: $ GO TO 160
306: *
307: * Perform QR iterations on rows and columns ILO to I until a
308: * submatrix of order 1 or 2 splits off at the bottom because a
309: * subdiagonal element has become negligible.
310: *
311: DO 140 ITS = 0, ITMAX
312: *
313: * Look for a single small subdiagonal element.
314: *
315: DO 30 K = I, L + 1, -1
316: IF( ABS( H( K, K-1 ) ).LE.SMLNUM )
317: $ GO TO 40
318: TST = ABS( H( K-1, K-1 ) ) + ABS( H( K, K ) )
319: IF( TST.EQ.ZERO ) THEN
320: IF( K-2.GE.ILO )
321: $ TST = TST + ABS( H( K-1, K-2 ) )
322: IF( K+1.LE.IHI )
323: $ TST = TST + ABS( H( K+1, K ) )
324: END IF
325: * ==== The following is a conservative small subdiagonal
326: * . deflation criterion due to Ahues & Tisseur (LAWN 122,
327: * . 1997). It has better mathematical foundation and
328: * . improves accuracy in some cases. ====
329: IF( ABS( H( K, K-1 ) ).LE.ULP*TST ) THEN
330: AB = MAX( ABS( H( K, K-1 ) ), ABS( H( K-1, K ) ) )
331: BA = MIN( ABS( H( K, K-1 ) ), ABS( H( K-1, K ) ) )
332: AA = MAX( ABS( H( K, K ) ),
333: $ ABS( H( K-1, K-1 )-H( K, K ) ) )
334: BB = MIN( ABS( H( K, K ) ),
335: $ ABS( H( K-1, K-1 )-H( K, K ) ) )
336: S = AA + AB
337: IF( BA*( AB / S ).LE.MAX( SMLNUM,
338: $ ULP*( BB*( AA / S ) ) ) )GO TO 40
339: END IF
340: 30 CONTINUE
341: 40 CONTINUE
342: L = K
343: IF( L.GT.ILO ) THEN
344: *
345: * H(L,L-1) is negligible
346: *
347: H( L, L-1 ) = ZERO
348: END IF
349: *
350: * Exit from loop if a submatrix of order 1 or 2 has split off.
351: *
352: IF( L.GE.I-1 )
353: $ GO TO 150
354: *
355: * Now the active submatrix is in rows and columns L to I. If
356: * eigenvalues only are being computed, only the active submatrix
357: * need be transformed.
358: *
359: IF( .NOT.WANTT ) THEN
360: I1 = L
361: I2 = I
362: END IF
363: *
364: IF( ITS.EQ.10 ) THEN
365: *
366: * Exceptional shift.
367: *
368: S = ABS( H( L+1, L ) ) + ABS( H( L+2, L+1 ) )
369: H11 = DAT1*S + H( L, L )
370: H12 = DAT2*S
371: H21 = S
372: H22 = H11
373: ELSE IF( ITS.EQ.20 ) THEN
374: *
375: * Exceptional shift.
376: *
377: S = ABS( H( I, I-1 ) ) + ABS( H( I-1, I-2 ) )
378: H11 = DAT1*S + H( I, I )
379: H12 = DAT2*S
380: H21 = S
381: H22 = H11
382: ELSE
383: *
384: * Prepare to use Francis' double shift
385: * (i.e. 2nd degree generalized Rayleigh quotient)
386: *
387: H11 = H( I-1, I-1 )
388: H21 = H( I, I-1 )
389: H12 = H( I-1, I )
390: H22 = H( I, I )
391: END IF
392: S = ABS( H11 ) + ABS( H12 ) + ABS( H21 ) + ABS( H22 )
393: IF( S.EQ.ZERO ) THEN
394: RT1R = ZERO
395: RT1I = ZERO
396: RT2R = ZERO
397: RT2I = ZERO
398: ELSE
399: H11 = H11 / S
400: H21 = H21 / S
401: H12 = H12 / S
402: H22 = H22 / S
403: TR = ( H11+H22 ) / TWO
404: DET = ( H11-TR )*( H22-TR ) - H12*H21
405: RTDISC = SQRT( ABS( DET ) )
406: IF( DET.GE.ZERO ) THEN
407: *
408: * ==== complex conjugate shifts ====
409: *
410: RT1R = TR*S
411: RT2R = RT1R
412: RT1I = RTDISC*S
413: RT2I = -RT1I
414: ELSE
415: *
416: * ==== real shifts (use only one of them) ====
417: *
418: RT1R = TR + RTDISC
419: RT2R = TR - RTDISC
420: IF( ABS( RT1R-H22 ).LE.ABS( RT2R-H22 ) ) THEN
421: RT1R = RT1R*S
422: RT2R = RT1R
423: ELSE
424: RT2R = RT2R*S
425: RT1R = RT2R
426: END IF
427: RT1I = ZERO
428: RT2I = ZERO
429: END IF
430: END IF
431: *
432: * Look for two consecutive small subdiagonal elements.
433: *
434: DO 50 M = I - 2, L, -1
435: * Determine the effect of starting the double-shift QR
436: * iteration at row M, and see if this would make H(M,M-1)
437: * negligible. (The following uses scaling to avoid
438: * overflows and most underflows.)
439: *
440: H21S = H( M+1, M )
441: S = ABS( H( M, M )-RT2R ) + ABS( RT2I ) + ABS( H21S )
442: H21S = H( M+1, M ) / S
443: V( 1 ) = H21S*H( M, M+1 ) + ( H( M, M )-RT1R )*
444: $ ( ( H( M, M )-RT2R ) / S ) - RT1I*( RT2I / S )
445: V( 2 ) = H21S*( H( M, M )+H( M+1, M+1 )-RT1R-RT2R )
446: V( 3 ) = H21S*H( M+2, M+1 )
447: S = ABS( V( 1 ) ) + ABS( V( 2 ) ) + ABS( V( 3 ) )
448: V( 1 ) = V( 1 ) / S
449: V( 2 ) = V( 2 ) / S
450: V( 3 ) = V( 3 ) / S
451: IF( M.EQ.L )
452: $ GO TO 60
453: IF( ABS( H( M, M-1 ) )*( ABS( V( 2 ) )+ABS( V( 3 ) ) ).LE.
454: $ ULP*ABS( V( 1 ) )*( ABS( H( M-1, M-1 ) )+ABS( H( M,
455: $ M ) )+ABS( H( M+1, M+1 ) ) ) )GO TO 60
456: 50 CONTINUE
457: 60 CONTINUE
458: *
459: * Double-shift QR step
460: *
461: DO 130 K = M, I - 1
462: *
463: * The first iteration of this loop determines a reflection G
464: * from the vector V and applies it from left and right to H,
465: * thus creating a nonzero bulge below the subdiagonal.
466: *
467: * Each subsequent iteration determines a reflection G to
468: * restore the Hessenberg form in the (K-1)th column, and thus
469: * chases the bulge one step toward the bottom of the active
470: * submatrix. NR is the order of G.
471: *
472: NR = MIN( 3, I-K+1 )
473: IF( K.GT.M )
474: $ CALL DCOPY( NR, H( K, K-1 ), 1, V, 1 )
475: CALL DLARFG( NR, V( 1 ), V( 2 ), 1, T1 )
476: IF( K.GT.M ) THEN
477: H( K, K-1 ) = V( 1 )
478: H( K+1, K-1 ) = ZERO
479: IF( K.LT.I-1 )
480: $ H( K+2, K-1 ) = ZERO
481: ELSE IF( M.GT.L ) THEN
482: * ==== Use the following instead of
483: * . H( K, K-1 ) = -H( K, K-1 ) to
484: * . avoid a bug when v(2) and v(3)
485: * . underflow. ====
486: H( K, K-1 ) = H( K, K-1 )*( ONE-T1 )
487: END IF
488: V2 = V( 2 )
489: T2 = T1*V2
490: IF( NR.EQ.3 ) THEN
491: V3 = V( 3 )
492: T3 = T1*V3
493: *
494: * Apply G from the left to transform the rows of the matrix
495: * in columns K to I2.
496: *
497: DO 70 J = K, I2
498: SUM = H( K, J ) + V2*H( K+1, J ) + V3*H( K+2, J )
499: H( K, J ) = H( K, J ) - SUM*T1
500: H( K+1, J ) = H( K+1, J ) - SUM*T2
501: H( K+2, J ) = H( K+2, J ) - SUM*T3
502: 70 CONTINUE
503: *
504: * Apply G from the right to transform the columns of the
505: * matrix in rows I1 to min(K+3,I).
506: *
507: DO 80 J = I1, MIN( K+3, I )
508: SUM = H( J, K ) + V2*H( J, K+1 ) + V3*H( J, K+2 )
509: H( J, K ) = H( J, K ) - SUM*T1
510: H( J, K+1 ) = H( J, K+1 ) - SUM*T2
511: H( J, K+2 ) = H( J, K+2 ) - SUM*T3
512: 80 CONTINUE
513: *
514: IF( WANTZ ) THEN
515: *
516: * Accumulate transformations in the matrix Z
517: *
518: DO 90 J = ILOZ, IHIZ
519: SUM = Z( J, K ) + V2*Z( J, K+1 ) + V3*Z( J, K+2 )
520: Z( J, K ) = Z( J, K ) - SUM*T1
521: Z( J, K+1 ) = Z( J, K+1 ) - SUM*T2
522: Z( J, K+2 ) = Z( J, K+2 ) - SUM*T3
523: 90 CONTINUE
524: END IF
525: ELSE IF( NR.EQ.2 ) THEN
526: *
527: * Apply G from the left to transform the rows of the matrix
528: * in columns K to I2.
529: *
530: DO 100 J = K, I2
531: SUM = H( K, J ) + V2*H( K+1, J )
532: H( K, J ) = H( K, J ) - SUM*T1
533: H( K+1, J ) = H( K+1, J ) - SUM*T2
534: 100 CONTINUE
535: *
536: * Apply G from the right to transform the columns of the
537: * matrix in rows I1 to min(K+3,I).
538: *
539: DO 110 J = I1, I
540: SUM = H( J, K ) + V2*H( J, K+1 )
541: H( J, K ) = H( J, K ) - SUM*T1
542: H( J, K+1 ) = H( J, K+1 ) - SUM*T2
543: 110 CONTINUE
544: *
545: IF( WANTZ ) THEN
546: *
547: * Accumulate transformations in the matrix Z
548: *
549: DO 120 J = ILOZ, IHIZ
550: SUM = Z( J, K ) + V2*Z( J, K+1 )
551: Z( J, K ) = Z( J, K ) - SUM*T1
552: Z( J, K+1 ) = Z( J, K+1 ) - SUM*T2
553: 120 CONTINUE
554: END IF
555: END IF
556: 130 CONTINUE
557: *
558: 140 CONTINUE
559: *
560: * Failure to converge in remaining number of iterations
561: *
562: INFO = I
563: RETURN
564: *
565: 150 CONTINUE
566: *
567: IF( L.EQ.I ) THEN
568: *
569: * H(I,I-1) is negligible: one eigenvalue has converged.
570: *
571: WR( I ) = H( I, I )
572: WI( I ) = ZERO
573: ELSE IF( L.EQ.I-1 ) THEN
574: *
575: * H(I-1,I-2) is negligible: a pair of eigenvalues have converged.
576: *
577: * Transform the 2-by-2 submatrix to standard Schur form,
578: * and compute and store the eigenvalues.
579: *
580: CALL DLANV2( H( I-1, I-1 ), H( I-1, I ), H( I, I-1 ),
581: $ H( I, I ), WR( I-1 ), WI( I-1 ), WR( I ), WI( I ),
582: $ CS, SN )
583: *
584: IF( WANTT ) THEN
585: *
586: * Apply the transformation to the rest of H.
587: *
588: IF( I2.GT.I )
589: $ CALL DROT( I2-I, H( I-1, I+1 ), LDH, H( I, I+1 ), LDH,
590: $ CS, SN )
591: CALL DROT( I-I1-1, H( I1, I-1 ), 1, H( I1, I ), 1, CS, SN )
592: END IF
593: IF( WANTZ ) THEN
594: *
595: * Apply the transformation to Z.
596: *
597: CALL DROT( NZ, Z( ILOZ, I-1 ), 1, Z( ILOZ, I ), 1, CS, SN )
598: END IF
599: END IF
600: *
601: * return to start of the main loop with new value of I.
602: *
603: I = L - 1
604: GO TO 20
605: *
606: 160 CONTINUE
607: RETURN
608: *
609: * End of DLAHQR
610: *
611: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dlahqr.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack