LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
dlasd5.f
Go to the documentation of this file.
1 *> \brief \b DLASD5 computes the square root of the i-th eigenvalue of a positive symmetric rank-one modification of a 2-by-2 diagonal matrix. Used by sbdsdc.
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download DLASD5 + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasd5.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasd5.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasd5.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE DLASD5( I, D, Z, DELTA, RHO, DSIGMA, WORK )
22 *
23 * .. Scalar Arguments ..
24 * INTEGER I
25 * DOUBLE PRECISION DSIGMA, RHO
26 * ..
27 * .. Array Arguments ..
28 * DOUBLE PRECISION D( 2 ), DELTA( 2 ), WORK( 2 ), Z( 2 )
29 * ..
30 *
31 *
32 *> \par Purpose:
33 * =============
34 *>
35 *> \verbatim
36 *>
37 *> This subroutine computes the square root of the I-th eigenvalue
38 *> of a positive symmetric rank-one modification of a 2-by-2 diagonal
39 *> matrix
40 *>
41 *> diag( D ) * diag( D ) + RHO * Z * transpose(Z) .
42 *>
43 *> The diagonal entries in the array D are assumed to satisfy
44 *>
45 *> 0 <= D(i) < D(j) for i < j .
46 *>
47 *> We also assume RHO > 0 and that the Euclidean norm of the vector
48 *> Z is one.
49 *> \endverbatim
50 *
51 * Arguments:
52 * ==========
53 *
54 *> \param[in] I
55 *> \verbatim
56 *> I is INTEGER
57 *> The index of the eigenvalue to be computed. I = 1 or I = 2.
58 *> \endverbatim
59 *>
60 *> \param[in] D
61 *> \verbatim
62 *> D is DOUBLE PRECISION array, dimension ( 2 )
63 *> The original eigenvalues. We assume 0 <= D(1) < D(2).
64 *> \endverbatim
65 *>
66 *> \param[in] Z
67 *> \verbatim
68 *> Z is DOUBLE PRECISION array, dimension ( 2 )
69 *> The components of the updating vector.
70 *> \endverbatim
71 *>
72 *> \param[out] DELTA
73 *> \verbatim
74 *> DELTA is DOUBLE PRECISION array, dimension ( 2 )
75 *> Contains (D(j) - sigma_I) in its j-th component.
76 *> The vector DELTA contains the information necessary
77 *> to construct the eigenvectors.
78 *> \endverbatim
79 *>
80 *> \param[in] RHO
81 *> \verbatim
82 *> RHO is DOUBLE PRECISION
83 *> The scalar in the symmetric updating formula.
84 *> \endverbatim
85 *>
86 *> \param[out] DSIGMA
87 *> \verbatim
88 *> DSIGMA is DOUBLE PRECISION
89 *> The computed sigma_I, the I-th updated eigenvalue.
90 *> \endverbatim
91 *>
92 *> \param[out] WORK
93 *> \verbatim
94 *> WORK is DOUBLE PRECISION array, dimension ( 2 )
95 *> WORK contains (D(j) + sigma_I) in its j-th component.
96 *> \endverbatim
97 *
98 * Authors:
99 * ========
100 *
101 *> \author Univ. of Tennessee
102 *> \author Univ. of California Berkeley
103 *> \author Univ. of Colorado Denver
104 *> \author NAG Ltd.
105 *
106 *> \ingroup OTHERauxiliary
107 *
108 *> \par Contributors:
109 * ==================
110 *>
111 *> Ren-Cang Li, Computer Science Division, University of California
112 *> at Berkeley, USA
113 *>
114 * =====================================================================
115  SUBROUTINE dlasd5( I, D, Z, DELTA, RHO, DSIGMA, WORK )
116 *
117 * -- LAPACK auxiliary routine --
118 * -- LAPACK is a software package provided by Univ. of Tennessee, --
119 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
120 *
121 * .. Scalar Arguments ..
122  INTEGER I
123  DOUBLE PRECISION DSIGMA, RHO
124 * ..
125 * .. Array Arguments ..
126  DOUBLE PRECISION D( 2 ), DELTA( 2 ), WORK( 2 ), Z( 2 )
127 * ..
128 *
129 * =====================================================================
130 *
131 * .. Parameters ..
132  DOUBLE PRECISION ZERO, ONE, TWO, THREE, FOUR
133  parameter( zero = 0.0d+0, one = 1.0d+0, two = 2.0d+0,
134  $ three = 3.0d+0, four = 4.0d+0 )
135 * ..
136 * .. Local Scalars ..
137  DOUBLE PRECISION B, C, DEL, DELSQ, TAU, W
138 * ..
139 * .. Intrinsic Functions ..
140  INTRINSIC abs, sqrt
141 * ..
142 * .. Executable Statements ..
143 *
144  del = d( 2 ) - d( 1 )
145  delsq = del*( d( 2 )+d( 1 ) )
146  IF( i.EQ.1 ) THEN
147  w = one + four*rho*( z( 2 )*z( 2 ) / ( d( 1 )+three*d( 2 ) )-
148  $ z( 1 )*z( 1 ) / ( three*d( 1 )+d( 2 ) ) ) / del
149  IF( w.GT.zero ) THEN
150  b = delsq + rho*( z( 1 )*z( 1 )+z( 2 )*z( 2 ) )
151  c = rho*z( 1 )*z( 1 )*delsq
152 *
153 * B > ZERO, always
154 *
155 * The following TAU is DSIGMA * DSIGMA - D( 1 ) * D( 1 )
156 *
157  tau = two*c / ( b+sqrt( abs( b*b-four*c ) ) )
158 *
159 * The following TAU is DSIGMA - D( 1 )
160 *
161  tau = tau / ( d( 1 )+sqrt( d( 1 )*d( 1 )+tau ) )
162  dsigma = d( 1 ) + tau
163  delta( 1 ) = -tau
164  delta( 2 ) = del - tau
165  work( 1 ) = two*d( 1 ) + tau
166  work( 2 ) = ( d( 1 )+tau ) + d( 2 )
167 * DELTA( 1 ) = -Z( 1 ) / TAU
168 * DELTA( 2 ) = Z( 2 ) / ( DEL-TAU )
169  ELSE
170  b = -delsq + rho*( z( 1 )*z( 1 )+z( 2 )*z( 2 ) )
171  c = rho*z( 2 )*z( 2 )*delsq
172 *
173 * The following TAU is DSIGMA * DSIGMA - D( 2 ) * D( 2 )
174 *
175  IF( b.GT.zero ) THEN
176  tau = -two*c / ( b+sqrt( b*b+four*c ) )
177  ELSE
178  tau = ( b-sqrt( b*b+four*c ) ) / two
179  END IF
180 *
181 * The following TAU is DSIGMA - D( 2 )
182 *
183  tau = tau / ( d( 2 )+sqrt( abs( d( 2 )*d( 2 )+tau ) ) )
184  dsigma = d( 2 ) + tau
185  delta( 1 ) = -( del+tau )
186  delta( 2 ) = -tau
187  work( 1 ) = d( 1 ) + tau + d( 2 )
188  work( 2 ) = two*d( 2 ) + tau
189 * DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU )
190 * DELTA( 2 ) = -Z( 2 ) / TAU
191  END IF
192 * TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) )
193 * DELTA( 1 ) = DELTA( 1 ) / TEMP
194 * DELTA( 2 ) = DELTA( 2 ) / TEMP
195  ELSE
196 *
197 * Now I=2
198 *
199  b = -delsq + rho*( z( 1 )*z( 1 )+z( 2 )*z( 2 ) )
200  c = rho*z( 2 )*z( 2 )*delsq
201 *
202 * The following TAU is DSIGMA * DSIGMA - D( 2 ) * D( 2 )
203 *
204  IF( b.GT.zero ) THEN
205  tau = ( b+sqrt( b*b+four*c ) ) / two
206  ELSE
207  tau = two*c / ( -b+sqrt( b*b+four*c ) )
208  END IF
209 *
210 * The following TAU is DSIGMA - D( 2 )
211 *
212  tau = tau / ( d( 2 )+sqrt( d( 2 )*d( 2 )+tau ) )
213  dsigma = d( 2 ) + tau
214  delta( 1 ) = -( del+tau )
215  delta( 2 ) = -tau
216  work( 1 ) = d( 1 ) + tau + d( 2 )
217  work( 2 ) = two*d( 2 ) + tau
218 * DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU )
219 * DELTA( 2 ) = -Z( 2 ) / TAU
220 * TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) )
221 * DELTA( 1 ) = DELTA( 1 ) / TEMP
222 * DELTA( 2 ) = DELTA( 2 ) / TEMP
223  END IF
224  RETURN
225 *
226 * End of DLASD5
227 *
228  END
subroutine dlasd5(I, D, Z, DELTA, RHO, DSIGMA, WORK)
DLASD5 computes the square root of the i-th eigenvalue of a positive symmetric rank-one modification ...
Definition: dlasd5.f:116