1*> \brief \b CCHKTZ
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 CCHKTZ( DOTYPE, NM, MVAL, NN, NVAL, THRESH, TSTERR, A,
12*                          COPYA, S, TAU, WORK, RWORK, NOUT )
13*
14*       .. Scalar Arguments ..
15*       LOGICAL            TSTERR
16*       INTEGER            NM, NN, NOUT
17*       REAL               THRESH
18*       ..
19*       .. Array Arguments ..
20*       LOGICAL            DOTYPE( * )
21*       INTEGER            MVAL( * ), NVAL( * )
22*       REAL               S( * ), RWORK( * )
23*       COMPLEX            A( * ), COPYA( * ), TAU( * ), WORK( * )
24*       ..
25*
26*
27*> \par Purpose:
28*  =============
29*>
30*> \verbatim
31*>
32*> CCHKTZ tests CTZRQF and CTZRZF.
33*> \endverbatim
34*
35*  Arguments:
36*  ==========
37*
38*> \param[in] DOTYPE
39*> \verbatim
40*>          DOTYPE is LOGICAL array, dimension (NTYPES)
41*>          The matrix types to be used for testing.  Matrices of type j
42*>          (for 1 <= j <= NTYPES) are used for testing if DOTYPE(j) =
43*>          .TRUE.; if DOTYPE(j) = .FALSE., then type j is not used.
44*> \endverbatim
45*>
46*> \param[in] NM
47*> \verbatim
48*>          NM is INTEGER
49*>          The number of values of M contained in the vector MVAL.
50*> \endverbatim
51*>
52*> \param[in] MVAL
53*> \verbatim
54*>          MVAL is INTEGER array, dimension (NM)
55*>          The values of the matrix row dimension M.
56*> \endverbatim
57*>
58*> \param[in] NN
59*> \verbatim
60*>          NN is INTEGER
61*>          The number of values of N contained in the vector NVAL.
62*> \endverbatim
63*>
64*> \param[in] NVAL
65*> \verbatim
66*>          NVAL is INTEGER array, dimension (NN)
67*>          The values of the matrix column dimension N.
68*> \endverbatim
69*>
70*> \param[in] THRESH
71*> \verbatim
72*>          THRESH is REAL
73*>          The threshold value for the test ratios.  A result is
74*>          included in the output file if RESULT >= THRESH.  To have
75*>          every test ratio printed, use THRESH = 0.
76*> \endverbatim
77*>
78*> \param[in] TSTERR
79*> \verbatim
80*>          TSTERR is LOGICAL
81*>          Flag that indicates whether error exits are to be tested.
82*> \endverbatim
83*>
84*> \param[out] A
85*> \verbatim
86*>          A is COMPLEX array, dimension (MMAX*NMAX)
87*>          where MMAX is the maximum value of M in MVAL and NMAX is the
88*>          maximum value of N in NVAL.
89*> \endverbatim
90*>
91*> \param[out] COPYA
92*> \verbatim
93*>          COPYA is COMPLEX array, dimension (MMAX*NMAX)
94*> \endverbatim
95*>
96*> \param[out] S
97*> \verbatim
98*>          S is REAL array, dimension
99*>                      (min(MMAX,NMAX))
100*> \endverbatim
101*>
102*> \param[out] TAU
103*> \verbatim
104*>          TAU is COMPLEX array, dimension (MMAX)
105*> \endverbatim
106*>
107*> \param[out] WORK
108*> \verbatim
109*>          WORK is COMPLEX array, dimension
110*>                      (MMAX*NMAX + 4*NMAX + MMAX)
111*> \endverbatim
112*>
113*> \param[out] RWORK
114*> \verbatim
115*>          RWORK is REAL array, dimension (2*NMAX)
116*> \endverbatim
117*>
118*> \param[in] NOUT
119*> \verbatim
120*>          NOUT is INTEGER
121*>          The unit number for output.
122*> \endverbatim
123*
124*  Authors:
125*  ========
126*
127*> \author Univ. of Tennessee
128*> \author Univ. of California Berkeley
129*> \author Univ. of Colorado Denver
130*> \author NAG Ltd.
131*
132*> \date November 2011
133*
134*> \ingroup complex_lin
135*
136*  =====================================================================
137      SUBROUTINE CCHKTZ( DOTYPE, NM, MVAL, NN, NVAL, THRESH, TSTERR, A,
138     $                   COPYA, S, TAU, WORK, RWORK, NOUT )
139*
140*  -- LAPACK test routine (version 3.4.0) --
141*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
142*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
143*     November 2011
144*
145*     .. Scalar Arguments ..
146      LOGICAL            TSTERR
147      INTEGER            NM, NN, NOUT
148      REAL               THRESH
149*     ..
150*     .. Array Arguments ..
151      LOGICAL            DOTYPE( * )
152      INTEGER            MVAL( * ), NVAL( * )
153      REAL               S( * ), RWORK( * )
154      COMPLEX            A( * ), COPYA( * ), TAU( * ), WORK( * )
155*     ..
156*
157*  =====================================================================
158*
159*     .. Parameters ..
160      INTEGER            NTYPES
161      PARAMETER          ( NTYPES = 3 )
162      INTEGER            NTESTS
163      PARAMETER          ( NTESTS = 6 )
164      REAL               ONE, ZERO
165      PARAMETER          ( ONE = 1.0E0, ZERO = 0.0E0 )
166*     ..
167*     .. Local Scalars ..
168      CHARACTER*3        PATH
169      INTEGER            I, IM, IMODE, IN, INFO, K, LDA, LWORK, M,
170     $                   MNMIN, MODE, N, NERRS, NFAIL, NRUN
171      REAL               EPS
172*     ..
173*     .. Local Arrays ..
174      INTEGER            ISEED( 4 ), ISEEDY( 4 )
175      REAL               RESULT( NTESTS )
176*     ..
177*     .. External Functions ..
178      REAL               CQRT12, CRZT01, CRZT02, CTZT01, CTZT02, SLAMCH
179      EXTERNAL           CQRT12, CRZT01, CRZT02, CTZT01, CTZT02, SLAMCH
180*     ..
181*     .. External Subroutines ..
182      EXTERNAL           ALAHD, ALASUM, CERRTZ, CGEQR2, CLACPY, CLASET,
183     $                   CLATMS, CTZRQF, CTZRZF, SLAORD
184*     ..
185*     .. Intrinsic Functions ..
186      INTRINSIC          CMPLX, MAX, MIN
187*     ..
188*     .. Scalars in Common ..
189      LOGICAL            LERR, OK
190      CHARACTER*32       SRNAMT
191      INTEGER            INFOT, IOUNIT
192*     ..
193*     .. Common blocks ..
194      COMMON             / INFOC / INFOT, IOUNIT, OK, LERR
195      COMMON             / SRNAMC / SRNAMT
196*     ..
197*     .. Data statements ..
198      DATA               ISEEDY / 1988, 1989, 1990, 1991 /
199*     ..
200*     .. Executable Statements ..
201*
202*     Initialize constants and the random number seed.
203*
204      PATH( 1: 1 ) = 'Complex precision'
205      PATH( 2: 3 ) = 'TZ'
206      NRUN = 0
207      NFAIL = 0
208      NERRS = 0
209      DO 10 I = 1, 4
210         ISEED( I ) = ISEEDY( I )
211   10 CONTINUE
212      EPS = SLAMCH( 'Epsilon' )
213*
214*     Test the error exits
215*
216      IF( TSTERR )
217     $   CALL CERRTZ( PATH, NOUT )
218      INFOT = 0
219*
220      DO 70 IM = 1, NM
221*
222*        Do for each value of M in MVAL.
223*
224         M = MVAL( IM )
225         LDA = MAX( 1, M )
226*
227         DO 60 IN = 1, NN
228*
229*           Do for each value of N in NVAL for which M .LE. N.
230*
231            N = NVAL( IN )
232            MNMIN = MIN( M, N )
233            LWORK = MAX( 1, N*N+4*M+N )
234*
235            IF( M.LE.N ) THEN
236               DO 50 IMODE = 1, NTYPES
237                  IF( .NOT.DOTYPE( IMODE ) )
238     $               GO TO 50
239*
240*                 Do for each type of singular value distribution.
241*                    0:  zero matrix
242*                    1:  one small singular value
243*                    2:  exponential distribution
244*
245                  MODE = IMODE - 1
246*
247*                 Test CTZRQF
248*
249*                 Generate test matrix of size m by n using
250*                 singular value distribution indicated by `mode'.
251*
252                  IF( MODE.EQ.0 ) THEN
253                     CALL CLASET( 'Full', M, N, CMPLX( ZERO ),
254     $                            CMPLX( ZERO ), A, LDA )
255                     DO 20 I = 1, MNMIN
256                        S( I ) = ZERO
257   20                CONTINUE
258                  ELSE
259                     CALL CLATMS( M, N, 'Uniform', ISEED,
260     $                            'Nonsymmetric', S, IMODE,
261     $                            ONE / EPS, ONE, M, N, 'No packing', A,
262     $                            LDA, WORK, INFO )
263                     CALL CGEQR2( M, N, A, LDA, WORK, WORK( MNMIN+1 ),
264     $                            INFO )
265                     CALL CLASET( 'Lower', M-1, N, CMPLX( ZERO ),
266     $                            CMPLX( ZERO ), A( 2 ), LDA )
267                     CALL SLAORD( 'Decreasing', MNMIN, S, 1 )
268                  END IF
269*
270*                 Save A and its singular values
271*
272                  CALL CLACPY( 'All', M, N, A, LDA, COPYA, LDA )
273*
274*                 Call CTZRQF to reduce the upper trapezoidal matrix to
275*                 upper triangular form.
276*
277                  SRNAMT = 'CTZRQF'
278                  CALL CTZRQF( M, N, A, LDA, TAU, INFO )
279*
280*                 Compute norm(svd(a) - svd(r))
281*
282                  RESULT( 1 ) = CQRT12( M, M, A, LDA, S, WORK,
283     $                          LWORK, RWORK )
284*
285*                 Compute norm( A - R*Q )
286*
287                  RESULT( 2 ) = CTZT01( M, N, COPYA, A, LDA, TAU, WORK,
288     $                          LWORK )
289*
290*                 Compute norm(Q'*Q - I).
291*
292                  RESULT( 3 ) = CTZT02( M, N, A, LDA, TAU, WORK, LWORK )
293*
294*                 Test CTZRZF
295*
296*                 Generate test matrix of size m by n using
297*                 singular value distribution indicated by `mode'.
298*
299                  IF( MODE.EQ.0 ) THEN
300                     CALL CLASET( 'Full', M, N, CMPLX( ZERO ),
301     $                            CMPLX( ZERO ), A, LDA )
302                     DO 30 I = 1, MNMIN
303                        S( I ) = ZERO
304   30                CONTINUE
305                  ELSE
306                     CALL CLATMS( M, N, 'Uniform', ISEED,
307     $                            'Nonsymmetric', S, IMODE,
308     $                            ONE / EPS, ONE, M, N, 'No packing', A,
309     $                            LDA, WORK, INFO )
310                     CALL CGEQR2( M, N, A, LDA, WORK, WORK( MNMIN+1 ),
311     $                            INFO )
312                     CALL CLASET( 'Lower', M-1, N, CMPLX( ZERO ),
313     $                            CMPLX( ZERO ), A( 2 ), LDA )
314                     CALL SLAORD( 'Decreasing', MNMIN, S, 1 )
315                  END IF
316*
317*                 Save A and its singular values
318*
319                  CALL CLACPY( 'All', M, N, A, LDA, COPYA, LDA )
320*
321*                 Call CTZRZF to reduce the upper trapezoidal matrix to
322*                 upper triangular form.
323*
324                  SRNAMT = 'CTZRZF'
325                  CALL CTZRZF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
326*
327*                 Compute norm(svd(a) - svd(r))
328*
329                  RESULT( 4 ) = CQRT12( M, M, A, LDA, S, WORK,
330     $                          LWORK, RWORK )
331*
332*                 Compute norm( A - R*Q )
333*
334                  RESULT( 5 ) = CRZT01( M, N, COPYA, A, LDA, TAU, WORK,
335     $                          LWORK )
336*
337*                 Compute norm(Q'*Q - I).
338*
339                  RESULT( 6 ) = CRZT02( M, N, A, LDA, TAU, WORK, LWORK )
340*
341*                 Print information about the tests that did not pass
342*                 the threshold.
343*
344                  DO 40 K = 1, 6
345                     IF( RESULT( K ).GE.THRESH ) THEN
346                        IF( NFAIL.EQ.0 .AND. NERRS.EQ.0 )
347     $                     CALL ALAHD( NOUT, PATH )
348                        WRITE( NOUT, FMT = 9999 )M, N, IMODE, K,
349     $                     RESULT( K )
350                        NFAIL = NFAIL + 1
351                     END IF
352   40             CONTINUE
353                  NRUN = NRUN + 6
354   50          CONTINUE
355            END IF
356   60    CONTINUE
357   70 CONTINUE
358*
359*     Print a summary of the results.
360*
361      CALL ALASUM( PATH, NOUT, NFAIL, NRUN, NERRS )
362*
363 9999 FORMAT( ' M =', I5, ', N =', I5, ', type ', I2, ', test ', I2,
364     $      ', ratio =', G12.5 )
365*
366*     End if CCHKTZ
367*
368      END
369