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

◆ dgeqpf()

subroutine dgeqpf ( integer m,
integer n,
double precision, dimension( lda, * ) a,
integer lda,
integer, dimension( * ) jpvt,
double precision, dimension( * ) tau,
double precision, dimension( * ) work,
integer info )

DGEQPF

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

Purpose:
!>
!> This routine is deprecated and has been replaced by routine DGEQP3.
!>
!> DGEQPF computes a QR factorization with column pivoting of a
!> real M-by-N matrix A: A*P = Q*R.
!>
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,out]A
!>          A is DOUBLE PRECISION array, dimension (LDA,N)
!>          On entry, the M-by-N matrix A.
!>          On exit, the upper triangle of the array contains the
!>          min(M,N)-by-N upper triangular matrix R; the elements
!>          below the diagonal, together with the array TAU,
!>          represent the orthogonal matrix Q as a product of
!>          min(m,n) elementary reflectors.
!>
[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.
!>
[out]WORK
!>          WORK is DOUBLE PRECISION array, dimension (3*N)
!>
[out]INFO
!>          INFO is INTEGER
!>          = 0:  successful exit
!>          < 0:  if INFO = -i, the i-th argument had an illegal value
!>
Author
Univ. of Tennessee
Univ. of California Berkeley
Univ. of Colorado Denver
NAG Ltd.
Further Details:
!>
!>  The matrix Q is represented as a product of elementary reflectors
!>
!>     Q = H(1) H(2) . . . H(n)
!>
!>  Each H(i) has the form
!>
!>     H = I - tau * v * v**T
!>
!>  where tau is a real scalar, and v is a real vector with
!>  v(1:i-1) = 0 and v(i) = 1; v(i+1:m) is stored on exit in A(i+1:m,i).
!>
!>  The matrix P is represented in jpvt as follows: If
!>     jpvt(j) = i
!>  then the jth column of P is the ith canonical unit vector.
!>
!>  Partial column norm updating strategy modified by
!>    Z. Drmac and Z. Bujanovic, Dept. of Mathematics,
!>    University of Zagreb, Croatia.
!>  -- April 2011                                                      --
!>  For more details see LAPACK Working Note 176.
!>

Definition at line 139 of file dgeqpf.f.

140*
141* -- LAPACK computational routine --
142* -- LAPACK is a software package provided by Univ. of Tennessee, --
143* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
144*
145* .. Scalar Arguments ..
146 INTEGER INFO, LDA, M, N
147* ..
148* .. Array Arguments ..
149 INTEGER JPVT( * )
150 DOUBLE PRECISION A( LDA, * ), TAU( * ), WORK( * )
151* ..
152*
153* =====================================================================
154*
155* .. Parameters ..
156 DOUBLE PRECISION ZERO, ONE
157 parameter( zero = 0.0d+0, one = 1.0d+0 )
158* ..
159* .. Local Scalars ..
160 INTEGER I, ITEMP, J, MA, MN, PVT
161 DOUBLE PRECISION AII, TEMP, TEMP2, TOL3Z
162* ..
163* .. External Subroutines ..
164 EXTERNAL dgeqr2, dlarf, dlarfg, dorm2r, dswap, xerbla
165* ..
166* .. Intrinsic Functions ..
167 INTRINSIC abs, max, min, sqrt
168* ..
169* .. External Functions ..
170 INTEGER IDAMAX
171 DOUBLE PRECISION DLAMCH, DNRM2
172 EXTERNAL idamax, dlamch, dnrm2
173* ..
174* .. Executable Statements ..
175*
176* Test the input arguments
177*
178 info = 0
179 IF( m.LT.0 ) THEN
180 info = -1
181 ELSE IF( n.LT.0 ) THEN
182 info = -2
183 ELSE IF( lda.LT.max( 1, m ) ) THEN
184 info = -4
185 END IF
186 IF( info.NE.0 ) THEN
187 CALL xerbla( 'DGEQPF', -info )
188 RETURN
189 END IF
190*
191 mn = min( m, n )
192 tol3z = sqrt(dlamch('Epsilon'))
193*
194* Move initial columns up front
195*
196 itemp = 1
197 DO 10 i = 1, n
198 IF( jpvt( i ).NE.0 ) THEN
199 IF( i.NE.itemp ) THEN
200 CALL dswap( m, a( 1, i ), 1, a( 1, itemp ), 1 )
201 jpvt( i ) = jpvt( itemp )
202 jpvt( itemp ) = i
203 ELSE
204 jpvt( i ) = i
205 END IF
206 itemp = itemp + 1
207 ELSE
208 jpvt( i ) = i
209 END IF
210 10 CONTINUE
211 itemp = itemp - 1
212*
213* Compute the QR factorization and update remaining columns
214*
215 IF( itemp.GT.0 ) THEN
216 ma = min( itemp, m )
217 CALL dgeqr2( m, ma, a, lda, tau, work, info )
218 IF( ma.LT.n ) THEN
219 CALL dorm2r( 'Left', 'Transpose', m, n-ma, ma, a, lda,
220 $ tau, a( 1, ma+1 ), lda, work, info )
221 END IF
222 END IF
223*
224 IF( itemp.LT.mn ) THEN
225*
226* Initialize partial column norms. The first n elements of
227* work store the exact column norms.
228*
229 DO 20 i = itemp + 1, n
230 work( i ) = dnrm2( m-itemp, a( itemp+1, i ), 1 )
231 work( n+i ) = work( i )
232 20 CONTINUE
233*
234* Compute factorization
235*
236 DO 40 i = itemp + 1, mn
237*
238* Determine ith pivot column and swap if necessary
239*
240 pvt = ( i-1 ) + idamax( n-i+1, work( i ), 1 )
241*
242 IF( pvt.NE.i ) THEN
243 CALL dswap( m, a( 1, pvt ), 1, a( 1, i ), 1 )
244 itemp = jpvt( pvt )
245 jpvt( pvt ) = jpvt( i )
246 jpvt( i ) = itemp
247 work( pvt ) = work( i )
248 work( n+pvt ) = work( n+i )
249 END IF
250*
251* Generate elementary reflector H(i)
252*
253 IF( i.LT.m ) THEN
254 CALL dlarfg( m-i+1, a( i, i ), a( i+1, i ), 1,
255 $ tau( i ) )
256 ELSE
257 CALL dlarfg( 1, a( m, m ), a( m, m ), 1, tau( m ) )
258 END IF
259*
260 IF( i.LT.n ) THEN
261*
262* Apply H(i) to A(i:m,i+1:n) from the left
263*
264 aii = a( i, i )
265 a( i, i ) = one
266 CALL dlarf( 'LEFT', m-i+1, n-i, a( i, i ), 1,
267 $ tau( i ), a( i, i+1 ), lda, work( 2*n+1 ) )
268 a( i, i ) = aii
269 END IF
270*
271* Update partial column norms
272*
273 DO 30 j = i + 1, n
274 IF( work( j ).NE.zero ) THEN
275*
276* NOTE: The following 4 lines follow from the analysis in
277* Lapack Working Note 176.
278*
279 temp = abs( a( i, j ) ) / work( j )
280 temp = max( zero, ( one+temp )*( one-temp ) )
281 temp2 = temp*( work( j ) / work( n+j ) )**2
282 IF( temp2 .LE. tol3z ) THEN
283 IF( m-i.GT.0 ) THEN
284 work( j ) = dnrm2( m-i, a( i+1, j ), 1 )
285 work( n+j ) = work( j )
286 ELSE
287 work( j ) = zero
288 work( n+j ) = zero
289 END IF
290 ELSE
291 work( j ) = work( j )*sqrt( temp )
292 END IF
293 END IF
294 30 CONTINUE
295*
296 40 CONTINUE
297 END IF
298 RETURN
299*
300* End of DGEQPF
301*
xerbla
subroutine xerbla(srname, info)
Definition cblat2.f:3285
dgeqr2
subroutine dgeqr2(m, n, a, lda, tau, work, info)
DGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm.
Definition dgeqr2.f:128
idamax
integer function idamax(n, dx, incx)
IDAMAX
Definition idamax.f:71
dlamch
double precision function dlamch(cmach)
DLAMCH
Definition dlamch.f:69
dlarf
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
dlarfg
subroutine dlarfg(n, alpha, x, incx, tau)
DLARFG generates an elementary reflector (Householder matrix).
Definition dlarfg.f:104
dnrm2
real(wp) function dnrm2(n, x, incx)
DNRM2
Definition dnrm2.f90:89
dswap
subroutine dswap(n, dx, incx, dy, incy)
DSWAP
Definition dswap.f:82
dorm2r
subroutine dorm2r(side, trans, m, n, k, a, lda, tau, c, ldc, work, info)
DORM2R multiplies a general matrix by the orthogonal matrix from a QR factorization determined by sge...
Definition dorm2r.f:156
Here is the call graph for this function:
Here is the caller graph for this function: