1 SUBROUTINE ZUNMBR( VECT, SIDE, TRANS, M, N, K, A, LDA, TAU, C,
2 $ LDC, WORK, LWORK, INFO )
3 *
4 * -- LAPACK routine (version 3.2) --
5 * -- LAPACK is a software package provided by Univ. of Tennessee, --
6 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7 * November 2006
8 *
9 * .. Scalar Arguments ..
10 CHARACTER SIDE, TRANS, VECT
11 INTEGER INFO, K, LDA, LDC, LWORK, M, N
12 * ..
13 * .. Array Arguments ..
14 COMPLEX*16 A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * If VECT = 'Q', ZUNMBR overwrites the general complex M-by-N matrix C
21 * with
22 * SIDE = 'L' SIDE = 'R'
23 * TRANS = 'N': Q * C C * Q
24 * TRANS = 'C': Q**H * C C * Q**H
25 *
26 * If VECT = 'P', ZUNMBR overwrites the general complex M-by-N matrix C
27 * with
28 * SIDE = 'L' SIDE = 'R'
29 * TRANS = 'N': P * C C * P
30 * TRANS = 'C': P**H * C C * P**H
31 *
32 * Here Q and P**H are the unitary matrices determined by ZGEBRD when
33 * reducing a complex matrix A to bidiagonal form: A = Q * B * P**H. Q
34 * and P**H are defined as products of elementary reflectors H(i) and
35 * G(i) respectively.
36 *
37 * Let nq = m if SIDE = 'L' and nq = n if SIDE = 'R'. Thus nq is the
38 * order of the unitary matrix Q or P**H that is applied.
39 *
40 * If VECT = 'Q', A is assumed to have been an NQ-by-K matrix:
41 * if nq >= k, Q = H(1) H(2) . . . H(k);
42 * if nq < k, Q = H(1) H(2) . . . H(nq-1).
43 *
44 * If VECT = 'P', A is assumed to have been a K-by-NQ matrix:
45 * if k < nq, P = G(1) G(2) . . . G(k);
46 * if k >= nq, P = G(1) G(2) . . . G(nq-1).
47 *
48 * Arguments
49 * =========
50 *
51 * VECT (input) CHARACTER*1
52 * = 'Q': apply Q or Q**H;
53 * = 'P': apply P or P**H.
54 *
55 * SIDE (input) CHARACTER*1
56 * = 'L': apply Q, Q**H, P or P**H from the Left;
57 * = 'R': apply Q, Q**H, P or P**H from the Right.
58 *
59 * TRANS (input) CHARACTER*1
60 * = 'N': No transpose, apply Q or P;
61 * = 'C': Conjugate transpose, apply Q**H or P**H.
62 *
63 * M (input) INTEGER
64 * The number of rows of the matrix C. M >= 0.
65 *
66 * N (input) INTEGER
67 * The number of columns of the matrix C. N >= 0.
68 *
69 * K (input) INTEGER
70 * If VECT = 'Q', the number of columns in the original
71 * matrix reduced by ZGEBRD.
72 * If VECT = 'P', the number of rows in the original
73 * matrix reduced by ZGEBRD.
74 * K >= 0.
75 *
76 * A (input) COMPLEX*16 array, dimension
77 * (LDA,min(nq,K)) if VECT = 'Q'
78 * (LDA,nq) if VECT = 'P'
79 * The vectors which define the elementary reflectors H(i) and
80 * G(i), whose products determine the matrices Q and P, as
81 * returned by ZGEBRD.
82 *
83 * LDA (input) INTEGER
84 * The leading dimension of the array A.
85 * If VECT = 'Q', LDA >= max(1,nq);
86 * if VECT = 'P', LDA >= max(1,min(nq,K)).
87 *
88 * TAU (input) COMPLEX*16 array, dimension (min(nq,K))
89 * TAU(i) must contain the scalar factor of the elementary
90 * reflector H(i) or G(i) which determines Q or P, as returned
91 * by ZGEBRD in the array argument TAUQ or TAUP.
92 *
93 * C (input/output) COMPLEX*16 array, dimension (LDC,N)
94 * On entry, the M-by-N matrix C.
95 * On exit, C is overwritten by Q*C or Q**H*C or C*Q**H or C*Q
96 * or P*C or P**H*C or C*P or C*P**H.
97 *
98 * LDC (input) INTEGER
99 * The leading dimension of the array C. LDC >= max(1,M).
100 *
101 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
102 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
103 *
104 * LWORK (input) INTEGER
105 * The dimension of the array WORK.
106 * If SIDE = 'L', LWORK >= max(1,N);
107 * if SIDE = 'R', LWORK >= max(1,M);
108 * if N = 0 or M = 0, LWORK >= 1.
109 * For optimum performance LWORK >= max(1,N*NB) if SIDE = 'L',
110 * and LWORK >= max(1,M*NB) if SIDE = 'R', where NB is the
111 * optimal blocksize. (NB = 0 if M = 0 or N = 0.)
112 *
113 * If LWORK = -1, then a workspace query is assumed; the routine
114 * only calculates the optimal size of the WORK array, returns
115 * this value as the first entry of the WORK array, and no error
116 * message related to LWORK is issued by XERBLA.
117 *
118 * INFO (output) INTEGER
119 * = 0: successful exit
120 * < 0: if INFO = -i, the i-th argument had an illegal value
121 *
122 * =====================================================================
123 *
124 * .. Local Scalars ..
125 LOGICAL APPLYQ, LEFT, LQUERY, NOTRAN
126 CHARACTER TRANST
127 INTEGER I1, I2, IINFO, LWKOPT, MI, NB, NI, NQ, NW
128 * ..
129 * .. External Functions ..
130 LOGICAL LSAME
131 INTEGER ILAENV
132 EXTERNAL LSAME, ILAENV
133 * ..
134 * .. External Subroutines ..
135 EXTERNAL XERBLA, ZUNMLQ, ZUNMQR
136 * ..
137 * .. Intrinsic Functions ..
138 INTRINSIC MAX, MIN
139 * ..
140 * .. Executable Statements ..
141 *
142 * Test the input arguments
143 *
144 INFO = 0
145 APPLYQ = LSAME( VECT, 'Q' )
146 LEFT = LSAME( SIDE, 'L' )
147 NOTRAN = LSAME( TRANS, 'N' )
148 LQUERY = ( LWORK.EQ.-1 )
149 *
150 * NQ is the order of Q or P and NW is the minimum dimension of WORK
151 *
152 IF( LEFT ) THEN
153 NQ = M
154 NW = N
155 ELSE
156 NQ = N
157 NW = M
158 END IF
159 IF( M.EQ.0 .OR. N.EQ.0 ) THEN
160 NW = 0
161 END IF
162 IF( .NOT.APPLYQ .AND. .NOT.LSAME( VECT, 'P' ) ) THEN
163 INFO = -1
164 ELSE IF( .NOT.LEFT .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
165 INFO = -2
166 ELSE IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN
167 INFO = -3
168 ELSE IF( M.LT.0 ) THEN
169 INFO = -4
170 ELSE IF( N.LT.0 ) THEN
171 INFO = -5
172 ELSE IF( K.LT.0 ) THEN
173 INFO = -6
174 ELSE IF( ( APPLYQ .AND. LDA.LT.MAX( 1, NQ ) ) .OR.
175 $ ( .NOT.APPLYQ .AND. LDA.LT.MAX( 1, MIN( NQ, K ) ) ) )
176 $ THEN
177 INFO = -8
178 ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
179 INFO = -11
180 ELSE IF( LWORK.LT.MAX( 1, NW ) .AND. .NOT.LQUERY ) THEN
181 INFO = -13
182 END IF
183 *
184 IF( INFO.EQ.0 ) THEN
185 IF( NW.GT.0 ) THEN
186 IF( APPLYQ ) THEN
187 IF( LEFT ) THEN
188 NB = ILAENV( 1, 'ZUNMQR', SIDE // TRANS, M-1, N, M-1,
189 $ -1 )
190 ELSE
191 NB = ILAENV( 1, 'ZUNMQR', SIDE // TRANS, M, N-1, N-1,
192 $ -1 )
193 END IF
194 ELSE
195 IF( LEFT ) THEN
196 NB = ILAENV( 1, 'ZUNMLQ', SIDE // TRANS, M-1, N, M-1,
197 $ -1 )
198 ELSE
199 NB = ILAENV( 1, 'ZUNMLQ', SIDE // TRANS, M, N-1, N-1,
200 $ -1 )
201 END IF
202 END IF
203 LWKOPT = MAX( 1, NW*NB )
204 ELSE
205 LWKOPT = 1
206 END IF
207 WORK( 1 ) = LWKOPT
208 END IF
209 *
210 IF( INFO.NE.0 ) THEN
211 CALL XERBLA( 'ZUNMBR', -INFO )
212 RETURN
213 ELSE IF( LQUERY ) THEN
214 RETURN
215 END IF
216 *
217 * Quick return if possible
218 *
219 IF( M.EQ.0 .OR. N.EQ.0 )
220 $ RETURN
221 *
222 IF( APPLYQ ) THEN
223 *
224 * Apply Q
225 *
226 IF( NQ.GE.K ) THEN
227 *
228 * Q was determined by a call to ZGEBRD with nq >= k
229 *
230 CALL ZUNMQR( SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC,
231 $ WORK, LWORK, IINFO )
232 ELSE IF( NQ.GT.1 ) THEN
233 *
234 * Q was determined by a call to ZGEBRD with nq < k
235 *
236 IF( LEFT ) THEN
237 MI = M - 1
238 NI = N
239 I1 = 2
240 I2 = 1
241 ELSE
242 MI = M
243 NI = N - 1
244 I1 = 1
245 I2 = 2
246 END IF
247 CALL ZUNMQR( SIDE, TRANS, MI, NI, NQ-1, A( 2, 1 ), LDA, TAU,
248 $ C( I1, I2 ), LDC, WORK, LWORK, IINFO )
249 END IF
250 ELSE
251 *
252 * Apply P
253 *
254 IF( NOTRAN ) THEN
255 TRANST = 'C'
256 ELSE
257 TRANST = 'N'
258 END IF
259 IF( NQ.GT.K ) THEN
260 *
261 * P was determined by a call to ZGEBRD with nq > k
262 *
263 CALL ZUNMLQ( SIDE, TRANST, M, N, K, A, LDA, TAU, C, LDC,
264 $ WORK, LWORK, IINFO )
265 ELSE IF( NQ.GT.1 ) THEN
266 *
267 * P was determined by a call to ZGEBRD with nq <= k
268 *
269 IF( LEFT ) THEN
270 MI = M - 1
271 NI = N
272 I1 = 2
273 I2 = 1
274 ELSE
275 MI = M
276 NI = N - 1
277 I1 = 1
278 I2 = 2
279 END IF
280 CALL ZUNMLQ( SIDE, TRANST, MI, NI, NQ-1, A( 1, 2 ), LDA,
281 $ TAU, C( I1, I2 ), LDC, WORK, LWORK, IINFO )
282 END IF
283 END IF
284 WORK( 1 ) = LWKOPT
285 RETURN
286 *
287 * End of ZUNMBR
288 *
289 END
2 $ LDC, WORK, LWORK, INFO )
3 *
4 * -- LAPACK routine (version 3.2) --
5 * -- LAPACK is a software package provided by Univ. of Tennessee, --
6 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
7 * November 2006
8 *
9 * .. Scalar Arguments ..
10 CHARACTER SIDE, TRANS, VECT
11 INTEGER INFO, K, LDA, LDC, LWORK, M, N
12 * ..
13 * .. Array Arguments ..
14 COMPLEX*16 A( LDA, * ), C( LDC, * ), TAU( * ), WORK( * )
15 * ..
16 *
17 * Purpose
18 * =======
19 *
20 * If VECT = 'Q', ZUNMBR overwrites the general complex M-by-N matrix C
21 * with
22 * SIDE = 'L' SIDE = 'R'
23 * TRANS = 'N': Q * C C * Q
24 * TRANS = 'C': Q**H * C C * Q**H
25 *
26 * If VECT = 'P', ZUNMBR overwrites the general complex M-by-N matrix C
27 * with
28 * SIDE = 'L' SIDE = 'R'
29 * TRANS = 'N': P * C C * P
30 * TRANS = 'C': P**H * C C * P**H
31 *
32 * Here Q and P**H are the unitary matrices determined by ZGEBRD when
33 * reducing a complex matrix A to bidiagonal form: A = Q * B * P**H. Q
34 * and P**H are defined as products of elementary reflectors H(i) and
35 * G(i) respectively.
36 *
37 * Let nq = m if SIDE = 'L' and nq = n if SIDE = 'R'. Thus nq is the
38 * order of the unitary matrix Q or P**H that is applied.
39 *
40 * If VECT = 'Q', A is assumed to have been an NQ-by-K matrix:
41 * if nq >= k, Q = H(1) H(2) . . . H(k);
42 * if nq < k, Q = H(1) H(2) . . . H(nq-1).
43 *
44 * If VECT = 'P', A is assumed to have been a K-by-NQ matrix:
45 * if k < nq, P = G(1) G(2) . . . G(k);
46 * if k >= nq, P = G(1) G(2) . . . G(nq-1).
47 *
48 * Arguments
49 * =========
50 *
51 * VECT (input) CHARACTER*1
52 * = 'Q': apply Q or Q**H;
53 * = 'P': apply P or P**H.
54 *
55 * SIDE (input) CHARACTER*1
56 * = 'L': apply Q, Q**H, P or P**H from the Left;
57 * = 'R': apply Q, Q**H, P or P**H from the Right.
58 *
59 * TRANS (input) CHARACTER*1
60 * = 'N': No transpose, apply Q or P;
61 * = 'C': Conjugate transpose, apply Q**H or P**H.
62 *
63 * M (input) INTEGER
64 * The number of rows of the matrix C. M >= 0.
65 *
66 * N (input) INTEGER
67 * The number of columns of the matrix C. N >= 0.
68 *
69 * K (input) INTEGER
70 * If VECT = 'Q', the number of columns in the original
71 * matrix reduced by ZGEBRD.
72 * If VECT = 'P', the number of rows in the original
73 * matrix reduced by ZGEBRD.
74 * K >= 0.
75 *
76 * A (input) COMPLEX*16 array, dimension
77 * (LDA,min(nq,K)) if VECT = 'Q'
78 * (LDA,nq) if VECT = 'P'
79 * The vectors which define the elementary reflectors H(i) and
80 * G(i), whose products determine the matrices Q and P, as
81 * returned by ZGEBRD.
82 *
83 * LDA (input) INTEGER
84 * The leading dimension of the array A.
85 * If VECT = 'Q', LDA >= max(1,nq);
86 * if VECT = 'P', LDA >= max(1,min(nq,K)).
87 *
88 * TAU (input) COMPLEX*16 array, dimension (min(nq,K))
89 * TAU(i) must contain the scalar factor of the elementary
90 * reflector H(i) or G(i) which determines Q or P, as returned
91 * by ZGEBRD in the array argument TAUQ or TAUP.
92 *
93 * C (input/output) COMPLEX*16 array, dimension (LDC,N)
94 * On entry, the M-by-N matrix C.
95 * On exit, C is overwritten by Q*C or Q**H*C or C*Q**H or C*Q
96 * or P*C or P**H*C or C*P or C*P**H.
97 *
98 * LDC (input) INTEGER
99 * The leading dimension of the array C. LDC >= max(1,M).
100 *
101 * WORK (workspace/output) COMPLEX*16 array, dimension (MAX(1,LWORK))
102 * On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
103 *
104 * LWORK (input) INTEGER
105 * The dimension of the array WORK.
106 * If SIDE = 'L', LWORK >= max(1,N);
107 * if SIDE = 'R', LWORK >= max(1,M);
108 * if N = 0 or M = 0, LWORK >= 1.
109 * For optimum performance LWORK >= max(1,N*NB) if SIDE = 'L',
110 * and LWORK >= max(1,M*NB) if SIDE = 'R', where NB is the
111 * optimal blocksize. (NB = 0 if M = 0 or N = 0.)
112 *
113 * If LWORK = -1, then a workspace query is assumed; the routine
114 * only calculates the optimal size of the WORK array, returns
115 * this value as the first entry of the WORK array, and no error
116 * message related to LWORK is issued by XERBLA.
117 *
118 * INFO (output) INTEGER
119 * = 0: successful exit
120 * < 0: if INFO = -i, the i-th argument had an illegal value
121 *
122 * =====================================================================
123 *
124 * .. Local Scalars ..
125 LOGICAL APPLYQ, LEFT, LQUERY, NOTRAN
126 CHARACTER TRANST
127 INTEGER I1, I2, IINFO, LWKOPT, MI, NB, NI, NQ, NW
128 * ..
129 * .. External Functions ..
130 LOGICAL LSAME
131 INTEGER ILAENV
132 EXTERNAL LSAME, ILAENV
133 * ..
134 * .. External Subroutines ..
135 EXTERNAL XERBLA, ZUNMLQ, ZUNMQR
136 * ..
137 * .. Intrinsic Functions ..
138 INTRINSIC MAX, MIN
139 * ..
140 * .. Executable Statements ..
141 *
142 * Test the input arguments
143 *
144 INFO = 0
145 APPLYQ = LSAME( VECT, 'Q' )
146 LEFT = LSAME( SIDE, 'L' )
147 NOTRAN = LSAME( TRANS, 'N' )
148 LQUERY = ( LWORK.EQ.-1 )
149 *
150 * NQ is the order of Q or P and NW is the minimum dimension of WORK
151 *
152 IF( LEFT ) THEN
153 NQ = M
154 NW = N
155 ELSE
156 NQ = N
157 NW = M
158 END IF
159 IF( M.EQ.0 .OR. N.EQ.0 ) THEN
160 NW = 0
161 END IF
162 IF( .NOT.APPLYQ .AND. .NOT.LSAME( VECT, 'P' ) ) THEN
163 INFO = -1
164 ELSE IF( .NOT.LEFT .AND. .NOT.LSAME( SIDE, 'R' ) ) THEN
165 INFO = -2
166 ELSE IF( .NOT.NOTRAN .AND. .NOT.LSAME( TRANS, 'C' ) ) THEN
167 INFO = -3
168 ELSE IF( M.LT.0 ) THEN
169 INFO = -4
170 ELSE IF( N.LT.0 ) THEN
171 INFO = -5
172 ELSE IF( K.LT.0 ) THEN
173 INFO = -6
174 ELSE IF( ( APPLYQ .AND. LDA.LT.MAX( 1, NQ ) ) .OR.
175 $ ( .NOT.APPLYQ .AND. LDA.LT.MAX( 1, MIN( NQ, K ) ) ) )
176 $ THEN
177 INFO = -8
178 ELSE IF( LDC.LT.MAX( 1, M ) ) THEN
179 INFO = -11
180 ELSE IF( LWORK.LT.MAX( 1, NW ) .AND. .NOT.LQUERY ) THEN
181 INFO = -13
182 END IF
183 *
184 IF( INFO.EQ.0 ) THEN
185 IF( NW.GT.0 ) THEN
186 IF( APPLYQ ) THEN
187 IF( LEFT ) THEN
188 NB = ILAENV( 1, 'ZUNMQR', SIDE // TRANS, M-1, N, M-1,
189 $ -1 )
190 ELSE
191 NB = ILAENV( 1, 'ZUNMQR', SIDE // TRANS, M, N-1, N-1,
192 $ -1 )
193 END IF
194 ELSE
195 IF( LEFT ) THEN
196 NB = ILAENV( 1, 'ZUNMLQ', SIDE // TRANS, M-1, N, M-1,
197 $ -1 )
198 ELSE
199 NB = ILAENV( 1, 'ZUNMLQ', SIDE // TRANS, M, N-1, N-1,
200 $ -1 )
201 END IF
202 END IF
203 LWKOPT = MAX( 1, NW*NB )
204 ELSE
205 LWKOPT = 1
206 END IF
207 WORK( 1 ) = LWKOPT
208 END IF
209 *
210 IF( INFO.NE.0 ) THEN
211 CALL XERBLA( 'ZUNMBR', -INFO )
212 RETURN
213 ELSE IF( LQUERY ) THEN
214 RETURN
215 END IF
216 *
217 * Quick return if possible
218 *
219 IF( M.EQ.0 .OR. N.EQ.0 )
220 $ RETURN
221 *
222 IF( APPLYQ ) THEN
223 *
224 * Apply Q
225 *
226 IF( NQ.GE.K ) THEN
227 *
228 * Q was determined by a call to ZGEBRD with nq >= k
229 *
230 CALL ZUNMQR( SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC,
231 $ WORK, LWORK, IINFO )
232 ELSE IF( NQ.GT.1 ) THEN
233 *
234 * Q was determined by a call to ZGEBRD with nq < k
235 *
236 IF( LEFT ) THEN
237 MI = M - 1
238 NI = N
239 I1 = 2
240 I2 = 1
241 ELSE
242 MI = M
243 NI = N - 1
244 I1 = 1
245 I2 = 2
246 END IF
247 CALL ZUNMQR( SIDE, TRANS, MI, NI, NQ-1, A( 2, 1 ), LDA, TAU,
248 $ C( I1, I2 ), LDC, WORK, LWORK, IINFO )
249 END IF
250 ELSE
251 *
252 * Apply P
253 *
254 IF( NOTRAN ) THEN
255 TRANST = 'C'
256 ELSE
257 TRANST = 'N'
258 END IF
259 IF( NQ.GT.K ) THEN
260 *
261 * P was determined by a call to ZGEBRD with nq > k
262 *
263 CALL ZUNMLQ( SIDE, TRANST, M, N, K, A, LDA, TAU, C, LDC,
264 $ WORK, LWORK, IINFO )
265 ELSE IF( NQ.GT.1 ) THEN
266 *
267 * P was determined by a call to ZGEBRD with nq <= k
268 *
269 IF( LEFT ) THEN
270 MI = M - 1
271 NI = N
272 I1 = 2
273 I2 = 1
274 ELSE
275 MI = M
276 NI = N - 1
277 I1 = 1
278 I2 = 2
279 END IF
280 CALL ZUNMLQ( SIDE, TRANST, MI, NI, NQ-1, A( 1, 2 ), LDA,
281 $ TAU, C( I1, I2 ), LDC, WORK, LWORK, IINFO )
282 END IF
283 END IF
284 WORK( 1 ) = LWKOPT
285 RETURN
286 *
287 * End of ZUNMBR
288 *
289 END