LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
cget22.f
Go to the documentation of this file.
1*> \brief \b CGET22
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 CGET22( TRANSA, TRANSE, TRANSW, N, A, LDA, E, LDE, W,
12* WORK, RWORK, RESULT )
13*
14* .. Scalar Arguments ..
15* CHARACTER TRANSA, TRANSE, TRANSW
16* INTEGER LDA, LDE, N
17* ..
18* .. Array Arguments ..
19* REAL RESULT( 2 ), RWORK( * )
20* COMPLEX A( LDA, * ), E( LDE, * ), W( * ), WORK( * )
21* ..
22*
23*
24*> \par Purpose:
25* =============
26*>
27*> \verbatim
28*>
29*> CGET22 does an eigenvector check.
30*>
31*> The basic test is:
32*>
33*> RESULT(1) = | A E - E W | / ( |A| |E| ulp )
34*>
35*> using the 1-norm. It also tests the normalization of E:
36*>
37*> RESULT(2) = max | m-norm(E(j)) - 1 | / ( n ulp )
38*> j
39*>
40*> where E(j) is the j-th eigenvector, and m-norm is the max-norm of a
41*> vector. The max-norm of a complex n-vector x in this case is the
42*> maximum of |re(x(i)| + |im(x(i)| over i = 1, ..., n.
43*> \endverbatim
44*
45* Arguments:
46* ==========
47*
48*> \param[in] TRANSA
49*> \verbatim
50*> TRANSA is CHARACTER*1
51*> Specifies whether or not A is transposed.
52*> = 'N': No transpose
53*> = 'T': Transpose
54*> = 'C': Conjugate transpose
55*> \endverbatim
56*>
57*> \param[in] TRANSE
58*> \verbatim
59*> TRANSE is CHARACTER*1
60*> Specifies whether or not E is transposed.
61*> = 'N': No transpose, eigenvectors are in columns of E
62*> = 'T': Transpose, eigenvectors are in rows of E
63*> = 'C': Conjugate transpose, eigenvectors are in rows of E
64*> \endverbatim
65*>
66*> \param[in] TRANSW
67*> \verbatim
68*> TRANSW is CHARACTER*1
69*> Specifies whether or not W is transposed.
70*> = 'N': No transpose
71*> = 'T': Transpose, same as TRANSW = 'N'
72*> = 'C': Conjugate transpose, use -WI(j) instead of WI(j)
73*> \endverbatim
74*>
75*> \param[in] N
76*> \verbatim
77*> N is INTEGER
78*> The order of the matrix A. N >= 0.
79*> \endverbatim
80*>
81*> \param[in] A
82*> \verbatim
83*> A is COMPLEX array, dimension (LDA,N)
84*> The matrix whose eigenvectors are in E.
85*> \endverbatim
86*>
87*> \param[in] LDA
88*> \verbatim
89*> LDA is INTEGER
90*> The leading dimension of the array A. LDA >= max(1,N).
91*> \endverbatim
92*>
93*> \param[in] E
94*> \verbatim
95*> E is COMPLEX array, dimension (LDE,N)
96*> The matrix of eigenvectors. If TRANSE = 'N', the eigenvectors
97*> are stored in the columns of E, if TRANSE = 'T' or 'C', the
98*> eigenvectors are stored in the rows of E.
99*> \endverbatim
100*>
101*> \param[in] LDE
102*> \verbatim
103*> LDE is INTEGER
104*> The leading dimension of the array E. LDE >= max(1,N).
105*> \endverbatim
106*>
107*> \param[in] W
108*> \verbatim
109*> W is COMPLEX array, dimension (N)
110*> The eigenvalues of A.
111*> \endverbatim
112*>
113*> \param[out] WORK
114*> \verbatim
115*> WORK is COMPLEX array, dimension (N*N)
116*> \endverbatim
117*>
118*> \param[out] RWORK
119*> \verbatim
120*> RWORK is REAL array, dimension (N)
121*> \endverbatim
122*>
123*> \param[out] RESULT
124*> \verbatim
125*> RESULT is REAL array, dimension (2)
126*> RESULT(1) = | A E - E W | / ( |A| |E| ulp )
127*> RESULT(2) = max | m-norm(E(j)) - 1 | / ( n ulp )
128*> j
129*> \endverbatim
130*
131* Authors:
132* ========
133*
134*> \author Univ. of Tennessee
135*> \author Univ. of California Berkeley
136*> \author Univ. of Colorado Denver
137*> \author NAG Ltd.
138*
139*> \ingroup complex_eig
140*
141* =====================================================================
142 SUBROUTINE cget22( TRANSA, TRANSE, TRANSW, N, A, LDA, E, LDE, W,
143 $ WORK, RWORK, RESULT )
144*
145* -- LAPACK test routine --
146* -- LAPACK is a software package provided by Univ. of Tennessee, --
147* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
148*
149* .. Scalar Arguments ..
150 CHARACTER TRANSA, TRANSE, TRANSW
151 INTEGER LDA, LDE, N
152* ..
153* .. Array Arguments ..
154 REAL RESULT( 2 ), RWORK( * )
155 COMPLEX A( LDA, * ), E( LDE, * ), W( * ), WORK( * )
156* ..
157*
158* =====================================================================
159*
160* .. Parameters ..
161 REAL ZERO, ONE
162 parameter( zero = 0.0e+0, one = 1.0e+0 )
163 COMPLEX CZERO, CONE
164 parameter( czero = ( 0.0e+0, 0.0e+0 ),
165 $ cone = ( 1.0e+0, 0.0e+0 ) )
166* ..
167* .. Local Scalars ..
168 CHARACTER NORMA, NORME
169 INTEGER ITRNSE, ITRNSW, J, JCOL, JOFF, JROW, JVEC
170 REAL ANORM, ENORM, ENRMAX, ENRMIN, ERRNRM, TEMP1,
171 $ ulp, unfl
172 COMPLEX WTEMP
173* ..
174* .. External Functions ..
175 LOGICAL LSAME
176 REAL CLANGE, SLAMCH
177 EXTERNAL lsame, clange, slamch
178* ..
179* .. External Subroutines ..
180 EXTERNAL cgemm, claset
181* ..
182* .. Intrinsic Functions ..
183 INTRINSIC abs, aimag, conjg, max, min, real
184* ..
185* .. Executable Statements ..
186*
187* Initialize RESULT (in case N=0)
188*
189 result( 1 ) = zero
190 result( 2 ) = zero
191 IF( n.LE.0 )
192 $ RETURN
193*
194 unfl = slamch( 'Safe minimum' )
195 ulp = slamch( 'Precision' )
196*
197 itrnse = 0
198 itrnsw = 0
199 norma = 'O'
200 norme = 'O'
201*
202 IF( lsame( transa, 'T' ) .OR. lsame( transa, 'C' ) ) THEN
203 norma = 'I'
204 END IF
205*
206 IF( lsame( transe, 'T' ) ) THEN
207 itrnse = 1
208 norme = 'I'
209 ELSE IF( lsame( transe, 'C' ) ) THEN
210 itrnse = 2
211 norme = 'I'
212 END IF
213*
214 IF( lsame( transw, 'C' ) ) THEN
215 itrnsw = 1
216 END IF
217*
218* Normalization of E:
219*
220 enrmin = one / ulp
221 enrmax = zero
222 IF( itrnse.EQ.0 ) THEN
223 DO 20 jvec = 1, n
224 temp1 = zero
225 DO 10 j = 1, n
226 temp1 = max( temp1, abs( real( e( j, jvec ) ) )+
227 $ abs( aimag( e( j, jvec ) ) ) )
228 10 CONTINUE
229 enrmin = min( enrmin, temp1 )
230 enrmax = max( enrmax, temp1 )
231 20 CONTINUE
232 ELSE
233 DO 30 jvec = 1, n
234 rwork( jvec ) = zero
235 30 CONTINUE
236*
237 DO 50 j = 1, n
238 DO 40 jvec = 1, n
239 rwork( jvec ) = max( rwork( jvec ),
240 $ abs( real( e( jvec, j ) ) )+
241 $ abs( aimag( e( jvec, j ) ) ) )
242 40 CONTINUE
243 50 CONTINUE
244*
245 DO 60 jvec = 1, n
246 enrmin = min( enrmin, rwork( jvec ) )
247 enrmax = max( enrmax, rwork( jvec ) )
248 60 CONTINUE
249 END IF
250*
251* Norm of A:
252*
253 anorm = max( clange( norma, n, n, a, lda, rwork ), unfl )
254*
255* Norm of E:
256*
257 enorm = max( clange( norme, n, n, e, lde, rwork ), ulp )
258*
259* Norm of error:
260*
261* Error = AE - EW
262*
263 CALL claset( 'Full', n, n, czero, czero, work, n )
264*
265 joff = 0
266 DO 100 jcol = 1, n
267 IF( itrnsw.EQ.0 ) THEN
268 wtemp = w( jcol )
269 ELSE
270 wtemp = conjg( w( jcol ) )
271 END IF
272*
273 IF( itrnse.EQ.0 ) THEN
274 DO 70 jrow = 1, n
275 work( joff+jrow ) = e( jrow, jcol )*wtemp
276 70 CONTINUE
277 ELSE IF( itrnse.EQ.1 ) THEN
278 DO 80 jrow = 1, n
279 work( joff+jrow ) = e( jcol, jrow )*wtemp
280 80 CONTINUE
281 ELSE
282 DO 90 jrow = 1, n
283 work( joff+jrow ) = conjg( e( jcol, jrow ) )*wtemp
284 90 CONTINUE
285 END IF
286 joff = joff + n
287 100 CONTINUE
288*
289 CALL cgemm( transa, transe, n, n, n, cone, a, lda, e, lde, -cone,
290 $ work, n )
291*
292 errnrm = clange( 'One', n, n, work, n, rwork ) / enorm
293*
294* Compute RESULT(1) (avoiding under/overflow)
295*
296 IF( anorm.GT.errnrm ) THEN
297 result( 1 ) = ( errnrm / anorm ) / ulp
298 ELSE
299 IF( anorm.LT.one ) THEN
300 result( 1 ) = one / ulp
301 ELSE
302 result( 1 ) = min( errnrm / anorm, one ) / ulp
303 END IF
304 END IF
305*
306* Compute RESULT(2) : the normalization error in E.
307*
308 result( 2 ) = max( abs( enrmax-one ), abs( enrmin-one ) ) /
309 $ ( real( n )*ulp )
310*
311 RETURN
312*
313* End of CGET22
314*
315 END
subroutine cget22(transa, transe, transw, n, a, lda, e, lde, w, work, rwork, result)
CGET22
Definition cget22.f:144
subroutine cgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
CGEMM
Definition cgemm.f:188
subroutine claset(uplo, m, n, alpha, beta, a, lda)
CLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition claset.f:106