Return to dlas2.f CVS log

Up to [local] / rpl / lapack / lapack

Annotation of rpl/lapack/lapack/dlas2.f, revision 1.1

1.1     ! bertrand    1:       SUBROUTINE DLAS2( F, G, H, SSMIN, SSMAX )
        !             2: *
        !             3: *  -- LAPACK auxiliary routine (version 3.2) --
        !             4: *  -- LAPACK is a software package provided by Univ. of Tennessee,    --
        !             5: *  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
        !             6: *     November 2006
        !             7: *
        !             8: *     .. Scalar Arguments ..
        !             9:       DOUBLE PRECISION   F, G, H, SSMAX, SSMIN
        !            10: *     ..
        !            11: *
        !            12: *  Purpose
        !            13: *  =======
        !            14: *
        !            15: *  DLAS2  computes the singular values of the 2-by-2 matrix
        !            16: *     [  F   G  ]
        !            17: *     [  0   H  ].
        !            18: *  On return, SSMIN is the smaller singular value and SSMAX is the
        !            19: *  larger singular value.
        !            20: *
        !            21: *  Arguments
        !            22: *  =========
        !            23: *
        !            24: *  F       (input) DOUBLE PRECISION
        !            25: *          The (1,1) element of the 2-by-2 matrix.
        !            26: *
        !            27: *  G       (input) DOUBLE PRECISION
        !            28: *          The (1,2) element of the 2-by-2 matrix.
        !            29: *
        !            30: *  H       (input) DOUBLE PRECISION
        !            31: *          The (2,2) element of the 2-by-2 matrix.
        !            32: *
        !            33: *  SSMIN   (output) DOUBLE PRECISION
        !            34: *          The smaller singular value.
        !            35: *
        !            36: *  SSMAX   (output) DOUBLE PRECISION
        !            37: *          The larger singular value.
        !            38: *
        !            39: *  Further Details
        !            40: *  ===============
        !            41: *
        !            42: *  Barring over/underflow, all output quantities are correct to within
        !            43: *  a few units in the last place (ulps), even in the absence of a guard
        !            44: *  digit in addition/subtraction.
        !            45: *
        !            46: *  In IEEE arithmetic, the code works correctly if one matrix element is
        !            47: *  infinite.
        !            48: *
        !            49: *  Overflow will not occur unless the largest singular value itself
        !            50: *  overflows, or is within a few ulps of overflow. (On machines with
        !            51: *  partial overflow, like the Cray, overflow may occur if the largest
        !            52: *  singular value is within a factor of 2 of overflow.)
        !            53: *
        !            54: *  Underflow is harmless if underflow is gradual. Otherwise, results
        !            55: *  may correspond to a matrix modified by perturbations of size near
        !            56: *  the underflow threshold.
        !            57: *
        !            58: *  ====================================================================
        !            59: *
        !            60: *     .. Parameters ..
        !            61:       DOUBLE PRECISION   ZERO
        !            62:       PARAMETER          ( ZERO = 0.0D0 )
        !            63:       DOUBLE PRECISION   ONE
        !            64:       PARAMETER          ( ONE = 1.0D0 )
        !            65:       DOUBLE PRECISION   TWO
        !            66:       PARAMETER          ( TWO = 2.0D0 )
        !            67: *     ..
        !            68: *     .. Local Scalars ..
        !            69:       DOUBLE PRECISION   AS, AT, AU, C, FA, FHMN, FHMX, GA, HA
        !            70: *     ..
        !            71: *     .. Intrinsic Functions ..
        !            72:       INTRINSIC          ABS, MAX, MIN, SQRT
        !            73: *     ..
        !            74: *     .. Executable Statements ..
        !            75: *
        !            76:       FA = ABS( F )
        !            77:       GA = ABS( G )
        !            78:       HA = ABS( H )
        !            79:       FHMN = MIN( FA, HA )
        !            80:       FHMX = MAX( FA, HA )
        !            81:       IF( FHMN.EQ.ZERO ) THEN
        !            82:          SSMIN = ZERO
        !            83:          IF( FHMX.EQ.ZERO ) THEN
        !            84:             SSMAX = GA
        !            85:          ELSE
        !            86:             SSMAX = MAX( FHMX, GA )*SQRT( ONE+
        !            87:      $              ( MIN( FHMX, GA ) / MAX( FHMX, GA ) )**2 )
        !            88:          END IF
        !            89:       ELSE
        !            90:          IF( GA.LT.FHMX ) THEN
        !            91:             AS = ONE + FHMN / FHMX
        !            92:             AT = ( FHMX-FHMN ) / FHMX
        !            93:             AU = ( GA / FHMX )**2
        !            94:             C = TWO / ( SQRT( AS*AS+AU )+SQRT( AT*AT+AU ) )
        !            95:             SSMIN = FHMN*C
        !            96:             SSMAX = FHMX / C
        !            97:          ELSE
        !            98:             AU = FHMX / GA
        !            99:             IF( AU.EQ.ZERO ) THEN
        !           100: *
        !           101: *              Avoid possible harmful underflow if exponent range
        !           102: *              asymmetric (true SSMIN may not underflow even if
        !           103: *              AU underflows)
        !           104: *
        !           105:                SSMIN = ( FHMN*FHMX ) / GA
        !           106:                SSMAX = GA
        !           107:             ELSE
        !           108:                AS = ONE + FHMN / FHMX
        !           109:                AT = ( FHMX-FHMN ) / FHMX
        !           110:                C = ONE / ( SQRT( ONE+( AS*AU )**2 )+
        !           111:      $             SQRT( ONE+( AT*AU )**2 ) )
        !           112:                SSMIN = ( FHMN*C )*AU
        !           113:                SSMIN = SSMIN + SSMIN
        !           114:                SSMAX = GA / ( C+C )
        !           115:             END IF
        !           116:          END IF
        !           117:       END IF
        !           118:       RETURN
        !           119: *
        !           120: *     End of DLAS2
        !           121: *
        !           122:       END