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

## ◆ zlaqps()

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

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

Purpose:
``` ZLAQPS 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*16 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*16 array, dimension (KB) The scalar factors of the elementary reflectors.``` [in,out] VN1 ``` VN1 is DOUBLE PRECISION array, dimension (N) The vector with the partial column norms.``` [in,out] VN2 ``` VN2 is DOUBLE PRECISION array, dimension (N) The vector with the exact column norms.``` [in,out] AUXV ``` AUXV is COMPLEX*16 array, dimension (NB) Auxiliary vector.``` [in,out] F ``` F is COMPLEX*16 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).```
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 175 of file zlaqps.f.

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