Annotation of rpl/lapack/lapack/dlasd5.f, revision 1.8
1.8 ! bertrand 1: *> \brief \b DLASD5
! 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: *> \date November 2011
! 107: *
! 108: *> \ingroup auxOTHERauxiliary
! 109: *
! 110: *> \par Contributors:
! 111: * ==================
! 112: *>
! 113: *> Ren-Cang Li, Computer Science Division, University of California
! 114: *> at Berkeley, USA
! 115: *>
! 116: * =====================================================================
1.1 bertrand 117: SUBROUTINE DLASD5( I, D, Z, DELTA, RHO, DSIGMA, WORK )
118: *
1.8 ! bertrand 119: * -- LAPACK auxiliary routine (version 3.4.0) --
1.1 bertrand 120: * -- LAPACK is a software package provided by Univ. of Tennessee, --
121: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8 ! bertrand 122: * November 2011
1.1 bertrand 123: *
124: * .. Scalar Arguments ..
125: INTEGER I
126: DOUBLE PRECISION DSIGMA, RHO
127: * ..
128: * .. Array Arguments ..
129: DOUBLE PRECISION D( 2 ), DELTA( 2 ), WORK( 2 ), Z( 2 )
130: * ..
131: *
132: * =====================================================================
133: *
134: * .. Parameters ..
135: DOUBLE PRECISION ZERO, ONE, TWO, THREE, FOUR
136: PARAMETER ( ZERO = 0.0D+0, ONE = 1.0D+0, TWO = 2.0D+0,
137: $ THREE = 3.0D+0, FOUR = 4.0D+0 )
138: * ..
139: * .. Local Scalars ..
140: DOUBLE PRECISION B, C, DEL, DELSQ, TAU, W
141: * ..
142: * .. Intrinsic Functions ..
143: INTRINSIC ABS, SQRT
144: * ..
145: * .. Executable Statements ..
146: *
147: DEL = D( 2 ) - D( 1 )
148: DELSQ = DEL*( D( 2 )+D( 1 ) )
149: IF( I.EQ.1 ) THEN
150: W = ONE + FOUR*RHO*( Z( 2 )*Z( 2 ) / ( D( 1 )+THREE*D( 2 ) )-
151: $ Z( 1 )*Z( 1 ) / ( THREE*D( 1 )+D( 2 ) ) ) / DEL
152: IF( W.GT.ZERO ) THEN
153: B = DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
154: C = RHO*Z( 1 )*Z( 1 )*DELSQ
155: *
156: * B > ZERO, always
157: *
158: * The following TAU is DSIGMA * DSIGMA - D( 1 ) * D( 1 )
159: *
160: TAU = TWO*C / ( B+SQRT( ABS( B*B-FOUR*C ) ) )
161: *
162: * The following TAU is DSIGMA - D( 1 )
163: *
164: TAU = TAU / ( D( 1 )+SQRT( D( 1 )*D( 1 )+TAU ) )
165: DSIGMA = D( 1 ) + TAU
166: DELTA( 1 ) = -TAU
167: DELTA( 2 ) = DEL - TAU
168: WORK( 1 ) = TWO*D( 1 ) + TAU
169: WORK( 2 ) = ( D( 1 )+TAU ) + D( 2 )
170: * DELTA( 1 ) = -Z( 1 ) / TAU
171: * DELTA( 2 ) = Z( 2 ) / ( DEL-TAU )
172: ELSE
173: B = -DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
174: C = RHO*Z( 2 )*Z( 2 )*DELSQ
175: *
176: * The following TAU is DSIGMA * DSIGMA - D( 2 ) * D( 2 )
177: *
178: IF( B.GT.ZERO ) THEN
179: TAU = -TWO*C / ( B+SQRT( B*B+FOUR*C ) )
180: ELSE
181: TAU = ( B-SQRT( B*B+FOUR*C ) ) / TWO
182: END IF
183: *
184: * The following TAU is DSIGMA - D( 2 )
185: *
186: TAU = TAU / ( D( 2 )+SQRT( ABS( D( 2 )*D( 2 )+TAU ) ) )
187: DSIGMA = D( 2 ) + TAU
188: DELTA( 1 ) = -( DEL+TAU )
189: DELTA( 2 ) = -TAU
190: WORK( 1 ) = D( 1 ) + TAU + D( 2 )
191: WORK( 2 ) = TWO*D( 2 ) + TAU
192: * DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU )
193: * DELTA( 2 ) = -Z( 2 ) / TAU
194: END IF
195: * TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) )
196: * DELTA( 1 ) = DELTA( 1 ) / TEMP
197: * DELTA( 2 ) = DELTA( 2 ) / TEMP
198: ELSE
199: *
200: * Now I=2
201: *
202: B = -DELSQ + RHO*( Z( 1 )*Z( 1 )+Z( 2 )*Z( 2 ) )
203: C = RHO*Z( 2 )*Z( 2 )*DELSQ
204: *
205: * The following TAU is DSIGMA * DSIGMA - D( 2 ) * D( 2 )
206: *
207: IF( B.GT.ZERO ) THEN
208: TAU = ( B+SQRT( B*B+FOUR*C ) ) / TWO
209: ELSE
210: TAU = TWO*C / ( -B+SQRT( B*B+FOUR*C ) )
211: END IF
212: *
213: * The following TAU is DSIGMA - D( 2 )
214: *
215: TAU = TAU / ( D( 2 )+SQRT( D( 2 )*D( 2 )+TAU ) )
216: DSIGMA = D( 2 ) + TAU
217: DELTA( 1 ) = -( DEL+TAU )
218: DELTA( 2 ) = -TAU
219: WORK( 1 ) = D( 1 ) + TAU + D( 2 )
220: WORK( 2 ) = TWO*D( 2 ) + TAU
221: * DELTA( 1 ) = -Z( 1 ) / ( DEL+TAU )
222: * DELTA( 2 ) = -Z( 2 ) / TAU
223: * TEMP = SQRT( DELTA( 1 )*DELTA( 1 )+DELTA( 2 )*DELTA( 2 ) )
224: * DELTA( 1 ) = DELTA( 1 ) / TEMP
225: * DELTA( 2 ) = DELTA( 2 ) / TEMP
226: END IF
227: RETURN
228: *
229: * End of DLASD5
230: *
231: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dlasd5.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack