LAPACK  3.9.0 LAPACK: Linear Algebra PACKage

## ◆ claqps()

 subroutine claqps ( integer M, integer N, integer OFFSET, integer NB, integer KB, complex, dimension( lda, * ) A, integer LDA, integer, dimension( * ) JPVT, complex, dimension( * ) TAU, real, dimension( * ) VN1, real, dimension( * ) VN2, complex, dimension( * ) AUXV, complex, dimension( ldf, * ) F, integer LDF )

CLAQPS computes a step of QR factorization with column pivoting of a real m-by-n matrix A by using BLAS level 3.

Purpose:
``` CLAQPS computes a step of QR factorization with column pivoting
of a complex M-by-N matrix A by using Blas-3.  It tries to factorize
NB columns from A starting from the row OFFSET+1, and updates all
of the matrix with Blas-3 xGEMM.

In some cases, due to catastrophic cancellations, it cannot
factorize NB columns.  Hence, the actual number of factorized
columns is returned in KB.

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 A that have been factorized in previous steps.``` [in] NB ``` NB is INTEGER The number of columns to factorize.``` [out] KB ``` KB is INTEGER The number of columns actually factorized.``` [in,out] A ``` A is COMPLEX array, dimension (LDA,N) On entry, the M-by-N matrix A. On exit, block A(OFFSET+1:M,1:KB) is the triangular factor obtained and block A(1:OFFSET,1:N) has been accordingly pivoted, but no factorized. The rest of the matrix, block A(OFFSET+1:M,KB+1:N) has been updated.``` [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) JPVT(I) = K <==> Column K of the full matrix A has been permuted into position I in AP.``` [out] TAU ``` TAU is COMPLEX array, dimension (KB) The scalar factors of the elementary reflectors.``` [in,out] VN1 ``` VN1 is REAL array, dimension (N) The vector with the partial column norms.``` [in,out] VN2 ``` VN2 is REAL array, dimension (N) The vector with the exact column norms.``` [in,out] AUXV ``` AUXV is COMPLEX array, dimension (NB) Auxiliary vector.``` [in,out] F ``` F is COMPLEX array, dimension (LDF,NB) Matrix F**H = L * Y**H * A.``` [in] LDF ``` LDF is INTEGER The leading dimension of the array F. LDF >= max(1,N).```
Date
December 2016
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 180 of file claqps.f.

180 *
181 * -- LAPACK auxiliary routine (version 3.7.0) --
182 * -- LAPACK is a software package provided by Univ. of Tennessee, --
183 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
184 * December 2016
185 *
186 * .. Scalar Arguments ..
187  INTEGER KB, LDA, LDF, M, N, NB, OFFSET
188 * ..
189 * .. Array Arguments ..
190  INTEGER JPVT( * )
191  REAL VN1( * ), VN2( * )
192  COMPLEX A( LDA, * ), AUXV( * ), F( LDF, * ), TAU( * )
193 * ..
194 *
195 * =====================================================================
196 *
197 * .. Parameters ..
198  REAL ZERO, ONE
199  COMPLEX CZERO, CONE
200  parameter( zero = 0.0e+0, one = 1.0e+0,
201  \$ czero = ( 0.0e+0, 0.0e+0 ),
202  \$ cone = ( 1.0e+0, 0.0e+0 ) )
203 * ..
204 * .. Local Scalars ..
205  INTEGER ITEMP, J, K, LASTRK, LSTICC, PVT, RK
206  REAL TEMP, TEMP2, TOL3Z
207  COMPLEX AKK
208 * ..
209 * .. External Subroutines ..
210  EXTERNAL cgemm, cgemv, clarfg, cswap
211 * ..
212 * .. Intrinsic Functions ..
213  INTRINSIC abs, conjg, max, min, nint, real, sqrt
214 * ..
215 * .. External Functions ..
216  INTEGER ISAMAX
217  REAL SCNRM2, SLAMCH
218  EXTERNAL isamax, scnrm2, slamch
219 * ..
220 * .. Executable Statements ..
221 *
222  lastrk = min( m, n+offset )
223  lsticc = 0
224  k = 0
225  tol3z = sqrt(slamch('Epsilon'))
226 *
227 * Beginning of while loop.
228 *
229  10 CONTINUE
230  IF( ( k.LT.nb ) .AND. ( lsticc.EQ.0 ) ) THEN
231  k = k + 1
232  rk = offset + k
233 *
234 * Determine ith pivot column and swap if necessary
235 *
236  pvt = ( k-1 ) + isamax( n-k+1, vn1( k ), 1 )
237  IF( pvt.NE.k ) THEN
238  CALL cswap( m, a( 1, pvt ), 1, a( 1, k ), 1 )
239  CALL cswap( k-1, f( pvt, 1 ), ldf, f( k, 1 ), ldf )
240  itemp = jpvt( pvt )
241  jpvt( pvt ) = jpvt( k )
242  jpvt( k ) = itemp
243  vn1( pvt ) = vn1( k )
244  vn2( pvt ) = vn2( k )
245  END IF
246 *
247 * Apply previous Householder reflectors to column K:
248 * A(RK:M,K) := A(RK:M,K) - A(RK:M,1:K-1)*F(K,1:K-1)**H.
249 *
250  IF( k.GT.1 ) THEN
251  DO 20 j = 1, k - 1
252  f( k, j ) = conjg( f( k, j ) )
253  20 CONTINUE
254  CALL cgemv( 'No transpose', m-rk+1, k-1, -cone, a( rk, 1 ),
255  \$ lda, f( k, 1 ), ldf, cone, a( rk, k ), 1 )
256  DO 30 j = 1, k - 1
257  f( k, j ) = conjg( f( k, j ) )
258  30 CONTINUE
259  END IF
260 *
261 * Generate elementary reflector H(k).
262 *
263  IF( rk.LT.m ) THEN
264  CALL clarfg( m-rk+1, a( rk, k ), a( rk+1, k ), 1, tau( k ) )
265  ELSE
266  CALL clarfg( 1, a( rk, k ), a( rk, k ), 1, tau( k ) )
267  END IF
268 *
269  akk = a( rk, k )
270  a( rk, k ) = cone
271 *
272 * Compute Kth column of F:
273 *
274 * Compute F(K+1:N,K) := tau(K)*A(RK:M,K+1:N)**H*A(RK:M,K).
275 *
276  IF( k.LT.n ) THEN
277  CALL cgemv( 'Conjugate transpose', m-rk+1, n-k, tau( k ),
278  \$ a( rk, k+1 ), lda, a( rk, k ), 1, czero,
279  \$ f( k+1, k ), 1 )
280  END IF
281 *
282 * Padding F(1:K,K) with zeros.
283 *
284  DO 40 j = 1, k
285  f( j, k ) = czero
286  40 CONTINUE
287 *
288 * Incremental updating of F:
289 * F(1:N,K) := F(1:N,K) - tau(K)*F(1:N,1:K-1)*A(RK:M,1:K-1)**H
290 * *A(RK:M,K).
291 *
292  IF( k.GT.1 ) THEN
293  CALL cgemv( 'Conjugate transpose', m-rk+1, k-1, -tau( k ),
294  \$ a( rk, 1 ), lda, a( rk, k ), 1, czero,
295  \$ auxv( 1 ), 1 )
296 *
297  CALL cgemv( 'No transpose', n, k-1, cone, f( 1, 1 ), ldf,
298  \$ auxv( 1 ), 1, cone, f( 1, k ), 1 )
299  END IF
300 *
301 * Update the current row of A:
302 * A(RK,K+1:N) := A(RK,K+1:N) - A(RK,1:K)*F(K+1:N,1:K)**H.
303 *
304  IF( k.LT.n ) THEN
305  CALL cgemm( 'No transpose', 'Conjugate transpose', 1, n-k,
306  \$ k, -cone, a( rk, 1 ), lda, f( k+1, 1 ), ldf,
307  \$ cone, a( rk, k+1 ), lda )
308  END IF
309 *
310 * Update partial column norms.
311 *
312  IF( rk.LT.lastrk ) THEN
313  DO 50 j = k + 1, n
314  IF( vn1( j ).NE.zero ) THEN
315 *
316 * NOTE: The following 4 lines follow from the analysis in
317 * Lapack Working Note 176.
318 *
319  temp = abs( a( rk, j ) ) / vn1( j )
320  temp = max( zero, ( one+temp )*( one-temp ) )
321  temp2 = temp*( vn1( j ) / vn2( j ) )**2
322  IF( temp2 .LE. tol3z ) THEN
323  vn2( j ) = real( lsticc )
324  lsticc = j
325  ELSE
326  vn1( j ) = vn1( j )*sqrt( temp )
327  END IF
328  END IF
329  50 CONTINUE
330  END IF
331 *
332  a( rk, k ) = akk
333 *
334 * End of while loop.
335 *
336  GO TO 10
337  END IF
338  kb = k
339  rk = offset + kb
340 *
341 * Apply the block reflector to the rest of the matrix:
342 * A(OFFSET+KB+1:M,KB+1:N) := A(OFFSET+KB+1:M,KB+1:N) -
343 * A(OFFSET+KB+1:M,1:KB)*F(KB+1:N,1:KB)**H.
344 *
345  IF( kb.LT.min( n, m-offset ) ) THEN
346  CALL cgemm( 'No transpose', 'Conjugate transpose', m-rk, n-kb,
347  \$ kb, -cone, a( rk+1, 1 ), lda, f( kb+1, 1 ), ldf,
348  \$ cone, a( rk+1, kb+1 ), lda )
349  END IF
350 *
351 * Recomputation of difficult columns.
352 *
353  60 CONTINUE
354  IF( lsticc.GT.0 ) THEN
355  itemp = nint( vn2( lsticc ) )
356  vn1( lsticc ) = scnrm2( m-rk, a( rk+1, lsticc ), 1 )
357 *
358 * NOTE: The computation of VN1( LSTICC ) relies on the fact that
359 * SNRM2 does not fail on vectors with norm below the value of
360 * SQRT(DLAMCH('S'))
361 *
362  vn2( lsticc ) = vn1( lsticc )
363  lsticc = itemp
364  GO TO 60
365  END IF
366 *
367  RETURN
368 *
369 * End of CLAQPS
370 *
Here is the call graph for this function:
Here is the caller graph for this function:
clarfg
subroutine clarfg(N, ALPHA, X, INCX, TAU)
CLARFG generates an elementary reflector (Householder matrix).
Definition: clarfg.f:108
cgemm
subroutine cgemm(TRANSA, TRANSB, M, N, K, ALPHA, A, LDA, B, LDB, BETA, C, LDC)
CGEMM
Definition: cgemm.f:189
scnrm2
real function scnrm2(N, X, INCX)
SCNRM2
Definition: scnrm2.f:77
cgemv
subroutine cgemv(TRANS, M, N, ALPHA, A, LDA, X, INCX, BETA, Y, INCY)
CGEMV
Definition: cgemv.f:160
isamax
integer function isamax(N, SX, INCX)
ISAMAX
Definition: isamax.f:73
cswap
subroutine cswap(N, CX, INCX, CY, INCY)
CSWAP
Definition: cswap.f:83
slamch
real function slamch(CMACH)
SLAMCH
Definition: slamch.f:70