dlas2.f

Go to the documentation of this file.

*> \brief \b DLAS2 computes singular values of a 2-by-2 triangular matrix.
*
*  =========== DOCUMENTATION ===========
*
* Online html documentation available at
*            http://www.netlib.org/lapack/explore-html/
*
*> \htmlonly
*> Download DLAS2 + dependencies
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dlas2.f">
*> [TGZ]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dlas2.f">
*> [ZIP]</a>
*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dlas2.f">
*> [TXT]</a>
*> \endhtmlonly
*
*  Definition:
*  ===========
*
*       SUBROUTINE DLAS2( F, G, H, SSMIN, SSMAX )
*
*       .. Scalar Arguments ..
*       DOUBLE PRECISION   F, G, H, SSMAX, SSMIN
*       ..
*
*
*> \par Purpose:
*  =============
*>
*> \verbatim
*>
*> DLAS2  computes the singular values of the 2-by-2 matrix
*>    [  F   G  ]
*>    [  0   H  ].
*> On return, SSMIN is the smaller singular value and SSMAX is the
*> larger singular value.
*> \endverbatim
*
*  Arguments:
*  ==========
*
*> \param[in] F
*> \verbatim
*>          F is DOUBLE PRECISION
*>          The (1,1) element of the 2-by-2 matrix.
*> \endverbatim
*>
*> \param[in] G
*> \verbatim
*>          G is DOUBLE PRECISION
*>          The (1,2) element of the 2-by-2 matrix.
*> \endverbatim
*>
*> \param[in] H
*> \verbatim
*>          H is DOUBLE PRECISION
*>          The (2,2) element of the 2-by-2 matrix.
*> \endverbatim
*>
*> \param[out] SSMIN
*> \verbatim
*>          SSMIN is DOUBLE PRECISION
*>          The smaller singular value.
*> \endverbatim
*>
*> \param[out] SSMAX
*> \verbatim
*>          SSMAX is DOUBLE PRECISION
*>          The larger singular value.
*> \endverbatim
*
*  Authors:
*  ========
*
*> \author Univ. of Tennessee
*> \author Univ. of California Berkeley
*> \author Univ. of Colorado Denver
*> \author NAG Ltd.
*
*> \ingroup las2
*
*> \par Further Details:
*  =====================
*>
*> \verbatim
*>
*>  Barring over/underflow, all output quantities are correct to within
*>  a few units in the last place (ulps), even in the absence of a guard
*>  digit in addition/subtraction.
*>
*>  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.
*>
*>  Underflow is harmless if underflow is gradual. Otherwise, results
*>  may correspond to a matrix modified by perturbations of size near
*>  the underflow threshold.
*> \endverbatim
*>
*  =====================================================================
      SUBROUTINE dlas2( F, G, H, SSMIN, SSMAX )
*
*  -- LAPACK auxiliary routine --
*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
*
*     .. Scalar Arguments ..
      DOUBLE PRECISION   F, G, H, SSMAX, SSMIN
*     ..
*
*  ====================================================================
*
*     .. Parameters ..
      DOUBLE PRECISION   ZERO
      parameter( zero = 0.0d0 )
      DOUBLE PRECISION   ONE
      parameter( one = 1.0d0 )
      DOUBLE PRECISION   TWO
      parameter( two = 2.0d0 )
*     ..
*     .. Local Scalars ..
      DOUBLE PRECISION   AS, AT, AU, C, FA, FHMN, FHMX, GA, HA
*     ..
*     .. Intrinsic Functions ..
      INTRINSIC          abs, max, min, sqrt
*     ..
*     .. Executable Statements ..
*
      fa = abs( f )
      ga = abs( g )
      ha = abs( h )
      fhmn = min( fa, ha )
      fhmx = max( fa, ha )
      IF( fhmn.EQ.zero ) THEN
         ssmin = zero
         IF( fhmx.EQ.zero ) THEN
            ssmax = ga
         ELSE
            ssmax = max( fhmx, ga )*sqrt( one+
     $              ( min( fhmx, ga ) / max( fhmx, ga ) )**2 )
         END IF
      ELSE
         IF( ga.LT.fhmx ) THEN
            as = one + fhmn / fhmx
            at = ( fhmx-fhmn ) / fhmx
            au = ( ga / fhmx )**2
            c = two / ( sqrt( as*as+au )+sqrt( at*at+au ) )
            ssmin = fhmn*c
            ssmax = fhmx / c
         ELSE
            au = fhmx / ga
            IF( au.EQ.zero ) THEN
*
*              Avoid possible harmful underflow if exponent range
*              asymmetric (true SSMIN may not underflow even if
*              AU underflows)
*
               ssmin = ( fhmn*fhmx ) / ga
               ssmax = ga
            ELSE
               as = one + fhmn / fhmx
               at = ( fhmx-fhmn ) / fhmx
               c = one / ( sqrt( one+( as*au )**2 )+
     $             sqrt( one+( at*au )**2 ) )
               ssmin = ( fhmn*c )*au
               ssmin = ssmin + ssmin
               ssmax = ga / ( c+c )
            END IF
         END IF
      END IF
      RETURN
*
*     End of DLAS2
*
      END