1*> \brief \b CUNGTSQR
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download CUNGTSQR + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cuntsqr.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cungtsqr.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cungtsqr.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*>
18*  Definition:
19*  ===========
20*
21*       SUBROUTINE CUNGTSQR( M, N, MB, NB, A, LDA, T, LDT, WORK, LWORK,
22*      $                     INFO )
23*
24*       .. Scalar Arguments ..
25*       INTEGER           INFO, LDA, LDT, LWORK, M, N, MB, NB
26*       ..
27*       .. Array Arguments ..
28*       COMPLEX           A( LDA, * ), T( LDT, * ), WORK( * )
29*       ..
30*
31*> \par Purpose:
32*  =============
33*>
34*> \verbatim
35*>
36*> CUNGTSQR generates an M-by-N complex matrix Q_out with orthonormal
37*> columns, which are the first N columns of a product of comlpex unitary
38*> matrices of order M which are returned by CLATSQR
39*>
40*>      Q_out = first_N_columns_of( Q(1)_in * Q(2)_in * ... * Q(k)_in ).
41*>
42*> See the documentation for CLATSQR.
43*> \endverbatim
44*
45*  Arguments:
46*  ==========
47*
48*> \param[in] M
49*> \verbatim
50*>          M is INTEGER
51*>          The number of rows of the matrix A.  M >= 0.
52*> \endverbatim
53*>
54*> \param[in] N
55*> \verbatim
56*>          N is INTEGER
57*>          The number of columns of the matrix A. M >= N >= 0.
58*> \endverbatim
59*>
60*> \param[in] MB
61*> \verbatim
62*>          MB is INTEGER
63*>          The row block size used by CLATSQR to return
64*>          arrays A and T. MB > N.
65*>          (Note that if MB > M, then M is used instead of MB
66*>          as the row block size).
67*> \endverbatim
68*>
69*> \param[in] NB
70*> \verbatim
71*>          NB is INTEGER
72*>          The column block size used by CLATSQR to return
73*>          arrays A and T. NB >= 1.
74*>          (Note that if NB > N, then N is used instead of NB
75*>          as the column block size).
76*> \endverbatim
77*>
78*> \param[in,out] A
79*> \verbatim
80*>          A is COMPLEX array, dimension (LDA,N)
81*>
82*>          On entry:
83*>
84*>             The elements on and above the diagonal are not accessed.
85*>             The elements below the diagonal represent the unit
86*>             lower-trapezoidal blocked matrix V computed by CLATSQR
87*>             that defines the input matrices Q_in(k) (ones on the
88*>             diagonal are not stored) (same format as the output A
89*>             below the diagonal in CLATSQR).
90*>
91*>          On exit:
92*>
93*>             The array A contains an M-by-N orthonormal matrix Q_out,
94*>             i.e the columns of A are orthogonal unit vectors.
95*> \endverbatim
96*>
97*> \param[in] LDA
98*> \verbatim
99*>          LDA is INTEGER
100*>          The leading dimension of the array A.  LDA >= max(1,M).
101*> \endverbatim
102*>
103*> \param[in] T
104*> \verbatim
105*>          T is COMPLEX array,
106*>          dimension (LDT, N * NIRB)
107*>          where NIRB = Number_of_input_row_blocks
108*>                     = MAX( 1, CEIL((M-N)/(MB-N)) )
109*>          Let NICB = Number_of_input_col_blocks
110*>                   = CEIL(N/NB)
111*>
112*>          The upper-triangular block reflectors used to define the
113*>          input matrices Q_in(k), k=(1:NIRB*NICB). The block
114*>          reflectors are stored in compact form in NIRB block
115*>          reflector sequences. Each of NIRB block reflector sequences
116*>          is stored in a larger NB-by-N column block of T and consists
117*>          of NICB smaller NB-by-NB upper-triangular column blocks.
118*>          (same format as the output T in CLATSQR).
119*> \endverbatim
120*>
121*> \param[in] LDT
122*> \verbatim
123*>          LDT is INTEGER
124*>          The leading dimension of the array T.
125*>          LDT >= max(1,min(NB1,N)).
126*> \endverbatim
127*>
128*> \param[out] WORK
129*> \verbatim
130*>          (workspace) COMPLEX array, dimension (MAX(2,LWORK))
131*>          On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
132*> \endverbatim
133*>
134*> \param[in] LWORK
135*> \verbatim
136*>          The dimension of the array WORK.  LWORK >= (M+NB)*N.
137*>          If LWORK = -1, then a workspace query is assumed.
138*>          The routine only calculates the optimal size of the WORK
139*>          array, returns this value as the first entry of the WORK
140*>          array, and no error message related to LWORK is issued
141*>          by XERBLA.
142*> \endverbatim
143*>
144*> \param[out] INFO
145*> \verbatim
146*>          INFO is INTEGER
147*>          = 0:  successful exit
148*>          < 0:  if INFO = -i, the i-th argument had an illegal value
149*> \endverbatim
150*>
151*  Authors:
152*  ========
153*
154*> \author Univ. of Tennessee
155*> \author Univ. of California Berkeley
156*> \author Univ. of Colorado Denver
157*> \author NAG Ltd.
158*
159*> \ingroup complexOTHERcomputational
160*
161*> \par Contributors:
162*  ==================
163*>
164*> \verbatim
165*>
166*> November 2019, Igor Kozachenko,
167*>                Computer Science Division,
168*>                University of California, Berkeley
169*>
170*> \endverbatim
171*
172*  =====================================================================
173      SUBROUTINE CUNGTSQR( M, N, MB, NB, A, LDA, T, LDT, WORK, LWORK,
174     $                     INFO )
175      IMPLICIT NONE
176*
177*  -- LAPACK computational routine --
178*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
179*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
180*
181*     .. Scalar Arguments ..
182      INTEGER           INFO, LDA, LDT, LWORK, M, N, MB, NB
183*     ..
184*     .. Array Arguments ..
185      COMPLEX        A( LDA, * ), T( LDT, * ), WORK( * )
186*     ..
187*
188*  =====================================================================
189*
190*     .. Parameters ..
191      COMPLEX            CONE, CZERO
192      PARAMETER          ( CONE = ( 1.0E+0, 0.0E+0 ),
193     $                     CZERO = ( 0.0E+0, 0.0E+0 ) )
194*     ..
195*     .. Local Scalars ..
196      LOGICAL            LQUERY
197      INTEGER            IINFO, LDC, LWORKOPT, LC, LW, NBLOCAL, J
198*     ..
199*     .. External Subroutines ..
200      EXTERNAL           CCOPY, CLAMTSQR, CLASET, XERBLA
201*     ..
202*     .. Intrinsic Functions ..
203      INTRINSIC          CMPLX, MAX, MIN
204*     ..
205*     .. Executable Statements ..
206*
207*     Test the input parameters
208*
209      LQUERY  = LWORK.EQ.-1
210      INFO = 0
211      IF( M.LT.0 ) THEN
212         INFO = -1
213      ELSE IF( N.LT.0 .OR. M.LT.N ) THEN
214         INFO = -2
215      ELSE IF( MB.LE.N ) THEN
216         INFO = -3
217      ELSE IF( NB.LT.1 ) THEN
218         INFO = -4
219      ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
220         INFO = -6
221      ELSE IF( LDT.LT.MAX( 1, MIN( NB, N ) ) ) THEN
222         INFO = -8
223      ELSE
224*
225*        Test the input LWORK for the dimension of the array WORK.
226*        This workspace is used to store array C(LDC, N) and WORK(LWORK)
227*        in the call to CLAMTSQR. See the documentation for CLAMTSQR.
228*
229         IF( LWORK.LT.2 .AND. (.NOT.LQUERY) ) THEN
230            INFO = -10
231         ELSE
232*
233*           Set block size for column blocks
234*
235            NBLOCAL = MIN( NB, N )
236*
237*           LWORK = -1, then set the size for the array C(LDC,N)
238*           in CLAMTSQR call and set the optimal size of the work array
239*           WORK(LWORK) in CLAMTSQR call.
240*
241            LDC = M
242            LC = LDC*N
243            LW = N * NBLOCAL
244*
245            LWORKOPT = LC+LW
246*
247            IF( ( LWORK.LT.MAX( 1, LWORKOPT ) ).AND.(.NOT.LQUERY) ) THEN
248               INFO = -10
249            END IF
250         END IF
251*
252      END IF
253*
254*     Handle error in the input parameters and return workspace query.
255*
256      IF( INFO.NE.0 ) THEN
257         CALL XERBLA( 'CUNGTSQR', -INFO )
258         RETURN
259      ELSE IF ( LQUERY ) THEN
260         WORK( 1 ) = CMPLX( LWORKOPT )
261         RETURN
262      END IF
263*
264*     Quick return if possible
265*
266      IF( MIN( M, N ).EQ.0 ) THEN
267         WORK( 1 ) = CMPLX( LWORKOPT )
268         RETURN
269      END IF
270*
271*     (1) Form explicitly the tall-skinny M-by-N left submatrix Q1_in
272*     of M-by-M orthogonal matrix Q_in, which is implicitly stored in
273*     the subdiagonal part of input array A and in the input array T.
274*     Perform by the following operation using the routine CLAMTSQR.
275*
276*         Q1_in = Q_in * ( I ), where I is a N-by-N identity matrix,
277*                        ( 0 )        0 is a (M-N)-by-N zero matrix.
278*
279*     (1a) Form M-by-N matrix in the array WORK(1:LDC*N) with ones
280*     on the diagonal and zeros elsewhere.
281*
282      CALL CLASET( 'F', M, N, CZERO, CONE, WORK, LDC )
283*
284*     (1b)  On input, WORK(1:LDC*N) stores ( I );
285*                                          ( 0 )
286*
287*           On output, WORK(1:LDC*N) stores Q1_in.
288*
289      CALL CLAMTSQR( 'L', 'N', M, N, N, MB, NBLOCAL, A, LDA, T, LDT,
290     $               WORK, LDC, WORK( LC+1 ), LW, IINFO )
291*
292*     (2) Copy the result from the part of the work array (1:M,1:N)
293*     with the leading dimension LDC that starts at WORK(1) into
294*     the output array A(1:M,1:N) column-by-column.
295*
296      DO J = 1, N
297         CALL CCOPY( M, WORK( (J-1)*LDC + 1 ), 1, A( 1, J ), 1 )
298      END DO
299*
300      WORK( 1 ) = CMPLX( LWORKOPT )
301      RETURN
302*
303*     End of CUNGTSQR
304*
305      END
306