Annotation of rpl/lapack/lapack/dlasv2.f, revision 1.8
1.8 ! bertrand 1: *> \brief \b DLASV2
! 2: *
! 3: * =========== DOCUMENTATION ===========
! 4: *
! 5: * Online html documentation available at
! 6: * http://www.netlib.org/lapack/explore-html/
! 7: *
! 8: *> \htmlonly
! 9: *> Download DLASV2 + dependencies
! 10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlasv2.f">
! 11: *> [TGZ]</a>
! 12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlasv2.f">
! 13: *> [ZIP]</a>
! 14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlasv2.f">
! 15: *> [TXT]</a>
! 16: *> \endhtmlonly
! 17: *
! 18: * Definition:
! 19: * ===========
! 20: *
! 21: * SUBROUTINE DLASV2( F, G, H, SSMIN, SSMAX, SNR, CSR, SNL, CSL )
! 22: *
! 23: * .. Scalar Arguments ..
! 24: * DOUBLE PRECISION CSL, CSR, F, G, H, SNL, SNR, SSMAX, SSMIN
! 25: * ..
! 26: *
! 27: *
! 28: *> \par Purpose:
! 29: * =============
! 30: *>
! 31: *> \verbatim
! 32: *>
! 33: *> DLASV2 computes the singular value decomposition of a 2-by-2
! 34: *> triangular matrix
! 35: *> [ F G ]
! 36: *> [ 0 H ].
! 37: *> On return, abs(SSMAX) is the larger singular value, abs(SSMIN) is the
! 38: *> smaller singular value, and (CSL,SNL) and (CSR,SNR) are the left and
! 39: *> right singular vectors for abs(SSMAX), giving the decomposition
! 40: *>
! 41: *> [ CSL SNL ] [ F G ] [ CSR -SNR ] = [ SSMAX 0 ]
! 42: *> [-SNL CSL ] [ 0 H ] [ SNR CSR ] [ 0 SSMIN ].
! 43: *> \endverbatim
! 44: *
! 45: * Arguments:
! 46: * ==========
! 47: *
! 48: *> \param[in] F
! 49: *> \verbatim
! 50: *> F is DOUBLE PRECISION
! 51: *> The (1,1) element of the 2-by-2 matrix.
! 52: *> \endverbatim
! 53: *>
! 54: *> \param[in] G
! 55: *> \verbatim
! 56: *> G is DOUBLE PRECISION
! 57: *> The (1,2) element of the 2-by-2 matrix.
! 58: *> \endverbatim
! 59: *>
! 60: *> \param[in] H
! 61: *> \verbatim
! 62: *> H is DOUBLE PRECISION
! 63: *> The (2,2) element of the 2-by-2 matrix.
! 64: *> \endverbatim
! 65: *>
! 66: *> \param[out] SSMIN
! 67: *> \verbatim
! 68: *> SSMIN is DOUBLE PRECISION
! 69: *> abs(SSMIN) is the smaller singular value.
! 70: *> \endverbatim
! 71: *>
! 72: *> \param[out] SSMAX
! 73: *> \verbatim
! 74: *> SSMAX is DOUBLE PRECISION
! 75: *> abs(SSMAX) is the larger singular value.
! 76: *> \endverbatim
! 77: *>
! 78: *> \param[out] SNL
! 79: *> \verbatim
! 80: *> SNL is DOUBLE PRECISION
! 81: *> \endverbatim
! 82: *>
! 83: *> \param[out] CSL
! 84: *> \verbatim
! 85: *> CSL is DOUBLE PRECISION
! 86: *> The vector (CSL, SNL) is a unit left singular vector for the
! 87: *> singular value abs(SSMAX).
! 88: *> \endverbatim
! 89: *>
! 90: *> \param[out] SNR
! 91: *> \verbatim
! 92: *> SNR is DOUBLE PRECISION
! 93: *> \endverbatim
! 94: *>
! 95: *> \param[out] CSR
! 96: *> \verbatim
! 97: *> CSR is DOUBLE PRECISION
! 98: *> The vector (CSR, SNR) is a unit right singular vector for the
! 99: *> singular value abs(SSMAX).
! 100: *> \endverbatim
! 101: *
! 102: * Authors:
! 103: * ========
! 104: *
! 105: *> \author Univ. of Tennessee
! 106: *> \author Univ. of California Berkeley
! 107: *> \author Univ. of Colorado Denver
! 108: *> \author NAG Ltd.
! 109: *
! 110: *> \date November 2011
! 111: *
! 112: *> \ingroup auxOTHERauxiliary
! 113: *
! 114: *> \par Further Details:
! 115: * =====================
! 116: *>
! 117: *> \verbatim
! 118: *>
! 119: *> Any input parameter may be aliased with any output parameter.
! 120: *>
! 121: *> Barring over/underflow and assuming a guard digit in subtraction, all
! 122: *> output quantities are correct to within a few units in the last
! 123: *> place (ulps).
! 124: *>
! 125: *> In IEEE arithmetic, the code works correctly if one matrix element is
! 126: *> infinite.
! 127: *>
! 128: *> Overflow will not occur unless the largest singular value itself
! 129: *> overflows or is within a few ulps of overflow. (On machines with
! 130: *> partial overflow, like the Cray, overflow may occur if the largest
! 131: *> singular value is within a factor of 2 of overflow.)
! 132: *>
! 133: *> Underflow is harmless if underflow is gradual. Otherwise, results
! 134: *> may correspond to a matrix modified by perturbations of size near
! 135: *> the underflow threshold.
! 136: *> \endverbatim
! 137: *>
! 138: * =====================================================================
1.1 bertrand 139: SUBROUTINE DLASV2( F, G, H, SSMIN, SSMAX, SNR, CSR, SNL, CSL )
140: *
1.8 ! bertrand 141: * -- LAPACK auxiliary routine (version 3.4.0) --
1.1 bertrand 142: * -- LAPACK is a software package provided by Univ. of Tennessee, --
143: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
1.8 ! bertrand 144: * November 2011
1.1 bertrand 145: *
146: * .. Scalar Arguments ..
147: DOUBLE PRECISION CSL, CSR, F, G, H, SNL, SNR, SSMAX, SSMIN
148: * ..
149: *
150: * =====================================================================
151: *
152: * .. Parameters ..
153: DOUBLE PRECISION ZERO
154: PARAMETER ( ZERO = 0.0D0 )
155: DOUBLE PRECISION HALF
156: PARAMETER ( HALF = 0.5D0 )
157: DOUBLE PRECISION ONE
158: PARAMETER ( ONE = 1.0D0 )
159: DOUBLE PRECISION TWO
160: PARAMETER ( TWO = 2.0D0 )
161: DOUBLE PRECISION FOUR
162: PARAMETER ( FOUR = 4.0D0 )
163: * ..
164: * .. Local Scalars ..
165: LOGICAL GASMAL, SWAP
166: INTEGER PMAX
167: DOUBLE PRECISION A, CLT, CRT, D, FA, FT, GA, GT, HA, HT, L, M,
168: $ MM, R, S, SLT, SRT, T, TEMP, TSIGN, TT
169: * ..
170: * .. Intrinsic Functions ..
171: INTRINSIC ABS, SIGN, SQRT
172: * ..
173: * .. External Functions ..
174: DOUBLE PRECISION DLAMCH
175: EXTERNAL DLAMCH
176: * ..
177: * .. Executable Statements ..
178: *
179: FT = F
180: FA = ABS( FT )
181: HT = H
182: HA = ABS( H )
183: *
184: * PMAX points to the maximum absolute element of matrix
185: * PMAX = 1 if F largest in absolute values
186: * PMAX = 2 if G largest in absolute values
187: * PMAX = 3 if H largest in absolute values
188: *
189: PMAX = 1
190: SWAP = ( HA.GT.FA )
191: IF( SWAP ) THEN
192: PMAX = 3
193: TEMP = FT
194: FT = HT
195: HT = TEMP
196: TEMP = FA
197: FA = HA
198: HA = TEMP
199: *
200: * Now FA .ge. HA
201: *
202: END IF
203: GT = G
204: GA = ABS( GT )
205: IF( GA.EQ.ZERO ) THEN
206: *
207: * Diagonal matrix
208: *
209: SSMIN = HA
210: SSMAX = FA
211: CLT = ONE
212: CRT = ONE
213: SLT = ZERO
214: SRT = ZERO
215: ELSE
216: GASMAL = .TRUE.
217: IF( GA.GT.FA ) THEN
218: PMAX = 2
219: IF( ( FA / GA ).LT.DLAMCH( 'EPS' ) ) THEN
220: *
221: * Case of very large GA
222: *
223: GASMAL = .FALSE.
224: SSMAX = GA
225: IF( HA.GT.ONE ) THEN
226: SSMIN = FA / ( GA / HA )
227: ELSE
228: SSMIN = ( FA / GA )*HA
229: END IF
230: CLT = ONE
231: SLT = HT / GT
232: SRT = ONE
233: CRT = FT / GT
234: END IF
235: END IF
236: IF( GASMAL ) THEN
237: *
238: * Normal case
239: *
240: D = FA - HA
241: IF( D.EQ.FA ) THEN
242: *
243: * Copes with infinite F or H
244: *
245: L = ONE
246: ELSE
247: L = D / FA
248: END IF
249: *
250: * Note that 0 .le. L .le. 1
251: *
252: M = GT / FT
253: *
254: * Note that abs(M) .le. 1/macheps
255: *
256: T = TWO - L
257: *
258: * Note that T .ge. 1
259: *
260: MM = M*M
261: TT = T*T
262: S = SQRT( TT+MM )
263: *
264: * Note that 1 .le. S .le. 1 + 1/macheps
265: *
266: IF( L.EQ.ZERO ) THEN
267: R = ABS( M )
268: ELSE
269: R = SQRT( L*L+MM )
270: END IF
271: *
272: * Note that 0 .le. R .le. 1 + 1/macheps
273: *
274: A = HALF*( S+R )
275: *
276: * Note that 1 .le. A .le. 1 + abs(M)
277: *
278: SSMIN = HA / A
279: SSMAX = FA*A
280: IF( MM.EQ.ZERO ) THEN
281: *
282: * Note that M is very tiny
283: *
284: IF( L.EQ.ZERO ) THEN
285: T = SIGN( TWO, FT )*SIGN( ONE, GT )
286: ELSE
287: T = GT / SIGN( D, FT ) + M / T
288: END IF
289: ELSE
290: T = ( M / ( S+T )+M / ( R+L ) )*( ONE+A )
291: END IF
292: L = SQRT( T*T+FOUR )
293: CRT = TWO / L
294: SRT = T / L
295: CLT = ( CRT+SRT*M ) / A
296: SLT = ( HT / FT )*SRT / A
297: END IF
298: END IF
299: IF( SWAP ) THEN
300: CSL = SRT
301: SNL = CRT
302: CSR = SLT
303: SNR = CLT
304: ELSE
305: CSL = CLT
306: SNL = SLT
307: CSR = CRT
308: SNR = SRT
309: END IF
310: *
311: * Correct signs of SSMAX and SSMIN
312: *
313: IF( PMAX.EQ.1 )
314: $ TSIGN = SIGN( ONE, CSR )*SIGN( ONE, CSL )*SIGN( ONE, F )
315: IF( PMAX.EQ.2 )
316: $ TSIGN = SIGN( ONE, SNR )*SIGN( ONE, CSL )*SIGN( ONE, G )
317: IF( PMAX.EQ.3 )
318: $ TSIGN = SIGN( ONE, SNR )*SIGN( ONE, SNL )*SIGN( ONE, H )
319: SSMAX = SIGN( SSMAX, TSIGN )
320: SSMIN = SIGN( SSMIN, TSIGN*SIGN( ONE, F )*SIGN( ONE, H ) )
321: RETURN
322: *
323: * End of DLASV2
324: *
325: END
CVSweb interface <joel.bertrand@systella.fr>
![[BACK]](/cvsweb/icons/back.gif){kind=icon}
Return to dlasv2.f CVS log ![[TXT]](/cvsweb/icons/text.gif){kind=icon}
![[DIR]](/cvsweb/icons/dir.gif){kind=icon}
Up to [local] / rpl / lapack / lapack