SUBROUTINE SLASD5( I, D, Z, DELTA, RHO, DSIGMA, WORK ) * * -- LAPACK auxiliary routine (instrumented to count ops, version 3.0) -- * Univ. of Tennessee, Oak Ridge National Lab, Argonne National Lab, * Courant Institute, NAG Ltd., and Rice University * June 30, 1999 * * .. Scalar Arguments .. INTEGER I REAL DSIGMA, RHO * .. * .. Array Arguments .. REAL D( 2 ), DELTA( 2 ), WORK( 2 ), Z( 2 ) * .. * .. Common block to return operation count .. COMMON / LATIME / OPS, ITCNT * .. * .. Scalars in Common .. REAL ITCNT, OPS * .. * * Purpose * ======= * * This subroutine computes the square root of the I-th eigenvalue * of a positive symmetric rank-one modification of a 2-by-2 diagonal * matrix * * diag( D ) * diag( D ) + RHO * Z * transpose(Z) . * * The diagonal entries in the array D are assumed to satisfy * * 0 <= D(i) < D(j) for i < j . * * We also assume RHO > 0 and that the Euclidean norm of the vector * Z is one. * * Arguments * ========= * * I (input) INTEGER * The index of the eigenvalue to be computed. I = 1 or I = 2. * * D (input) REAL array, dimension ( 2 ) * The original eigenvalues. We assume 0 <= D(1) < D(2). * * Z (input) REAL array, dimension ( 2 ) * The components of the updating vector. * * DELTA (output) REAL array, dimension ( 2 ) * Contains (D(j) - lambda_I) in its j-th component. * The vector DELTA contains the information necessary * to construct the eigenvectors. * * RHO (input) REAL * The scalar in the symmetric updating formula. * * DSIGMA (output) REAL * The computed lambda_I, the I-th updated eigenvalue. * * WORK (workspace) REAL array, dimension ( 2 ) * WORK contains (D(j) + sigma_I) in its j-th component. * * Further Details * =============== * * Based on contributions by * Ren-Cang Li, Computer Science Division, University of California * at Berkeley, USA * * ===================================================================== * * .. Parameters .. REAL ZERO, ONE, TWO, THREE, FOUR PARAMETER ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0, $ THREE = 3.0E0, FOUR = 4.0E0 ) * .. * .. Local Scalars .. REAL B, C, DEL, DELSQ, TAU, W * .. * .. Intrinsic Functions .. INTRINSIC REAL, ABS, SQRT * .. * .. Executable Statements .. * OPS = OPS + REAL( 3 ) DEL = D( 2 ) - D( 1 ) DELSQ = DEL*( D( 2 )+D( 1 ) ) IF( I.EQ.1 ) THEN OPS = OPS + REAL( 13 ) W = ONE + FOUR*RHO*( Z( 2 )*Z( 2 ) / ( D( 1 )+THREE*D( 2 ) )- $ Z( 1 )*Z( 1 ) / ( THREE*D( 1 )+D( 2 ) ) ) / DEL IF( W.GT.ZERO ) THEN OPS = OPS + REAL( 8 ) B = DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) ) C = RHO*Z( 1 )*Z( 1 )*DELSQ * * B > ZERO, always * * The following TAU is DSIGMA * DSIGMA - D( 1 ) * D( 1 ) * OPS = OPS + REAL( 7 ) TAU = TWO*C / ( B+SQRT( ABS( B*B-FOUR*C ) ) ) * * The following TAU is DSIGMA - D( 1 ) * OPS = OPS + REAL( 14 ) TAU = TAU / ( D( 1 )+SQRT( D( 1 )*D( 1 )+TAU ) ) DSIGMA = D( 1 ) + TAU DELTA( 1 ) = -TAU DELTA( 2 ) = DEL - TAU WORK( 1 ) = TWO*D( 1 ) + TAU WORK( 2 ) = ( D( 1 )+TAU ) + D( 2 ) * DELTA( 1 ) = -Z( 1 ) / TAU * DELTA( 2 ) = Z( 2 ) / ( DEL-TAU ) ELSE OPS = OPS + REAL( 8 ) B = -DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) ) C = RHO*Z( 2 )*Z( 2 )*DELSQ * * The following TAU is DSIGMA * DSIGMA - D( 2 ) * D( 2 ) * IF( B.GT.ZERO ) THEN OPS = OPS + REAL( 7 ) TAU = -TWO*C / ( B+SQRT( B*B+FOUR*C ) ) ELSE OPS = OPS + REAL( 6 ) TAU = ( B-SQRT( B*B+FOUR*C ) ) / TWO END IF * * The following TAU is DSIGMA - D( 2 ) * OPS = OPS + REAL( 14 ) TAU = TAU / ( D( 2 )+SQRT( ABS( D( 2 )*D( 2 )+TAU ) ) ) DSIGMA = D( 2 ) + TAU DELTA( 1 ) = -( DEL+TAU ) DELTA( 2 ) = -TAU WORK( 1 ) = D( 1 ) + TAU + D( 2 ) WORK( 2 ) = TWO*D( 2 ) + TAU * DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU ) * DELTA( 2 ) = -Z( 2 ) / TAU END IF OPS = OPS + REAL( 6 ) * TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) ) * DELTA( 1 ) = DELTA( 1 ) / TEMP * DELTA( 2 ) = DELTA( 2 ) / TEMP ELSE * * Now I=2 * OPS = OPS + REAL( 8 ) B = -DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) ) C = RHO*Z( 2 )*Z( 2 )*DELSQ * * The following TAU is DSIGMA * DSIGMA - D( 2 ) * D( 2 ) * IF( B.GT.ZERO ) THEN OPS = OPS + REAL( 6 ) TAU = ( B+SQRT( B*B+FOUR*C ) ) / TWO ELSE OPS = OPS + REAL( 7 ) TAU = TWO*C / ( -B+SQRT( B*B+FOUR*C ) ) END IF * * The following TAU is DSIGMA - D( 2 ) * OPS = OPS + REAL( 20 ) TAU = TAU / ( D( 2 )+SQRT( D( 2 )*D( 2 )+TAU ) ) DSIGMA = D( 2 ) + TAU DELTA( 1 ) = -( DEL+TAU ) DELTA( 2 ) = -TAU WORK( 1 ) = D( 1 ) + TAU + D( 2 ) WORK( 2 ) = TWO*D( 2 ) + TAU * DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU ) * DELTA( 2 ) = -Z( 2 ) / TAU * TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) ) * DELTA( 1 ) = DELTA( 1 ) / TEMP * DELTA( 2 ) = DELTA( 2 ) / TEMP END IF RETURN * * End of SLASD5 * END