LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
 All Classes Files Functions Variables Typedefs Macros
cget52.f
Go to the documentation of this file.
1 *> \brief \b CGET52
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 * Definition:
9 * ===========
10 *
11 * SUBROUTINE CGET52( LEFT, N, A, LDA, B, LDB, E, LDE, ALPHA, BETA,
12 * WORK, RWORK, RESULT )
13 *
14 * .. Scalar Arguments ..
15 * LOGICAL LEFT
16 * INTEGER LDA, LDB, LDE, N
17 * ..
18 * .. Array Arguments ..
19 * REAL RESULT( 2 ), RWORK( * )
20 * COMPLEX A( LDA, * ), ALPHA( * ), B( LDB, * ),
21 * $ BETA( * ), E( LDE, * ), WORK( * )
22 * ..
23 *
24 *
25 *> \par Purpose:
26 * =============
27 *>
28 *> \verbatim
29 *>
30 *> CGET52 does an eigenvector check for the generalized eigenvalue
31 *> problem.
32 *>
33 *> The basic test for right eigenvectors is:
34 *>
35 *> | b(i) A E(i) - a(i) B E(i) |
36 *> RESULT(1) = max -------------------------------
37 *> i n ulp max( |b(i) A|, |a(i) B| )
38 *>
39 *> using the 1-norm. Here, a(i)/b(i) = w is the i-th generalized
40 *> eigenvalue of A - w B, or, equivalently, b(i)/a(i) = m is the i-th
41 *> generalized eigenvalue of m A - B.
42 *>
43 *> H H _ _
44 *> For left eigenvectors, A , B , a, and b are used.
45 *>
46 *> CGET52 also tests the normalization of E. Each eigenvector is
47 *> supposed to be normalized so that the maximum "absolute value"
48 *> of its elements is 1, where in this case, "absolute value"
49 *> of a complex value x is |Re(x)| + |Im(x)| ; let us call this
50 *> maximum "absolute value" norm of a vector v M(v).
51 *> if a(i)=b(i)=0, then the eigenvector is set to be the jth coordinate
52 *> vector. The normalization test is:
53 *>
54 *> RESULT(2) = max | M(v(i)) - 1 | / ( n ulp )
55 *> eigenvectors v(i)
56 *> \endverbatim
57 *
58 * Arguments:
59 * ==========
60 *
61 *> \param[in] LEFT
62 *> \verbatim
63 *> LEFT is LOGICAL
64 *> =.TRUE.: The eigenvectors in the columns of E are assumed
65 *> to be *left* eigenvectors.
66 *> =.FALSE.: The eigenvectors in the columns of E are assumed
67 *> to be *right* eigenvectors.
68 *> \endverbatim
69 *>
70 *> \param[in] N
71 *> \verbatim
72 *> N is INTEGER
73 *> The size of the matrices. If it is zero, CGET52 does
74 *> nothing. It must be at least zero.
75 *> \endverbatim
76 *>
77 *> \param[in] A
78 *> \verbatim
79 *> A is COMPLEX array, dimension (LDA, N)
80 *> The matrix A.
81 *> \endverbatim
82 *>
83 *> \param[in] LDA
84 *> \verbatim
85 *> LDA is INTEGER
86 *> The leading dimension of A. It must be at least 1
87 *> and at least N.
88 *> \endverbatim
89 *>
90 *> \param[in] B
91 *> \verbatim
92 *> B is COMPLEX array, dimension (LDB, N)
93 *> The matrix B.
94 *> \endverbatim
95 *>
96 *> \param[in] LDB
97 *> \verbatim
98 *> LDB is INTEGER
99 *> The leading dimension of B. It must be at least 1
100 *> and at least N.
101 *> \endverbatim
102 *>
103 *> \param[in] E
104 *> \verbatim
105 *> E is COMPLEX array, dimension (LDE, N)
106 *> The matrix of eigenvectors. It must be O( 1 ).
107 *> \endverbatim
108 *>
109 *> \param[in] LDE
110 *> \verbatim
111 *> LDE is INTEGER
112 *> The leading dimension of E. It must be at least 1 and at
113 *> least N.
114 *> \endverbatim
115 *>
116 *> \param[in] ALPHA
117 *> \verbatim
118 *> ALPHA is COMPLEX array, dimension (N)
119 *> The values a(i) as described above, which, along with b(i),
120 *> define the generalized eigenvalues.
121 *> \endverbatim
122 *>
123 *> \param[in] BETA
124 *> \verbatim
125 *> BETA is COMPLEX array, dimension (N)
126 *> The values b(i) as described above, which, along with a(i),
127 *> define the generalized eigenvalues.
128 *> \endverbatim
129 *>
130 *> \param[out] WORK
131 *> \verbatim
132 *> WORK is COMPLEX array, dimension (N**2)
133 *> \endverbatim
134 *>
135 *> \param[out] RWORK
136 *> \verbatim
137 *> RWORK is REAL array, dimension (N)
138 *> \endverbatim
139 *>
140 *> \param[out] RESULT
141 *> \verbatim
142 *> RESULT is REAL array, dimension (2)
143 *> The values computed by the test described above. If A E or
144 *> B E is likely to overflow, then RESULT(1:2) is set to
145 *> 10 / ulp.
146 *> \endverbatim
147 *
148 * Authors:
149 * ========
150 *
151 *> \author Univ. of Tennessee
152 *> \author Univ. of California Berkeley
153 *> \author Univ. of Colorado Denver
154 *> \author NAG Ltd.
155 *
156 *> \date November 2011
157 *
158 *> \ingroup complex_eig
159 *
160 * =====================================================================
161  SUBROUTINE cget52( LEFT, N, A, LDA, B, LDB, E, LDE, ALPHA, BETA,
162  $ work, rwork, result )
163 *
164 * -- LAPACK test routine (version 3.4.0) --
165 * -- LAPACK is a software package provided by Univ. of Tennessee, --
166 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
167 * November 2011
168 *
169 * .. Scalar Arguments ..
170  LOGICAL left
171  INTEGER lda, ldb, lde, n
172 * ..
173 * .. Array Arguments ..
174  REAL result( 2 ), rwork( * )
175  COMPLEX a( lda, * ), alpha( * ), b( ldb, * ),
176  $ beta( * ), e( lde, * ), work( * )
177 * ..
178 *
179 * =====================================================================
180 *
181 * .. Parameters ..
182  REAL zero, one
183  parameter( zero = 0.0e+0, one = 1.0e+0 )
184  COMPLEX czero, cone
185  parameter( czero = ( 0.0e+0, 0.0e+0 ),
186  $ cone = ( 1.0e+0, 0.0e+0 ) )
187 * ..
188 * .. Local Scalars ..
189  CHARACTER normab, trans
190  INTEGER j, jvec
191  REAL abmax, alfmax, anorm, betmax, bnorm, enorm,
192  $ enrmer, errnrm, safmax, safmin, scale, temp1,
193  $ ulp
194  COMPLEX acoeff, alphai, bcoeff, betai, x
195 * ..
196 * .. External Functions ..
197  REAL clange, slamch
198  EXTERNAL clange, slamch
199 * ..
200 * .. External Subroutines ..
201  EXTERNAL cgemv
202 * ..
203 * .. Intrinsic Functions ..
204  INTRINSIC abs, aimag, conjg, max, real
205 * ..
206 * .. Statement Functions ..
207  REAL abs1
208 * ..
209 * .. Statement Function definitions ..
210  abs1( x ) = abs( REAL( X ) ) + abs( aimag( x ) )
211 * ..
212 * .. Executable Statements ..
213 *
214  result( 1 ) = zero
215  result( 2 ) = zero
216  IF( n.LE.0 )
217  $ RETURN
218 *
219  safmin = slamch( 'Safe minimum' )
220  safmax = one / safmin
221  ulp = slamch( 'Epsilon' )*slamch( 'Base' )
222 *
223  IF( left ) THEN
224  trans = 'C'
225  normab = 'I'
226  ELSE
227  trans = 'N'
228  normab = 'O'
229  END IF
230 *
231 * Norm of A, B, and E:
232 *
233  anorm = max( clange( normab, n, n, a, lda, rwork ), safmin )
234  bnorm = max( clange( normab, n, n, b, ldb, rwork ), safmin )
235  enorm = max( clange( 'O', n, n, e, lde, rwork ), ulp )
236  alfmax = safmax / max( one, bnorm )
237  betmax = safmax / max( one, anorm )
238 *
239 * Compute error matrix.
240 * Column i = ( b(i) A - a(i) B ) E(i) / max( |a(i) B| |b(i) A| )
241 *
242  DO 10 jvec = 1, n
243  alphai = alpha( jvec )
244  betai = beta( jvec )
245  abmax = max( abs1( alphai ), abs1( betai ) )
246  IF( abs1( alphai ).GT.alfmax .OR. abs1( betai ).GT.betmax .OR.
247  $ abmax.LT.one ) THEN
248  scale = one / max( abmax, safmin )
249  alphai = scale*alphai
250  betai = scale*betai
251  END IF
252  scale = one / max( abs1( alphai )*bnorm, abs1( betai )*anorm,
253  $ safmin )
254  acoeff = scale*betai
255  bcoeff = scale*alphai
256  IF( left ) THEN
257  acoeff = conjg( acoeff )
258  bcoeff = conjg( bcoeff )
259  END IF
260  CALL cgemv( trans, n, n, acoeff, a, lda, e( 1, jvec ), 1,
261  $ czero, work( n*( jvec-1 )+1 ), 1 )
262  CALL cgemv( trans, n, n, -bcoeff, b, lda, e( 1, jvec ), 1,
263  $ cone, work( n*( jvec-1 )+1 ), 1 )
264  10 CONTINUE
265 *
266  errnrm = clange( 'One', n, n, work, n, rwork ) / enorm
267 *
268 * Compute RESULT(1)
269 *
270  result( 1 ) = errnrm / ulp
271 *
272 * Normalization of E:
273 *
274  enrmer = zero
275  DO 30 jvec = 1, n
276  temp1 = zero
277  DO 20 j = 1, n
278  temp1 = max( temp1, abs1( e( j, jvec ) ) )
279  20 CONTINUE
280  enrmer = max( enrmer, temp1-one )
281  30 CONTINUE
282 *
283 * Compute RESULT(2) : the normalization error in E.
284 *
285  result( 2 ) = enrmer / ( REAL( n )*ulp )
286 *
287  RETURN
288 *
289 * End of CGET52
290 *
291  END