LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
claesy.f
Go to the documentation of this file.
1 *> \brief \b CLAESY computes the eigenvalues and eigenvectors of a 2-by-2 complex symmetric matrix.
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download CLAESY + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/claesy.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/claesy.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/claesy.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE CLAESY( A, B, C, RT1, RT2, EVSCAL, CS1, SN1 )
22 *
23 * .. Scalar Arguments ..
24 * COMPLEX A, B, C, CS1, EVSCAL, RT1, RT2, SN1
25 * ..
26 *
27 *
28 *> \par Purpose:
29 * =============
30 *>
31 *> \verbatim
32 *>
33 *> CLAESY computes the eigendecomposition of a 2-by-2 symmetric matrix
34 *> ( ( A, B );( B, C ) )
35 *> provided the norm of the matrix of eigenvectors is larger than
36 *> some threshold value.
37 *>
38 *> RT1 is the eigenvalue of larger absolute value, and RT2 of
39 *> smaller absolute value. If the eigenvectors are computed, then
40 *> on return ( CS1, SN1 ) is the unit eigenvector for RT1, hence
41 *>
42 *> [ CS1 SN1 ] . [ A B ] . [ CS1 -SN1 ] = [ RT1 0 ]
43 *> [ -SN1 CS1 ] [ B C ] [ SN1 CS1 ] [ 0 RT2 ]
44 *> \endverbatim
45 *
46 * Arguments:
47 * ==========
48 *
49 *> \param[in] A
50 *> \verbatim
51 *> A is COMPLEX
52 *> The ( 1, 1 ) element of input matrix.
53 *> \endverbatim
54 *>
55 *> \param[in] B
56 *> \verbatim
57 *> B is COMPLEX
58 *> The ( 1, 2 ) element of input matrix. The ( 2, 1 ) element
59 *> is also given by B, since the 2-by-2 matrix is symmetric.
60 *> \endverbatim
61 *>
62 *> \param[in] C
63 *> \verbatim
64 *> C is COMPLEX
65 *> The ( 2, 2 ) element of input matrix.
66 *> \endverbatim
67 *>
68 *> \param[out] RT1
69 *> \verbatim
70 *> RT1 is COMPLEX
71 *> The eigenvalue of larger modulus.
72 *> \endverbatim
73 *>
74 *> \param[out] RT2
75 *> \verbatim
76 *> RT2 is COMPLEX
77 *> The eigenvalue of smaller modulus.
78 *> \endverbatim
79 *>
80 *> \param[out] EVSCAL
81 *> \verbatim
82 *> EVSCAL is COMPLEX
83 *> The complex value by which the eigenvector matrix was scaled
84 *> to make it orthonormal. If EVSCAL is zero, the eigenvectors
85 *> were not computed. This means one of two things: the 2-by-2
86 *> matrix could not be diagonalized, or the norm of the matrix
87 *> of eigenvectors before scaling was larger than the threshold
88 *> value THRESH (set below).
89 *> \endverbatim
90 *>
91 *> \param[out] CS1
92 *> \verbatim
93 *> CS1 is COMPLEX
94 *> \endverbatim
95 *>
96 *> \param[out] SN1
97 *> \verbatim
98 *> SN1 is COMPLEX
99 *> If EVSCAL .NE. 0, ( CS1, SN1 ) is the unit right eigenvector
100 *> for RT1.
101 *> \endverbatim
102 *
103 * Authors:
104 * ========
105 *
106 *> \author Univ. of Tennessee
107 *> \author Univ. of California Berkeley
108 *> \author Univ. of Colorado Denver
109 *> \author NAG Ltd.
110 *
111 *> \ingroup complexSYauxiliary
112 *
113 * =====================================================================
114  SUBROUTINE claesy( A, B, C, RT1, RT2, EVSCAL, CS1, SN1 )
115 *
116 * -- LAPACK auxiliary routine --
117 * -- LAPACK is a software package provided by Univ. of Tennessee, --
118 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
119 *
120 * .. Scalar Arguments ..
121  COMPLEX A, B, C, CS1, EVSCAL, RT1, RT2, SN1
122 * ..
123 *
124 * =====================================================================
125 *
126 * .. Parameters ..
127  REAL ZERO
128  parameter( zero = 0.0e0 )
129  REAL ONE
130  parameter( one = 1.0e0 )
131  COMPLEX CONE
132  parameter( cone = ( 1.0e0, 0.0e0 ) )
133  REAL HALF
134  parameter( half = 0.5e0 )
135  REAL THRESH
136  parameter( thresh = 0.1e0 )
137 * ..
138 * .. Local Scalars ..
139  REAL BABS, EVNORM, TABS, Z
140  COMPLEX S, T, TMP
141 * ..
142 * .. Intrinsic Functions ..
143  INTRINSIC abs, max, sqrt
144 * ..
145 * .. Executable Statements ..
146 *
147 *
148 * Special case: The matrix is actually diagonal.
149 * To avoid divide by zero later, we treat this case separately.
150 *
151  IF( abs( b ).EQ.zero ) THEN
152  rt1 = a
153  rt2 = c
154  IF( abs( rt1 ).LT.abs( rt2 ) ) THEN
155  tmp = rt1
156  rt1 = rt2
157  rt2 = tmp
158  cs1 = zero
159  sn1 = one
160  ELSE
161  cs1 = one
162  sn1 = zero
163  END IF
164  ELSE
165 *
166 * Compute the eigenvalues and eigenvectors.
167 * The characteristic equation is
168 * lambda **2 - (A+C) lambda + (A*C - B*B)
169 * and we solve it using the quadratic formula.
170 *
171  s = ( a+c )*half
172  t = ( a-c )*half
173 *
174 * Take the square root carefully to avoid over/under flow.
175 *
176  babs = abs( b )
177  tabs = abs( t )
178  z = max( babs, tabs )
179  IF( z.GT.zero )
180  $ t = z*sqrt( ( t / z )**2+( b / z )**2 )
181 *
182 * Compute the two eigenvalues. RT1 and RT2 are exchanged
183 * if necessary so that RT1 will have the greater magnitude.
184 *
185  rt1 = s + t
186  rt2 = s - t
187  IF( abs( rt1 ).LT.abs( rt2 ) ) THEN
188  tmp = rt1
189  rt1 = rt2
190  rt2 = tmp
191  END IF
192 *
193 * Choose CS1 = 1 and SN1 to satisfy the first equation, then
194 * scale the components of this eigenvector so that the matrix
195 * of eigenvectors X satisfies X * X**T = I . (No scaling is
196 * done if the norm of the eigenvalue matrix is less than THRESH.)
197 *
198  sn1 = ( rt1-a ) / b
199  tabs = abs( sn1 )
200  IF( tabs.GT.one ) THEN
201  t = tabs*sqrt( ( one / tabs )**2+( sn1 / tabs )**2 )
202  ELSE
203  t = sqrt( cone+sn1*sn1 )
204  END IF
205  evnorm = abs( t )
206  IF( evnorm.GE.thresh ) THEN
207  evscal = cone / t
208  cs1 = evscal
209  sn1 = sn1*evscal
210  ELSE
211  evscal = zero
212  END IF
213  END IF
214  RETURN
215 *
216 * End of CLAESY
217 *
218  END
subroutine claesy(A, B, C, RT1, RT2, EVSCAL, CS1, SN1)
CLAESY computes the eigenvalues and eigenvectors of a 2-by-2 complex symmetric matrix.
Definition: claesy.f:115