LAPACK  3.6.1
LAPACK: Linear Algebra PACKage
subroutine zgeqp3 ( integer  M,
integer  N,
complex*16, dimension( lda, * )  A,
integer  LDA,
integer, dimension( * )  JPVT,
complex*16, dimension( * )  TAU,
complex*16, dimension( * )  WORK,
integer  LWORK,
double precision, dimension( * )  RWORK,
integer  INFO 
)

ZGEQP3

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

Purpose:
 ZGEQP3 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*16 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*16 array, dimension (min(M,N))
          The scalar factors of the elementary reflectors.
[out]WORK
          WORK is COMPLEX*16 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 DOUBLE PRECISION 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
November 2015
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 zgeqp3.f.

161 *
162 * -- LAPACK computational routine (version 3.6.0) --
163 * -- LAPACK is a software package provided by Univ. of Tennessee, --
164 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
165 * November 2015
166 *
167 * .. Scalar Arguments ..
168  INTEGER info, lda, lwork, m, n
169 * ..
170 * .. Array Arguments ..
171  INTEGER jpvt( * )
172  DOUBLE PRECISION rwork( * )
173  COMPLEX*16 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 xerbla, zgeqrf, zlaqp2, zlaqps, zswap, zunmqr
189 * ..
190 * .. External Functions ..
191  INTEGER ilaenv
192  DOUBLE PRECISION dznrm2
193  EXTERNAL ilaenv, dznrm2
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, 'ZGEQRF', ' ', m, n, -1, -1 )
221  lwkopt = ( n + 1 )*nb
222  END IF
223  work( 1 ) = dcmplx( 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( 'ZGEQP3', -info )
232  RETURN
233  ELSE IF( lquery ) THEN
234  RETURN
235  END IF
236 *
237 * Move initial columns up front.
238 *
239  nfxd = 1
240  DO 10 j = 1, n
241  IF( jpvt( j ).NE.0 ) THEN
242  IF( j.NE.nfxd ) THEN
243  CALL zswap( m, a( 1, j ), 1, a( 1, nfxd ), 1 )
244  jpvt( j ) = jpvt( nfxd )
245  jpvt( nfxd ) = j
246  ELSE
247  jpvt( j ) = j
248  END IF
249  nfxd = nfxd + 1
250  ELSE
251  jpvt( j ) = j
252  END IF
253  10 CONTINUE
254  nfxd = nfxd - 1
255 *
256 * Factorize fixed columns
257 * =======================
258 *
259 * Compute the QR factorization of fixed columns and update
260 * remaining columns.
261 *
262  IF( nfxd.GT.0 ) THEN
263  na = min( m, nfxd )
264 *CC CALL ZGEQR2( M, NA, A, LDA, TAU, WORK, INFO )
265  CALL zgeqrf( m, na, a, lda, tau, work, lwork, info )
266  iws = max( iws, int( work( 1 ) ) )
267  IF( na.LT.n ) THEN
268 *CC CALL ZUNM2R( 'Left', 'Conjugate Transpose', M, N-NA,
269 *CC $ NA, A, LDA, TAU, A( 1, NA+1 ), LDA, WORK,
270 *CC $ INFO )
271  CALL zunmqr( 'Left', 'Conjugate Transpose', m, n-na, na, a,
272  $ lda, tau, a( 1, na+1 ), lda, work, lwork,
273  $ info )
274  iws = max( iws, int( work( 1 ) ) )
275  END IF
276  END IF
277 *
278 * Factorize free columns
279 * ======================
280 *
281  IF( nfxd.LT.minmn ) THEN
282 *
283  sm = m - nfxd
284  sn = n - nfxd
285  sminmn = minmn - nfxd
286 *
287 * Determine the block size.
288 *
289  nb = ilaenv( inb, 'ZGEQRF', ' ', sm, sn, -1, -1 )
290  nbmin = 2
291  nx = 0
292 *
293  IF( ( nb.GT.1 ) .AND. ( nb.LT.sminmn ) ) THEN
294 *
295 * Determine when to cross over from blocked to unblocked code.
296 *
297  nx = max( 0, ilaenv( ixover, 'ZGEQRF', ' ', sm, sn, -1,
298  $ -1 ) )
299 *
300 *
301  IF( nx.LT.sminmn ) THEN
302 *
303 * Determine if workspace is large enough for blocked code.
304 *
305  minws = ( sn+1 )*nb
306  iws = max( iws, minws )
307  IF( lwork.LT.minws ) THEN
308 *
309 * Not enough workspace to use optimal NB: Reduce NB and
310 * determine the minimum value of NB.
311 *
312  nb = lwork / ( sn+1 )
313  nbmin = max( 2, ilaenv( inbmin, 'ZGEQRF', ' ', sm, sn,
314  $ -1, -1 ) )
315 *
316 *
317  END IF
318  END IF
319  END IF
320 *
321 * Initialize partial column norms. The first N elements of work
322 * store the exact column norms.
323 *
324  DO 20 j = nfxd + 1, n
325  rwork( j ) = dznrm2( sm, a( nfxd+1, j ), 1 )
326  rwork( n+j ) = rwork( j )
327  20 CONTINUE
328 *
329  IF( ( nb.GE.nbmin ) .AND. ( nb.LT.sminmn ) .AND.
330  $ ( nx.LT.sminmn ) ) THEN
331 *
332 * Use blocked code initially.
333 *
334  j = nfxd + 1
335 *
336 * Compute factorization: while loop.
337 *
338 *
339  topbmn = minmn - nx
340  30 CONTINUE
341  IF( j.LE.topbmn ) THEN
342  jb = min( nb, topbmn-j+1 )
343 *
344 * Factorize JB columns among columns J:N.
345 *
346  CALL zlaqps( m, n-j+1, j-1, jb, fjb, a( 1, j ), lda,
347  $ jpvt( j ), tau( j ), rwork( j ),
348  $ rwork( n+j ), work( 1 ), work( jb+1 ),
349  $ n-j+1 )
350 *
351  j = j + fjb
352  GO TO 30
353  END IF
354  ELSE
355  j = nfxd + 1
356  END IF
357 *
358 * Use unblocked code to factor the last or only block.
359 *
360 *
361  IF( j.LE.minmn )
362  $ CALL zlaqp2( m, n-j+1, j-1, a( 1, j ), lda, jpvt( j ),
363  $ tau( j ), rwork( j ), rwork( n+j ), work( 1 ) )
364 *
365  END IF
366 *
367  work( 1 ) = dcmplx( lwkopt )
368  RETURN
369 *
370 * End of ZGEQP3
371 *
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:151
subroutine zgeqrf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
ZGEQRF VARIANT: left-looking Level 3 BLAS of the algorithm.
Definition: zgeqrf.f:151
subroutine zswap(N, ZX, INCX, ZY, INCY)
ZSWAP
Definition: zswap.f:52
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine zunmqr(SIDE, TRANS, M, N, K, A, LDA, TAU, C, LDC, WORK, LWORK, INFO)
ZUNMQR
Definition: zunmqr.f:169
double precision function dznrm2(N, X, INCX)
DZNRM2
Definition: dznrm2.f:56
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine zlaqps(M, N, OFFSET, NB, KB, A, LDA, JPVT, TAU, VN1, VN2, AUXV, F, LDF)
ZLAQPS computes a step of QR factorization with column pivoting of a real m-by-n matrix A by using BL...
Definition: zlaqps.f:179

Here is the call graph for this function:

Here is the caller graph for this function: