LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
zlaqp2.f
Go to the documentation of this file.
1*> \brief \b ZLAQP2 computes a QR factorization with column pivoting of the matrix block.
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> Download ZLAQP2 + dependencies
9*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/zlaqp2.f">
10*> [TGZ]</a>
11*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/zlaqp2.f">
12*> [ZIP]</a>
13*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/zlaqp2.f">
14*> [TXT]</a>
15*
16* Definition:
17* ===========
18*
19* SUBROUTINE ZLAQP2( M, N, OFFSET, A, LDA, JPVT, TAU, VN1, VN2,
20* WORK )
21*
22* .. Scalar Arguments ..
23* INTEGER LDA, M, N, OFFSET
24* ..
25* .. Array Arguments ..
26* INTEGER JPVT( * )
27* DOUBLE PRECISION VN1( * ), VN2( * )
28* COMPLEX*16 A( LDA, * ), TAU( * ), WORK( * )
29* ..
30*
31*
32*> \par Purpose:
33* =============
34*>
35*> \verbatim
36*>
37*> ZLAQP2 computes a QR factorization with column pivoting of
38*> the block A(OFFSET+1:M,1:N).
39*> The block A(1:OFFSET,1:N) is accordingly pivoted, but not factorized.
40*> \endverbatim
41*
42* Arguments:
43* ==========
44*
45*> \param[in] M
46*> \verbatim
47*> M is INTEGER
48*> The number of rows of the matrix A. M >= 0.
49*> \endverbatim
50*>
51*> \param[in] N
52*> \verbatim
53*> N is INTEGER
54*> The number of columns of the matrix A. N >= 0.
55*> \endverbatim
56*>
57*> \param[in] OFFSET
58*> \verbatim
59*> OFFSET is INTEGER
60*> The number of rows of the matrix A that must be pivoted
61*> but no factorized. OFFSET >= 0.
62*> \endverbatim
63*>
64*> \param[in,out] A
65*> \verbatim
66*> A is COMPLEX*16 array, dimension (LDA,N)
67*> On entry, the M-by-N matrix A.
68*> On exit, the upper triangle of block A(OFFSET+1:M,1:N) is
69*> the triangular factor obtained; the elements in block
70*> A(OFFSET+1:M,1:N) below the diagonal, together with the
71*> array TAU, represent the orthogonal matrix Q as a product of
72*> elementary reflectors. Block A(1:OFFSET,1:N) has been
73*> accordingly pivoted, but no factorized.
74*> \endverbatim
75*>
76*> \param[in] LDA
77*> \verbatim
78*> LDA is INTEGER
79*> The leading dimension of the array A. LDA >= max(1,M).
80*> \endverbatim
81*>
82*> \param[in,out] JPVT
83*> \verbatim
84*> JPVT is INTEGER array, dimension (N)
85*> On entry, if JPVT(i) .ne. 0, the i-th column of A is permuted
86*> to the front of A*P (a leading column); if JPVT(i) = 0,
87*> the i-th column of A is a free column.
88*> On exit, if JPVT(i) = k, then the i-th column of A*P
89*> was the k-th column of A.
90*> \endverbatim
91*>
92*> \param[out] TAU
93*> \verbatim
94*> TAU is COMPLEX*16 array, dimension (min(M,N))
95*> The scalar factors of the elementary reflectors.
96*> \endverbatim
97*>
98*> \param[in,out] VN1
99*> \verbatim
100*> VN1 is DOUBLE PRECISION array, dimension (N)
101*> The vector with the partial column norms.
102*> \endverbatim
103*>
104*> \param[in,out] VN2
105*> \verbatim
106*> VN2 is DOUBLE PRECISION array, dimension (N)
107*> The vector with the exact column norms.
108*> \endverbatim
109*>
110*> \param[out] WORK
111*> \verbatim
112*> WORK is COMPLEX*16 array, dimension (N)
113*> \endverbatim
114*
115* Authors:
116* ========
117*
118*> \author Univ. of Tennessee
119*> \author Univ. of California Berkeley
120*> \author Univ. of Colorado Denver
121*> \author NAG Ltd.
122*
123*> \ingroup laqp2
124*
125*> \par Contributors:
126* ==================
127*>
128*> G. Quintana-Orti, Depto. de Informatica, Universidad Jaime I, Spain
129*> X. Sun, Computer Science Dept., Duke University, USA
130*> \n
131*> Partial column norm updating strategy modified on April 2011
132*> Z. Drmac and Z. Bujanovic, Dept. of Mathematics,
133*> University of Zagreb, Croatia.
134*
135*> \par References:
136* ================
137*>
138*> LAPACK Working Note 176
139*
140*> <a href="http://www.netlib.org/lapack/lawnspdf/lawn176.pdf">[PDF]</a>
141*
142* =====================================================================
143 SUBROUTINE zlaqp2( M, N, OFFSET, A, LDA, JPVT, TAU, VN1, VN2,
144 $ WORK )
145*
146* -- LAPACK auxiliary routine --
147* -- LAPACK is a software package provided by Univ. of Tennessee, --
148* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
149*
150* .. Scalar Arguments ..
151 INTEGER LDA, M, N, OFFSET
152* ..
153* .. Array Arguments ..
154 INTEGER JPVT( * )
155 DOUBLE PRECISION VN1( * ), VN2( * )
156 COMPLEX*16 A( LDA, * ), TAU( * ), WORK( * )
157* ..
158*
159* =====================================================================
160*
161* .. Parameters ..
162 DOUBLE PRECISION ZERO, ONE
163 COMPLEX*16 CONE
164 parameter( zero = 0.0d+0, one = 1.0d+0,
165 $ cone = ( 1.0d+0, 0.0d+0 ) )
166* ..
167* .. Local Scalars ..
168 INTEGER I, ITEMP, J, MN, OFFPI, PVT
169 DOUBLE PRECISION TEMP, TEMP2, TOL3Z
170* ..
171* .. External Subroutines ..
172 EXTERNAL zlarf1f, zlarfg, zswap
173* ..
174* .. Intrinsic Functions ..
175 INTRINSIC abs, dconjg, max, min, sqrt
176* ..
177* .. External Functions ..
178 INTEGER IDAMAX
179 DOUBLE PRECISION DLAMCH, DZNRM2
180 EXTERNAL idamax, dlamch, dznrm2
181* ..
182* .. Executable Statements ..
183*
184 mn = min( m-offset, n )
185 tol3z = sqrt(dlamch('Epsilon'))
186*
187* Compute factorization.
188*
189 DO 20 i = 1, mn
190*
191 offpi = offset + i
192*
193* Determine ith pivot column and swap if necessary.
194*
195 pvt = ( i-1 ) + idamax( n-i+1, vn1( i ), 1 )
196*
197 IF( pvt.NE.i ) THEN
198 CALL zswap( m, a( 1, pvt ), 1, a( 1, i ), 1 )
199 itemp = jpvt( pvt )
200 jpvt( pvt ) = jpvt( i )
201 jpvt( i ) = itemp
202 vn1( pvt ) = vn1( i )
203 vn2( pvt ) = vn2( i )
204 END IF
205*
206* Generate elementary reflector H(i).
207*
208 IF( offpi.LT.m ) THEN
209 CALL zlarfg( m-offpi+1, a( offpi, i ), a( offpi+1, i ),
210 $ 1,
211 $ tau( i ) )
212 ELSE
213 CALL zlarfg( 1, a( m, i ), a( m, i ), 1, tau( i ) )
214 END IF
215*
216 IF( i.LT.n ) THEN
217*
218* Apply H(i)**H to A(offset+i:m,i+1:n) from the left.
219*
220 CALL zlarf1f( 'Left', m-offpi+1, n-i, a( offpi, i ), 1,
221 $ conjg( tau( i ) ), a( offpi, i+1 ), lda,
222 $ work( 1 ) )
223 END IF
224*
225* Update partial column norms.
226*
227 DO 10 j = i + 1, n
228 IF( vn1( j ).NE.zero ) THEN
229*
230* NOTE: The following 4 lines follow from the analysis in
231* Lapack Working Note 176.
232*
233 temp = one - ( abs( a( offpi, j ) ) / vn1( j ) )**2
234 temp = max( temp, zero )
235 temp2 = temp*( vn1( j ) / vn2( j ) )**2
236 IF( temp2 .LE. tol3z ) THEN
237 IF( offpi.LT.m ) THEN
238 vn1( j ) = dznrm2( m-offpi, a( offpi+1, j ), 1 )
239 vn2( j ) = vn1( j )
240 ELSE
241 vn1( j ) = zero
242 vn2( j ) = zero
243 END IF
244 ELSE
245 vn1( j ) = vn1( j )*sqrt( temp )
246 END IF
247 END IF
248 10 CONTINUE
249*
250 20 CONTINUE
251*
252 RETURN
253*
254* End of ZLAQP2
255*
256 END
subroutine zlaqp2(m, n, offset, a, lda, jpvt, tau, vn1, vn2, work)
ZLAQP2 computes a QR factorization with column pivoting of the matrix block.
Definition zlaqp2.f:145
subroutine zlarf1f(side, m, n, v, incv, tau, c, ldc, work)
ZLARF1F applies an elementary reflector to a general rectangular
Definition zlarf1f.f:157
subroutine zlarfg(n, alpha, x, incx, tau)
ZLARFG generates an elementary reflector (Householder matrix).
Definition zlarfg.f:104
subroutine zswap(n, zx, incx, zy, incy)
ZSWAP
Definition zswap.f:81