LAPACK  3.5.0
LAPACK: Linear Algebra PACKage
cgeqp3.f File Reference

Go to the source code of this file.

Functions/Subroutines

subroutine cgeqp3 (M, N, A, LDA, JPVT, TAU, WORK, LWORK, RWORK, INFO)
 CGEQP3 More...
 

Function/Subroutine Documentation

subroutine cgeqp3 ( integer  M,
integer  N,
complex, dimension( lda, * )  A,
integer  LDA,
integer, dimension( * )  JPVT,
complex, dimension( * )  TAU,
complex, dimension( * )  WORK,
integer  LWORK,
real, dimension( * )  RWORK,
integer  INFO 
)

CGEQP3

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

Purpose:
 CGEQP3 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 COMPLEX 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
          unitary 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 COMPLEX array, dimension (min(M,N))
          The scalar factors of the elementary reflectors.
[out]WORK
          WORK is COMPLEX 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 >= N+1.
          For optimal performance LWORK >= ( 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]RWORK
          RWORK is REAL array, dimension (2*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.
Date
September 2012
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**H

  where tau is a complex 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 161 of file cgeqp3.f.

161 *
162 * -- LAPACK computational routine (version 3.4.2) --
163 * -- LAPACK is a software package provided by Univ. of Tennessee, --
164 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
165 * September 2012
166 *
167 * .. Scalar Arguments ..
168  INTEGER info, lda, lwork, m, n
169 * ..
170 * .. Array Arguments ..
171  INTEGER jpvt( * )
172  REAL rwork( * )
173  COMPLEX a( lda, * ), tau( * ), work( * )
174 * ..
175 *
176 * =====================================================================
177 *
178 * .. Parameters ..
179  INTEGER inb, inbmin, ixover
180  parameter( inb = 1, inbmin = 2, ixover = 3 )
181 * ..
182 * .. Local Scalars ..
183  LOGICAL lquery
184  INTEGER fjb, iws, j, jb, lwkopt, minmn, minws, na, nb,
185  $ nbmin, nfxd, nx, sm, sminmn, sn, topbmn
186 * ..
187 * .. External Subroutines ..
188  EXTERNAL cgeqrf, claqp2, claqps, cswap, cunmqr, xerbla
189 * ..
190 * .. External Functions ..
191  INTEGER ilaenv
192  REAL scnrm2
193  EXTERNAL ilaenv, scnrm2
194 * ..
195 * .. Intrinsic Functions ..
196  INTRINSIC int, max, min
197 * ..
198 * .. Executable Statements ..
199 *
200 * Test input arguments
201 * ====================
202 *
203  info = 0
204  lquery = ( lwork.EQ.-1 )
205  IF( m.LT.0 ) THEN
206  info = -1
207  ELSE IF( n.LT.0 ) THEN
208  info = -2
209  ELSE IF( lda.LT.max( 1, m ) ) THEN
210  info = -4
211  END IF
212 *
213  IF( info.EQ.0 ) THEN
214  minmn = min( m, n )
215  IF( minmn.EQ.0 ) THEN
216  iws = 1
217  lwkopt = 1
218  ELSE
219  iws = n + 1
220  nb = ilaenv( inb, 'CGEQRF', ' ', m, n, -1, -1 )
221  lwkopt = ( n + 1 )*nb
222  END IF
223  work( 1 ) = lwkopt
224 *
225  IF( ( lwork.LT.iws ) .AND. .NOT.lquery ) THEN
226  info = -8
227  END IF
228  END IF
229 *
230  IF( info.NE.0 ) THEN
231  CALL xerbla( 'CGEQP3', -info )
232  RETURN
233  ELSE IF( lquery ) THEN
234  RETURN
235  END IF
236 *
237 * Quick return if possible.
238 *
239  IF( minmn.EQ.0 ) THEN
240  RETURN
241  END IF
242 *
243 * Move initial columns up front.
244 *
245  nfxd = 1
246  DO 10 j = 1, n
247  IF( jpvt( j ).NE.0 ) THEN
248  IF( j.NE.nfxd ) THEN
249  CALL cswap( m, a( 1, j ), 1, a( 1, nfxd ), 1 )
250  jpvt( j ) = jpvt( nfxd )
251  jpvt( nfxd ) = j
252  ELSE
253  jpvt( j ) = j
254  END IF
255  nfxd = nfxd + 1
256  ELSE
257  jpvt( j ) = j
258  END IF
259  10 CONTINUE
260  nfxd = nfxd - 1
261 *
262 * Factorize fixed columns
263 * =======================
264 *
265 * Compute the QR factorization of fixed columns and update
266 * remaining columns.
267 *
268  IF( nfxd.GT.0 ) THEN
269  na = min( m, nfxd )
270 *CC CALL CGEQR2( M, NA, A, LDA, TAU, WORK, INFO )
271  CALL cgeqrf( m, na, a, lda, tau, work, lwork, info )
272  iws = max( iws, int( work( 1 ) ) )
273  IF( na.LT.n ) THEN
274 *CC CALL CUNM2R( 'Left', 'Conjugate Transpose', M, N-NA,
275 *CC $ NA, A, LDA, TAU, A( 1, NA+1 ), LDA, WORK,
276 *CC $ INFO )
277  CALL cunmqr( 'Left', 'Conjugate Transpose', m, n-na, na, a,
278  $ lda, tau, a( 1, na+1 ), lda, work, lwork,
279  $ info )
280  iws = max( iws, int( work( 1 ) ) )
281  END IF
282  END IF
283 *
284 * Factorize free columns
285 * ======================
286 *
287  IF( nfxd.LT.minmn ) THEN
288 *
289  sm = m - nfxd
290  sn = n - nfxd
291  sminmn = minmn - nfxd
292 *
293 * Determine the block size.
294 *
295  nb = ilaenv( inb, 'CGEQRF', ' ', sm, sn, -1, -1 )
296  nbmin = 2
297  nx = 0
298 *
299  IF( ( nb.GT.1 ) .AND. ( nb.LT.sminmn ) ) THEN
300 *
301 * Determine when to cross over from blocked to unblocked code.
302 *
303  nx = max( 0, ilaenv( ixover, 'CGEQRF', ' ', sm, sn, -1,
304  $ -1 ) )
305 *
306 *
307  IF( nx.LT.sminmn ) THEN
308 *
309 * Determine if workspace is large enough for blocked code.
310 *
311  minws = ( sn+1 )*nb
312  iws = max( iws, minws )
313  IF( lwork.LT.minws ) THEN
314 *
315 * Not enough workspace to use optimal NB: Reduce NB and
316 * determine the minimum value of NB.
317 *
318  nb = lwork / ( sn+1 )
319  nbmin = max( 2, ilaenv( inbmin, 'CGEQRF', ' ', sm, sn,
320  $ -1, -1 ) )
321 *
322 *
323  END IF
324  END IF
325  END IF
326 *
327 * Initialize partial column norms. The first N elements of work
328 * store the exact column norms.
329 *
330  DO 20 j = nfxd + 1, n
331  rwork( j ) = scnrm2( sm, a( nfxd+1, j ), 1 )
332  rwork( n+j ) = rwork( j )
333  20 CONTINUE
334 *
335  IF( ( nb.GE.nbmin ) .AND. ( nb.LT.sminmn ) .AND.
336  $ ( nx.LT.sminmn ) ) THEN
337 *
338 * Use blocked code initially.
339 *
340  j = nfxd + 1
341 *
342 * Compute factorization: while loop.
343 *
344 *
345  topbmn = minmn - nx
346  30 CONTINUE
347  IF( j.LE.topbmn ) THEN
348  jb = min( nb, topbmn-j+1 )
349 *
350 * Factorize JB columns among columns J:N.
351 *
352  CALL claqps( m, n-j+1, j-1, jb, fjb, a( 1, j ), lda,
353  $ jpvt( j ), tau( j ), rwork( j ),
354  $ rwork( n+j ), work( 1 ), work( jb+1 ),
355  $ n-j+1 )
356 *
357  j = j + fjb
358  GO TO 30
359  END IF
360  ELSE
361  j = nfxd + 1
362  END IF
363 *
364 * Use unblocked code to factor the last or only block.
365 *
366 *
367  IF( j.LE.minmn )
368  $ CALL claqp2( m, n-j+1, j-1, a( 1, j ), lda, jpvt( j ),
369  $ tau( j ), rwork( j ), rwork( n+j ), work( 1 ) )
370 *
371  END IF
372 *
373  work( 1 ) = iws
374  RETURN
375 *
376 * End of CGEQP3
377 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine claqp2(M, N, OFFSET, A, LDA, JPVT, TAU, VN1, VN2, WORK)
CLAQP2 computes a QR factorization with column pivoting of the matrix block.
Definition: claqp2.f:151
subroutine cunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
CUNMQR
Definition: cunmqr.f:172
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
Definition: cswap.f:52
subroutine claqps(M, N, OFFSET, NB, KB, A, LDA, JPVT, TAU, VN1, VN2, AUXV, F, LDF)
CLAQPS computes a step of QR factorization with column pivoting of a real m-by-n matrix A by using BL...
Definition: claqps.f:180
real function scnrm2(N, X, INCX)
SCNRM2
Definition: scnrm2.f:56
subroutine cgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
CGEQRF
Definition: cgeqrf.f:138
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83

Here is the call graph for this function:

Here is the caller graph for this function: