SRC/slaed5.f

   1 *> \brief \b SLAED5 used by sstedc. Solves the 2-by-2 secular equation.
   2 *
   3 *  =========== DOCUMENTATION ===========
   4 *
   5 * Online html documentation available at
   6 *            http://www.netlib.org/lapack/explore-html/
   7 *
   8 *> \htmlonly
   9 *> Download SLAED5 + dependencies
  10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slaed5.f">
  11 *> [TGZ]</a>
  12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slaed5.f">
  13 *> [ZIP]</a>
  14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slaed5.f">
  15 *> [TXT]</a>
  16 *> \endhtmlonly
  17 *
  18 *  Definition:
  19 *  ===========
  20 *
  21 *       SUBROUTINE SLAED5( I, D, Z, DELTA, RHO, DLAM )
  22 *
  23 *       .. Scalar Arguments ..
  24 *       INTEGER            I
  25 *       REAL               DLAM, RHO
  26 *       ..
  27 *       .. Array Arguments ..
  28 *       REAL               D( 2 ), DELTA( 2 ), Z( 2 )
  29 *       ..
  30 *
  31 *
  32 *> \par Purpose:
  33 *  =============
  34 *>
  35 *> \verbatim
  36 *>
  37 *> This subroutine computes the I-th eigenvalue of a symmetric rank-one
  38 *> modification of a 2-by-2 diagonal matrix
  39 *>
  40 *>            diag( D )  +  RHO * Z * transpose(Z) .
  41 *>
  42 *> The diagonal elements in the array D are assumed to satisfy
  43 *>
  44 *>            D(i) < D(j)  for  i < j .
  45 *>
  46 *> We also assume RHO > 0 and that the Euclidean norm of the vector
  47 *> Z is one.
  48 *> \endverbatim
  49 *
  50 *  Arguments:
  51 *  ==========
  52 *
  53 *> \param[in] I
  54 *> \verbatim
  55 *>          I is INTEGER
  56 *>         The index of the eigenvalue to be computed.  I = 1 or I = 2.
  57 *> \endverbatim
  58 *>
  59 *> \param[in] D
  60 *> \verbatim
  61 *>          D is REAL array, dimension (2)
  62 *>         The original eigenvalues.  We assume D(1) < D(2).
  63 *> \endverbatim
  64 *>
  65 *> \param[in] Z
  66 *> \verbatim
  67 *>          Z is REAL array, dimension (2)
  68 *>         The components of the updating vector.
  69 *> \endverbatim
  70 *>
  71 *> \param[out] DELTA
  72 *> \verbatim
  73 *>          DELTA is REAL array, dimension (2)
  74 *>         The vector DELTA contains the information necessary
  75 *>         to construct the eigenvectors.
  76 *> \endverbatim
  77 *>
  78 *> \param[in] RHO
  79 *> \verbatim
  80 *>          RHO is REAL
  81 *>         The scalar in the symmetric updating formula.
  82 *> \endverbatim
  83 *>
  84 *> \param[out] DLAM
  85 *> \verbatim
  86 *>          DLAM is REAL
  87 *>         The computed lambda_I, the I-th updated eigenvalue.
  88 *> \endverbatim
  89 *
  90 *  Authors:
  91 *  ========
  92 *
  93 *> \author Univ. of Tennessee
  94 *> \author Univ. of California Berkeley
  95 *> \author Univ. of Colorado Denver
  96 *> \author NAG Ltd.
  97 *
  98 *> \date September 2012
  99 *
 100 *> \ingroup auxOTHERcomputational
 101 *
 102 *> \par Contributors:
 103 *  ==================
 104 *>
 105 *>     Ren-Cang Li, Computer Science Division, University of California
 106 *>     at Berkeley, USA
 107 *>
 108 *  =====================================================================
 109       SUBROUTINE SLAED5( I, D, Z, DELTA, RHO, DLAM )
 110 *
 111 *  -- LAPACK computational routine (version 3.4.2) --
 112 *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
 113 *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
 114 *     September 2012
 115 *
 116 *     .. Scalar Arguments ..
 117       INTEGER            I
 118       REAL               DLAM, RHO
 119 *     ..
 120 *     .. Array Arguments ..
 121       REAL               D( 2 ), DELTA( 2 ), Z( 2 )
 122 *     ..
 123 *
 124 *  =====================================================================
 125 *
 126 *     .. Parameters ..
 127       REAL               ZERO, ONE, TWO, FOUR
 128       PARAMETER          ( ZERO = 0.0E0, ONE = 1.0E0, TWO = 2.0E0,
 129      $                   FOUR = 4.0E0 )
 130 *     ..
 131 *     .. Local Scalars ..
 132       REAL               B, C, DEL, TAU, TEMP, W
 133 *     ..
 134 *     .. Intrinsic Functions ..
 135       INTRINSIC          ABS, SQRT
 136 *     ..
 137 *     .. Executable Statements ..
 138 *
 139       DEL = D( 2 ) - D( 1 )
 140       IF( I.EQ.1 ) THEN
 141          W = ONE + TWO*RHO*( Z( 2 )*Z( 2 )-Z( 1 )*Z( 1 ) ) / DEL
 142          IF( W.GT.ZERO ) THEN
 143             B = DEL + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
 144             C = RHO*Z( 1 )*Z( 1 )*DEL
 145 *
 146 *           B > ZERO, always
 147 *
 148             TAU = TWO*C / ( B+SQRT( ABS( B*B-FOUR*C ) ) )
 149             DLAM = D( 1 ) + TAU
 150             DELTA( 1 ) = -Z( 1 ) / TAU
 151             DELTA( 2 ) = Z( 2 ) / ( DEL-TAU )
 152          ELSE
 153             B = -DEL + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
 154             C = RHO*Z( 2 )*Z( 2 )*DEL
 155             IF( B.GT.ZERO ) THEN
 156                TAU = -TWO*C / ( B+SQRT( B*B+FOUR*C ) )
 157             ELSE
 158                TAU = ( B-SQRT( B*B+FOUR*C ) ) / TWO
 159             END IF
 160             DLAM = D( 2 ) + TAU
 161             DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU )
 162             DELTA( 2 ) = -Z( 2 ) / TAU
 163          END IF
 164          TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) )
 165          DELTA( 1 ) = DELTA( 1 ) / TEMP
 166          DELTA( 2 ) = DELTA( 2 ) / TEMP
 167       ELSE
 168 *
 169 *     Now I=2
 170 *
 171          B = -DEL + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
 172          C = RHO*Z( 2 )*Z( 2 )*DEL
 173          IF( B.GT.ZERO ) THEN
 174             TAU = ( B+SQRT( B*B+FOUR*C ) ) / TWO
 175          ELSE
 176             TAU = TWO*C / ( -B+SQRT( B*B+FOUR*C ) )
 177          END IF
 178          DLAM = D( 2 ) + TAU
 179          DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU )
 180          DELTA( 2 ) = -Z( 2 ) / TAU
 181          TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) )
 182          DELTA( 1 ) = DELTA( 1 ) / TEMP
 183          DELTA( 2 ) = DELTA( 2 ) / TEMP
 184       END IF
 185       RETURN
 186 *
 187 *     End OF SLAED5
 188 *
 189       END