1*> \brief \b SLAEV2 computes the eigenvalues and eigenvectors of a 2-by-2 symmetric/Hermitian matrix.
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download SLAEV2 + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/slaev2.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/slaev2.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/slaev2.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18*  Definition:
19*  ===========
20*
21*       SUBROUTINE SLAEV2( A, B, C, RT1, RT2, CS1, SN1 )
22*
23*       .. Scalar Arguments ..
24*       REAL               A, B, C, CS1, RT1, RT2, SN1
25*       ..
26*
27*
28*> \par Purpose:
29*  =============
30*>
31*> \verbatim
32*>
33*> SLAEV2 computes the eigendecomposition of a 2-by-2 symmetric matrix
34*>    [  A   B  ]
35*>    [  B   C  ].
36*> On return, RT1 is the eigenvalue of larger absolute value, RT2 is the
37*> eigenvalue of smaller absolute value, and (CS1,SN1) is the unit right
38*> eigenvector for RT1, giving the decomposition
39*>
40*>    [ CS1  SN1 ] [  A   B  ] [ CS1 -SN1 ]  =  [ RT1  0  ]
41*>    [-SN1  CS1 ] [  B   C  ] [ SN1  CS1 ]     [  0  RT2 ].
42*> \endverbatim
43*
44*  Arguments:
45*  ==========
46*
47*> \param[in] A
48*> \verbatim
49*>          A is REAL
50*>          The (1,1) element of the 2-by-2 matrix.
51*> \endverbatim
52*>
53*> \param[in] B
54*> \verbatim
55*>          B is REAL
56*>          The (1,2) element and the conjugate of the (2,1) element of
57*>          the 2-by-2 matrix.
58*> \endverbatim
59*>
60*> \param[in] C
61*> \verbatim
62*>          C is REAL
63*>          The (2,2) element of the 2-by-2 matrix.
64*> \endverbatim
65*>
66*> \param[out] RT1
67*> \verbatim
68*>          RT1 is REAL
69*>          The eigenvalue of larger absolute value.
70*> \endverbatim
71*>
72*> \param[out] RT2
73*> \verbatim
74*>          RT2 is REAL
75*>          The eigenvalue of smaller absolute value.
76*> \endverbatim
77*>
78*> \param[out] CS1
79*> \verbatim
80*>          CS1 is REAL
81*> \endverbatim
82*>
83*> \param[out] SN1
84*> \verbatim
85*>          SN1 is REAL
86*>          The vector (CS1, SN1) is a unit right eigenvector for RT1.
87*> \endverbatim
88*
89*  Authors:
90*  ========
91*
92*> \author Univ. of Tennessee
93*> \author Univ. of California Berkeley
94*> \author Univ. of Colorado Denver
95*> \author NAG Ltd.
96*
97*> \ingroup OTHERauxiliary
98*
99*> \par Further Details:
100*  =====================
101*>
102*> \verbatim
103*>
104*>  RT1 is accurate to a few ulps barring over/underflow.
105*>
106*>  RT2 may be inaccurate if there is massive cancellation in the
107*>  determinant A*C-B*B; higher precision or correctly rounded or
108*>  correctly truncated arithmetic would be needed to compute RT2
109*>  accurately in all cases.
110*>
111*>  CS1 and SN1 are accurate to a few ulps barring over/underflow.
112*>
113*>  Overflow is possible only if RT1 is within a factor of 5 of overflow.
114*>  Underflow is harmless if the input data is 0 or exceeds
115*>     underflow_threshold / macheps.
116*> \endverbatim
117*>
118*  =====================================================================
119      SUBROUTINE SLAEV2( A, B, C, RT1, RT2, CS1, SN1 )
120*
121*  -- LAPACK auxiliary routine --
122*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
123*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
124*
125*     .. Scalar Arguments ..
126      REAL               A, B, C, CS1, RT1, RT2, SN1
127*     ..
128*
129* =====================================================================
130*
131*     .. Parameters ..
132      REAL               ONE
133      PARAMETER          ( ONE = 1.0E0 )
134      REAL               TWO
135      PARAMETER          ( TWO = 2.0E0 )
136      REAL               ZERO
137      PARAMETER          ( ZERO = 0.0E0 )
138      REAL               HALF
139      PARAMETER          ( HALF = 0.5E0 )
140*     ..
141*     .. Local Scalars ..
142      INTEGER            SGN1, SGN2
143      REAL               AB, ACMN, ACMX, ACS, ADF, CS, CT, DF, RT, SM,
144     $                   TB, TN
145*     ..
146*     .. Intrinsic Functions ..
147      INTRINSIC          ABS, SQRT
148*     ..
149*     .. Executable Statements ..
150*
151*     Compute the eigenvalues
152*
153      SM = A + C
154      DF = A - C
155      ADF = ABS( DF )
156      TB = B + B
157      AB = ABS( TB )
158      IF( ABS( A ).GT.ABS( C ) ) THEN
159         ACMX = A
160         ACMN = C
161      ELSE
162         ACMX = C
163         ACMN = A
164      END IF
165      IF( ADF.GT.AB ) THEN
166         RT = ADF*SQRT( ONE+( AB / ADF )**2 )
167      ELSE IF( ADF.LT.AB ) THEN
168         RT = AB*SQRT( ONE+( ADF / AB )**2 )
169      ELSE
170*
171*        Includes case AB=ADF=0
172*
173         RT = AB*SQRT( TWO )
174      END IF
175      IF( SM.LT.ZERO ) THEN
176         RT1 = HALF*( SM-RT )
177         SGN1 = -1
178*
179*        Order of execution important.
180*        To get fully accurate smaller eigenvalue,
181*        next line needs to be executed in higher precision.
182*
183         RT2 = ( ACMX / RT1 )*ACMN - ( B / RT1 )*B
184      ELSE IF( SM.GT.ZERO ) THEN
185         RT1 = HALF*( SM+RT )
186         SGN1 = 1
187*
188*        Order of execution important.
189*        To get fully accurate smaller eigenvalue,
190*        next line needs to be executed in higher precision.
191*
192         RT2 = ( ACMX / RT1 )*ACMN - ( B / RT1 )*B
193      ELSE
194*
195*        Includes case RT1 = RT2 = 0
196*
197         RT1 = HALF*RT
198         RT2 = -HALF*RT
199         SGN1 = 1
200      END IF
201*
202*     Compute the eigenvector
203*
204      IF( DF.GE.ZERO ) THEN
205         CS = DF + RT
206         SGN2 = 1
207      ELSE
208         CS = DF - RT
209         SGN2 = -1
210      END IF
211      ACS = ABS( CS )
212      IF( ACS.GT.AB ) THEN
213         CT = -TB / CS
214         SN1 = ONE / SQRT( ONE+CT*CT )
215         CS1 = CT*SN1
216      ELSE
217         IF( AB.EQ.ZERO ) THEN
218            CS1 = ONE
219            SN1 = ZERO
220         ELSE
221            TN = -CS / TB
222            CS1 = ONE / SQRT( ONE+TN*TN )
223            SN1 = TN*CS1
224         END IF
225      END IF
226      IF( SGN1.EQ.SGN2 ) THEN
227         TN = CS1
228         CS1 = -SN1
229         SN1 = TN
230      END IF
231      RETURN
232*
233*     End of SLAEV2
234*
235      END
236