LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
cgeqpf.f
Go to the documentation of this file.
1*> \brief \b CGEQPF
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download CGEQPF + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/cgeqpf.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/cgeqpf.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/cgeqpf.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE CGEQPF( M, N, A, LDA, JPVT, TAU, WORK, RWORK, INFO )
22*
23* .. Scalar Arguments ..
24* INTEGER INFO, LDA, M, N
25* ..
26* .. Array Arguments ..
27* INTEGER JPVT( * )
28* REAL RWORK( * )
29* COMPLEX A( LDA, * ), TAU( * ), WORK( * )
30* ..
31*
32*
33*> \par Purpose:
34* =============
35*>
36*> \verbatim
37*>
38*> This routine is deprecated and has been replaced by routine CGEQP3.
39*>
40*> CGEQPF computes a QR factorization with column pivoting of a
41*> complex M-by-N matrix A: A*P = Q*R.
42*> \endverbatim
43*
44* Arguments:
45* ==========
46*
47*> \param[in] M
48*> \verbatim
49*> M is INTEGER
50*> The number of rows of the matrix A. M >= 0.
51*> \endverbatim
52*>
53*> \param[in] N
54*> \verbatim
55*> N is INTEGER
56*> The number of columns of the matrix A. N >= 0
57*> \endverbatim
58*>
59*> \param[in,out] A
60*> \verbatim
61*> A is COMPLEX array, dimension (LDA,N)
62*> On entry, the M-by-N matrix A.
63*> On exit, the upper triangle of the array contains the
64*> min(M,N)-by-N upper triangular matrix R; the elements
65*> below the diagonal, together with the array TAU,
66*> represent the unitary matrix Q as a product of
67*> min(m,n) elementary reflectors.
68*> \endverbatim
69*>
70*> \param[in] LDA
71*> \verbatim
72*> LDA is INTEGER
73*> The leading dimension of the array A. LDA >= max(1,M).
74*> \endverbatim
75*>
76*> \param[in,out] JPVT
77*> \verbatim
78*> JPVT is INTEGER array, dimension (N)
79*> On entry, if JPVT(i) .ne. 0, the i-th column of A is permuted
80*> to the front of A*P (a leading column); if JPVT(i) = 0,
81*> the i-th column of A is a free column.
82*> On exit, if JPVT(i) = k, then the i-th column of A*P
83*> was the k-th column of A.
84*> \endverbatim
85*>
86*> \param[out] TAU
87*> \verbatim
88*> TAU is COMPLEX array, dimension (min(M,N))
89*> The scalar factors of the elementary reflectors.
90*> \endverbatim
91*>
92*> \param[out] WORK
93*> \verbatim
94*> WORK is COMPLEX array, dimension (N)
95*> \endverbatim
96*>
97*> \param[out] RWORK
98*> \verbatim
99*> RWORK is REAL array, dimension (2*N)
100*> \endverbatim
101*>
102*> \param[out] INFO
103*> \verbatim
104*> INFO is INTEGER
105*> = 0: successful exit
106*> < 0: if INFO = -i, the i-th argument had an illegal value
107*> \endverbatim
108*
109* Authors:
110* ========
111*
112*> \author Univ. of Tennessee
113*> \author Univ. of California Berkeley
114*> \author Univ. of Colorado Denver
115*> \author NAG Ltd.
116*
117*> \ingroup complexGEcomputational
118*
119*> \par Further Details:
120* =====================
121*>
122*> \verbatim
123*>
124*> The matrix Q is represented as a product of elementary reflectors
125*>
126*> Q = H(1) H(2) . . . H(n)
127*>
128*> Each H(i) has the form
129*>
130*> H = I - tau * v * v**H
131*>
132*> where tau is a complex scalar, and v is a complex vector with
133*> v(1:i-1) = 0 and v(i) = 1; v(i+1:m) is stored on exit in A(i+1:m,i).
134*>
135*> The matrix P is represented in jpvt as follows: If
136*> jpvt(j) = i
137*> then the jth column of P is the ith canonical unit vector.
138*>
139*> Partial column norm updating strategy modified by
140*> Z. Drmac and Z. Bujanovic, Dept. of Mathematics,
141*> University of Zagreb, Croatia.
142*> -- April 2011 --
143*> For more details see LAPACK Working Note 176.
144*> \endverbatim
145*>
146* =====================================================================
147 SUBROUTINE cgeqpf( M, N, A, LDA, JPVT, TAU, WORK, RWORK, INFO )
148*
149* -- LAPACK computational routine --
150* -- LAPACK is a software package provided by Univ. of Tennessee, --
151* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
152*
153* .. Scalar Arguments ..
154 INTEGER INFO, LDA, M, N
155* ..
156* .. Array Arguments ..
157 INTEGER JPVT( * )
158 REAL RWORK( * )
159 COMPLEX A( LDA, * ), TAU( * ), WORK( * )
160* ..
161*
162* =====================================================================
163*
164* .. Parameters ..
165 REAL ZERO, ONE
166 parameter( zero = 0.0e+0, one = 1.0e+0 )
167* ..
168* .. Local Scalars ..
169 INTEGER I, ITEMP, J, MA, MN, PVT
170 REAL TEMP, TEMP2, TOL3Z
171 COMPLEX AII
172* ..
173* .. External Subroutines ..
174 EXTERNAL cgeqr2, clarf, clarfg, cswap, cunm2r, xerbla
175* ..
176* .. Intrinsic Functions ..
177 INTRINSIC abs, cmplx, conjg, max, min, sqrt
178* ..
179* .. External Functions ..
180 INTEGER ISAMAX
181 REAL SCNRM2, SLAMCH
182 EXTERNAL isamax, scnrm2, slamch
183* ..
184* .. Executable Statements ..
185*
186* Test the input arguments
187*
188 info = 0
189 IF( m.LT.0 ) THEN
190 info = -1
191 ELSE IF( n.LT.0 ) THEN
192 info = -2
193 ELSE IF( lda.LT.max( 1, m ) ) THEN
194 info = -4
195 END IF
196 IF( info.NE.0 ) THEN
197 CALL xerbla( 'CGEQPF', -info )
198 RETURN
199 END IF
200*
201 mn = min( m, n )
202 tol3z = sqrt(slamch('Epsilon'))
203*
204* Move initial columns up front
205*
206 itemp = 1
207 DO 10 i = 1, n
208 IF( jpvt( i ).NE.0 ) THEN
209 IF( i.NE.itemp ) THEN
210 CALL cswap( m, a( 1, i ), 1, a( 1, itemp ), 1 )
211 jpvt( i ) = jpvt( itemp )
212 jpvt( itemp ) = i
213 ELSE
214 jpvt( i ) = i
215 END IF
216 itemp = itemp + 1
217 ELSE
218 jpvt( i ) = i
219 END IF
220 10 CONTINUE
221 itemp = itemp - 1
222*
223* Compute the QR factorization and update remaining columns
224*
225 IF( itemp.GT.0 ) THEN
226 ma = min( itemp, m )
227 CALL cgeqr2( m, ma, a, lda, tau, work, info )
228 IF( ma.LT.n ) THEN
229 CALL cunm2r( 'Left', 'Conjugate transpose', m, n-ma, ma, a,
230 $ lda, tau, a( 1, ma+1 ), lda, work, info )
231 END IF
232 END IF
233*
234 IF( itemp.LT.mn ) THEN
235*
236* Initialize partial column norms. The first n elements of
237* work store the exact column norms.
238*
239 DO 20 i = itemp + 1, n
240 rwork( i ) = scnrm2( m-itemp, a( itemp+1, i ), 1 )
241 rwork( n+i ) = rwork( i )
242 20 CONTINUE
243*
244* Compute factorization
245*
246 DO 40 i = itemp + 1, mn
247*
248* Determine ith pivot column and swap if necessary
249*
250 pvt = ( i-1 ) + isamax( n-i+1, rwork( i ), 1 )
251*
252 IF( pvt.NE.i ) THEN
253 CALL cswap( m, a( 1, pvt ), 1, a( 1, i ), 1 )
254 itemp = jpvt( pvt )
255 jpvt( pvt ) = jpvt( i )
256 jpvt( i ) = itemp
257 rwork( pvt ) = rwork( i )
258 rwork( n+pvt ) = rwork( n+i )
259 END IF
260*
261* Generate elementary reflector H(i)
262*
263 aii = a( i, i )
264 CALL clarfg( m-i+1, aii, a( min( i+1, m ), i ), 1,
265 $ tau( i ) )
266 a( i, i ) = aii
267*
268 IF( i.LT.n ) THEN
269*
270* Apply H(i) to A(i:m,i+1:n) from the left
271*
272 aii = a( i, i )
273 a( i, i ) = cmplx( one )
274 CALL clarf( 'Left', m-i+1, n-i, a( i, i ), 1,
275 $ conjg( tau( i ) ), a( i, i+1 ), lda, work )
276 a( i, i ) = aii
277 END IF
278*
279* Update partial column norms
280*
281 DO 30 j = i + 1, n
282 IF( rwork( j ).NE.zero ) THEN
283*
284* NOTE: The following 4 lines follow from the analysis in
285* Lapack Working Note 176.
286*
287 temp = abs( a( i, j ) ) / rwork( j )
288 temp = max( zero, ( one+temp )*( one-temp ) )
289 temp2 = temp*( rwork( j ) / rwork( n+j ) )**2
290 IF( temp2 .LE. tol3z ) THEN
291 IF( m-i.GT.0 ) THEN
292 rwork( j ) = scnrm2( m-i, a( i+1, j ), 1 )
293 rwork( n+j ) = rwork( j )
294 ELSE
295 rwork( j ) = zero
296 rwork( n+j ) = zero
297 END IF
298 ELSE
299 rwork( j ) = rwork( j )*sqrt( temp )
300 END IF
301 END IF
302 30 CONTINUE
303*
304 40 CONTINUE
305 END IF
306 RETURN
307*
308* End of CGEQPF
309*
310 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine cgeqpf(m, n, a, lda, jpvt, tau, work, rwork, info)
CGEQPF
Definition cgeqpf.f:148
subroutine cgeqr2(m, n, a, lda, tau, work, info)
CGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm.
Definition cgeqr2.f:130
subroutine clarf(side, m, n, v, incv, tau, c, ldc, work)
CLARF applies an elementary reflector to a general rectangular matrix.
Definition clarf.f:128
subroutine clarfg(n, alpha, x, incx, tau)
CLARFG generates an elementary reflector (Householder matrix).
Definition clarfg.f:106
subroutine cswap(n, cx, incx, cy, incy)
CSWAP
Definition cswap.f:81
subroutine cunm2r(side, trans, m, n, k, a, lda, tau, c, ldc, work, info)
CUNM2R multiplies a general matrix by the unitary matrix from a QR factorization determined by cgeqrf...
Definition cunm2r.f:159