1*> \brief \b ZTPLQT2 computes a LQ factorization of a real or complex "triangular-pentagonal" matrix, which is composed of a triangular block and a pentagonal block, using the compact WY representation for Q.
2*
3*  =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6*            http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download ZTPLQT2 + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/ztplqt2.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/ztplqt2.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/ztplqt2.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18*  Definition:
19*  ===========
20*
21*       SUBROUTINE ZTPLQT2( M, N, L, A, LDA, B, LDB, T, LDT, INFO )
22*
23*       .. Scalar Arguments ..
24*       INTEGER   INFO, LDA, LDB, LDT, N, M, L
25*       ..
26*       .. Array Arguments ..
27*       COMPLEX*16   A( LDA, * ), B( LDB, * ), T( LDT, * )
28*       ..
29*
30*
31*> \par Purpose:
32*  =============
33*>
34*> \verbatim
35*>
36*> ZTPLQT2 computes a LQ a factorization of a complex "triangular-pentagonal"
37*> matrix C, which is composed of a triangular block A and pentagonal block B,
38*> using the compact WY representation for Q.
39*> \endverbatim
40*
41*  Arguments:
42*  ==========
43*
44*> \param[in] M
45*> \verbatim
46*>          M is INTEGER
47*>          The total number of rows of the matrix B.
48*>          M >= 0.
49*> \endverbatim
50*>
51*> \param[in] N
52*> \verbatim
53*>          N is INTEGER
54*>          The number of columns of the matrix B, and the order of
55*>          the triangular matrix A.
56*>          N >= 0.
57*> \endverbatim
58*>
59*> \param[in] L
60*> \verbatim
61*>          L is INTEGER
62*>          The number of rows of the lower trapezoidal part of B.
63*>          MIN(M,N) >= L >= 0.  See Further Details.
64*> \endverbatim
65*>
66*> \param[in,out] A
67*> \verbatim
68*>          A is COMPLEX*16 array, dimension (LDA,M)
69*>          On entry, the lower triangular M-by-M matrix A.
70*>          On exit, the elements on and below the diagonal of the array
71*>          contain the lower triangular matrix L.
72*> \endverbatim
73*>
74*> \param[in] LDA
75*> \verbatim
76*>          LDA is INTEGER
77*>          The leading dimension of the array A.  LDA >= max(1,M).
78*> \endverbatim
79*>
80*> \param[in,out] B
81*> \verbatim
82*>          B is COMPLEX*16 array, dimension (LDB,N)
83*>          On entry, the pentagonal M-by-N matrix B.  The first N-L columns
84*>          are rectangular, and the last L columns are lower trapezoidal.
85*>          On exit, B contains the pentagonal matrix V.  See Further Details.
86*> \endverbatim
87*>
88*> \param[in] LDB
89*> \verbatim
90*>          LDB is INTEGER
91*>          The leading dimension of the array B.  LDB >= max(1,M).
92*> \endverbatim
93*>
94*> \param[out] T
95*> \verbatim
96*>          T is COMPLEX*16 array, dimension (LDT,M)
97*>          The N-by-N upper triangular factor T of the block reflector.
98*>          See Further Details.
99*> \endverbatim
100*>
101*> \param[in] LDT
102*> \verbatim
103*>          LDT is INTEGER
104*>          The leading dimension of the array T.  LDT >= max(1,M)
105*> \endverbatim
106*>
107*> \param[out] INFO
108*> \verbatim
109*>          INFO is INTEGER
110*>          = 0: successful exit
111*>          < 0: if INFO = -i, the i-th argument had an illegal value
112*> \endverbatim
113*
114*  Authors:
115*  ========
116*
117*> \author Univ. of Tennessee
118*> \author Univ. of California Berkeley
119*> \author Univ. of Colorado Denver
120*> \author NAG Ltd.
121*
122*> \ingroup doubleOTHERcomputational
123*
124*> \par Further Details:
125*  =====================
126*>
127*> \verbatim
128*>
129*>  The input matrix C is a M-by-(M+N) matrix
130*>
131*>               C = [ A ][ B ]
132*>
133*>
134*>  where A is an lower triangular M-by-M matrix, and B is M-by-N pentagonal
135*>  matrix consisting of a M-by-(N-L) rectangular matrix B1 left of a M-by-L
136*>  upper trapezoidal matrix B2:
137*>
138*>               B = [ B1 ][ B2 ]
139*>                   [ B1 ]  <-     M-by-(N-L) rectangular
140*>                   [ B2 ]  <-     M-by-L lower trapezoidal.
141*>
142*>  The lower trapezoidal matrix B2 consists of the first L columns of a
143*>  N-by-N lower triangular matrix, where 0 <= L <= MIN(M,N).  If L=0,
144*>  B is rectangular M-by-N; if M=L=N, B is lower triangular.
145*>
146*>  The matrix W stores the elementary reflectors H(i) in the i-th row
147*>  above the diagonal (of A) in the M-by-(M+N) input matrix C
148*>
149*>               C = [ A ][ B ]
150*>                   [ A ]  <- lower triangular M-by-M
151*>                   [ B ]  <- M-by-N pentagonal
152*>
153*>  so that W can be represented as
154*>
155*>               W = [ I ][ V ]
156*>                   [ I ]  <- identity, M-by-M
157*>                   [ V ]  <- M-by-N, same form as B.
158*>
159*>  Thus, all of information needed for W is contained on exit in B, which
160*>  we call V above.  Note that V has the same form as B; that is,
161*>
162*>               W = [ V1 ][ V2 ]
163*>                   [ V1 ] <-     M-by-(N-L) rectangular
164*>                   [ V2 ] <-     M-by-L lower trapezoidal.
165*>
166*>  The rows of V represent the vectors which define the H(i)'s.
167*>  The (M+N)-by-(M+N) block reflector H is then given by
168*>
169*>               H = I - W**T * T * W
170*>
171*>  where W^H is the conjugate transpose of W and T is the upper triangular
172*>  factor of the block reflector.
173*> \endverbatim
174*>
175*  =====================================================================
176      SUBROUTINE ZTPLQT2( M, N, L, A, LDA, B, LDB, T, LDT, INFO )
177*
178*  -- LAPACK computational routine --
179*  -- LAPACK is a software package provided by Univ. of Tennessee,    --
180*  -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
181*
182*     .. Scalar Arguments ..
183      INTEGER        INFO, LDA, LDB, LDT, N, M, L
184*     ..
185*     .. Array Arguments ..
186      COMPLEX*16     A( LDA, * ), B( LDB, * ), T( LDT, * )
187*     ..
188*
189*  =====================================================================
190*
191*     .. Parameters ..
192      COMPLEX*16  ONE, ZERO
193      PARAMETER( ZERO = ( 0.0D+0, 0.0D+0 ),ONE  = ( 1.0D+0, 0.0D+0 ) )
194*     ..
195*     .. Local Scalars ..
196      INTEGER   I, J, P, MP, NP
197      COMPLEX*16   ALPHA
198*     ..
199*     .. External Subroutines ..
200      EXTERNAL  ZLARFG, ZGEMV, ZGERC, ZTRMV, XERBLA
201*     ..
202*     .. Intrinsic Functions ..
203      INTRINSIC MAX, MIN
204*     ..
205*     .. Executable Statements ..
206*
207*     Test the input arguments
208*
209      INFO = 0
210      IF( M.LT.0 ) THEN
211         INFO = -1
212      ELSE IF( N.LT.0 ) THEN
213         INFO = -2
214      ELSE IF( L.LT.0 .OR. L.GT.MIN(M,N) ) THEN
215         INFO = -3
216      ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
217         INFO = -5
218      ELSE IF( LDB.LT.MAX( 1, M ) ) THEN
219         INFO = -7
220      ELSE IF( LDT.LT.MAX( 1, M ) ) THEN
221         INFO = -9
222      END IF
223      IF( INFO.NE.0 ) THEN
224         CALL XERBLA( 'ZTPLQT2', -INFO )
225         RETURN
226      END IF
227*
228*     Quick return if possible
229*
230      IF( N.EQ.0 .OR. M.EQ.0 ) RETURN
231*
232      DO I = 1, M
233*
234*        Generate elementary reflector H(I) to annihilate B(I,:)
235*
236         P = N-L+MIN( L, I )
237         CALL ZLARFG( P+1, A( I, I ), B( I, 1 ), LDB, T( 1, I ) )
238         T(1,I)=CONJG(T(1,I))
239         IF( I.LT.M ) THEN
240            DO J = 1, P
241               B( I, J ) = CONJG(B(I,J))
242            END DO
243*
244*           W(M-I:1) := C(I+1:M,I:N) * C(I,I:N) [use W = T(M,:)]
245*
246            DO J = 1, M-I
247               T( M, J ) = (A( I+J, I ))
248            END DO
249            CALL ZGEMV( 'N', M-I, P, ONE, B( I+1, 1 ), LDB,
250     $                  B( I, 1 ), LDB, ONE, T( M, 1 ), LDT )
251*
252*           C(I+1:M,I:N) = C(I+1:M,I:N) + alpha * C(I,I:N)*W(M-1:1)^H
253*
254            ALPHA = -(T( 1, I ))
255            DO J = 1, M-I
256               A( I+J, I ) = A( I+J, I ) + ALPHA*(T( M, J ))
257            END DO
258            CALL ZGERC( M-I, P, (ALPHA),  T( M, 1 ), LDT,
259     $          B( I, 1 ), LDB, B( I+1, 1 ), LDB )
260            DO J = 1, P
261               B( I, J ) = CONJG(B(I,J))
262            END DO
263         END IF
264      END DO
265*
266      DO I = 2, M
267*
268*        T(I,1:I-1) := C(I:I-1,1:N)**H * (alpha * C(I,I:N))
269*
270         ALPHA = -(T( 1, I ))
271         DO J = 1, I-1
272            T( I, J ) = ZERO
273         END DO
274         P = MIN( I-1, L )
275         NP = MIN( N-L+1, N )
276         MP = MIN( P+1, M )
277         DO J = 1, N-L+P
278           B(I,J)=CONJG(B(I,J))
279         END DO
280*
281*        Triangular part of B2
282*
283         DO J = 1, P
284            T( I, J ) = (ALPHA*B( I, N-L+J ))
285         END DO
286         CALL ZTRMV( 'L', 'N', 'N', P, B( 1, NP ), LDB,
287     $               T( I, 1 ), LDT )
288*
289*        Rectangular part of B2
290*
291         CALL ZGEMV( 'N', I-1-P, L,  ALPHA, B( MP, NP ), LDB,
292     $               B( I, NP ), LDB, ZERO, T( I,MP ), LDT )
293*
294*        B1
295
296*
297         CALL ZGEMV( 'N', I-1, N-L, ALPHA, B, LDB, B( I, 1 ), LDB,
298     $               ONE, T( I, 1 ), LDT )
299*
300
301*
302*        T(1:I-1,I) := T(1:I-1,1:I-1) * T(I,1:I-1)
303*
304         DO J = 1, I-1
305            T(I,J)=CONJG(T(I,J))
306         END DO
307         CALL ZTRMV( 'L', 'C', 'N', I-1, T, LDT, T( I, 1 ), LDT )
308         DO J = 1, I-1
309            T(I,J)=CONJG(T(I,J))
310         END DO
311         DO J = 1, N-L+P
312            B(I,J)=CONJG(B(I,J))
313         END DO
314*
315*        T(I,I) = tau(I)
316*
317         T( I, I ) = T( 1, I )
318         T( 1, I ) = ZERO
319      END DO
320      DO I=1,M
321         DO J= I+1,M
322            T(I,J)=(T(J,I))
323            T(J,I)=ZERO
324         END DO
325      END DO
326
327*
328*     End of ZTPLQT2
329*
330      END
331