LAPACK 3.12.1
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches

◆ dlaqp2()

subroutine dlaqp2 ( integer m,
integer n,
integer offset,
double precision, dimension( lda, * ) a,
integer lda,
integer, dimension( * ) jpvt,
double precision, dimension( * ) tau,
double precision, dimension( * ) vn1,
double precision, dimension( * ) vn2,
double precision, dimension( * ) work )

DLAQP2 computes a QR factorization with column pivoting of the matrix block.

Download DLAQP2 + dependencies [TGZ] [ZIP] [TXT]

Purpose:
!>
!> DLAQP2 computes a QR factorization with column pivoting of
!> the block A(OFFSET+1:M,1:N).
!> The block A(1:OFFSET,1:N) is accordingly pivoted, but not factorized.
!> 
Parameters
[in]M
!>          M is INTEGER
!>          The number of rows of the matrix A. M >= 0.
!> 
[in]N
!>          N is INTEGER
!>          The number of columns of the matrix A. N >= 0.
!> 
[in]OFFSET
!>          OFFSET is INTEGER
!>          The number of rows of the matrix A that must be pivoted
!>          but no factorized. OFFSET >= 0.
!> 
[in,out]A
!>          A is DOUBLE PRECISION array, dimension (LDA,N)
!>          On entry, the M-by-N matrix A.
!>          On exit, the upper triangle of block A(OFFSET+1:M,1:N) is
!>          the triangular factor obtained; the elements in block
!>          A(OFFSET+1:M,1:N) below the diagonal, together with the
!>          array TAU, represent the orthogonal matrix Q as a product of
!>          elementary reflectors. Block A(1:OFFSET,1:N) has been
!>          accordingly pivoted, but no factorized.
!> 
[in]LDA
!>          LDA is INTEGER
!>          The leading dimension of the array A. LDA >= max(1,M).
!> 
[in,out]JPVT
!>          JPVT is INTEGER array, dimension (N)
!>          On entry, if JPVT(i) .ne. 0, the i-th column of A is permuted
!>          to the front of A*P (a leading column); if JPVT(i) = 0,
!>          the i-th column of A is a free column.
!>          On exit, if JPVT(i) = k, then the i-th column of A*P
!>          was the k-th column of A.
!> 
[out]TAU
!>          TAU is DOUBLE PRECISION array, dimension (min(M,N))
!>          The scalar factors of the elementary reflectors.
!> 
[in,out]VN1
!>          VN1 is DOUBLE PRECISION array, dimension (N)
!>          The vector with the partial column norms.
!> 
[in,out]VN2
!>          VN2 is DOUBLE PRECISION array, dimension (N)
!>          The vector with the exact column norms.
!> 
[out]WORK
!>          WORK is DOUBLE PRECISION array, dimension (N)
!> 
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Contributors:
G. Quintana-Orti, Depto. de Informatica, Universidad Jaime I, Spain X. Sun, Computer Science Dept., Duke University, USA
Partial column norm updating strategy modified on April 2011 Z. Drmac and Z. Bujanovic, Dept. of Mathematics, University of Zagreb, Croatia.
References:
LAPACK Working Note 176 [PDF]

Definition at line 143 of file dlaqp2.f.

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 A( LDA, * ), TAU( * ), VN1( * ), VN2( * ),
156 $ WORK( * )
157* ..
158*
159* =====================================================================
160*
161* .. Parameters ..
162 DOUBLE PRECISION ZERO, ONE
163 parameter( zero = 0.0d+0, one = 1.0d+0 )
164* ..
165* .. Local Scalars ..
166 INTEGER I, ITEMP, J, MN, OFFPI, PVT
167 DOUBLE PRECISION TEMP, TEMP2, TOL3Z
168* ..
169* .. External Subroutines ..
170 EXTERNAL dlarf, dlarfg, dswap
171* ..
172* .. Intrinsic Functions ..
173 INTRINSIC abs, max, min, sqrt
174* ..
175* .. External Functions ..
176 INTEGER IDAMAX
177 DOUBLE PRECISION DLAMCH, DNRM2
178 EXTERNAL idamax, dlamch, dnrm2
179* ..
180* .. Executable Statements ..
181*
182 mn = min( m-offset, n )
183 tol3z = sqrt(dlamch('Epsilon'))
184*
185* Compute factorization.
186*
187 DO 20 i = 1, mn
188*
189 offpi = offset + i
190*
191* Determine ith pivot column and swap if necessary.
192*
193 pvt = ( i-1 ) + idamax( n-i+1, vn1( i ), 1 )
194*
195 IF( pvt.NE.i ) THEN
196 CALL dswap( m, a( 1, pvt ), 1, a( 1, i ), 1 )
197 itemp = jpvt( pvt )
198 jpvt( pvt ) = jpvt( i )
199 jpvt( i ) = itemp
200 vn1( pvt ) = vn1( i )
201 vn2( pvt ) = vn2( i )
202 END IF
203*
204* Generate elementary reflector H(i).
205*
206 IF( offpi.LT.m ) THEN
207 CALL dlarfg( m-offpi+1, a( offpi, i ), a( offpi+1, i ),
208 $ 1,
209 $ tau( i ) )
210 ELSE
211 CALL dlarfg( 1, a( m, i ), a( m, i ), 1, tau( i ) )
212 END IF
213*
214 IF( i.LT.n ) THEN
215*
216* Apply H(i)**T to A(offset+i:m,i+1:n) from the left.
217*
218 CALL dlarf1f( 'Left', m-offpi+1, n-i, a( offpi, i ), 1,
219 $ tau( i ), a( offpi, i+1 ), lda, work( 1 ) )
220 END IF
221*
222* Update partial column norms.
223*
224 DO 10 j = i + 1, n
225 IF( vn1( j ).NE.zero ) THEN
226*
227* NOTE: The following 4 lines follow from the analysis in
228* Lapack Working Note 176.
229*
230 temp = one - ( abs( a( offpi, j ) ) / vn1( j ) )**2
231 temp = max( temp, zero )
232 temp2 = temp*( vn1( j ) / vn2( j ) )**2
233 IF( temp2 .LE. tol3z ) THEN
234 IF( offpi.LT.m ) THEN
235 vn1( j ) = dnrm2( m-offpi, a( offpi+1, j ), 1 )
236 vn2( j ) = vn1( j )
237 ELSE
238 vn1( j ) = zero
239 vn2( j ) = zero
240 END IF
241 ELSE
242 vn1( j ) = vn1( j )*sqrt( temp )
243 END IF
244 END IF
245 10 CONTINUE
246*
247 20 CONTINUE
248*
249 RETURN
250*
251* End of DLAQP2
252*
integer function idamax(n, dx, incx)
IDAMAX
Definition idamax.f:71
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
subroutine dlarf1f(side, m, n, v, incv, tau, c, ldc, work)
DLARF1F applies an elementary reflector to a general rectangular
Definition dlarf1f.f:157
subroutine dlarf(side, m, n, v, incv, tau, c, ldc, work)
DLARF applies an elementary reflector to a general rectangular matrix.
Definition dlarf.f:122
subroutine dlarfg(n, alpha, x, incx, tau)
DLARFG generates an elementary reflector (Householder matrix).
Definition dlarfg.f:104
real(wp) function dnrm2(n, x, incx)
DNRM2
Definition dnrm2.f90:89
subroutine dswap(n, dx, incx, dy, incy)
DSWAP
Definition dswap.f:82
Here is the call graph for this function:
Here is the caller graph for this function: