LAPACK  3.4.2
LAPACK: Linear Algebra PACKage
 All Files Functions Groups
cggqrf.f
Go to the documentation of this file.
1 *> \brief \b CGGQRF
2 *
3 * =========== DOCUMENTATION ===========
4 *
5 * Online html documentation available at
6 * http://www.netlib.org/lapack/explore-html/
7 *
8 *> \htmlonly
9 *> Download CGGQRF + dependencies
10 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cggqrf.f">
11 *> [TGZ]</a>
12 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cggqrf.f">
13 *> [ZIP]</a>
14 *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cggqrf.f">
15 *> [TXT]</a>
16 *> \endhtmlonly
17 *
18 * Definition:
19 * ===========
20 *
21 * SUBROUTINE CGGQRF( N, M, P, A, LDA, TAUA, B, LDB, TAUB, WORK,
22 * LWORK, INFO )
23 *
24 * .. Scalar Arguments ..
25 * INTEGER INFO, LDA, LDB, LWORK, M, N, P
26 * ..
27 * .. Array Arguments ..
28 * COMPLEX A( LDA, * ), B( LDB, * ), TAUA( * ), TAUB( * ),
29 * $ WORK( * )
30 * ..
31 *
32 *
33 *> \par Purpose:
34 * =============
35 *>
36 *> \verbatim
37 *>
38 *> CGGQRF computes a generalized QR factorization of an N-by-M matrix A
39 *> and an N-by-P matrix B:
40 *>
41 *> A = Q*R, B = Q*T*Z,
42 *>
43 *> where Q is an N-by-N unitary matrix, Z is a P-by-P unitary matrix,
44 *> and R and T assume one of the forms:
45 *>
46 *> if N >= M, R = ( R11 ) M , or if N < M, R = ( R11 R12 ) N,
47 *> ( 0 ) N-M N M-N
48 *> M
49 *>
50 *> where R11 is upper triangular, and
51 *>
52 *> if N <= P, T = ( 0 T12 ) N, or if N > P, T = ( T11 ) N-P,
53 *> P-N N ( T21 ) P
54 *> P
55 *>
56 *> where T12 or T21 is upper triangular.
57 *>
58 *> In particular, if B is square and nonsingular, the GQR factorization
59 *> of A and B implicitly gives the QR factorization of inv(B)*A:
60 *>
61 *> inv(B)*A = Z**H * (inv(T)*R)
62 *>
63 *> where inv(B) denotes the inverse of the matrix B, and Z' denotes the
64 *> conjugate transpose of matrix Z.
65 *> \endverbatim
66 *
67 * Arguments:
68 * ==========
69 *
70 *> \param[in] N
71 *> \verbatim
72 *> N is INTEGER
73 *> The number of rows of the matrices A and B. N >= 0.
74 *> \endverbatim
75 *>
76 *> \param[in] M
77 *> \verbatim
78 *> M is INTEGER
79 *> The number of columns of the matrix A. M >= 0.
80 *> \endverbatim
81 *>
82 *> \param[in] P
83 *> \verbatim
84 *> P is INTEGER
85 *> The number of columns of the matrix B. P >= 0.
86 *> \endverbatim
87 *>
88 *> \param[in,out] A
89 *> \verbatim
90 *> A is COMPLEX array, dimension (LDA,M)
91 *> On entry, the N-by-M matrix A.
92 *> On exit, the elements on and above the diagonal of the array
93 *> contain the min(N,M)-by-M upper trapezoidal matrix R (R is
94 *> upper triangular if N >= M); the elements below the diagonal,
95 *> with the array TAUA, represent the unitary matrix Q as a
96 *> product of min(N,M) elementary reflectors (see Further
97 *> Details).
98 *> \endverbatim
99 *>
100 *> \param[in] LDA
101 *> \verbatim
102 *> LDA is INTEGER
103 *> The leading dimension of the array A. LDA >= max(1,N).
104 *> \endverbatim
105 *>
106 *> \param[out] TAUA
107 *> \verbatim
108 *> TAUA is COMPLEX array, dimension (min(N,M))
109 *> The scalar factors of the elementary reflectors which
110 *> represent the unitary matrix Q (see Further Details).
111 *> \endverbatim
112 *>
113 *> \param[in,out] B
114 *> \verbatim
115 *> B is COMPLEX array, dimension (LDB,P)
116 *> On entry, the N-by-P matrix B.
117 *> On exit, if N <= P, the upper triangle of the subarray
118 *> B(1:N,P-N+1:P) contains the N-by-N upper triangular matrix T;
119 *> if N > P, the elements on and above the (N-P)-th subdiagonal
120 *> contain the N-by-P upper trapezoidal matrix T; the remaining
121 *> elements, with the array TAUB, represent the unitary
122 *> matrix Z as a product of elementary reflectors (see Further
123 *> Details).
124 *> \endverbatim
125 *>
126 *> \param[in] LDB
127 *> \verbatim
128 *> LDB is INTEGER
129 *> The leading dimension of the array B. LDB >= max(1,N).
130 *> \endverbatim
131 *>
132 *> \param[out] TAUB
133 *> \verbatim
134 *> TAUB is COMPLEX array, dimension (min(N,P))
135 *> The scalar factors of the elementary reflectors which
136 *> represent the unitary matrix Z (see Further Details).
137 *> \endverbatim
138 *>
139 *> \param[out] WORK
140 *> \verbatim
141 *> WORK is COMPLEX array, dimension (MAX(1,LWORK))
142 *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
143 *> \endverbatim
144 *>
145 *> \param[in] LWORK
146 *> \verbatim
147 *> LWORK is INTEGER
148 *> The dimension of the array WORK. LWORK >= max(1,N,M,P).
149 *> For optimum performance LWORK >= max(N,M,P)*max(NB1,NB2,NB3),
150 *> where NB1 is the optimal blocksize for the QR factorization
151 *> of an N-by-M matrix, NB2 is the optimal blocksize for the
152 *> RQ factorization of an N-by-P matrix, and NB3 is the optimal
153 *> blocksize for a call of CUNMQR.
154 *>
155 *> If LWORK = -1, then a workspace query is assumed; the routine
156 *> only calculates the optimal size of the WORK array, returns
157 *> this value as the first entry of the WORK array, and no error
158 *> message related to LWORK is issued by XERBLA.
159 *> \endverbatim
160 *>
161 *> \param[out] INFO
162 *> \verbatim
163 *> INFO is INTEGER
164 *> = 0: successful exit
165 *> < 0: if INFO = -i, the i-th argument had an illegal value.
166 *> \endverbatim
167 *
168 * Authors:
169 * ========
170 *
171 *> \author Univ. of Tennessee
172 *> \author Univ. of California Berkeley
173 *> \author Univ. of Colorado Denver
174 *> \author NAG Ltd.
175 *
176 *> \date November 2011
177 *
178 *> \ingroup complexOTHERcomputational
179 *
180 *> \par Further Details:
181 * =====================
182 *>
183 *> \verbatim
184 *>
185 *> The matrix Q is represented as a product of elementary reflectors
186 *>
187 *> Q = H(1) H(2) . . . H(k), where k = min(n,m).
188 *>
189 *> Each H(i) has the form
190 *>
191 *> H(i) = I - taua * v * v**H
192 *>
193 *> where taua is a complex scalar, and v is a complex vector with
194 *> v(1:i-1) = 0 and v(i) = 1; v(i+1:n) is stored on exit in A(i+1:n,i),
195 *> and taua in TAUA(i).
196 *> To form Q explicitly, use LAPACK subroutine CUNGQR.
197 *> To use Q to update another matrix, use LAPACK subroutine CUNMQR.
198 *>
199 *> The matrix Z is represented as a product of elementary reflectors
200 *>
201 *> Z = H(1) H(2) . . . H(k), where k = min(n,p).
202 *>
203 *> Each H(i) has the form
204 *>
205 *> H(i) = I - taub * v * v**H
206 *>
207 *> where taub is a complex scalar, and v is a complex vector with
208 *> v(p-k+i+1:p) = 0 and v(p-k+i) = 1; v(1:p-k+i-1) is stored on exit in
209 *> B(n-k+i,1:p-k+i-1), and taub in TAUB(i).
210 *> To form Z explicitly, use LAPACK subroutine CUNGRQ.
211 *> To use Z to update another matrix, use LAPACK subroutine CUNMRQ.
212 *> \endverbatim
213 *>
214 * =====================================================================
215  SUBROUTINE cggqrf( N, M, P, A, LDA, TAUA, B, LDB, TAUB, WORK,
216  $ lwork, info )
217 *
218 * -- LAPACK computational routine (version 3.4.0) --
219 * -- LAPACK is a software package provided by Univ. of Tennessee, --
220 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
221 * November 2011
222 *
223 * .. Scalar Arguments ..
224  INTEGER info, lda, ldb, lwork, m, n, p
225 * ..
226 * .. Array Arguments ..
227  COMPLEX a( lda, * ), b( ldb, * ), taua( * ), taub( * ),
228  $ work( * )
229 * ..
230 *
231 * =====================================================================
232 *
233 * .. Local Scalars ..
234  LOGICAL lquery
235  INTEGER lopt, lwkopt, nb, nb1, nb2, nb3
236 * ..
237 * .. External Subroutines ..
238  EXTERNAL cgeqrf, cgerqf, cunmqr, xerbla
239 * ..
240 * .. External Functions ..
241  INTEGER ilaenv
242  EXTERNAL ilaenv
243 * ..
244 * .. Intrinsic Functions ..
245  INTRINSIC int, max, min
246 * ..
247 * .. Executable Statements ..
248 *
249 * Test the input parameters
250 *
251  info = 0
252  nb1 = ilaenv( 1, 'CGEQRF', ' ', n, m, -1, -1 )
253  nb2 = ilaenv( 1, 'CGERQF', ' ', n, p, -1, -1 )
254  nb3 = ilaenv( 1, 'CUNMQR', ' ', n, m, p, -1 )
255  nb = max( nb1, nb2, nb3 )
256  lwkopt = max( n, m, p)*nb
257  work( 1 ) = lwkopt
258  lquery = ( lwork.EQ.-1 )
259  IF( n.LT.0 ) THEN
260  info = -1
261  ELSE IF( m.LT.0 ) THEN
262  info = -2
263  ELSE IF( p.LT.0 ) THEN
264  info = -3
265  ELSE IF( lda.LT.max( 1, n ) ) THEN
266  info = -5
267  ELSE IF( ldb.LT.max( 1, n ) ) THEN
268  info = -8
269  ELSE IF( lwork.LT.max( 1, n, m, p ) .AND. .NOT.lquery ) THEN
270  info = -11
271  END IF
272  IF( info.NE.0 ) THEN
273  CALL xerbla( 'CGGQRF', -info )
274  return
275  ELSE IF( lquery ) THEN
276  return
277  END IF
278 *
279 * QR factorization of N-by-M matrix A: A = Q*R
280 *
281  CALL cgeqrf( n, m, a, lda, taua, work, lwork, info )
282  lopt = work( 1 )
283 *
284 * Update B := Q**H*B.
285 *
286  CALL cunmqr( 'Left', 'Conjugate Transpose', n, p, min( n, m ), a,
287  $ lda, taua, b, ldb, work, lwork, info )
288  lopt = max( lopt, int( work( 1 ) ) )
289 *
290 * RQ factorization of N-by-P matrix B: B = T*Z.
291 *
292  CALL cgerqf( n, p, b, ldb, taub, work, lwork, info )
293  work( 1 ) = max( lopt, int( work( 1 ) ) )
294 *
295  return
296 *
297 * End of CGGQRF
298 *
299  END