OpenRadioss 2025.1.11
OpenRadioss project
Loading...
Searching...
No Matches
ctsqr01.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine ctsqr01 (tssw, m, n, mb, nb, result)
 CTSQR01

Function/Subroutine Documentation

◆ ctsqr01()

subroutine ctsqr01 ( character tssw,
integer m,
integer n,
integer mb,
integer nb,
real, dimension(6) result )

CTSQR01

Purpose:
!>
!> DTSQR01 tests DGEQR , DGELQ, DGEMLQ and DGEMQR.
!>
Parameters
[in]TSSW
!>          TSSW is CHARACTER
!>          'TS' for testing tall skinny QR
!>               and anything else for testing short wide LQ
!>
[in]M
!>          M is INTEGER
!>          Number of rows in test matrix.
!>
[in]N
!>          N is INTEGER
!>          Number of columns in test matrix.
!>
[in]MB
!>          MB is INTEGER
!>          Number of row in row block in test matrix.
!>
[in]NB
!>          NB is INTEGER
!>          Number of columns in column block test matrix.
!>
[out]RESULT
!>          RESULT is REAL array, dimension (6)
!>          Results of each of the six tests below.
!>
!>          RESULT(1) = | A - Q R | or | A - L Q |
!>          RESULT(2) = | I - Q^H Q | or | I - Q Q^H |
!>          RESULT(3) = | Q C - Q C |
!>          RESULT(4) = | Q^H C - Q^H C |
!>          RESULT(5) = | C Q - C Q |
!>          RESULT(6) = | C Q^H - C Q^H |
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.

Definition at line 81 of file ctsqr01.f.

82 IMPLICIT NONE
83*
84* -- LAPACK test routine --
85* -- LAPACK is a software package provided by Univ. of Tennessee, --
86* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
87*
88* .. Scalar Arguments ..
89 CHARACTER TSSW
90 INTEGER M, N, MB, NB
91* .. Return values ..
92 REAL RESULT(6)
93*
94* =====================================================================
95*
96* ..
97* .. Local allocatable arrays
98 COMPLEX, ALLOCATABLE :: AF(:,:), Q(:,:),
99 $ R(:,:), RWORK(:), WORK( : ), T(:),
100 $ CF(:,:), DF(:,:), A(:,:), C(:,:), D(:,:), LQ(:,:)
101*
102* .. Parameters ..
103 REAL ZERO
104 COMPLEX ONE, CZERO
105 parameter( zero = 0.0, one = (1.0,0.0), czero=(0.0,0.0) )
106* ..
107* .. Local Scalars ..
108 LOGICAL TESTZEROS, TS
109 INTEGER INFO, J, K, L, LWORK, TSIZE, MNB
110 REAL ANORM, EPS, RESID, CNORM, DNORM
111* ..
112* .. Local Arrays ..
113 INTEGER ISEED( 4 )
114 COMPLEX TQUERY( 5 ), WORKQUERY( 1 )
115* ..
116* .. External Functions ..
117 REAL SLAMCH, CLANGE, CLANSY
118 LOGICAL LSAME
119 INTEGER ILAENV
120 EXTERNAL slamch, clange, clansy, lsame, ilaenv
121* ..
122* .. Intrinsic Functions ..
123 INTRINSIC max, min
124* .. Scalars in Common ..
125 CHARACTER*32 srnamt
126* ..
127* .. Common blocks ..
128 COMMON / srnamc / srnamt
129* ..
130* .. Data statements ..
131 DATA iseed / 1988, 1989, 1990, 1991 /
132*
133* TEST TALL SKINNY OR SHORT WIDE
134*
135 ts = lsame(tssw, 'TS')
136*
137* TEST MATRICES WITH HALF OF MATRIX BEING ZEROS
138*
139 testzeros = .false.
140*
141 eps = slamch( 'Epsilon' )
142 k = min(m,n)
143 l = max(m,n,1)
144 mnb = max( mb, nb)
145 lwork = max(3,l)*mnb
146*
147* Dynamically allocate local arrays
148*
149 ALLOCATE ( a(m,n), af(m,n), q(l,l), r(m,l), rwork(l),
150 $ c(m,n), cf(m,n),
151 $ d(n,m), df(n,m), lq(l,n) )
152*
153* Put random numbers into A and copy to AF
154*
155 DO j=1,n
156 CALL clarnv( 2, iseed, m, a( 1, j ) )
157 END DO
158 IF (testzeros) THEN
159 IF (m.GE.4) THEN
160 DO j=1,n
161 CALL clarnv( 2, iseed, m/2, a( m/4, j ) )
162 END DO
163 END IF
164 END IF
165 CALL clacpy( 'Full', m, n, a, m, af, m )
166*
167 IF (ts) THEN
168*
169* Factor the matrix A in the array AF.
170*
171 CALL cgeqr( m, n, af, m, tquery, -1, workquery, -1, info )
172 tsize = int( tquery( 1 ) )
173 lwork = int( workquery( 1 ) )
174 CALL cgemqr( 'L', 'N', m, m, k, af, m, tquery, tsize, cf, m,
175 $ workquery, -1, info)
176 lwork = max( lwork, int( workquery( 1 ) ) )
177 CALL cgemqr( 'L', 'N', m, n, k, af, m, tquery, tsize, cf, m,
178 $ workquery, -1, info)
179 lwork = max( lwork, int( workquery( 1 ) ) )
180 CALL cgemqr( 'L', 'C', m, n, k, af, m, tquery, tsize, cf, m,
181 $ workquery, -1, info)
182 lwork = max( lwork, int( workquery( 1 ) ) )
183 CALL cgemqr( 'R', 'N', n, m, k, af, m, tquery, tsize, df, n,
184 $ workquery, -1, info)
185 lwork = max( lwork, int( workquery( 1 ) ) )
186 CALL cgemqr( 'R', 'C', n, m, k, af, m, tquery, tsize, df, n,
187 $ workquery, -1, info)
188 lwork = max( lwork, int( workquery( 1 ) ) )
189 ALLOCATE ( t( tsize ) )
190 ALLOCATE ( work( lwork ) )
191 srnamt = 'CGEQR'
192 CALL cgeqr( m, n, af, m, t, tsize, work, lwork, info )
193*
194* Generate the m-by-m matrix Q
195*
196 CALL claset( 'Full', m, m, czero, one, q, m )
197 srnamt = 'CGEMQR'
198 CALL cgemqr( 'L', 'N', m, m, k, af, m, t, tsize, q, m,
199 $ work, lwork, info )
200*
201* Copy R
202*
203 CALL claset( 'Full', m, n, czero, czero, r, m )
204 CALL clacpy( 'Upper', m, n, af, m, r, m )
205*
206* Compute |R - Q'*A| / |A| and store in RESULT(1)
207*
208 CALL cgemm( 'C', 'N', m, n, m, -one, q, m, a, m, one, r, m )
209 anorm = clange( '1', m, n, a, m, rwork )
210 resid = clange( '1', m, n, r, m, rwork )
211 IF( anorm.GT.zero ) THEN
212 result( 1 ) = resid / (eps*max(1,m)*anorm)
213 ELSE
214 result( 1 ) = zero
215 END IF
216*
217* Compute |I - Q'*Q| and store in RESULT(2)
218*
219 CALL claset( 'Full', m, m, czero, one, r, m )
220 CALL cherk( 'U', 'C', m, m, real(-one), q, m, real(one), r, m )
221 resid = clansy( '1', 'Upper', m, r, m, rwork )
222 result( 2 ) = resid / (eps*max(1,m))
223*
224* Generate random m-by-n matrix C and a copy CF
225*
226 DO j=1,n
227 CALL clarnv( 2, iseed, m, c( 1, j ) )
228 END DO
229 cnorm = clange( '1', m, n, c, m, rwork)
230 CALL clacpy( 'Full', m, n, c, m, cf, m )
231*
232* Apply Q to C as Q*C
233*
234 srnamt = 'CGEMQR'
235 CALL cgemqr( 'L', 'N', m, n, k, af, m, t, tsize, cf, m,
236 $ work, lwork, info)
237*
238* Compute |Q*C - Q*C| / |C|
239*
240 CALL cgemm( 'N', 'N', m, n, m, -one, q, m, c, m, one, cf, m )
241 resid = clange( '1', m, n, cf, m, rwork )
242 IF( cnorm.GT.zero ) THEN
243 result( 3 ) = resid / (eps*max(1,m)*cnorm)
244 ELSE
245 result( 3 ) = zero
246 END IF
247*
248* Copy C into CF again
249*
250 CALL clacpy( 'Full', m, n, c, m, cf, m )
251*
252* Apply Q to C as QT*C
253*
254 srnamt = 'CGEMQR'
255 CALL cgemqr( 'L', 'C', m, n, k, af, m, t, tsize, cf, m,
256 $ work, lwork, info)
257*
258* Compute |QT*C - QT*C| / |C|
259*
260 CALL cgemm( 'C', 'N', m, n, m, -one, q, m, c, m, one, cf, m )
261 resid = clange( '1', m, n, cf, m, rwork )
262 IF( cnorm.GT.zero ) THEN
263 result( 4 ) = resid / (eps*max(1,m)*cnorm)
264 ELSE
265 result( 4 ) = zero
266 END IF
267*
268* Generate random n-by-m matrix D and a copy DF
269*
270 DO j=1,m
271 CALL clarnv( 2, iseed, n, d( 1, j ) )
272 END DO
273 dnorm = clange( '1', n, m, d, n, rwork)
274 CALL clacpy( 'Full', n, m, d, n, df, n )
275*
276* Apply Q to D as D*Q
277*
278 srnamt = 'CGEMQR'
279 CALL cgemqr( 'R', 'N', n, m, k, af, m, t, tsize, df, n,
280 $ work, lwork, info)
281*
282* Compute |D*Q - D*Q| / |D|
283*
284 CALL cgemm( 'N', 'N', n, m, m, -one, d, n, q, m, one, df, n )
285 resid = clange( '1', n, m, df, n, rwork )
286 IF( dnorm.GT.zero ) THEN
287 result( 5 ) = resid / (eps*max(1,m)*dnorm)
288 ELSE
289 result( 5 ) = zero
290 END IF
291*
292* Copy D into DF again
293*
294 CALL clacpy( 'Full', n, m, d, n, df, n )
295*
296* Apply Q to D as D*QT
297*
298 CALL cgemqr( 'R', 'C', n, m, k, af, m, t, tsize, df, n,
299 $ work, lwork, info)
300*
301* Compute |D*QT - D*QT| / |D|
302*
303 CALL cgemm( 'N', 'C', n, m, m, -one, d, n, q, m, one, df, n )
304 resid = clange( '1', n, m, df, n, rwork )
305 IF( cnorm.GT.zero ) THEN
306 result( 6 ) = resid / (eps*max(1,m)*dnorm)
307 ELSE
308 result( 6 ) = zero
309 END IF
310*
311* Short and wide
312*
313 ELSE
314 CALL cgelq( m, n, af, m, tquery, -1, workquery, -1, info )
315 tsize = int( tquery( 1 ) )
316 lwork = int( workquery( 1 ) )
317 CALL cgemlq( 'R', 'N', n, n, k, af, m, tquery, tsize, q, n,
318 $ workquery, -1, info )
319 lwork = max( lwork, int( workquery( 1 ) ) )
320 CALL cgemlq( 'L', 'N', n, m, k, af, m, tquery, tsize, df, n,
321 $ workquery, -1, info)
322 lwork = max( lwork, int( workquery( 1 ) ) )
323 CALL cgemlq( 'L', 'C', n, m, k, af, m, tquery, tsize, df, n,
324 $ workquery, -1, info)
325 lwork = max( lwork, int( workquery( 1 ) ) )
326 CALL cgemlq( 'R', 'N', m, n, k, af, m, tquery, tsize, cf, m,
327 $ workquery, -1, info)
328 lwork = max( lwork, int( workquery( 1 ) ) )
329 CALL cgemlq( 'R', 'C', m, n, k, af, m, tquery, tsize, cf, m,
330 $ workquery, -1, info)
331 lwork = max( lwork, int( workquery( 1 ) ) )
332 ALLOCATE ( t( tsize ) )
333 ALLOCATE ( work( lwork ) )
334 srnamt = 'CGELQ'
335 CALL cgelq( m, n, af, m, t, tsize, work, lwork, info )
336*
337*
338* Generate the n-by-n matrix Q
339*
340 CALL claset( 'Full', n, n, czero, one, q, n )
341 srnamt = 'CGEMLQ'
342 CALL cgemlq( 'R', 'N', n, n, k, af, m, t, tsize, q, n,
343 $ work, lwork, info )
344*
345* Copy R
346*
347 CALL claset( 'Full', m, n, czero, czero, lq, l )
348 CALL clacpy( 'Lower', m, n, af, m, lq, l )
349*
350* Compute |L - A*Q'| / |A| and store in RESULT(1)
351*
352 CALL cgemm( 'N', 'C', m, n, n, -one, a, m, q, n, one, lq, l )
353 anorm = clange( '1', m, n, a, m, rwork )
354 resid = clange( '1', m, n, lq, l, rwork )
355 IF( anorm.GT.zero ) THEN
356 result( 1 ) = resid / (eps*max(1,n)*anorm)
357 ELSE
358 result( 1 ) = zero
359 END IF
360*
361* Compute |I - Q'*Q| and store in RESULT(2)
362*
363 CALL claset( 'Full', n, n, czero, one, lq, l )
364 CALL cherk( 'U', 'C', n, n, real(-one), q, n, real(one), lq, l)
365 resid = clansy( '1', 'Upper', n, lq, l, rwork )
366 result( 2 ) = resid / (eps*max(1,n))
367*
368* Generate random m-by-n matrix C and a copy CF
369*
370 DO j=1,m
371 CALL clarnv( 2, iseed, n, d( 1, j ) )
372 END DO
373 dnorm = clange( '1', n, m, d, n, rwork)
374 CALL clacpy( 'Full', n, m, d, n, df, n )
375*
376* Apply Q to C as Q*C
377*
378 CALL cgemlq( 'L', 'N', n, m, k, af, m, t, tsize, df, n,
379 $ work, lwork, info)
380*
381* Compute |Q*D - Q*D| / |D|
382*
383 CALL cgemm( 'N', 'N', n, m, n, -one, q, n, d, n, one, df, n )
384 resid = clange( '1', n, m, df, n, rwork )
385 IF( dnorm.GT.zero ) THEN
386 result( 3 ) = resid / (eps*max(1,n)*dnorm)
387 ELSE
388 result( 3 ) = zero
389 END IF
390*
391* Copy D into DF again
392*
393 CALL clacpy( 'Full', n, m, d, n, df, n )
394*
395* Apply Q to D as QT*D
396*
397 CALL cgemlq( 'L', 'C', n, m, k, af, m, t, tsize, df, n,
398 $ work, lwork, info)
399*
400* Compute |QT*D - QT*D| / |D|
401*
402 CALL cgemm( 'C', 'N', n, m, n, -one, q, n, d, n, one, df, n )
403 resid = clange( '1', n, m, df, n, rwork )
404 IF( dnorm.GT.zero ) THEN
405 result( 4 ) = resid / (eps*max(1,n)*dnorm)
406 ELSE
407 result( 4 ) = zero
408 END IF
409*
410* Generate random n-by-m matrix D and a copy DF
411*
412 DO j=1,n
413 CALL clarnv( 2, iseed, m, c( 1, j ) )
414 END DO
415 cnorm = clange( '1', m, n, c, m, rwork)
416 CALL clacpy( 'Full', m, n, c, m, cf, m )
417*
418* Apply Q to C as C*Q
419*
420 CALL cgemlq( 'R', 'N', m, n, k, af, m, t, tsize, cf, m,
421 $ work, lwork, info)
422*
423* Compute |C*Q - C*Q| / |C|
424*
425 CALL cgemm( 'N', 'N', m, n, n, -one, c, m, q, n, one, cf, m )
426 resid = clange( '1', n, m, df, n, rwork )
427 IF( cnorm.GT.zero ) THEN
428 result( 5 ) = resid / (eps*max(1,n)*cnorm)
429 ELSE
430 result( 5 ) = zero
431 END IF
432*
433* Copy C into CF again
434*
435 CALL clacpy( 'Full', m, n, c, m, cf, m )
436*
437* Apply Q to D as D*QT
438*
439 CALL cgemlq( 'R', 'C', m, n, k, af, m, t, tsize, cf, m,
440 $ work, lwork, info)
441*
442* Compute |C*QT - C*QT| / |C|
443*
444 CALL cgemm( 'N', 'C', m, n, n, -one, c, m, q, n, one, cf, m )
445 resid = clange( '1', m, n, cf, m, rwork )
446 IF( cnorm.GT.zero ) THEN
447 result( 6 ) = resid / (eps*max(1,n)*cnorm)
448 ELSE
449 result( 6 ) = zero
450 END IF
451*
452 END IF
453*
454* Deallocate all arrays
455*
456 DEALLOCATE ( a, af, q, r, rwork, work, t, c, d, cf, df)
457*
458 RETURN
cgelq
subroutine cgelq(m, n, a, lda, t, tsize, work, lwork, info)
CGELQ
Definition cgelq.f:172
cgemlq
subroutine cgemlq(side, trans, m, n, k, a, lda, t, tsize, c, ldc, work, lwork, info)
CGEMLQ
Definition cgemlq.f:170
cgemqr
subroutine cgemqr(side, trans, m, n, k, a, lda, t, tsize, c, ldc, work, lwork, info)
CGEMQR
Definition cgemqr.f:172
cgeqr
subroutine cgeqr(m, n, a, lda, t, tsize, work, lwork, info)
CGEQR
Definition cgeqr.f:174
ilaenv
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
Definition ilaenv.f:162
lsame
logical function lsame(ca, cb)
LSAME
Definition lsame.f:53
clange
real function clange(norm, m, n, a, lda, work)
CLANGE returns the value of the 1-norm, Frobenius norm, infinity-norm, or the largest absolute value ...
Definition clange.f:115
clacpy
subroutine clacpy(uplo, m, n, a, lda, b, ldb)
CLACPY copies all or part of one two-dimensional array to another.
Definition clacpy.f:103
claset
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
clarnv
subroutine clarnv(idist, iseed, n, x)
CLARNV returns a vector of random numbers from a uniform or normal distribution.
Definition clarnv.f:99
clansy
real function clansy(norm, uplo, n, a, lda, work)
CLANSY returns the value of the 1-norm, or the Frobenius norm, or the infinity norm,...
Definition clansy.f:123
cherk
subroutine cherk(uplo, trans, n, k, alpha, a, lda, beta, c, ldc)
CHERK
Definition cherk.f:173
cgemm
subroutine cgemm(transa, transb, m, n, k, alpha, a, lda, b, ldb, beta, c, ldc)
CGEMM
Definition cgemm.f:187
slamch
real function slamch(cmach)
SLAMCH
Definition slamch.f:68
min
#define min(a, b)
Definition macros.h:20
max
#define max(a, b)
Definition macros.h:21