LAPACK 3.11.0 LAPACK: Linear Algebra PACKage
Searching...
No Matches
dtzrzf.f
Go to the documentation of this file.
1*> \brief \b DTZRZF
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dtzrzf.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dtzrzf.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dtzrzf.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DTZRZF( M, N, A, LDA, TAU, WORK, LWORK, INFO )
22*
23* .. Scalar Arguments ..
24* INTEGER INFO, LDA, LWORK, M, N
25* ..
26* .. Array Arguments ..
27* DOUBLE PRECISION A( LDA, * ), TAU( * ), WORK( * )
28* ..
29*
30*
31*> \par Purpose:
32* =============
33*>
34*> \verbatim
35*>
36*> DTZRZF reduces the M-by-N ( M<=N ) real upper trapezoidal matrix A
37*> to upper triangular form by means of orthogonal transformations.
38*>
39*> The upper trapezoidal matrix A is factored as
40*>
41*> A = ( R 0 ) * Z,
42*>
43*> where Z is an N-by-N orthogonal matrix and R is an M-by-M upper
44*> triangular matrix.
45*> \endverbatim
46*
47* Arguments:
48* ==========
49*
50*> \param[in] M
51*> \verbatim
52*> M is INTEGER
53*> The number of rows of the matrix A. M >= 0.
54*> \endverbatim
55*>
56*> \param[in] N
57*> \verbatim
58*> N is INTEGER
59*> The number of columns of the matrix A. N >= M.
60*> \endverbatim
61*>
62*> \param[in,out] A
63*> \verbatim
64*> A is DOUBLE PRECISION array, dimension (LDA,N)
65*> On entry, the leading M-by-N upper trapezoidal part of the
66*> array A must contain the matrix to be factorized.
67*> On exit, the leading M-by-M upper triangular part of A
68*> contains the upper triangular matrix R, and elements M+1 to
69*> N of the first M rows of A, with the array TAU, represent the
70*> orthogonal matrix Z as a product of M elementary reflectors.
71*> \endverbatim
72*>
73*> \param[in] LDA
74*> \verbatim
75*> LDA is INTEGER
76*> The leading dimension of the array A. LDA >= max(1,M).
77*> \endverbatim
78*>
79*> \param[out] TAU
80*> \verbatim
81*> TAU is DOUBLE PRECISION array, dimension (M)
82*> The scalar factors of the elementary reflectors.
83*> \endverbatim
84*>
85*> \param[out] WORK
86*> \verbatim
87*> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
88*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
89*> \endverbatim
90*>
91*> \param[in] LWORK
92*> \verbatim
93*> LWORK is INTEGER
94*> The dimension of the array WORK. LWORK >= max(1,M).
95*> For optimum performance LWORK >= M*NB, where NB is
96*> the optimal blocksize.
97*>
98*> If LWORK = -1, then a workspace query is assumed; the routine
99*> only calculates the optimal size of the WORK array, returns
100*> this value as the first entry of the WORK array, and no error
101*> message related to LWORK is issued by XERBLA.
102*> \endverbatim
103*>
104*> \param[out] INFO
105*> \verbatim
106*> INFO is INTEGER
107*> = 0: successful exit
108*> < 0: if INFO = -i, the i-th argument had an illegal value
109*> \endverbatim
110*
111* Authors:
112* ========
113*
114*> \author Univ. of Tennessee
115*> \author Univ. of California Berkeley
116*> \author Univ. of Colorado Denver
117*> \author NAG Ltd.
118*
119*> \ingroup doubleOTHERcomputational
120*
121*> \par Contributors:
122* ==================
123*>
124*> A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA
125*
126*> \par Further Details:
127* =====================
128*>
129*> \verbatim
130*>
131*> The N-by-N matrix Z can be computed by
132*>
133*> Z = Z(1)*Z(2)* ... *Z(M)
134*>
135*> where each N-by-N Z(k) is given by
136*>
137*> Z(k) = I - tau(k)*v(k)*v(k)**T
138*>
139*> with v(k) is the kth row vector of the M-by-N matrix
140*>
141*> V = ( I A(:,M+1:N) )
142*>
143*> I is the M-by-M identity matrix, A(:,M+1:N)
144*> is the output stored in A on exit from DTZRZF,
145*> and tau(k) is the kth element of the array TAU.
146*>
147*> \endverbatim
148*>
149* =====================================================================
150 SUBROUTINE dtzrzf( M, N, A, LDA, TAU, WORK, LWORK, INFO )
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 DOUBLE PRECISION A( LDA, * ), TAU( * ), WORK( * )
161* ..
162*
163* =====================================================================
164*
165* .. Parameters ..
166 DOUBLE PRECISION ZERO
167 parameter( zero = 0.0d+0 )
168* ..
169* .. Local Scalars ..
170 LOGICAL LQUERY
171 INTEGER I, IB, IWS, KI, KK, LDWORK, LWKMIN, LWKOPT,
172 \$ M1, MU, NB, NBMIN, NX
173* ..
174* .. External Subroutines ..
175 EXTERNAL xerbla, dlarzb, dlarzt, dlatrz
176* ..
177* .. Intrinsic Functions ..
178 INTRINSIC max, min
179* ..
180* .. External Functions ..
181 INTEGER ILAENV
182 EXTERNAL ilaenv
183* ..
184* .. Executable Statements ..
185*
186* Test the input arguments
187*
188 info = 0
189 lquery = ( lwork.EQ.-1 )
190 IF( m.LT.0 ) THEN
191 info = -1
192 ELSE IF( n.LT.m ) THEN
193 info = -2
194 ELSE IF( lda.LT.max( 1, m ) ) THEN
195 info = -4
196 END IF
197*
198 IF( info.EQ.0 ) THEN
199 IF( m.EQ.0 .OR. m.EQ.n ) THEN
200 lwkopt = 1
201 lwkmin = 1
202 ELSE
203*
204* Determine the block size.
205*
206 nb = ilaenv( 1, 'DGERQF', ' ', m, n, -1, -1 )
207 lwkopt = m*nb
208 lwkmin = max( 1, m )
209 END IF
210 work( 1 ) = lwkopt
211*
212 IF( lwork.LT.lwkmin .AND. .NOT.lquery ) THEN
213 info = -7
214 END IF
215 END IF
216*
217 IF( info.NE.0 ) THEN
218 CALL xerbla( 'DTZRZF', -info )
219 RETURN
220 ELSE IF( lquery ) THEN
221 RETURN
222 END IF
223*
224* Quick return if possible
225*
226 IF( m.EQ.0 ) THEN
227 RETURN
228 ELSE IF( m.EQ.n ) THEN
229 DO 10 i = 1, n
230 tau( i ) = zero
231 10 CONTINUE
232 RETURN
233 END IF
234*
235 nbmin = 2
236 nx = 1
237 iws = m
238 IF( nb.GT.1 .AND. nb.LT.m ) THEN
239*
240* Determine when to cross over from blocked to unblocked code.
241*
242 nx = max( 0, ilaenv( 3, 'DGERQF', ' ', m, n, -1, -1 ) )
243 IF( nx.LT.m ) THEN
244*
245* Determine if workspace is large enough for blocked code.
246*
247 ldwork = m
248 iws = ldwork*nb
249 IF( lwork.LT.iws ) THEN
250*
251* Not enough workspace to use optimal NB: reduce NB and
252* determine the minimum value of NB.
253*
254 nb = lwork / ldwork
255 nbmin = max( 2, ilaenv( 2, 'DGERQF', ' ', m, n, -1,
256 \$ -1 ) )
257 END IF
258 END IF
259 END IF
260*
261 IF( nb.GE.nbmin .AND. nb.LT.m .AND. nx.LT.m ) THEN
262*
263* Use blocked code initially.
264* The last kk rows are handled by the block method.
265*
266 m1 = min( m+1, n )
267 ki = ( ( m-nx-1 ) / nb )*nb
268 kk = min( m, ki+nb )
269*
270 DO 20 i = m - kk + ki + 1, m - kk + 1, -nb
271 ib = min( m-i+1, nb )
272*
273* Compute the TZ factorization of the current block
274* A(i:i+ib-1,i:n)
275*
276 CALL dlatrz( ib, n-i+1, n-m, a( i, i ), lda, tau( i ),
277 \$ work )
278 IF( i.GT.1 ) THEN
279*
280* Form the triangular factor of the block reflector
281* H = H(i+ib-1) . . . H(i+1) H(i)
282*
283 CALL dlarzt( 'Backward', 'Rowwise', n-m, ib, a( i, m1 ),
284 \$ lda, tau( i ), work, ldwork )
285*
286* Apply H to A(1:i-1,i:n) from the right
287*
288 CALL dlarzb( 'Right', 'No transpose', 'Backward',
289 \$ 'Rowwise', i-1, n-i+1, ib, n-m, a( i, m1 ),
290 \$ lda, work, ldwork, a( 1, i ), lda,
291 \$ work( ib+1 ), ldwork )
292 END IF
293 20 CONTINUE
294 mu = i + nb - 1
295 ELSE
296 mu = m
297 END IF
298*
299* Use unblocked code to factor the last or only block
300*
301 IF( mu.GT.0 )
302 \$ CALL dlatrz( mu, n, n-m, a, lda, tau, work )
303*
304 work( 1 ) = lwkopt
305*
306 RETURN
307*
308* End of DTZRZF
309*
310 END
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
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:185
subroutine dtzrzf(M, N, A, LDA, TAU, WORK, LWORK, INFO)
DTZRZF
Definition: dtzrzf.f:151
subroutine dlatrz(M, N, L, A, LDA, TAU, WORK)
DLATRZ factors an upper trapezoidal matrix by means of orthogonal transformations.
Definition: dlatrz.f:140
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:183