1 SUBROUTINE SLAE2( A, B, C, RT1, RT2 ) 2* 3* -- LAPACK auxiliary routine (version 3.0) -- 4* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., 5* Courant Institute, Argonne National Lab, and Rice University 6* October 31, 1992 7* 8* .. Scalar Arguments .. 9 REAL A, B, C, RT1, RT2 10* .. 11* 12* Purpose 13* ======= 14* 15* SLAE2 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) REAL 25* The (1,1) element of the 2-by-2 matrix. 26* 27* B (input) REAL 28* The (1,2) and (2,1) elements of the 2-by-2 matrix. 29* 30* C (input) REAL 31* The (2,2) element of the 2-by-2 matrix. 32* 33* RT1 (output) REAL 34* The eigenvalue of larger absolute value. 35* 36* RT2 (output) REAL 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 REAL ONE 57 PARAMETER ( ONE = 1.0E0 ) 58 REAL TWO 59 PARAMETER ( TWO = 2.0E0 ) 60 REAL ZERO 61 PARAMETER ( ZERO = 0.0E0 ) 62 REAL HALF 63 PARAMETER ( HALF = 0.5E0 ) 64* .. 65* .. Local Scalars .. 66 REAL 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 SLAE2 123* 124 END 125