LAPACK  3.10.0
LAPACK: Linear Algebra PACKage
zunm22.f
Go to the documentation of this file.
1 *> \brief \b ZUNM22 multiplies a general matrix by a banded unitary matrix.
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download ZUNM22 + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zunm22.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zunm22.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zunm22.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE ZUNM22( SIDE, TRANS, M, N, N1, N2, Q, LDQ, C, LDC,
22 * $ WORK, LWORK, INFO )
23 *
24 * .. Scalar Arguments ..
25 * CHARACTER SIDE, TRANS
26 * INTEGER M, N, N1, N2, LDQ, LDC, LWORK, INFO
27 * ..
28 * .. Array Arguments ..
29 * COMPLEX*16 Q( LDQ, * ), C( LDC, * ), WORK( * )
30 * ..
31 *
32 *> \par Purpose
33 * ============
34 *>
35 *> \verbatim
36 *>
37 *> ZUNM22 overwrites the general complex M-by-N matrix C with
38 *>
39 *> SIDE = 'L' SIDE = 'R'
40 *> TRANS = 'N': Q * C C * Q
41 *> TRANS = 'C': Q**H * C C * Q**H
42 *>
43 *> where Q is a complex unitary matrix of order NQ, with NQ = M if
44 *> SIDE = 'L' and NQ = N if SIDE = 'R'.
45 *> The unitary matrix Q processes a 2-by-2 block structure
46 *>
47 *> [ Q11 Q12 ]
48 *> Q = [ ]
49 *> [ Q21 Q22 ],
50 *>
51 *> where Q12 is an N1-by-N1 lower triangular matrix and Q21 is an
52 *> N2-by-N2 upper triangular matrix.
53 *> \endverbatim
54 *
55 * Arguments:
56 * ==========
57 *
58 *> \param[in] SIDE
59 *> \verbatim
60 *> SIDE is CHARACTER*1
61 *> = 'L': apply Q or Q**H from the Left;
62 *> = 'R': apply Q or Q**H from the Right.
63 *> \endverbatim
64 *>
65 *> \param[in] TRANS
66 *> \verbatim
67 *> TRANS is CHARACTER*1
68 *> = 'N': apply Q (No transpose);
69 *> = 'C': apply Q**H (Conjugate transpose).
70 *> \endverbatim
71 *>
72 *> \param[in] M
73 *> \verbatim
74 *> M is INTEGER
75 *> The number of rows of the matrix C. M >= 0.
76 *> \endverbatim
77 *>
78 *> \param[in] N
79 *> \verbatim
80 *> N is INTEGER
81 *> The number of columns of the matrix C. N >= 0.
82 *> \endverbatim
83 *>
84 *> \param[in] N1
85 *> \param[in] N2
86 *> \verbatim
87 *> N1 is INTEGER
88 *> N2 is INTEGER
89 *> The dimension of Q12 and Q21, respectively. N1, N2 >= 0.
90 *> The following requirement must be satisfied:
91 *> N1 + N2 = M if SIDE = 'L' and N1 + N2 = N if SIDE = 'R'.
92 *> \endverbatim
93 *>
94 *> \param[in] Q
95 *> \verbatim
96 *> Q is COMPLEX*16 array, dimension
97 *> (LDQ,M) if SIDE = 'L'
98 *> (LDQ,N) if SIDE = 'R'
99 *> \endverbatim
100 *>
101 *> \param[in] LDQ
102 *> \verbatim
103 *> LDQ is INTEGER
104 *> The leading dimension of the array Q.
105 *> LDQ >= max(1,M) if SIDE = 'L'; LDQ >= max(1,N) if SIDE = 'R'.
106 *> \endverbatim
107 *>
108 *> \param[in,out] C
109 *> \verbatim
110 *> C is COMPLEX*16 array, dimension (LDC,N)
111 *> On entry, the M-by-N matrix C.
112 *> On exit, C is overwritten by Q*C or Q**H*C or C*Q**H or C*Q.
113 *> \endverbatim
114 *>
115 *> \param[in] LDC
116 *> \verbatim
117 *> LDC is INTEGER
118 *> The leading dimension of the array C. LDC >= max(1,M).
119 *> \endverbatim
120 *>
121 *> \param[out] WORK
122 *> \verbatim
123 *> WORK is COMPLEX*16 array, dimension (MAX(1,LWORK))
124 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
125 *> \endverbatim
126 *>
127 *> \param[in] LWORK
128 *> \verbatim
129 *> LWORK is INTEGER
130 *> The dimension of the array WORK.
131 *> If SIDE = 'L', LWORK >= max(1,N);
132 *> if SIDE = 'R', LWORK >= max(1,M).
133 *> For optimum performance LWORK >= M*N.
134 *>
135 *> If LWORK = -1, then a workspace query is assumed; the routine
136 *> only calculates the optimal size of the WORK array, returns
137 *> this value as the first entry of the WORK array, and no error
138 *> message related to LWORK is issued by XERBLA.
139 *> \endverbatim
140 *>
141 *> \param[out] INFO
142 *> \verbatim
143 *> INFO is INTEGER
144 *> = 0: successful exit
145 *> < 0: if INFO = -i, the i-th argument had an illegal value
146 *> \endverbatim
147 *
148 *
149 * Authors:
150 * ========
151 *
152 *> \author Univ. of Tennessee
153 *> \author Univ. of California Berkeley
154 *> \author Univ. of Colorado Denver
155 *> \author NAG Ltd.
156 *
157 *> \ingroup complexOTHERcomputational
158 *
159 * =====================================================================
160  SUBROUTINE zunm22( SIDE, TRANS, M, N, N1, N2, Q, LDQ, C, LDC,
161  $ WORK, LWORK, INFO )
162 *
163 * -- LAPACK computational routine --
164 * -- LAPACK is a software package provided by Univ. of Tennessee, --
165 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
166 *
167  IMPLICIT NONE
168 *
169 * .. Scalar Arguments ..
170  CHARACTER SIDE, TRANS
171  INTEGER M, N, N1, N2, LDQ, LDC, LWORK, INFO
172 * ..
173 * .. Array Arguments ..
174  COMPLEX*16 Q( LDQ, * ), C( LDC, * ), WORK( * )
175 * ..
176 *
177 * =====================================================================
178 *
179 * .. Parameters ..
180  COMPLEX*16 ONE
181  parameter( one = ( 1.0d+0, 0.0d+0 ) )
182 *
183 * .. Local Scalars ..
184  LOGICAL LEFT, LQUERY, NOTRAN
185  INTEGER I, LDWORK, LEN, LWKOPT, NB, NQ, NW
186 * ..
187 * .. External Functions ..
188  LOGICAL LSAME
189  EXTERNAL lsame
190 * ..
191 * .. External Subroutines ..
192  EXTERNAL zgemm, zlacpy, ztrmm, xerbla
193 * ..
194 * .. Intrinsic Functions ..
195  INTRINSIC dcmplx, max, min
196 * ..
197 * .. Executable Statements ..
198 *
199 * Test the input arguments
200 *
201  info = 0
202  left = lsame( side, 'L' )
203  notran = lsame( trans, 'N' )
204  lquery = ( lwork.EQ.-1 )
205 *
206 * NQ is the order of Q;
207 * NW is the minimum dimension of WORK.
208 *
209  IF( left ) THEN
210  nq = m
211  ELSE
212  nq = n
213  END IF
214  nw = nq
215  IF( n1.EQ.0 .OR. n2.EQ.0 ) nw = 1
216  IF( .NOT.left .AND. .NOT.lsame( side, 'R' ) ) THEN
217  info = -1
218  ELSE IF( .NOT.lsame( trans, 'N' ) .AND. .NOT.lsame( trans, 'C' ) )
219  $ THEN
220  info = -2
221  ELSE IF( m.LT.0 ) THEN
222  info = -3
223  ELSE IF( n.LT.0 ) THEN
224  info = -4
225  ELSE IF( n1.LT.0 .OR. n1+n2.NE.nq ) THEN
226  info = -5
227  ELSE IF( n2.LT.0 ) THEN
228  info = -6
229  ELSE IF( ldq.LT.max( 1, nq ) ) THEN
230  info = -8
231  ELSE IF( ldc.LT.max( 1, m ) ) THEN
232  info = -10
233  ELSE IF( lwork.LT.nw .AND. .NOT.lquery ) THEN
234  info = -12
235  END IF
236 *
237  IF( info.EQ.0 ) THEN
238  lwkopt = m*n
239  work( 1 ) = dcmplx( lwkopt )
240  END IF
241 *
242  IF( info.NE.0 ) THEN
243  CALL xerbla( 'ZUNM22', -info )
244  RETURN
245  ELSE IF( lquery ) THEN
246  RETURN
247  END IF
248 *
249 * Quick return if possible
250 *
251  IF( m.EQ.0 .OR. n.EQ.0 ) THEN
252  work( 1 ) = 1
253  RETURN
254  END IF
255 *
256 * Degenerate cases (N1 = 0 or N2 = 0) are handled using ZTRMM.
257 *
258  IF( n1.EQ.0 ) THEN
259  CALL ztrmm( side, 'Upper', trans, 'Non-Unit', m, n, one,
260  $ q, ldq, c, ldc )
261  work( 1 ) = one
262  RETURN
263  ELSE IF( n2.EQ.0 ) THEN
264  CALL ztrmm( side, 'Lower', trans, 'Non-Unit', m, n, one,
265  $ q, ldq, c, ldc )
266  work( 1 ) = one
267  RETURN
268  END IF
269 *
270 * Compute the largest chunk size available from the workspace.
271 *
272  nb = max( 1, min( lwork, lwkopt ) / nq )
273 *
274  IF( left ) THEN
275  IF( notran ) THEN
276  DO i = 1, n, nb
277  len = min( nb, n-i+1 )
278  ldwork = m
279 *
280 * Multiply bottom part of C by Q12.
281 *
282  CALL zlacpy( 'All', n1, len, c( n2+1, i ), ldc, work,
283  $ ldwork )
284  CALL ztrmm( 'Left', 'Lower', 'No Transpose', 'Non-Unit',
285  $ n1, len, one, q( 1, n2+1 ), ldq, work,
286  $ ldwork )
287 *
288 * Multiply top part of C by Q11.
289 *
290  CALL zgemm( 'No Transpose', 'No Transpose', n1, len, n2,
291  $ one, q, ldq, c( 1, i ), ldc, one, work,
292  $ ldwork )
293 *
294 * Multiply top part of C by Q21.
295 *
296  CALL zlacpy( 'All', n2, len, c( 1, i ), ldc,
297  $ work( n1+1 ), ldwork )
298  CALL ztrmm( 'Left', 'Upper', 'No Transpose', 'Non-Unit',
299  $ n2, len, one, q( n1+1, 1 ), ldq,
300  $ work( n1+1 ), ldwork )
301 *
302 * Multiply bottom part of C by Q22.
303 *
304  CALL zgemm( 'No Transpose', 'No Transpose', n2, len, n1,
305  $ one, q( n1+1, n2+1 ), ldq, c( n2+1, i ), ldc,
306  $ one, work( n1+1 ), ldwork )
307 *
308 * Copy everything back.
309 *
310  CALL zlacpy( 'All', m, len, work, ldwork, c( 1, i ),
311  $ ldc )
312  END DO
313  ELSE
314  DO i = 1, n, nb
315  len = min( nb, n-i+1 )
316  ldwork = m
317 *
318 * Multiply bottom part of C by Q21**H.
319 *
320  CALL zlacpy( 'All', n2, len, c( n1+1, i ), ldc, work,
321  $ ldwork )
322  CALL ztrmm( 'Left', 'Upper', 'Conjugate', 'Non-Unit',
323  $ n2, len, one, q( n1+1, 1 ), ldq, work,
324  $ ldwork )
325 *
326 * Multiply top part of C by Q11**H.
327 *
328  CALL zgemm( 'Conjugate', 'No Transpose', n2, len, n1,
329  $ one, q, ldq, c( 1, i ), ldc, one, work,
330  $ ldwork )
331 *
332 * Multiply top part of C by Q12**H.
333 *
334  CALL zlacpy( 'All', n1, len, c( 1, i ), ldc,
335  $ work( n2+1 ), ldwork )
336  CALL ztrmm( 'Left', 'Lower', 'Conjugate', 'Non-Unit',
337  $ n1, len, one, q( 1, n2+1 ), ldq,
338  $ work( n2+1 ), ldwork )
339 *
340 * Multiply bottom part of C by Q22**H.
341 *
342  CALL zgemm( 'Conjugate', 'No Transpose', n1, len, n2,
343  $ one, q( n1+1, n2+1 ), ldq, c( n1+1, i ), ldc,
344  $ one, work( n2+1 ), ldwork )
345 *
346 * Copy everything back.
347 *
348  CALL zlacpy( 'All', m, len, work, ldwork, c( 1, i ),
349  $ ldc )
350  END DO
351  END IF
352  ELSE
353  IF( notran ) THEN
354  DO i = 1, m, nb
355  len = min( nb, m-i+1 )
356  ldwork = len
357 *
358 * Multiply right part of C by Q21.
359 *
360  CALL zlacpy( 'All', len, n2, c( i, n1+1 ), ldc, work,
361  $ ldwork )
362  CALL ztrmm( 'Right', 'Upper', 'No Transpose', 'Non-Unit',
363  $ len, n2, one, q( n1+1, 1 ), ldq, work,
364  $ ldwork )
365 *
366 * Multiply left part of C by Q11.
367 *
368  CALL zgemm( 'No Transpose', 'No Transpose', len, n2, n1,
369  $ one, c( i, 1 ), ldc, q, ldq, one, work,
370  $ ldwork )
371 *
372 * Multiply left part of C by Q12.
373 *
374  CALL zlacpy( 'All', len, n1, c( i, 1 ), ldc,
375  $ work( 1 + n2*ldwork ), ldwork )
376  CALL ztrmm( 'Right', 'Lower', 'No Transpose', 'Non-Unit',
377  $ len, n1, one, q( 1, n2+1 ), ldq,
378  $ work( 1 + n2*ldwork ), ldwork )
379 *
380 * Multiply right part of C by Q22.
381 *
382  CALL zgemm( 'No Transpose', 'No Transpose', len, n1, n2,
383  $ one, c( i, n1+1 ), ldc, q( n1+1, n2+1 ), ldq,
384  $ one, work( 1 + n2*ldwork ), ldwork )
385 *
386 * Copy everything back.
387 *
388  CALL zlacpy( 'All', len, n, work, ldwork, c( i, 1 ),
389  $ ldc )
390  END DO
391  ELSE
392  DO i = 1, m, nb
393  len = min( nb, m-i+1 )
394  ldwork = len
395 *
396 * Multiply right part of C by Q12**H.
397 *
398  CALL zlacpy( 'All', len, n1, c( i, n2+1 ), ldc, work,
399  $ ldwork )
400  CALL ztrmm( 'Right', 'Lower', 'Conjugate', 'Non-Unit',
401  $ len, n1, one, q( 1, n2+1 ), ldq, work,
402  $ ldwork )
403 *
404 * Multiply left part of C by Q11**H.
405 *
406  CALL zgemm( 'No Transpose', 'Conjugate', len, n1, n2,
407  $ one, c( i, 1 ), ldc, q, ldq, one, work,
408  $ ldwork )
409 *
410 * Multiply left part of C by Q21**H.
411 *
412  CALL zlacpy( 'All', len, n2, c( i, 1 ), ldc,
413  $ work( 1 + n1*ldwork ), ldwork )
414  CALL ztrmm( 'Right', 'Upper', 'Conjugate', 'Non-Unit',
415  $ len, n2, one, q( n1+1, 1 ), ldq,
416  $ work( 1 + n1*ldwork ), ldwork )
417 *
418 * Multiply right part of C by Q22**H.
419 *
420  CALL zgemm( 'No Transpose', 'Conjugate', len, n2, n1,
421  $ one, c( i, n2+1 ), ldc, q( n1+1, n2+1 ), ldq,
422  $ one, work( 1 + n1*ldwork ), ldwork )
423 *
424 * Copy everything back.
425 *
426  CALL zlacpy( 'All', len, n, work, ldwork, c( i, 1 ),
427  $ ldc )
428  END DO
429  END IF
430  END IF
431 *
432  work( 1 ) = dcmplx( lwkopt )
433  RETURN
434 *
435 * End of ZUNM22
436 *
437  END
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine zgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
ZGEMM
Definition: zgemm.f:187
subroutine ztrmm(SIDE, UPLO, TRANSA, DIAG, M, N, ALPHA, A, LDA, B, LDB)
ZTRMM
Definition: ztrmm.f:177
subroutine zlacpy(UPLO, M, N, A, LDA, B, LDB)
ZLACPY copies all or part of one two-dimensional array to another.
Definition: zlacpy.f:103
subroutine zunm22(SIDE, TRANS, M, N, N1, N2, Q, LDQ, C, LDC, WORK, LWORK, INFO)
ZUNM22 multiplies a general matrix by a banded unitary matrix.
Definition: zunm22.f:162