LAPACK 3.12.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dgeqrfp.f
Go to the documentation of this file.
1*> \brief \b DGEQRFP
2*
3* =========== DOCUMENTATION ===========
4*
5* Online html documentation available at
6* http://www.netlib.org/lapack/explore-html/
7*
8*> \htmlonly
9*> Download DGEQRFP + dependencies
10*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dgeqrfp.f">
11*> [TGZ]</a>
12*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dgeqrfp.f">
13*> [ZIP]</a>
14*> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dgeqrfp.f">
15*> [TXT]</a>
16*> \endhtmlonly
17*
18* Definition:
19* ===========
20*
21* SUBROUTINE DGEQRFP( 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*> DGEQR2P computes a QR factorization of a real M-by-N matrix A:
37*>
38*> A = Q * ( R ),
39*> ( 0 )
40*>
41*> where:
42*>
43*> Q is a M-by-M orthogonal matrix;
44*> R is an upper-triangular N-by-N matrix with nonnegative diagonal
45*> entries;
46*> 0 is a (M-N)-by-N zero matrix, if M > N.
47*>
48*> \endverbatim
49*
50* Arguments:
51* ==========
52*
53*> \param[in] M
54*> \verbatim
55*> M is INTEGER
56*> The number of rows of the matrix A. M >= 0.
57*> \endverbatim
58*>
59*> \param[in] N
60*> \verbatim
61*> N is INTEGER
62*> The number of columns of the matrix A. N >= 0.
63*> \endverbatim
64*>
65*> \param[in,out] A
66*> \verbatim
67*> A is DOUBLE PRECISION array, dimension (LDA,N)
68*> On entry, the M-by-N matrix A.
69*> On exit, the elements on and above the diagonal of the array
70*> contain the min(M,N)-by-N upper trapezoidal matrix R (R is
71*> upper triangular if m >= n). The diagonal entries of R
72*> are nonnegative; the elements below the diagonal,
73*> with the array TAU, represent the orthogonal matrix Q as a
74*> product of min(m,n) elementary reflectors (see Further
75*> Details).
76*> \endverbatim
77*>
78*> \param[in] LDA
79*> \verbatim
80*> LDA is INTEGER
81*> The leading dimension of the array A. LDA >= max(1,M).
82*> \endverbatim
83*>
84*> \param[out] TAU
85*> \verbatim
86*> TAU is DOUBLE PRECISION array, dimension (min(M,N))
87*> The scalar factors of the elementary reflectors (see Further
88*> Details).
89*> \endverbatim
90*>
91*> \param[out] WORK
92*> \verbatim
93*> WORK is DOUBLE PRECISION array, dimension (MAX(1,LWORK))
94*> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
95*> \endverbatim
96*>
97*> \param[in] LWORK
98*> \verbatim
99*> LWORK is INTEGER
100*> The dimension of the array WORK. LWORK >= max(1,N).
101*> For optimum performance LWORK >= N*NB, where NB is
102*> the optimal blocksize.
103*>
104*> If LWORK = -1, then a workspace query is assumed; the routine
105*> only calculates the optimal size of the WORK array, returns
106*> this value as the first entry of the WORK array, and no error
107*> message related to LWORK is issued by XERBLA.
108*> \endverbatim
109*>
110*> \param[out] INFO
111*> \verbatim
112*> INFO is INTEGER
113*> = 0: successful exit
114*> < 0: if INFO = -i, the i-th argument had an illegal value
115*> \endverbatim
116*
117* Authors:
118* ========
119*
120*> \author Univ. of Tennessee
121*> \author Univ. of California Berkeley
122*> \author Univ. of Colorado Denver
123*> \author NAG Ltd.
124*
125*> \ingroup geqrfp
126*
127*> \par Further Details:
128* =====================
129*>
130*> \verbatim
131*>
132*> The matrix Q is represented as a product of elementary reflectors
133*>
134*> Q = H(1) H(2) . . . H(k), where k = min(m,n).
135*>
136*> Each H(i) has the form
137*>
138*> H(i) = I - tau * v * v**T
139*>
140*> where tau is a real scalar, and v is a real vector with
141*> v(1:i-1) = 0 and v(i) = 1; v(i+1:m) is stored on exit in A(i+1:m,i),
142*> and tau in TAU(i).
143*>
144*> See Lapack Working Note 203 for details
145*> \endverbatim
146*>
147* =====================================================================
148 SUBROUTINE dgeqrfp( M, N, A, LDA, TAU, WORK, LWORK, INFO )
149*
150* -- LAPACK computational routine --
151* -- LAPACK is a software package provided by Univ. of Tennessee, --
152* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
153*
154* .. Scalar Arguments ..
155 INTEGER INFO, LDA, LWORK, M, N
156* ..
157* .. Array Arguments ..
158 DOUBLE PRECISION A( LDA, * ), TAU( * ), WORK( * )
159* ..
160*
161* =====================================================================
162*
163* .. Local Scalars ..
164 LOGICAL LQUERY
165 INTEGER I, IB, IINFO, IWS, K, LDWORK, LWKOPT, NB,
166 $ NBMIN, NX
167* ..
168* .. External Subroutines ..
169 EXTERNAL dgeqr2p, dlarfb, dlarft, xerbla
170* ..
171* .. Intrinsic Functions ..
172 INTRINSIC max, min
173* ..
174* .. External Functions ..
175 INTEGER ILAENV
176 EXTERNAL ilaenv
177* ..
178* .. Executable Statements ..
179*
180* Test the input arguments
181*
182 info = 0
183 nb = ilaenv( 1, 'DGEQRF', ' ', m, n, -1, -1 )
184 lwkopt = n*nb
185 work( 1 ) = lwkopt
186 lquery = ( lwork.EQ.-1 )
187 IF( m.LT.0 ) THEN
188 info = -1
189 ELSE IF( n.LT.0 ) THEN
190 info = -2
191 ELSE IF( lda.LT.max( 1, m ) ) THEN
192 info = -4
193 ELSE IF( lwork.LT.max( 1, n ) .AND. .NOT.lquery ) THEN
194 info = -7
195 END IF
196 IF( info.NE.0 ) THEN
197 CALL xerbla( 'DGEQRFP', -info )
198 RETURN
199 ELSE IF( lquery ) THEN
200 RETURN
201 END IF
202*
203* Quick return if possible
204*
205 k = min( m, n )
206 IF( k.EQ.0 ) THEN
207 work( 1 ) = 1
208 RETURN
209 END IF
210*
211 nbmin = 2
212 nx = 0
213 iws = n
214 IF( nb.GT.1 .AND. nb.LT.k ) THEN
215*
216* Determine when to cross over from blocked to unblocked code.
217*
218 nx = max( 0, ilaenv( 3, 'DGEQRF', ' ', m, n, -1, -1 ) )
219 IF( nx.LT.k ) THEN
220*
221* Determine if workspace is large enough for blocked code.
222*
223 ldwork = n
224 iws = ldwork*nb
225 IF( lwork.LT.iws ) THEN
226*
227* Not enough workspace to use optimal NB: reduce NB and
228* determine the minimum value of NB.
229*
230 nb = lwork / ldwork
231 nbmin = max( 2, ilaenv( 2, 'DGEQRF', ' ', m, n, -1,
232 $ -1 ) )
233 END IF
234 END IF
235 END IF
236*
237 IF( nb.GE.nbmin .AND. nb.LT.k .AND. nx.LT.k ) THEN
238*
239* Use blocked code initially
240*
241 DO 10 i = 1, k - nx, nb
242 ib = min( k-i+1, nb )
243*
244* Compute the QR factorization of the current block
245* A(i:m,i:i+ib-1)
246*
247 CALL dgeqr2p( m-i+1, ib, a( i, i ), lda, tau( i ), work,
248 $ iinfo )
249 IF( i+ib.LE.n ) THEN
250*
251* Form the triangular factor of the block reflector
252* H = H(i) H(i+1) . . . H(i+ib-1)
253*
254 CALL dlarft( 'Forward', 'Columnwise', m-i+1, ib,
255 $ a( i, i ), lda, tau( i ), work, ldwork )
256*
257* Apply H**T to A(i:m,i+ib:n) from the left
258*
259 CALL dlarfb( 'Left', 'Transpose', 'Forward',
260 $ 'Columnwise', m-i+1, n-i-ib+1, ib,
261 $ a( i, i ), lda, work, ldwork, a( i, i+ib ),
262 $ lda, work( ib+1 ), ldwork )
263 END IF
264 10 CONTINUE
265 ELSE
266 i = 1
267 END IF
268*
269* Use unblocked code to factor the last or only block.
270*
271 IF( i.LE.k )
272 $ CALL dgeqr2p( m-i+1, n-i+1, a( i, i ), lda, tau( i ), work,
273 $ iinfo )
274*
275 work( 1 ) = iws
276 RETURN
277*
278* End of DGEQRFP
279*
280 END
subroutine xerbla(srname, info)
Definition cblat2.f:3285
subroutine dgeqr2p(m, n, a, lda, tau, work, info)
DGEQR2P computes the QR factorization of a general rectangular matrix with non-negative diagonal elem...
Definition dgeqr2p.f:134
subroutine dgeqrfp(m, n, a, lda, tau, work, lwork, info)
DGEQRFP
Definition dgeqrfp.f:149
subroutine dlarfb(side, trans, direct, storev, m, n, k, v, ldv, t, ldt, c, ldc, work, ldwork)
DLARFB applies a block reflector or its transpose to a general rectangular matrix.
Definition dlarfb.f:197
subroutine dlarft(direct, storev, n, k, v, ldv, tau, t, ldt)
DLARFT forms the triangular factor T of a block reflector H = I - vtvH
Definition dlarft.f:163