1*> \brief \b CQRT04
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 CQRT04(M,N,NB,RESULT)
12*
13*       .. Scalar Arguments ..
14*       INTEGER M, N, NB, LDT
15*       .. Return values ..
16*       REAL RESULT(6)
17*
18*
19*> \par Purpose:
20*  =============
21*>
22*> \verbatim
23*>
24*> CQRT04 tests CGEQRT and CGEMQRT.
25*> \endverbatim
26*
27*  Arguments:
28*  ==========
29*
30*> \param[in] M
31*> \verbatim
32*>          M is INTEGER
33*>          Number of rows in 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] NB
43*> \verbatim
44*>          NB is INTEGER
45*>          Block size of test matrix.  NB <= Min(M,N).
46*> \endverbatim
47*>
48*> \param[out] RESULT
49*> \verbatim
50*>          RESULT is REAL array, dimension (6)
51*>          Results of each of the six tests below.
52*>
53*>          RESULT(1) = | A - Q R |
54*>          RESULT(2) = | I - Q^H Q |
55*>          RESULT(3) = | Q C - Q C |
56*>          RESULT(4) = | Q^H C - Q^H C |
57*>          RESULT(5) = | C Q - C Q |
58*>          RESULT(6) = | C Q^H - C Q^H |
59*> \endverbatim
60*
61*  Authors:
62*  ========
63*
64*> \author Univ. of Tennessee
65*> \author Univ. of California Berkeley
66*> \author Univ. of Colorado Denver
67*> \author NAG Ltd.
68*
69*> \date April 2012
70*
71*> \ingroup complex_lin
72*
73*  =====================================================================
74      SUBROUTINE CQRT04(M,N,NB,RESULT)
75      IMPLICIT NONE
76*
77*  -- LAPACK test routine (version 3.4.1) --
78*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
79*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
80*     April 2012
81*
82*     .. Scalar Arguments ..
83      INTEGER M, N, NB, LDT
84*     .. Return values ..
85      REAL RESULT(6)
86*
87*  =====================================================================
88*
89*     ..
90*     .. Local allocatable arrays
91      COMPLEX, ALLOCATABLE :: AF(:,:), Q(:,:),
92     $  R(:,:), RWORK(:), WORK( : ), T(:,:),
93     $  CF(:,:), DF(:,:), A(:,:), C(:,:), D(:,:)
94*
95*     .. Parameters ..
96      REAL ZERO
97      COMPLEX ONE, CZERO
98      PARAMETER( ZERO = 0.0, ONE = (1.0,0.0), CZERO=(0.0,0.0) )
99*     ..
100*     .. Local Scalars ..
101      INTEGER INFO, J, K, L, LWORK
102      REAL   ANORM, EPS, RESID, CNORM, DNORM
103*     ..
104*     .. Local Arrays ..
105      INTEGER            ISEED( 4 )
106*     ..
107*     .. External Functions ..
108      REAL SLAMCH
109      REAL CLANGE, CLANSY
110      LOGICAL  LSAME
111      EXTERNAL SLAMCH, CLANGE, CLANSY, LSAME
112*     ..
113*     .. Intrinsic Functions ..
114      INTRINSIC  MAX, MIN
115*     ..
116*     .. Data statements ..
117      DATA ISEED / 1988, 1989, 1990, 1991 /
118*
119      EPS = SLAMCH( 'Epsilon' )
120      K = MIN(M,N)
121      L = MAX(M,N)
122      LWORK = MAX(2,L)*MAX(2,L)*NB
123*
124*     Dynamically allocate local arrays
125*
126      ALLOCATE ( A(M,N), AF(M,N), Q(M,M), R(M,L), RWORK(L),
127     $           WORK(LWORK), T(NB,N), C(M,N), CF(M,N),
128     $           D(N,M), DF(N,M) )
129*
130*     Put random numbers into A and copy to AF
131*
132      LDT=NB
133      DO J=1,N
134         CALL CLARNV( 2, ISEED, M, A( 1, J ) )
135      END DO
136      CALL CLACPY( 'Full', M, N, A, M, AF, M )
137*
138*     Factor the matrix A in the array AF.
139*
140      CALL CGEQRT( M, N, NB, AF, M, T, LDT, WORK, INFO )
141*
142*     Generate the m-by-m matrix Q
143*
144      CALL CLASET( 'Full', M, M, CZERO, ONE, Q, M )
145      CALL CGEMQRT( 'R', 'N', M, M, K, NB, AF, M, T, LDT, Q, M,
146     $              WORK, INFO )
147*
148*     Copy R
149*
150      CALL CLASET( 'Full', M, N, CZERO, CZERO, R, M )
151      CALL CLACPY( 'Upper', M, N, AF, M, R, M )
152*
153*     Compute |R - Q'*A| / |A| and store in RESULT(1)
154*
155      CALL CGEMM( 'C', 'N', M, N, M, -ONE, Q, M, A, M, ONE, R, M )
156      ANORM = CLANGE( '1', M, N, A, M, RWORK )
157      RESID = CLANGE( '1', M, N, R, M, RWORK )
158      IF( ANORM.GT.ZERO ) THEN
159         RESULT( 1 ) = RESID / (EPS*MAX(1,M)*ANORM)
160      ELSE
161         RESULT( 1 ) = ZERO
162      END IF
163*
164*     Compute |I - Q'*Q| and store in RESULT(2)
165*
166      CALL CLASET( 'Full', M, M, CZERO, ONE, R, M )
167      CALL CHERK( 'U', 'C', M, M, REAL(-ONE), Q, M, REAL(ONE), R, M )
168      RESID = CLANSY( '1', 'Upper', M, R, M, RWORK )
169      RESULT( 2 ) = RESID / (EPS*MAX(1,M))
170*
171*     Generate random m-by-n matrix C and a copy CF
172*
173      DO J=1,N
174         CALL CLARNV( 2, ISEED, M, C( 1, J ) )
175      END DO
176      CNORM = CLANGE( '1', M, N, C, M, RWORK)
177      CALL CLACPY( 'Full', M, N, C, M, CF, M )
178*
179*     Apply Q to C as Q*C
180*
181      CALL CGEMQRT( 'L', 'N', M, N, K, NB, AF, M, T, NB, CF, M,
182     $             WORK, INFO)
183*
184*     Compute |Q*C - Q*C| / |C|
185*
186      CALL CGEMM( 'N', 'N', M, N, M, -ONE, Q, M, C, M, ONE, CF, M )
187      RESID = CLANGE( '1', M, N, CF, M, RWORK )
188      IF( CNORM.GT.ZERO ) THEN
189         RESULT( 3 ) = RESID / (EPS*MAX(1,M)*CNORM)
190      ELSE
191         RESULT( 3 ) = ZERO
192      END IF
193*
194*     Copy C into CF again
195*
196      CALL CLACPY( 'Full', M, N, C, M, CF, M )
197*
198*     Apply Q to C as QT*C
199*
200      CALL CGEMQRT( 'L', 'C', M, N, K, NB, AF, M, T, NB, CF, M,
201     $             WORK, INFO)
202*
203*     Compute |QT*C - QT*C| / |C|
204*
205      CALL CGEMM( 'C', 'N', M, N, M, -ONE, Q, M, C, M, ONE, CF, M )
206      RESID = CLANGE( '1', M, N, CF, M, RWORK )
207      IF( CNORM.GT.ZERO ) THEN
208         RESULT( 4 ) = RESID / (EPS*MAX(1,M)*CNORM)
209      ELSE
210         RESULT( 4 ) = ZERO
211      END IF
212*
213*     Generate random n-by-m matrix D and a copy DF
214*
215      DO J=1,M
216         CALL CLARNV( 2, ISEED, N, D( 1, J ) )
217      END DO
218      DNORM = CLANGE( '1', N, M, D, N, RWORK)
219      CALL CLACPY( 'Full', N, M, D, N, DF, N )
220*
221*     Apply Q to D as D*Q
222*
223      CALL CGEMQRT( 'R', 'N', N, M, K, NB, AF, M, T, NB, DF, N,
224     $             WORK, INFO)
225*
226*     Compute |D*Q - D*Q| / |D|
227*
228      CALL CGEMM( 'N', 'N', N, M, M, -ONE, D, N, Q, M, ONE, DF, N )
229      RESID = CLANGE( '1', N, M, DF, N, RWORK )
230      IF( CNORM.GT.ZERO ) THEN
231         RESULT( 5 ) = RESID / (EPS*MAX(1,M)*DNORM)
232      ELSE
233         RESULT( 5 ) = ZERO
234      END IF
235*
236*     Copy D into DF again
237*
238      CALL CLACPY( 'Full', N, M, D, N, DF, N )
239*
240*     Apply Q to D as D*QT
241*
242      CALL CGEMQRT( 'R', 'C', N, M, K, NB, AF, M, T, NB, DF, N,
243     $             WORK, INFO)
244*
245*     Compute |D*QT - D*QT| / |D|
246*
247      CALL CGEMM( 'N', 'C', N, M, M, -ONE, D, N, Q, M, ONE, DF, N )
248      RESID = CLANGE( '1', N, M, DF, N, RWORK )
249      IF( CNORM.GT.ZERO ) THEN
250         RESULT( 6 ) = RESID / (EPS*MAX(1,M)*DNORM)
251      ELSE
252         RESULT( 6 ) = ZERO
253      END IF
254*
255*     Deallocate all arrays
256*
257      DEALLOCATE ( A, AF, Q, R, RWORK, WORK, T, C, D, CF, DF)
258*
259      RETURN
260      END
261
262