Return to dlae2.f CVS log

Up to [local] / rpl / lapack / lapack

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

1.1     ! bertrand    1:       SUBROUTINE DLAE2( A, B, C, RT1, RT2 )
        !             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   A, B, C, RT1, RT2
        !            10: *     ..
        !            11: *
        !            12: *  Purpose
        !            13: *  =======
        !            14: *
        !            15: *  DLAE2  computes the eigenvalues of a 2-by-2 symmetric matrix
        !            16: *     [  A   B  ]
        !            17: *     [  B   C  ].
        !            18: *  On return, RT1 is the eigenvalue of larger absolute value, and RT2
        !            19: *  is the eigenvalue of smaller absolute value.
        !            20: *
        !            21: *  Arguments
        !            22: *  =========
        !            23: *
        !            24: *  A       (input) DOUBLE PRECISION
        !            25: *          The (1,1) element of the 2-by-2 matrix.
        !            26: *
        !            27: *  B       (input) DOUBLE PRECISION
        !            28: *          The (1,2) and (2,1) elements of the 2-by-2 matrix.
        !            29: *
        !            30: *  C       (input) DOUBLE PRECISION
        !            31: *          The (2,2) element of the 2-by-2 matrix.
        !            32: *
        !            33: *  RT1     (output) DOUBLE PRECISION
        !            34: *          The eigenvalue of larger absolute value.
        !            35: *
        !            36: *  RT2     (output) DOUBLE PRECISION
        !            37: *          The eigenvalue of smaller absolute value.
        !            38: *
        !            39: *  Further Details
        !            40: *  ===============
        !            41: *
        !            42: *  RT1 is accurate to a few ulps barring over/underflow.
        !            43: *
        !            44: *  RT2 may be inaccurate if there is massive cancellation in the
        !            45: *  determinant A*C-B*B; higher precision or correctly rounded or
        !            46: *  correctly truncated arithmetic would be needed to compute RT2
        !            47: *  accurately in all cases.
        !            48: *
        !            49: *  Overflow is possible only if RT1 is within a factor of 5 of overflow.
        !            50: *  Underflow is harmless if the input data is 0 or exceeds
        !            51: *     underflow_threshold / macheps.
        !            52: *
        !            53: * =====================================================================
        !            54: *
        !            55: *     .. Parameters ..
        !            56:       DOUBLE PRECISION   ONE
        !            57:       PARAMETER          ( ONE = 1.0D0 )
        !            58:       DOUBLE PRECISION   TWO
        !            59:       PARAMETER          ( TWO = 2.0D0 )
        !            60:       DOUBLE PRECISION   ZERO
        !            61:       PARAMETER          ( ZERO = 0.0D0 )
        !            62:       DOUBLE PRECISION   HALF
        !            63:       PARAMETER          ( HALF = 0.5D0 )
        !            64: *     ..
        !            65: *     .. Local Scalars ..
        !            66:       DOUBLE PRECISION   AB, ACMN, ACMX, ADF, DF, RT, SM, TB
        !            67: *     ..
        !            68: *     .. Intrinsic Functions ..
        !            69:       INTRINSIC          ABS, SQRT
        !            70: *     ..
        !            71: *     .. Executable Statements ..
        !            72: *
        !            73: *     Compute the eigenvalues
        !            74: *
        !            75:       SM = A + C
        !            76:       DF = A - C
        !            77:       ADF = ABS( DF )
        !            78:       TB = B + B
        !            79:       AB = ABS( TB )
        !            80:       IF( ABS( A ).GT.ABS( C ) ) THEN
        !            81:          ACMX = A
        !            82:          ACMN = C
        !            83:       ELSE
        !            84:          ACMX = C
        !            85:          ACMN = A
        !            86:       END IF
        !            87:       IF( ADF.GT.AB ) THEN
        !            88:          RT = ADF*SQRT( ONE+( AB / ADF )**2 )
        !            89:       ELSE IF( ADF.LT.AB ) THEN
        !            90:          RT = AB*SQRT( ONE+( ADF / AB )**2 )
        !            91:       ELSE
        !            92: *
        !            93: *        Includes case AB=ADF=0
        !            94: *
        !            95:          RT = AB*SQRT( TWO )
        !            96:       END IF
        !            97:       IF( SM.LT.ZERO ) THEN
        !            98:          RT1 = HALF*( SM-RT )
        !            99: *
        !           100: *        Order of execution important.
        !           101: *        To get fully accurate smaller eigenvalue,
        !           102: *        next line needs to be executed in higher precision.
        !           103: *
        !           104:          RT2 = ( ACMX / RT1 )*ACMN - ( B / RT1 )*B
        !           105:       ELSE IF( SM.GT.ZERO ) THEN
        !           106:          RT1 = HALF*( SM+RT )
        !           107: *
        !           108: *        Order of execution important.
        !           109: *        To get fully accurate smaller eigenvalue,
        !           110: *        next line needs to be executed in higher precision.
        !           111: *
        !           112:          RT2 = ( ACMX / RT1 )*ACMN - ( B / RT1 )*B
        !           113:       ELSE
        !           114: *
        !           115: *        Includes case RT1 = RT2 = 0
        !           116: *
        !           117:          RT1 = HALF*RT
        !           118:          RT2 = -HALF*RT
        !           119:       END IF
        !           120:       RETURN
        !           121: *
        !           122: *     End of DLAE2
        !           123: *
        !           124:       END