LAPACK 3.12.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches

## ◆ sgeqp3()

 subroutine sgeqp3 ( integer m, integer n, real, dimension( lda, * ) a, integer lda, integer, dimension( * ) jpvt, real, dimension( * ) tau, real, dimension( * ) work, integer lwork, integer info )

SGEQP3

Purpose:
``` SGEQP3 computes a QR factorization with column pivoting of a
matrix A:  A*P = Q*R  using Level 3 BLAS.```
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 REAL 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 trapezoidal 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(J).ne.0, the J-th column of A is permuted to the front of A*P (a leading column); if JPVT(J)=0, the J-th column of A is a free column. On exit, if JPVT(J)=K, then the J-th column of A*P was the the K-th column of A.``` [out] TAU ``` TAU is REAL array, dimension (min(M,N)) The scalar factors of the elementary reflectors.``` [out] WORK ``` WORK is REAL array, dimension (MAX(1,LWORK)) On exit, if INFO=0, WORK(1) returns the optimal LWORK.``` [in] LWORK ``` LWORK is INTEGER The dimension of the array WORK. LWORK >= 3*N+1. For optimal performance LWORK >= 2*N+( N+1 )*NB, where NB is the optimal blocksize. If LWORK = -1, then a workspace query is assumed; the routine only calculates the optimal size of the WORK array, returns this value as the first entry of the WORK array, and no error message related to LWORK is issued by XERBLA.``` [out] INFO ``` INFO is INTEGER = 0: successful exit. < 0: if INFO = -i, the i-th argument had an illegal value.```
Further Details:
```  The matrix Q is represented as a product of elementary reflectors

Q = H(1) H(2) . . . H(k), where k = min(m,n).

Each H(i) has the form

H(i) = I - tau * v * v**T

where tau is a real scalar, and v is a real/complex 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), and tau in TAU(i).```
Contributors:
G. Quintana-Orti, Depto. de Informatica, Universidad Jaime I, Spain X. Sun, Computer Science Dept., Duke University, USA

Definition at line 150 of file sgeqp3.f.

151*
152* -- LAPACK computational routine --
153* -- LAPACK is a software package provided by Univ. of Tennessee, --
154* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
155*
156* .. Scalar Arguments ..
157 INTEGER INFO, LDA, LWORK, M, N
158* ..
159* .. Array Arguments ..
160 INTEGER JPVT( * )
161 REAL A( LDA, * ), TAU( * ), WORK( * )
162* ..
163*
164* =====================================================================
165*
166* .. Parameters ..
167 INTEGER INB, INBMIN, IXOVER
168 parameter( inb = 1, inbmin = 2, ixover = 3 )
169* ..
170* .. Local Scalars ..
171 LOGICAL LQUERY
172 INTEGER FJB, IWS, J, JB, LWKOPT, MINMN, MINWS, NA, NB,
173 \$ NBMIN, NFXD, NX, SM, SMINMN, SN, TOPBMN
174* ..
175* .. External Subroutines ..
176 EXTERNAL sgeqrf, slaqp2, slaqps, sormqr, sswap, xerbla
177* ..
178* .. External Functions ..
179 INTEGER ILAENV
180 REAL SNRM2, SROUNDUP_LWORK
181 EXTERNAL ilaenv, snrm2, sroundup_lwork
182* ..
183* .. Intrinsic Functions ..
184 INTRINSIC int, max, min
185* Test input arguments
186* ====================
187*
188 info = 0
189 lquery = ( lwork.EQ.-1 )
190 IF( m.LT.0 ) THEN
191 info = -1
192 ELSE IF( n.LT.0 ) THEN
193 info = -2
194 ELSE IF( lda.LT.max( 1, m ) ) THEN
195 info = -4
196 END IF
197*
198 IF( info.EQ.0 ) THEN
199 minmn = min( m, n )
200 IF( minmn.EQ.0 ) THEN
201 iws = 1
202 lwkopt = 1
203 ELSE
204 iws = 3*n + 1
205 nb = ilaenv( inb, 'SGEQRF', ' ', m, n, -1, -1 )
206 lwkopt = 2*n + ( n + 1 )*nb
207 END IF
208 work( 1 ) = sroundup_lwork(lwkopt)
209*
210 IF( ( lwork.LT.iws ) .AND. .NOT.lquery ) THEN
211 info = -8
212 END IF
213 END IF
214*
215 IF( info.NE.0 ) THEN
216 CALL xerbla( 'SGEQP3', -info )
217 RETURN
218 ELSE IF( lquery ) THEN
219 RETURN
220 END IF
221*
222* Move initial columns up front.
223*
224 nfxd = 1
225 DO 10 j = 1, n
226 IF( jpvt( j ).NE.0 ) THEN
227 IF( j.NE.nfxd ) THEN
228 CALL sswap( m, a( 1, j ), 1, a( 1, nfxd ), 1 )
229 jpvt( j ) = jpvt( nfxd )
230 jpvt( nfxd ) = j
231 ELSE
232 jpvt( j ) = j
233 END IF
234 nfxd = nfxd + 1
235 ELSE
236 jpvt( j ) = j
237 END IF
238 10 CONTINUE
239 nfxd = nfxd - 1
240*
241* Factorize fixed columns
242* =======================
243*
244* Compute the QR factorization of fixed columns and update
245* remaining columns.
246*
247 IF( nfxd.GT.0 ) THEN
248 na = min( m, nfxd )
249*CC CALL SGEQR2( M, NA, A, LDA, TAU, WORK, INFO )
250 CALL sgeqrf( m, na, a, lda, tau, work, lwork, info )
251 iws = max( iws, int( work( 1 ) ) )
252 IF( na.LT.n ) THEN
253*CC CALL SORM2R( 'Left', 'Transpose', M, N-NA, NA, A, LDA,
254*CC \$ TAU, A( 1, NA+1 ), LDA, WORK, INFO )
255 CALL sormqr( 'Left', 'Transpose', m, n-na, na, a, lda, tau,
256 \$ a( 1, na+1 ), lda, work, lwork, info )
257 iws = max( iws, int( work( 1 ) ) )
258 END IF
259 END IF
260*
261* Factorize free columns
262* ======================
263*
264 IF( nfxd.LT.minmn ) THEN
265*
266 sm = m - nfxd
267 sn = n - nfxd
268 sminmn = minmn - nfxd
269*
270* Determine the block size.
271*
272 nb = ilaenv( inb, 'SGEQRF', ' ', sm, sn, -1, -1 )
273 nbmin = 2
274 nx = 0
275*
276 IF( ( nb.GT.1 ) .AND. ( nb.LT.sminmn ) ) THEN
277*
278* Determine when to cross over from blocked to unblocked code.
279*
280 nx = max( 0, ilaenv( ixover, 'SGEQRF', ' ', sm, sn, -1,
281 \$ -1 ) )
282*
283*
284 IF( nx.LT.sminmn ) THEN
285*
286* Determine if workspace is large enough for blocked code.
287*
288 minws = 2*sn + ( sn+1 )*nb
289 iws = max( iws, minws )
290 IF( lwork.LT.minws ) THEN
291*
292* Not enough workspace to use optimal NB: Reduce NB and
293* determine the minimum value of NB.
294*
295 nb = ( lwork-2*sn ) / ( sn+1 )
296 nbmin = max( 2, ilaenv( inbmin, 'SGEQRF', ' ', sm, sn,
297 \$ -1, -1 ) )
298*
299*
300 END IF
301 END IF
302 END IF
303*
304* Initialize partial column norms. The first N elements of work
305* store the exact column norms.
306*
307 DO 20 j = nfxd + 1, n
308 work( j ) = snrm2( sm, a( nfxd+1, j ), 1 )
309 work( n+j ) = work( j )
310 20 CONTINUE
311*
312 IF( ( nb.GE.nbmin ) .AND. ( nb.LT.sminmn ) .AND.
313 \$ ( nx.LT.sminmn ) ) THEN
314*
315* Use blocked code initially.
316*
317 j = nfxd + 1
318*
319* Compute factorization: while loop.
320*
321*
322 topbmn = minmn - nx
323 30 CONTINUE
324 IF( j.LE.topbmn ) THEN
325 jb = min( nb, topbmn-j+1 )
326*
327* Factorize JB columns among columns J:N.
328*
329 CALL slaqps( m, n-j+1, j-1, jb, fjb, a( 1, j ), lda,
330 \$ jpvt( j ), tau( j ), work( j ), work( n+j ),
331 \$ work( 2*n+1 ), work( 2*n+jb+1 ), n-j+1 )
332*
333 j = j + fjb
334 GO TO 30
335 END IF
336 ELSE
337 j = nfxd + 1
338 END IF
339*
340* Use unblocked code to factor the last or only block.
341*
342*
343 IF( j.LE.minmn )
344 \$ CALL slaqp2( m, n-j+1, j-1, a( 1, j ), lda, jpvt( j ),
345 \$ tau( j ), work( j ), work( n+j ),
346 \$ work( 2*n+1 ) )
347*
348 END IF
349*
350 work( 1 ) = sroundup_lwork(iws)
351 RETURN
352*
353* End of SGEQP3
354*
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine sgeqrf(m, n, a, lda, tau, work, lwork, info)
SGEQRF
Definition sgeqrf.f:146
integer function ilaenv(ispec, name, opts, n1, n2, n3, n4)
ILAENV
Definition ilaenv.f:162
subroutine slaqp2(m, n, offset, a, lda, jpvt, tau, vn1, vn2, work)
SLAQP2 computes a QR factorization with column pivoting of the matrix block.
Definition slaqp2.f:149
subroutine slaqps(m, n, offset, nb, kb, a, lda, jpvt, tau, vn1, vn2, auxv, f, ldf)
SLAQPS computes a step of QR factorization with column pivoting of a real m-by-n matrix A by using BL...
Definition slaqps.f:178
real(wp) function snrm2(n, x, incx)
SNRM2
Definition snrm2.f90:89
real function sroundup_lwork(lwork)
SROUNDUP_LWORK
subroutine sswap(n, sx, incx, sy, incy)
SSWAP
Definition sswap.f:82
subroutine sormqr(side, trans, m, n, k, a, lda, tau, c, ldc, work, lwork, info)
SORMQR
Definition sormqr.f:168
Here is the call graph for this function:
Here is the caller graph for this function: