Return to dlaed5.f CVS log

Up to [local] / rpl / lapack / lapack

Annotation of rpl/lapack/lapack/dlaed5.f, revision 1.6

1.1       bertrand    1:       SUBROUTINE DLAED5( I, D, Z, DELTA, RHO, DLAM )
                      2: *
                      3: *  -- LAPACK routine (version 3.2) --
                      4: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
                      5: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
                      6: *     November 2006
                      7: *
                      8: *     .. Scalar Arguments ..
                      9:       INTEGER            I
                     10:       DOUBLE PRECISION   DLAM, RHO
                     11: *     ..
                     12: *     .. Array Arguments ..
                     13:       DOUBLE PRECISION   D( 2 ), DELTA( 2 ), Z( 2 )
                     14: *     ..
                     15: *
                     16: *  Purpose
                     17: *  =======
                     18: *
                     19: *  This subroutine computes the I-th eigenvalue of a symmetric rank-one
                     20: *  modification of a 2-by-2 diagonal matrix
                     21: *
                     22: *             diag( D )  +  RHO *  Z * transpose(Z) .
                     23: *
                     24: *  The diagonal elements in the array D are assumed to satisfy
                     25: *
                     26: *             D(i) < D(j)  for  i < j .
                     27: *
                     28: *  We also assume RHO > 0 and that the Euclidean norm of the vector
                     29: *  Z is one.
                     30: *
                     31: *  Arguments
                     32: *  =========
                     33: *
                     34: *  I      (input) INTEGER
                     35: *         The index of the eigenvalue to be computed.  I = 1 or I = 2.
                     36: *
                     37: *  D      (input) DOUBLE PRECISION array, dimension (2)
                     38: *         The original eigenvalues.  We assume D(1) < D(2).
                     39: *
                     40: *  Z      (input) DOUBLE PRECISION array, dimension (2)
                     41: *         The components of the updating vector.
                     42: *
                     43: *  DELTA  (output) DOUBLE PRECISION array, dimension (2)
                     44: *         The vector DELTA contains the information necessary
                     45: *         to construct the eigenvectors.
                     46: *
                     47: *  RHO    (input) DOUBLE PRECISION
                     48: *         The scalar in the symmetric updating formula.
                     49: *
                     50: *  DLAM   (output) DOUBLE PRECISION
                     51: *         The computed lambda_I, the I-th updated eigenvalue.
                     52: *
                     53: *  Further Details
                     54: *  ===============
                     55: *
                     56: *  Based on contributions by
                     57: *     Ren-Cang Li, Computer Science Division, University of California
                     58: *     at Berkeley, USA
                     59: *
                     60: *  =====================================================================
                     61: *
                     62: *     .. Parameters ..
                     63:       DOUBLE PRECISION   ZERO, ONE, TWO, FOUR
                     64:       PARAMETER          ( ZERO = 0.0D0, ONE = 1.0D0, TWO = 2.0D0,
                     65:      $                   FOUR = 4.0D0 )
                     66: *     ..
                     67: *     .. Local Scalars ..
                     68:       DOUBLE PRECISION   B, C, DEL, TAU, TEMP, W
                     69: *     ..
                     70: *     .. Intrinsic Functions ..
                     71:       INTRINSIC          ABS, SQRT
                     72: *     ..
                     73: *     .. Executable Statements ..
                     74: *
                     75:       DEL = D( 2 ) - D( 1 )
                     76:       IF( I.EQ.1 ) THEN
                     77:          W = ONE + TWO*RHO*( Z( 2 )*Z( 2 )-Z( 1 )*Z( 1 ) ) / DEL
                     78:          IF( W.GT.ZERO ) THEN
                     79:             B = DEL + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
                     80:             C = RHO*Z( 1 )*Z( 1 )*DEL
                     81: *
                     82: *           B > ZERO, always
                     83: *
                     84:             TAU = TWO*C / ( B+SQRT( ABS( B*B-FOUR*C ) ) )
                     85:             DLAM = D( 1 ) + TAU
                     86:             DELTA( 1 ) = -Z( 1 ) / TAU
                     87:             DELTA( 2 ) = Z( 2 ) / ( DEL-TAU )
                     88:          ELSE
                     89:             B = -DEL + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
                     90:             C = RHO*Z( 2 )*Z( 2 )*DEL
                     91:             IF( B.GT.ZERO ) THEN
                     92:                TAU = -TWO*C / ( B+SQRT( B*B+FOUR*C ) )
                     93:             ELSE
                     94:                TAU = ( B-SQRT( B*B+FOUR*C ) ) / TWO
                     95:             END IF
                     96:             DLAM = D( 2 ) + TAU
                     97:             DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU )
                     98:             DELTA( 2 ) = -Z( 2 ) / TAU
                     99:          END IF
                    100:          TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) )
                    101:          DELTA( 1 ) = DELTA( 1 ) / TEMP
                    102:          DELTA( 2 ) = DELTA( 2 ) / TEMP
                    103:       ELSE
                    104: *
                    105: *     Now I=2
                    106: *
                    107:          B = -DEL + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
                    108:          C = RHO*Z( 2 )*Z( 2 )*DEL
                    109:          IF( B.GT.ZERO ) THEN
                    110:             TAU = ( B+SQRT( B*B+FOUR*C ) ) / TWO
                    111:          ELSE
                    112:             TAU = TWO*C / ( -B+SQRT( B*B+FOUR*C ) )
                    113:          END IF
                    114:          DLAM = D( 2 ) + TAU
                    115:          DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU )
                    116:          DELTA( 2 ) = -Z( 2 ) / TAU
                    117:          TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) )
                    118:          DELTA( 1 ) = DELTA( 1 ) / TEMP
                    119:          DELTA( 2 ) = DELTA( 2 ) / TEMP
                    120:       END IF
                    121:       RETURN
                    122: *
                    123: *     End OF DLAED5
                    124: *
                    125:       END