LAPACK  3.10.0
LAPACK: Linear Algebra PACKage

◆ dlasv2()

subroutine dlasv2 ( double precision  F,
double precision  G,
double precision  H,
double precision  SSMIN,
double precision  SSMAX,
double precision  SNR,
double precision  CSR,
double precision  SNL,
double precision  CSL 
)

DLASV2 computes the singular value decomposition of a 2-by-2 triangular matrix.

Download DLASV2 + dependencies [TGZ] [ZIP] [TXT]

Purpose:
 DLASV2 computes the singular value decomposition of a 2-by-2
 triangular matrix
    [  F   G  ]
    [  0   H  ].
 On return, abs(SSMAX) is the larger singular value, abs(SSMIN) is the
 smaller singular value, and (CSL,SNL) and (CSR,SNR) are the left and
 right singular vectors for abs(SSMAX), giving the decomposition

    [ CSL  SNL ] [  F   G  ] [ CSR -SNR ]  =  [ SSMAX   0   ]
    [-SNL  CSL ] [  0   H  ] [ SNR  CSR ]     [  0    SSMIN ].
Parameters
[in]F
          F is DOUBLE PRECISION
          The (1,1) element of the 2-by-2 matrix.
[in]G
          G is DOUBLE PRECISION
          The (1,2) element of the 2-by-2 matrix.
[in]H
          H is DOUBLE PRECISION
          The (2,2) element of the 2-by-2 matrix.
[out]SSMIN
          SSMIN is DOUBLE PRECISION
          abs(SSMIN) is the smaller singular value.
[out]SSMAX
          SSMAX is DOUBLE PRECISION
          abs(SSMAX) is the larger singular value.
[out]SNL
          SNL is DOUBLE PRECISION
[out]CSL
          CSL is DOUBLE PRECISION
          The vector (CSL, SNL) is a unit left singular vector for the
          singular value abs(SSMAX).
[out]SNR
          SNR is DOUBLE PRECISION
[out]CSR
          CSR is DOUBLE PRECISION
          The vector (CSR, SNR) is a unit right singular vector for the
          singular value abs(SSMAX).
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
  Any input parameter may be aliased with any output parameter.

  Barring over/underflow and assuming a guard digit in subtraction, all
  output quantities are correct to within a few units in the last
  place (ulps).

  In IEEE arithmetic, the code works correctly if one matrix element is
  infinite.

  Overflow will not occur unless the largest singular value itself
  overflows or is within a few ulps of overflow. (On machines with
  partial overflow, like the Cray, overflow may occur if the largest
  singular value is within a factor of 2 of overflow.)

  Underflow is harmless if underflow is gradual. Otherwise, results
  may correspond to a matrix modified by perturbations of size near
  the underflow threshold.

Definition at line 137 of file dlasv2.f.

138 *
139 * -- LAPACK auxiliary routine --
140 * -- LAPACK is a software package provided by Univ. of Tennessee, --
141 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
142 *
143 * .. Scalar Arguments ..
144  DOUBLE PRECISION CSL, CSR, F, G, H, SNL, SNR, SSMAX, SSMIN
145 * ..
146 *
147 * =====================================================================
148 *
149 * .. Parameters ..
150  DOUBLE PRECISION ZERO
151  parameter( zero = 0.0d0 )
152  DOUBLE PRECISION HALF
153  parameter( half = 0.5d0 )
154  DOUBLE PRECISION ONE
155  parameter( one = 1.0d0 )
156  DOUBLE PRECISION TWO
157  parameter( two = 2.0d0 )
158  DOUBLE PRECISION FOUR
159  parameter( four = 4.0d0 )
160 * ..
161 * .. Local Scalars ..
162  LOGICAL GASMAL, SWAP
163  INTEGER PMAX
164  DOUBLE PRECISION A, CLT, CRT, D, FA, FT, GA, GT, HA, HT, L, M,
165  $ MM, R, S, SLT, SRT, T, TEMP, TSIGN, TT
166 * ..
167 * .. Intrinsic Functions ..
168  INTRINSIC abs, sign, sqrt
169 * ..
170 * .. External Functions ..
171  DOUBLE PRECISION DLAMCH
172  EXTERNAL dlamch
173 * ..
174 * .. Executable Statements ..
175 *
176  ft = f
177  fa = abs( ft )
178  ht = h
179  ha = abs( h )
180 *
181 * PMAX points to the maximum absolute element of matrix
182 * PMAX = 1 if F largest in absolute values
183 * PMAX = 2 if G largest in absolute values
184 * PMAX = 3 if H largest in absolute values
185 *
186  pmax = 1
187  swap = ( ha.GT.fa )
188  IF( swap ) THEN
189  pmax = 3
190  temp = ft
191  ft = ht
192  ht = temp
193  temp = fa
194  fa = ha
195  ha = temp
196 *
197 * Now FA .ge. HA
198 *
199  END IF
200  gt = g
201  ga = abs( gt )
202  IF( ga.EQ.zero ) THEN
203 *
204 * Diagonal matrix
205 *
206  ssmin = ha
207  ssmax = fa
208  clt = one
209  crt = one
210  slt = zero
211  srt = zero
212  ELSE
213  gasmal = .true.
214  IF( ga.GT.fa ) THEN
215  pmax = 2
216  IF( ( fa / ga ).LT.dlamch( 'EPS' ) ) THEN
217 *
218 * Case of very large GA
219 *
220  gasmal = .false.
221  ssmax = ga
222  IF( ha.GT.one ) THEN
223  ssmin = fa / ( ga / ha )
224  ELSE
225  ssmin = ( fa / ga )*ha
226  END IF
227  clt = one
228  slt = ht / gt
229  srt = one
230  crt = ft / gt
231  END IF
232  END IF
233  IF( gasmal ) THEN
234 *
235 * Normal case
236 *
237  d = fa - ha
238  IF( d.EQ.fa ) THEN
239 *
240 * Copes with infinite F or H
241 *
242  l = one
243  ELSE
244  l = d / fa
245  END IF
246 *
247 * Note that 0 .le. L .le. 1
248 *
249  m = gt / ft
250 *
251 * Note that abs(M) .le. 1/macheps
252 *
253  t = two - l
254 *
255 * Note that T .ge. 1
256 *
257  mm = m*m
258  tt = t*t
259  s = sqrt( tt+mm )
260 *
261 * Note that 1 .le. S .le. 1 + 1/macheps
262 *
263  IF( l.EQ.zero ) THEN
264  r = abs( m )
265  ELSE
266  r = sqrt( l*l+mm )
267  END IF
268 *
269 * Note that 0 .le. R .le. 1 + 1/macheps
270 *
271  a = half*( s+r )
272 *
273 * Note that 1 .le. A .le. 1 + abs(M)
274 *
275  ssmin = ha / a
276  ssmax = fa*a
277  IF( mm.EQ.zero ) THEN
278 *
279 * Note that M is very tiny
280 *
281  IF( l.EQ.zero ) THEN
282  t = sign( two, ft )*sign( one, gt )
283  ELSE
284  t = gt / sign( d, ft ) + m / t
285  END IF
286  ELSE
287  t = ( m / ( s+t )+m / ( r+l ) )*( one+a )
288  END IF
289  l = sqrt( t*t+four )
290  crt = two / l
291  srt = t / l
292  clt = ( crt+srt*m ) / a
293  slt = ( ht / ft )*srt / a
294  END IF
295  END IF
296  IF( swap ) THEN
297  csl = srt
298  snl = crt
299  csr = slt
300  snr = clt
301  ELSE
302  csl = clt
303  snl = slt
304  csr = crt
305  snr = srt
306  END IF
307 *
308 * Correct signs of SSMAX and SSMIN
309 *
310  IF( pmax.EQ.1 )
311  $ tsign = sign( one, csr )*sign( one, csl )*sign( one, f )
312  IF( pmax.EQ.2 )
313  $ tsign = sign( one, snr )*sign( one, csl )*sign( one, g )
314  IF( pmax.EQ.3 )
315  $ tsign = sign( one, snr )*sign( one, snl )*sign( one, h )
316  ssmax = sign( ssmax, tsign )
317  ssmin = sign( ssmin, tsign*sign( one, f )*sign( one, h ) )
318  RETURN
319 *
320 * End of DLASV2
321 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:69
Here is the caller graph for this function: