 LAPACK  3.10.1 LAPACK: Linear Algebra PACKage

## ◆ cgeqrf()

 subroutine cgeqrf ( integer M, integer N, complex, dimension( lda, * ) A, integer LDA, complex, dimension( * ) TAU, complex, dimension( * ) WORK, integer LWORK, integer INFO )

CGEQRF VARIANT: left-looking Level 3 BLAS version of the algorithm.

Purpose:

``` CGEQRF computes a QR factorization of a real M-by-N matrix A:
A = Q * R.

This is the left-looking Level 3 BLAS version of the algorithm.```
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 elements on and above the diagonal of the array contain the min(M,N)-by-N upper trapezoidal matrix R (R is upper triangular if m >= n); the elements below the diagonal, with the array TAU, represent the orthogonal matrix Q as a product of min(m,n) elementary reflectors (see Further Details).``` [in] LDA ``` LDA is INTEGER The leading dimension of the array A. LDA >= max(1,M).``` [out] TAU ``` TAU is COMPLEX array, dimension (min(M,N)) The scalar factors of the elementary reflectors (see Further Details).``` [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. The dimension can be divided into three parts.` ``` 1) The part for the triangular factor T. If the very last T is not bigger than any of the rest, then this part is NB x ceiling(K/NB), otherwise, NB x (K-NT), where K = min(M,N) and NT is the dimension of the very last T``` ``` 2) The part for the very last T when T is bigger than any of the rest T. The size of this part is NT x NT, where NT = K - ceiling ((K-NX)/NB) x NB, where K = min(M,N), NX is calculated by NX = MAX( 0, ILAENV( 3, 'CGEQRF', ' ', M, N, -1, -1 ) )``` ` 3) The part for dlarfb is of size max((N-M)*K, (N-M)*NB, K*NB, NB*NB)` ` So LWORK = part1 + part2 + part3` ``` 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```
Date
December 2016

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'

where tau is a real scalar, and v is a real 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).```

Definition at line 150 of file cgeqrf.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  COMPLEX A( LDA, * ), TAU( * ), WORK( * )
161 * ..
162 *
163 * =====================================================================
164 *
165 * .. Local Scalars ..
166  LOGICAL LQUERY
167  INTEGER I, IB, IINFO, IWS, J, K, LWKOPT, NB,
168  \$ NBMIN, NX, LBWORK, NT, LLWORK
169 * ..
170 * .. External Subroutines ..
171  EXTERNAL cgeqr2, clarfb, clarft, xerbla
172 * ..
173 * .. Intrinsic Functions ..
174  INTRINSIC max, min
175 * ..
176 * .. External Functions ..
177  INTEGER ILAENV
178  REAL SCEIL
179  EXTERNAL ilaenv, sceil
180 * ..
181 * .. Executable Statements ..
182
183  info = 0
184  nbmin = 2
185  nx = 0
186  iws = n
187  k = min( m, n )
188  nb = ilaenv( 1, 'CGEQRF', ' ', m, n, -1, -1 )
189
190  IF( nb.GT.1 .AND. nb.LT.k ) THEN
191 *
192 * Determine when to cross over from blocked to unblocked code.
193 *
194  nx = max( 0, ilaenv( 3, 'CGEQRF', ' ', m, n, -1, -1 ) )
195  END IF
196 *
197 * Get NT, the size of the very last T, which is the left-over from in-between K-NX and K to K, eg.:
198 *
199 * NB=3 2NB=6 K=10
200 * | | |
201 * 1--2--3--4--5--6--7--8--9--10
202 * | \________/
203 * K-NX=5 NT=4
204 *
205 * So here 4 x 4 is the last T stored in the workspace
206 *
207  nt = k-sceil(real(k-nx)/real(nb))*nb
208
209 *
210 * optimal workspace = space for dlarfb + space for normal T's + space for the last T
211 *
212  llwork = max(max((n-m)*k, (n-m)*nb), max(k*nb, nb*nb))
213  llwork = sceil(real(llwork)/real(nb))
214
215  IF ( nt.GT.nb ) THEN
216
217  lbwork = k-nt
218 *
219 * Optimal workspace for dlarfb = MAX(1,N)*NT
220 *
221  lwkopt = (lbwork+llwork)*nb
222  work( 1 ) = (lwkopt+nt*nt)
223
224  ELSE
225
226  lbwork = sceil(real(k)/real(nb))*nb
227  lwkopt = (lbwork+llwork-nb)*nb
228  work( 1 ) = lwkopt
229
230  END IF
231
232 *
233 * Test the input arguments
234 *
235  lquery = ( lwork.EQ.-1 )
236  IF( m.LT.0 ) THEN
237  info = -1
238  ELSE IF( n.LT.0 ) THEN
239  info = -2
240  ELSE IF( lda.LT.max( 1, m ) ) THEN
241  info = -4
242  ELSE IF( lwork.LT.max( 1, n ) .AND. .NOT.lquery ) THEN
243  info = -7
244  END IF
245  IF( info.NE.0 ) THEN
246  CALL xerbla( 'CGEQRF', -info )
247  RETURN
248  ELSE IF( lquery ) THEN
249  RETURN
250  END IF
251 *
252 * Quick return if possible
253 *
254  IF( k.EQ.0 ) THEN
255  work( 1 ) = 1
256  RETURN
257  END IF
258 *
259  IF( nb.GT.1 .AND. nb.LT.k ) THEN
260
261  IF( nx.LT.k ) THEN
262 *
263 * Determine if workspace is large enough for blocked code.
264 *
265  IF ( nt.LE.nb ) THEN
266  iws = (lbwork+llwork-nb)*nb
267  ELSE
268  iws = (lbwork+llwork)*nb+nt*nt
269  END IF
270
271  IF( lwork.LT.iws ) THEN
272 *
273 * Not enough workspace to use optimal NB: reduce NB and
274 * determine the minimum value of NB.
275 *
276  IF ( nt.LE.nb ) THEN
277  nb = lwork / (llwork+(lbwork-nb))
278  ELSE
279  nb = (lwork-nt*nt)/(lbwork+llwork)
280  END IF
281
282  nbmin = max( 2, ilaenv( 2, 'CGEQRF', ' ', m, n, -1,
283  \$ -1 ) )
284  END IF
285  END IF
286  END IF
287 *
288  IF( nb.GE.nbmin .AND. nb.LT.k .AND. nx.LT.k ) THEN
289 *
290 * Use blocked code initially
291 *
292  DO 10 i = 1, k - nx, nb
293  ib = min( k-i+1, nb )
294 *
295 * Update the current column using old T's
296 *
297  DO 20 j = 1, i - nb, nb
298 *
299 * Apply H' to A(J:M,I:I+IB-1) from the left
300 *
301  CALL clarfb( 'Left', 'Transpose', 'Forward',
302  \$ 'Columnwise', m-j+1, ib, nb,
303  \$ a( j, j ), lda, work(j), lbwork,
304  \$ a( j, i ), lda, work(lbwork*nb+nt*nt+1),
305  \$ ib)
306
307 20 CONTINUE
308 *
309 * Compute the QR factorization of the current block
310 * A(I:M,I:I+IB-1)
311 *
312  CALL cgeqr2( m-i+1, ib, a( i, i ), lda, tau( i ),
313  \$ work(lbwork*nb+nt*nt+1), iinfo )
314
315  IF( i+ib.LE.n ) THEN
316 *
317 * Form the triangular factor of the block reflector
318 * H = H(i) H(i+1) . . . H(i+ib-1)
319 *
320  CALL clarft( 'Forward', 'Columnwise', m-i+1, ib,
321  \$ a( i, i ), lda, tau( i ),
322  \$ work(i), lbwork )
323 *
324  END IF
325  10 CONTINUE
326  ELSE
327  i = 1
328  END IF
329 *
330 * Use unblocked code to factor the last or only block.
331 *
332  IF( i.LE.k ) THEN
333
334  IF ( i .NE. 1 ) THEN
335
336  DO 30 j = 1, i - nb, nb
337 *
338 * Apply H' to A(J:M,I:K) from the left
339 *
340  CALL clarfb( 'Left', 'Transpose', 'Forward',
341  \$ 'Columnwise', m-j+1, k-i+1, nb,
342  \$ a( j, j ), lda, work(j), lbwork,
343  \$ a( j, i ), lda, work(lbwork*nb+nt*nt+1),
344  \$ k-i+1)
345 30 CONTINUE
346
347  CALL cgeqr2( m-i+1, k-i+1, a( i, i ), lda, tau( i ),
348  \$ work(lbwork*nb+nt*nt+1),iinfo )
349
350  ELSE
351 *
352 * Use unblocked code to factor the last or only block.
353 *
354  CALL cgeqr2( m-i+1, n-i+1, a( i, i ), lda, tau( i ),
355  \$ work,iinfo )
356
357  END IF
358  END IF
359
360
361 *
362 * Apply update to the column M+1:N when N > M
363 *
364  IF ( m.LT.n .AND. i.NE.1) THEN
365 *
366 * Form the last triangular factor of the block reflector
367 * H = H(i) H(i+1) . . . H(i+ib-1)
368 *
369  IF ( nt .LE. nb ) THEN
370  CALL clarft( 'Forward', 'Columnwise', m-i+1, k-i+1,
371  \$ a( i, i ), lda, tau( i ), work(i), lbwork )
372  ELSE
373  CALL clarft( 'Forward', 'Columnwise', m-i+1, k-i+1,
374  \$ a( i, i ), lda, tau( i ),
375  \$ work(lbwork*nb+1), nt )
376  END IF
377
378 *
379 * Apply H' to A(1:M,M+1:N) from the left
380 *
381  DO 40 j = 1, k-nx, nb
382
383  ib = min( k-j+1, nb )
384
385  CALL clarfb( 'Left', 'Transpose', 'Forward',
386  \$ 'Columnwise', m-j+1, n-m, ib,
387  \$ a( j, j ), lda, work(j), lbwork,
388  \$ a( j, m+1 ), lda, work(lbwork*nb+nt*nt+1),
389  \$ n-m)
390
391 40 CONTINUE
392
393  IF ( nt.LE.nb ) THEN
394  CALL clarfb( 'Left', 'Transpose', 'Forward',
395  \$ 'Columnwise', m-j+1, n-m, k-j+1,
396  \$ a( j, j ), lda, work(j), lbwork,
397  \$ a( j, m+1 ), lda, work(lbwork*nb+nt*nt+1),
398  \$ n-m)
399  ELSE
400  CALL clarfb( 'Left', 'Transpose', 'Forward',
401  \$ 'Columnwise', m-j+1, n-m, k-j+1,
402  \$ a( j, j ), lda,
403  \$ work(lbwork*nb+1),
404  \$ nt, a( j, m+1 ), lda, work(lbwork*nb+nt*nt+1),
405  \$ n-m)
406  END IF
407
408  END IF
409
410  work( 1 ) = iws
411  RETURN
412 *
413 * End of CGEQRF
414 *
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
ILAENV
Definition: ilaenv.f:162
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine cgeqr2(M, N, A, LDA, TAU, WORK, INFO)
CGEQR2 computes the QR factorization of a general rectangular matrix using an unblocked algorithm.
Definition: cgeqr2.f:130
subroutine clarfb(SIDE, TRANS, DIRECT, STOREV, M, N, K, V, LDV, T, LDT, C, LDC, WORK, LDWORK)
CLARFB applies a block reflector or its conjugate-transpose to a general rectangular matrix.
Definition: clarfb.f:197
subroutine clarft(DIRECT, STOREV, N, K, V, LDV, TAU, T, LDT)
CLARFT forms the triangular factor T of a block reflector H = I - vtvH
Definition: clarft.f:163
real function sceil(A)
SCEIL
Definition: sceil.f:60
Here is the call graph for this function: