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

DTZRZF

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

Purpose:
 DTZRZF reduces the M-by-N ( M<=N ) real upper trapezoidal matrix A
 to upper triangular form by means of orthogonal transformations.

 The upper trapezoidal matrix A is factored as

    A = ( R  0 ) * Z,

 where Z is an N-by-N orthogonal matrix and R is an M-by-M upper
 triangular matrix.
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 >= M.
[in,out]A
          A is DOUBLE PRECISION array, dimension (LDA,N)
          On entry, the leading M-by-N upper trapezoidal part of the
          array A must contain the matrix to be factorized.
          On exit, the leading M-by-M upper triangular part of A
          contains the upper triangular matrix R, and elements M+1 to
          N of the first M rows of A, with the array TAU, represent the
          orthogonal matrix Z as a product of M elementary reflectors.
[in]LDA
          LDA is INTEGER
          The leading dimension of the array A.  LDA >= max(1,M).
[out]TAU
          TAU is DOUBLE PRECISION array, dimension (M)
          The scalar factors of the elementary reflectors.
[out]WORK
          WORK is DOUBLE PRECISION 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 >= max(1,M).
          For optimum performance LWORK >= M*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]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
April 2012
Contributors:
A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA
Further Details:
  The N-by-N matrix Z can be computed by

     Z =  Z(1)*Z(2)* ... *Z(M)

  where each N-by-N Z(k) is given by

     Z(k) = I - tau(k)*v(k)*v(k)**T

  with v(k) is the kth row vector of the M-by-N matrix

     V = ( I   A(:,M+1:N) )

  I is the M-by-M identity matrix, A(:,M+1:N) 
  is the output stored in A on exit from DTZRZF,
  and tau(k) is the kth element of the array TAU.

Definition at line 153 of file dtzrzf.f.

153 *
154 * -- LAPACK computational routine (version 3.4.1) --
155 * -- LAPACK is a software package provided by Univ. of Tennessee, --
156 * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
157 * April 2012
158 *
159 * .. Scalar Arguments ..
160  INTEGER info, lda, lwork, m, n
161 * ..
162 * .. Array Arguments ..
163  DOUBLE PRECISION a( lda, * ), tau( * ), work( * )
164 * ..
165 *
166 * =====================================================================
167 *
168 * .. Parameters ..
169  DOUBLE PRECISION zero
170  parameter ( zero = 0.0d+0 )
171 * ..
172 * .. Local Scalars ..
173  LOGICAL lquery
174  INTEGER i, ib, iws, ki, kk, ldwork, lwkmin, lwkopt,
175  $ m1, mu, nb, nbmin, nx
176 * ..
177 * .. External Subroutines ..
178  EXTERNAL xerbla, dlarzb, dlarzt, dlatrz
179 * ..
180 * .. Intrinsic Functions ..
181  INTRINSIC max, min
182 * ..
183 * .. External Functions ..
184  INTEGER ilaenv
185  EXTERNAL ilaenv
186 * ..
187 * .. Executable Statements ..
188 *
189 * Test the input arguments
190 *
191  info = 0
192  lquery = ( lwork.EQ.-1 )
193  IF( m.LT.0 ) THEN
194  info = -1
195  ELSE IF( n.LT.m ) THEN
196  info = -2
197  ELSE IF( lda.LT.max( 1, m ) ) THEN
198  info = -4
199  END IF
200 *
201  IF( info.EQ.0 ) THEN
202  IF( m.EQ.0 .OR. m.EQ.n ) THEN
203  lwkopt = 1
204  lwkmin = 1
205  ELSE
206 *
207 * Determine the block size.
208 *
209  nb = ilaenv( 1, 'DGERQF', ' ', m, n, -1, -1 )
210  lwkopt = m*nb
211  lwkmin = max( 1, m )
212  END IF
213  work( 1 ) = lwkopt
214 *
215  IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
216  info = -7
217  END IF
218  END IF
219 *
220  IF( info.NE.0 ) THEN
221  CALL xerbla( 'DTZRZF', -info )
222  RETURN
223  ELSE IF( lquery ) THEN
224  RETURN
225  END IF
226 *
227 * Quick return if possible
228 *
229  IF( m.EQ.0 ) THEN
230  RETURN
231  ELSE IF( m.EQ.n ) THEN
232  DO 10 i = 1, n
233  tau( i ) = zero
234  10 CONTINUE
235  RETURN
236  END IF
237 *
238  nbmin = 2
239  nx = 1
240  iws = m
241  IF( nb.GT.1 .AND. nb.LT.m ) THEN
242 *
243 * Determine when to cross over from blocked to unblocked code.
244 *
245  nx = max( 0, ilaenv( 3, 'DGERQF', ' ', m, n, -1, -1 ) )
246  IF( nx.LT.m ) THEN
247 *
248 * Determine if workspace is large enough for blocked code.
249 *
250  ldwork = m
251  iws = ldwork*nb
252  IF( lwork.LT.iws ) THEN
253 *
254 * Not enough workspace to use optimal NB: reduce NB and
255 * determine the minimum value of NB.
256 *
257  nb = lwork / ldwork
258  nbmin = max( 2, ilaenv( 2, 'DGERQF', ' ', m, n, -1,
259  $ -1 ) )
260  END IF
261  END IF
262  END IF
263 *
264  IF( nb.GE.nbmin .AND. nb.LT.m .AND. nx.LT.m ) THEN
265 *
266 * Use blocked code initially.
267 * The last kk rows are handled by the block method.
268 *
269  m1 = min( m+1, n )
270  ki = ( ( m-nx-1 ) / nb )*nb
271  kk = min( m, ki+nb )
272 *
273  DO 20 i = m - kk + ki + 1, m - kk + 1, -nb
274  ib = min( m-i+1, nb )
275 *
276 * Compute the TZ factorization of the current block
277 * A(i:i+ib-1,i:n)
278 *
279  CALL dlatrz( ib, n-i+1, n-m, a( i, i ), lda, tau( i ),
280  $ work )
281  IF( i.GT.1 ) THEN
282 *
283 * Form the triangular factor of the block reflector
284 * H = H(i+ib-1) . . . H(i+1) H(i)
285 *
286  CALL dlarzt( 'Backward', 'Rowwise', n-m, ib, a( i, m1 ),
287  $ lda, tau( i ), work, ldwork )
288 *
289 * Apply H to A(1:i-1,i:n) from the right
290 *
291  CALL dlarzb( 'Right', 'No transpose', 'Backward',
292  $ 'Rowwise', i-1, n-i+1, ib, n-m, a( i, m1 ),
293  $ lda, work, ldwork, a( 1, i ), lda,
294  $ work( ib+1 ), ldwork )
295  END IF
296  20 CONTINUE
297  mu = i + nb - 1
298  ELSE
299  mu = m
300  END IF
301 *
302 * Use unblocked code to factor the last or only block
303 *
304  IF( mu.GT.0 )
305  $ CALL dlatrz( mu, n, n-m, a, lda, tau, work )
306 *
307  work( 1 ) = lwkopt
308 *
309  RETURN
310 *
311 * End of DTZRZF
312 *
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:62
subroutine dlatrz(M, N, L, A, LDA, TAU, WORK)
DLATRZ factors an upper trapezoidal matrix by means of orthogonal transformations.
Definition: dlatrz.f:142
subroutine dlarzt(DIRECT, STOREV, N, K, V, LDV, TAU, T, LDT)
DLARZT forms the triangular factor T of a block reflector H = I - vtvH.
Definition: dlarzt.f:187
integer function ilaenv(ISPEC, NAME, OPTS, N1, N2, N3, N4)
Definition: tstiee.f:83
subroutine dlarzb(SIDE, TRANS, DIRECT, STOREV, M, N, K, L, V, LDV, T, LDT, C, LDC, WORK, LDWORK)
DLARZB applies a block reflector or its transpose to a general matrix.
Definition: dlarzb.f:185

Here is the call graph for this function:

Here is the caller graph for this function: