LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dlasv2.f
Go to the documentation of this file.
1*> \brief \b DLASV2 computes the singular value decomposition of a 2-by-2 triangular matrix.
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*> \ingroup lasv2
111*
112*> \par Further Details:
113* =====================
114*>
115*> \verbatim
116*>
117*> Any input parameter may be aliased with any output parameter.
118*>
119*> Barring over/underflow and assuming a guard digit in subtraction, all
120*> output quantities are correct to within a few units in the last
121*> place (ulps).
122*>
123*> In IEEE arithmetic, the code works correctly if one matrix element is
124*> infinite.
125*>
126*> Overflow will not occur unless the largest singular value itself
127*> overflows or is within a few ulps of overflow.
128*>
129*> Underflow is harmless if underflow is gradual. Otherwise, results
130*> may correspond to a matrix modified by perturbations of size near
131*> the underflow threshold.
132*> \endverbatim
133*>
134* =====================================================================
135 SUBROUTINE dlasv2( F, G, H, SSMIN, SSMAX, SNR, CSR, SNL, CSL )
136*
137* -- LAPACK auxiliary routine --
138* -- LAPACK is a software package provided by Univ. of Tennessee, --
139* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
140*
141* .. Scalar Arguments ..
142 DOUBLE PRECISION CSL, CSR, F, G, H, SNL, SNR, SSMAX, SSMIN
143* ..
144*
145* =====================================================================
146*
147* .. Parameters ..
148 DOUBLE PRECISION ZERO
149 parameter( zero = 0.0d0 )
150 DOUBLE PRECISION HALF
151 parameter( half = 0.5d0 )
152 DOUBLE PRECISION ONE
153 parameter( one = 1.0d0 )
154 DOUBLE PRECISION TWO
155 parameter( two = 2.0d0 )
156 DOUBLE PRECISION FOUR
157 parameter( four = 4.0d0 )
158* ..
159* .. Local Scalars ..
160 LOGICAL GASMAL, SWAP
161 INTEGER PMAX
162 DOUBLE PRECISION A, CLT, CRT, D, FA, FT, GA, GT, HA, HT, L, M,
163 $ MM, R, S, SLT, SRT, T, TEMP, TSIGN, TT
164* ..
165* .. Intrinsic Functions ..
166 INTRINSIC abs, sign, sqrt
167* ..
168* .. External Functions ..
169 DOUBLE PRECISION DLAMCH
170 EXTERNAL dlamch
171* ..
172* .. Executable Statements ..
173*
174 ft = f
175 fa = abs( ft )
176 ht = h
177 ha = abs( h )
178*
179* PMAX points to the maximum absolute element of matrix
180* PMAX = 1 if F largest in absolute values
181* PMAX = 2 if G largest in absolute values
182* PMAX = 3 if H largest in absolute values
183*
184 pmax = 1
185 swap = ( ha.GT.fa )
186 IF( swap ) THEN
187 pmax = 3
188 temp = ft
189 ft = ht
190 ht = temp
191 temp = fa
192 fa = ha
193 ha = temp
194*
195* Now FA .ge. HA
196*
197 END IF
198 gt = g
199 ga = abs( gt )
200 IF( ga.EQ.zero ) THEN
201*
202* Diagonal matrix
203*
204 ssmin = ha
205 ssmax = fa
206 clt = one
207 crt = one
208 slt = zero
209 srt = zero
210 ELSE
211 gasmal = .true.
212 IF( ga.GT.fa ) THEN
213 pmax = 2
214 IF( ( fa / ga ).LT.dlamch( 'EPS' ) ) THEN
215*
216* Case of very large GA
217*
218 gasmal = .false.
219 ssmax = ga
220 IF( ha.GT.one ) THEN
221 ssmin = fa / ( ga / ha )
222 ELSE
223 ssmin = ( fa / ga )*ha
224 END IF
225 clt = one
226 slt = ht / gt
227 srt = one
228 crt = ft / gt
229 END IF
230 END IF
231 IF( gasmal ) THEN
232*
233* Normal case
234*
235 d = fa - ha
236 IF( d.EQ.fa ) THEN
237*
238* Copes with infinite F or H
239*
240 l = one
241 ELSE
242 l = d / fa
243 END IF
244*
245* Note that 0 .le. L .le. 1
246*
247 m = gt / ft
248*
249* Note that abs(M) .le. 1/macheps
250*
251 t = two - l
252*
253* Note that T .ge. 1
254*
255 mm = m*m
256 tt = t*t
257 s = sqrt( tt+mm )
258*
259* Note that 1 .le. S .le. 1 + 1/macheps
260*
261 IF( l.EQ.zero ) THEN
262 r = abs( m )
263 ELSE
264 r = sqrt( l*l+mm )
265 END IF
266*
267* Note that 0 .le. R .le. 1 + 1/macheps
268*
269 a = half*( s+r )
270*
271* Note that 1 .le. A .le. 1 + abs(M)
272*
273 ssmin = ha / a
274 ssmax = fa*a
275 IF( mm.EQ.zero ) THEN
276*
277* Note that M is very tiny
278*
279 IF( l.EQ.zero ) THEN
280 t = sign( two, ft )*sign( one, gt )
281 ELSE
282 t = gt / sign( d, ft ) + m / t
283 END IF
284 ELSE
285 t = ( m / ( s+t )+m / ( r+l ) )*( one+a )
286 END IF
287 l = sqrt( t*t+four )
288 crt = two / l
289 srt = t / l
290 clt = ( crt+srt*m ) / a
291 slt = ( ht / ft )*srt / a
292 END IF
293 END IF
294 IF( swap ) THEN
295 csl = srt
296 snl = crt
297 csr = slt
298 snr = clt
299 ELSE
300 csl = clt
301 snl = slt
302 csr = crt
303 snr = srt
304 END IF
305*
306* Correct signs of SSMAX and SSMIN
307*
308 IF( pmax.EQ.1 )
309 $ tsign = sign( one, csr )*sign( one, csl )*sign( one, f )
310 IF( pmax.EQ.2 )
311 $ tsign = sign( one, snr )*sign( one, csl )*sign( one, g )
312 IF( pmax.EQ.3 )
313 $ tsign = sign( one, snr )*sign( one, snl )*sign( one, h )
314 ssmax = sign( ssmax, tsign )
315 ssmin = sign( ssmin, tsign*sign( one, f )*sign( one, h ) )
316 RETURN
317*
318* End of DLASV2
319*
320 END
subroutine dlasv2(f, g, h, ssmin, ssmax, snr, csr, snl, csl)
DLASV2 computes the singular value decomposition of a 2-by-2 triangular matrix.
Definition dlasv2.f:136