1*> \brief \b CTSQR01
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 CTSQR01(TSSW, M,N, MB, NB, RESULT)
12*
13*       .. Scalar Arguments ..
14*       INTEGER M, N, MB
15*       .. Return values ..
16*       REAL RESULT(6)
17*
18*
19*> \par Purpose:
20*  =============
21*>
22*> \verbatim
23*>
24*> DTSQR01 tests DGEQR , DGELQ, DGEMLQ and DGEMQR.
25*> \endverbatim
26*
27*  Arguments:
28*  ==========
29*
30*> \param[in] TSSW
31*> \verbatim
32*>          TSSW is CHARACTER
33*>          'TS' for testing tall skinny QR
34*>               and anything else for testing short wide LQ
35*> \endverbatim
36*> \param[in] M
37*> \verbatim
38*>          M is INTEGER
39*>          Number of rows in test matrix.
40*> \endverbatim
41*>
42*> \param[in] N
43*> \verbatim
44*>          N is INTEGER
45*>          Number of columns in test matrix.
46*> \endverbatim
47*> \param[in] MB
48*> \verbatim
49*>          MB is INTEGER
50*>          Number of row in row block in test matrix.
51*> \endverbatim
52*>
53*> \param[in] NB
54*> \verbatim
55*>          NB is INTEGER
56*>          Number of columns in column block test matrix.
57*> \endverbatim
58*>
59*> \param[out] RESULT
60*> \verbatim
61*>          RESULT is REAL array, dimension (6)
62*>          Results of each of the six tests below.
63*>
64*>          RESULT(1) = | A - Q R | or | A - L Q |
65*>          RESULT(2) = | I - Q^H Q | or | I - Q Q^H |
66*>          RESULT(3) = | Q C - Q C |
67*>          RESULT(4) = | Q^H C - Q^H C |
68*>          RESULT(5) = | C Q - C Q |
69*>          RESULT(6) = | C Q^H - C Q^H |
70*> \endverbatim
71*
72*  Authors:
73*  ========
74*
75*> \author Univ. of Tennessee
76*> \author Univ. of California Berkeley
77*> \author Univ. of Colorado Denver
78*> \author NAG Ltd.
79*
80*  =====================================================================
81      SUBROUTINE CTSQR01(TSSW, M, N, MB, NB, RESULT)
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
459      END
460