LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zhst01.f
Go to the documentation of this file.
1*> \brief \b ZHST01
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 ZHST01( N, ILO, IHI, A, LDA, H, LDH, Q, LDQ, WORK,
12* LWORK, RWORK, RESULT )
13*
14* .. Scalar Arguments ..
15* INTEGER IHI, ILO, LDA, LDH, LDQ, LWORK, N
16* ..
17* .. Array Arguments ..
18* DOUBLE PRECISION RESULT( 2 ), RWORK( * )
19* COMPLEX*16 A( LDA, * ), H( LDH, * ), Q( LDQ, * ),
20* $ WORK( LWORK )
21* ..
22*
23*
24*> \par Purpose:
25* =============
26*>
27*> \verbatim
28*>
29*> ZHST01 tests the reduction of a general matrix A to upper Hessenberg
30*> form: A = Q*H*Q'. Two test ratios are computed;
31*>
32*> RESULT(1) = norm( A - Q*H*Q' ) / ( norm(A) * N * EPS )
33*> RESULT(2) = norm( I - Q'*Q ) / ( N * EPS )
34*>
35*> The matrix Q is assumed to be given explicitly as it would be
36*> following ZGEHRD + ZUNGHR.
37*>
38*> In this version, ILO and IHI are not used, but they could be used
39*> to save some work if this is desired.
40*> \endverbatim
41*
42* Arguments:
43* ==========
44*
45*> \param[in] N
46*> \verbatim
47*> N is INTEGER
48*> The order of the matrix A. N >= 0.
49*> \endverbatim
50*>
51*> \param[in] ILO
52*> \verbatim
53*> ILO is INTEGER
54*> \endverbatim
55*>
56*> \param[in] IHI
57*> \verbatim
58*> IHI is INTEGER
59*>
60*> A is assumed to be upper triangular in rows and columns
61*> 1:ILO-1 and IHI+1:N, so Q differs from the identity only in
62*> rows and columns ILO+1:IHI.
63*> \endverbatim
64*>
65*> \param[in] A
66*> \verbatim
67*> A is COMPLEX*16 array, dimension (LDA,N)
68*> The original n by n matrix A.
69*> \endverbatim
70*>
71*> \param[in] LDA
72*> \verbatim
73*> LDA is INTEGER
74*> The leading dimension of the array A. LDA >= max(1,N).
75*> \endverbatim
76*>
77*> \param[in] H
78*> \verbatim
79*> H is COMPLEX*16 array, dimension (LDH,N)
80*> The upper Hessenberg matrix H from the reduction A = Q*H*Q'
81*> as computed by ZGEHRD. H is assumed to be zero below the
82*> first subdiagonal.
83*> \endverbatim
84*>
85*> \param[in] LDH
86*> \verbatim
87*> LDH is INTEGER
88*> The leading dimension of the array H. LDH >= max(1,N).
89*> \endverbatim
90*>
91*> \param[in] Q
92*> \verbatim
93*> Q is COMPLEX*16 array, dimension (LDQ,N)
94*> The orthogonal matrix Q from the reduction A = Q*H*Q' as
95*> computed by ZGEHRD + ZUNGHR.
96*> \endverbatim
97*>
98*> \param[in] LDQ
99*> \verbatim
100*> LDQ is INTEGER
101*> The leading dimension of the array Q. LDQ >= max(1,N).
102*> \endverbatim
103*>
104*> \param[out] WORK
105*> \verbatim
106*> WORK is COMPLEX*16 array, dimension (LWORK)
107*> \endverbatim
108*>
109*> \param[in] LWORK
110*> \verbatim
111*> LWORK is INTEGER
112*> The length of the array WORK. LWORK >= 2*N*N.
113*> \endverbatim
114*>
115*> \param[out] RWORK
116*> \verbatim
117*> RWORK is DOUBLE PRECISION array, dimension (N)
118*> \endverbatim
119*>
120*> \param[out] RESULT
121*> \verbatim
122*> RESULT is DOUBLE PRECISION array, dimension (2)
123*> RESULT(1) = norm( A - Q*H*Q' ) / ( norm(A) * N * EPS )
124*> RESULT(2) = norm( I - Q'*Q ) / ( N * EPS )
125*> \endverbatim
126*
127* Authors:
128* ========
129*
130*> \author Univ. of Tennessee
131*> \author Univ. of California Berkeley
132*> \author Univ. of Colorado Denver
133*> \author NAG Ltd.
134*
135*> \ingroup complex16_eig
136*
137* =====================================================================
138 SUBROUTINE zhst01( N, ILO, IHI, A, LDA, H, LDH, Q, LDQ, WORK,
139 $ LWORK, RWORK, RESULT )
140*
141* -- LAPACK test routine --
142* -- LAPACK is a software package provided by Univ. of Tennessee, --
143* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
144*
145* .. Scalar Arguments ..
146 INTEGER IHI, ILO, LDA, LDH, LDQ, LWORK, N
147* ..
148* .. Array Arguments ..
149 DOUBLE PRECISION RESULT( 2 ), RWORK( * )
150 COMPLEX*16 A( LDA, * ), H( LDH, * ), Q( LDQ, * ),
151 $ work( lwork )
152* ..
153*
154* =====================================================================
155*
156* .. Parameters ..
157 DOUBLE PRECISION ONE, ZERO
158 parameter( one = 1.0d+0, zero = 0.0d+0 )
159* ..
160* .. Local Scalars ..
161 INTEGER LDWORK
162 DOUBLE PRECISION ANORM, EPS, OVFL, SMLNUM, UNFL, WNORM
163* ..
164* .. External Functions ..
165 DOUBLE PRECISION DLAMCH, ZLANGE
166 EXTERNAL dlamch, zlange
167* ..
168* .. External Subroutines ..
169 EXTERNAL dlabad, zgemm, zlacpy, zunt01
170* ..
171* .. Intrinsic Functions ..
172 INTRINSIC dcmplx, max, min
173* ..
174* .. Executable Statements ..
175*
176* Quick return if possible
177*
178 IF( n.LE.0 ) THEN
179 result( 1 ) = zero
180 result( 2 ) = zero
181 RETURN
182 END IF
183*
184 unfl = dlamch( 'Safe minimum' )
185 eps = dlamch( 'Precision' )
186 ovfl = one / unfl
187 CALL dlabad( unfl, ovfl )
188 smlnum = unfl*n / eps
189*
190* Test 1: Compute norm( A - Q*H*Q' ) / ( norm(A) * N * EPS )
191*
192* Copy A to WORK
193*
194 ldwork = max( 1, n )
195 CALL zlacpy( ' ', n, n, a, lda, work, ldwork )
196*
197* Compute Q*H
198*
199 CALL zgemm( 'No transpose', 'No transpose', n, n, n,
200 $ dcmplx( one ), q, ldq, h, ldh, dcmplx( zero ),
201 $ work( ldwork*n+1 ), ldwork )
202*
203* Compute A - Q*H*Q'
204*
205 CALL zgemm( 'No transpose', 'Conjugate transpose', n, n, n,
206 $ dcmplx( -one ), work( ldwork*n+1 ), ldwork, q, ldq,
207 $ dcmplx( one ), work, ldwork )
208*
209 anorm = max( zlange( '1', n, n, a, lda, rwork ), unfl )
210 wnorm = zlange( '1', n, n, work, ldwork, rwork )
211*
212* Note that RESULT(1) cannot overflow and is bounded by 1/(N*EPS)
213*
214 result( 1 ) = min( wnorm, anorm ) / max( smlnum, anorm*eps ) / n
215*
216* Test 2: Compute norm( I - Q'*Q ) / ( N * EPS )
217*
218 CALL zunt01( 'Columns', n, n, q, ldq, work, lwork, rwork,
219 $ result( 2 ) )
220*
221 RETURN
222*
223* End of ZHST01
224*
225 END
subroutine dlabad(SMALL, LARGE)
DLABAD
Definition: dlabad.f:74
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:187
subroutine zhst01(N, ILO, IHI, A, LDA, H, LDH, Q, LDQ, WORK, LWORK, RWORK, RESULT)
ZHST01
Definition: zhst01.f:140
subroutine zunt01(ROWCOL, M, N, U, LDU, WORK, LWORK, RWORK, RESID)
ZUNT01
Definition: zunt01.f:126
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:103