LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zlqt05.f
Go to the documentation of this file.
1*> \brief \b ZLQT05
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 ZLQT05(M,N,L,NB,RESULT)
12*
13* .. Scalar Arguments ..
14* INTEGER LWORK, M, N, L, NB, LDT
15* .. Return values ..
16* DOUBLE PRECISION RESULT(6)
17*
18*
19*> \par Purpose:
20* =============
21*>
22*> \verbatim
23*>
24*> ZQRT05 tests ZTPLQT and ZTPMLQT.
25*> \endverbatim
26*
27* Arguments:
28* ==========
29*
30*> \param[in] M
31*> \verbatim
32*> M is INTEGER
33*> Number of rows in lower part of the test matrix.
34*> \endverbatim
35*>
36*> \param[in] N
37*> \verbatim
38*> N is INTEGER
39*> Number of columns in test matrix.
40*> \endverbatim
41*>
42*> \param[in] L
43*> \verbatim
44*> L is INTEGER
45*> The number of rows of the upper trapezoidal part the
46*> lower test matrix. 0 <= L <= M.
47*> \endverbatim
48*>
49*> \param[in] NB
50*> \verbatim
51*> NB is INTEGER
52*> Block size of test matrix. NB <= N.
53*> \endverbatim
54*>
55*> \param[out] RESULT
56*> \verbatim
57*> RESULT is DOUBLE PRECISION array, dimension (6)
58*> Results of each of the six tests below.
59*>
60*> RESULT(1) = | A - Q R |
61*> RESULT(2) = | I - Q^H Q |
62*> RESULT(3) = | Q C - Q C |
63*> RESULT(4) = | Q^H C - Q^H C |
64*> RESULT(5) = | C Q - C Q |
65*> RESULT(6) = | C Q^H - C Q^H |
66*> \endverbatim
67*
68* Authors:
69* ========
70*
71*> \author Univ. of Tennessee
72*> \author Univ. of California Berkeley
73*> \author Univ. of Colorado Denver
74*> \author NAG Ltd.
75*
76*> \ingroup double_lin
77*
78* =====================================================================
79 SUBROUTINE zlqt05(M,N,L,NB,RESULT)
80 IMPLICIT NONE
81*
82* -- LAPACK test routine --
83* -- LAPACK is a software package provided by Univ. of Tennessee, --
84* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
85*
86* .. Scalar Arguments ..
87 INTEGER LWORK, M, N, L, NB, LDT
88* .. Return values ..
89 DOUBLE PRECISION RESULT(6)
90*
91* =====================================================================
92*
93* ..
94* .. Local allocatable arrays
95 COMPLEX*16, ALLOCATABLE :: AF(:,:), Q(:,:),
96 $ R(:,:), WORK( : ), T(:,:),
97 $ CF(:,:), DF(:,:), A(:,:), C(:,:), D(:,:)
98 DOUBLE PRECISION, ALLOCATABLE :: RWORK(:)
99*
100* .. Parameters ..
101 DOUBLE PRECISION ZERO
102 COMPLEX*16 ONE, CZERO
103 parameter( zero = 0.0, one = (1.0,0.0), czero=(0.0,0.0) )
104* ..
105* .. Local Scalars ..
106 INTEGER INFO, J, K, N2, NP1,i
107 DOUBLE PRECISION ANORM, EPS, RESID, CNORM, DNORM
108* ..
109* .. Local Arrays ..
110 INTEGER ISEED( 4 )
111* ..
112* .. External Functions ..
113 DOUBLE PRECISION DLAMCH
114 DOUBLE PRECISION ZLANGE, ZLANSY
115 LOGICAL LSAME
116 EXTERNAL dlamch, zlange, zlansy, lsame
117* ..
118* .. Data statements ..
119 DATA iseed / 1988, 1989, 1990, 1991 /
120*
121 eps = dlamch( 'Epsilon' )
122 k = m
123 n2 = m+n
124 IF( n.GT.0 ) THEN
125 np1 = m+1
126 ELSE
127 np1 = 1
128 END IF
129 lwork = n2*n2*nb
130*
131* Dynamically allocate all arrays
132*
133 ALLOCATE(a(m,n2),af(m,n2),q(n2,n2),r(n2,n2),rwork(n2),
134 $ work(lwork),t(nb,m),c(n2,m),cf(n2,m),
135 $ d(m,n2),df(m,n2) )
136*
137* Put random stuff into A
138*
139 ldt=nb
140 CALL zlaset( 'Full', m, n2, czero, czero, a, m )
141 CALL zlaset( 'Full', nb, m, czero, czero, t, nb )
142 DO j=1,m
143 CALL zlarnv( 2, iseed, m-j+1, a( j, j ) )
144 END DO
145 IF( n.GT.0 ) THEN
146 DO j=1,n-l
147 CALL zlarnv( 2, iseed, m, a( 1, min(n+m,m+1) + j - 1 ) )
148 END DO
149 END IF
150 IF( l.GT.0 ) THEN
151 DO j=1,l
152 CALL zlarnv( 2, iseed, m-j+1, a( j, min(n+m,n+m-l+1)
153 $ + j - 1 ) )
154 END DO
155 END IF
156*
157* Copy the matrix A to the array AF.
158*
159 CALL zlacpy( 'Full', m, n2, a, m, af, m )
160*
161* Factor the matrix A in the array AF.
162*
163 CALL ztplqt( m,n,l,nb,af,m,af(1,np1),m,t,ldt,work,info)
164*
165* Generate the (M+N)-by-(M+N) matrix Q by applying H to I
166*
167 CALL zlaset( 'Full', n2, n2, czero, one, q, n2 )
168 CALL zgemlqt( 'L', 'N', n2, n2, k, nb, af, m, t, ldt, q, n2,
169 $ work, info )
170*
171* Copy L
172*
173 CALL zlaset( 'Full', n2, n2, czero, czero, r, n2 )
174 CALL zlacpy( 'Lower', m, n2, af, m, r, n2 )
175*
176* Compute |L - A*Q*C| / |A| and store in RESULT(1)
177*
178 CALL zgemm( 'N', 'C', m, n2, n2, -one, a, m, q, n2, one, r, n2)
179 anorm = zlange( '1', m, n2, a, m, rwork )
180 resid = zlange( '1', m, n2, r, n2, rwork )
181 IF( anorm.GT.zero ) THEN
182 result( 1 ) = resid / (eps*anorm*max(1,n2))
183 ELSE
184 result( 1 ) = zero
185 END IF
186*
187* Compute |I - Q*Q'| and store in RESULT(2)
188*
189 CALL zlaset( 'Full', n2, n2, czero, one, r, n2 )
190 CALL zherk( 'U', 'N', n2, n2, dreal(-one), q, n2, dreal(one),
191 $ r, n2 )
192 resid = zlansy( '1', 'Upper', n2, r, n2, rwork )
193 result( 2 ) = resid / (eps*max(1,n2))
194*
195* Generate random m-by-n matrix C and a copy CF
196*
197 CALL zlaset( 'Full', n2, m, czero, one, c, n2 )
198 DO j=1,m
199 CALL zlarnv( 2, iseed, n2, c( 1, j ) )
200 END DO
201 cnorm = zlange( '1', n2, m, c, n2, rwork)
202 CALL zlacpy( 'Full', n2, m, c, n2, cf, n2 )
203*
204* Apply Q to C as Q*C
205*
206 CALL ztpmlqt( 'L','N', n,m,k,l,nb,af(1, np1),m,t,ldt,cf,n2,
207 $ cf(np1,1),n2,work,info)
208*
209* Compute |Q*C - Q*C| / |C|
210*
211 CALL zgemm( 'N', 'N', n2, m, n2, -one, q, n2, c, n2, one, cf, n2 )
212 resid = zlange( '1', n2, m, cf, n2, rwork )
213 IF( cnorm.GT.zero ) THEN
214 result( 3 ) = resid / (eps*max(1,n2)*cnorm)
215 ELSE
216 result( 3 ) = zero
217 END IF
218
219*
220* Copy C into CF again
221*
222 CALL zlacpy( 'Full', n2, m, c, n2, cf, n2 )
223*
224* Apply Q to C as QT*C
225*
226 CALL ztpmlqt( 'L','C',n,m,k,l,nb,af(1,np1),m,t,ldt,cf,n2,
227 $ cf(np1,1),n2,work,info)
228*
229* Compute |QT*C - QT*C| / |C|
230*
231 CALL zgemm('C','N',n2,m,n2,-one,q,n2,c,n2,one,cf,n2)
232 resid = zlange( '1', n2, m, cf, n2, rwork )
233
234 IF( cnorm.GT.zero ) THEN
235 result( 4 ) = resid / (eps*max(1,n2)*cnorm)
236 ELSE
237 result( 4 ) = zero
238 END IF
239*
240* Generate random m-by-n matrix D and a copy DF
241*
242 DO j=1,n2
243 CALL zlarnv( 2, iseed, m, d( 1, j ) )
244 END DO
245 dnorm = zlange( '1', m, n2, d, m, rwork)
246 CALL zlacpy( 'Full', m, n2, d, m, df, m )
247*
248* Apply Q to D as D*Q
249*
250 CALL ztpmlqt('R','N',m,n,k,l,nb,af(1,np1),m,t,ldt,df,m,
251 $ df(1,np1),m,work,info)
252*
253* Compute |D*Q - D*Q| / |D|
254*
255 CALL zgemm('N','N',m,n2,n2,-one,d,m,q,n2,one,df,m)
256 resid = zlange('1',m, n2,df,m,rwork )
257 IF( cnorm.GT.zero ) THEN
258 result( 5 ) = resid / (eps*max(1,n2)*dnorm)
259 ELSE
260 result( 5 ) = zero
261 END IF
262*
263* Copy D into DF again
264*
265 CALL zlacpy('Full',m,n2,d,m,df,m )
266*
267* Apply Q to D as D*QT
268*
269 CALL ztpmlqt('R','C',m,n,k,l,nb,af(1,np1),m,t,ldt,df,m,
270 $ df(1,np1),m,work,info)
271
272*
273* Compute |D*QT - D*QT| / |D|
274*
275 CALL zgemm( 'N', 'C', m, n2, n2, -one, d, m, q, n2, one, df, m )
276 resid = zlange( '1', m, n2, df, m, rwork )
277 IF( cnorm.GT.zero ) THEN
278 result( 6 ) = resid / (eps*max(1,n2)*dnorm)
279 ELSE
280 result( 6 ) = zero
281 END IF
282*
283* Deallocate all arrays
284*
285 DEALLOCATE ( a, af, q, r, rwork, work, t, c, d, cf, df)
286 RETURN
287 END
subroutine zgemlqt(side, trans, m, n, k, mb, v, ldv, t, ldt, c, ldc, work, info)
ZGEMLQT
Definition zgemlqt.f:168
subroutine zgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
ZGEMM
Definition zgemm.f:188
subroutine zherk(uplo, trans, n, k, alpha, a, lda, beta, c, ldc)
ZHERK
Definition zherk.f:173
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
subroutine zlarnv(idist, iseed, n, x)
ZLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition zlarnv.f:99
subroutine zlaset(uplo, m, n, alpha, beta, a, lda)
ZLASET initializes the off-diagonal elements and the diagonal elements of a matrix to given values.
Definition zlaset.f:106
subroutine ztplqt(m, n, l, mb, a, lda, b, ldb, t, ldt, work, info)
ZTPLQT
Definition ztplqt.f:189
subroutine ztpmlqt(side, trans, m, n, k, l, mb, v, ldv, t, ldt, a, lda, b, ldb, work, info)
ZTPMLQT
Definition ztpmlqt.f:214
subroutine zlqt05(m, n, l, nb, result)
ZLQT05
Definition zlqt05.f:80