LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
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.
Date
September 2012
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 140 of file dlasv2.f.

140 *
141 * -- LAPACK auxiliary routine (version 3.4.2) --
142 * -- LAPACK is a software package provided by Univ. of Tennessee, --
143 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
144 * September 2012
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 *
double precision function dlamch(CMACH)
DLAMCH
Definition: dlamch.f:65

Here is the caller graph for this function: