LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
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*> \ingroup ggqrf
177*
178*> \par Further Details:
179* =====================
180*>
181*> \verbatim
182*>
183*> The matrix Q is represented as a product of elementary reflectors
184*>
185*> Q = H(1) H(2) . . . H(k), where k = min(n,m).
186*>
187*> Each H(i) has the form
188*>
189*> H(i) = I - taua * v * v**H
190*>
191*> where taua is a complex scalar, and v is a complex vector with
192*> v(1:i-1) = 0 and v(i) = 1; v(i+1:n) is stored on exit in A(i+1:n,i),
193*> and taua in TAUA(i).
194*> To form Q explicitly, use LAPACK subroutine CUNGQR.
195*> To use Q to update another matrix, use LAPACK subroutine CUNMQR.
196*>
197*> The matrix Z is represented as a product of elementary reflectors
198*>
199*> Z = H(1) H(2) . . . H(k), where k = min(n,p).
200*>
201*> Each H(i) has the form
202*>
203*> H(i) = I - taub * v * v**H
204*>
205*> where taub is a complex scalar, and v is a complex vector with
206*> v(p-k+i+1:p) = 0 and v(p-k+i) = 1; v(1:p-k+i-1) is stored on exit in
207*> B(n-k+i,1:p-k+i-1), and taub in TAUB(i).
208*> To form Z explicitly, use LAPACK subroutine CUNGRQ.
209*> To use Z to update another matrix, use LAPACK subroutine CUNMRQ.
210*> \endverbatim
211*>
212* =====================================================================
213 SUBROUTINE cggqrf( N, M, P, A, LDA, TAUA, B, LDB, TAUB, WORK,
214 $ LWORK, INFO )
215*
216* -- LAPACK computational routine --
217* -- LAPACK is a software package provided by Univ. of Tennessee, --
218* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
219*
220* .. Scalar Arguments ..
221 INTEGER INFO, LDA, LDB, LWORK, M, N, P
222* ..
223* .. Array Arguments ..
224 COMPLEX A( LDA, * ), B( LDB, * ), TAUA( * ), TAUB( * ),
225 $ work( * )
226* ..
227*
228* =====================================================================
229*
230* .. Local Scalars ..
231 LOGICAL LQUERY
232 INTEGER LOPT, LWKOPT, NB, NB1, NB2, NB3
233* ..
234* .. External Subroutines ..
235 EXTERNAL cgeqrf, cgerqf, cunmqr, xerbla
236* ..
237* .. External Functions ..
238 INTEGER ILAENV
239 REAL SROUNDUP_LWORK
240 EXTERNAL ilaenv, sroundup_lwork
241* ..
242* .. Intrinsic Functions ..
243 INTRINSIC int, max, min
244* ..
245* .. Executable Statements ..
246*
247* Test the input parameters
248*
249 info = 0
250 nb1 = ilaenv( 1, 'CGEQRF', ' ', n, m, -1, -1 )
251 nb2 = ilaenv( 1, 'CGERQF', ' ', n, p, -1, -1 )
252 nb3 = ilaenv( 1, 'CUNMQR', ' ', n, m, p, -1 )
253 nb = max( nb1, nb2, nb3 )
254 lwkopt = max( n, m, p)*nb
255 work( 1 ) = sroundup_lwork(lwkopt)
256 lquery = ( lwork.EQ.-1 )
257 IF( n.LT.0 ) THEN
258 info = -1
259 ELSE IF( m.LT.0 ) THEN
260 info = -2
261 ELSE IF( p.LT.0 ) THEN
262 info = -3
263 ELSE IF( lda.LT.max( 1, n ) ) THEN
264 info = -5
265 ELSE IF( ldb.LT.max( 1, n ) ) THEN
266 info = -8
267 ELSE IF( lwork.LT.max( 1, n, m, p ) .AND. .NOT.lquery ) THEN
268 info = -11
269 END IF
270 IF( info.NE.0 ) THEN
271 CALL xerbla( 'CGGQRF', -info )
272 RETURN
273 ELSE IF( lquery ) THEN
274 RETURN
275 END IF
276*
277* QR factorization of N-by-M matrix A: A = Q*R
278*
279 CALL cgeqrf( n, m, a, lda, taua, work, lwork, info )
280 lopt = int( work( 1 ) )
281*
282* Update B := Q**H*B.
283*
284 CALL cunmqr( 'Left', 'Conjugate Transpose', n, p, min( n, m ), a,
285 $ lda, taua, b, ldb, work, lwork, info )
286 lopt = max( lopt, int( work( 1 ) ) )
287*
288* RQ factorization of N-by-P matrix B: B = T*Z.
289*
290 CALL cgerqf( n, p, b, ldb, taub, work, lwork, info )
291 work( 1 ) = max( lopt, int( work( 1 ) ) )
292*
293 RETURN
294*
295* End of CGGQRF
296*
297 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine cgeqrf(m, n, a, lda, tau, work, lwork, info)
CGEQRF
Definition cgeqrf.f:146
subroutine cgerqf(m, n, a, lda, tau, work, lwork, info)
CGERQF
Definition cgerqf.f:139
subroutine cggqrf(n, m, p, a, lda, taua, b, ldb, taub, work, lwork, info)
CGGQRF
Definition cggqrf.f:215
subroutine cunmqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
CUNMQR
Definition cunmqr.f:168