LAPACK 3.11.0
LAPACK: Linear Algebra PACKage
Loading...
Searching...
No Matches
dgeqr.f
Go to the documentation of this file.
1*> \brief \b DGEQR
2*
3* Definition:
4* ===========
5*
6* SUBROUTINE DGEQR( M, N, A, LDA, T, TSIZE, WORK, LWORK,
7* INFO )
8*
9* .. Scalar Arguments ..
10* INTEGER INFO, LDA, M, N, TSIZE, LWORK
11* ..
12* .. Array Arguments ..
13* DOUBLE PRECISION A( LDA, * ), T( * ), WORK( * )
14* ..
15*
16*
17*> \par Purpose:
18* =============
19*>
20*> \verbatim
21*>
22*> DGEQR computes a QR factorization of a real M-by-N matrix A:
23*>
24*> A = Q * ( R ),
25*> ( 0 )
26*>
27*> where:
28*>
29*> Q is a M-by-M orthogonal matrix;
30*> R is an upper-triangular N-by-N matrix;
31*> 0 is a (M-N)-by-N zero matrix, if M > N.
32*>
33*> \endverbatim
34*
35* Arguments:
36* ==========
37*
38*> \param[in] M
39*> \verbatim
40*> M is INTEGER
41*> The number of rows of the matrix A. M >= 0.
42*> \endverbatim
43*>
44*> \param[in] N
45*> \verbatim
46*> N is INTEGER
47*> The number of columns of the matrix A. N >= 0.
48*> \endverbatim
49*>
50*> \param[in,out] A
51*> \verbatim
52*> A is DOUBLE PRECISION array, dimension (LDA,N)
53*> On entry, the M-by-N matrix A.
54*> On exit, the elements on and above the diagonal of the array
55*> contain the min(M,N)-by-N upper trapezoidal matrix R
56*> (R is upper triangular if M >= N);
57*> the elements below the diagonal are used to store part of the
58*> data structure to represent Q.
59*> \endverbatim
60*>
61*> \param[in] LDA
62*> \verbatim
63*> LDA is INTEGER
64*> The leading dimension of the array A. LDA >= max(1,M).
65*> \endverbatim
66*>
67*> \param[out] T
68*> \verbatim
69*> T is DOUBLE PRECISION array, dimension (MAX(5,TSIZE))
70*> On exit, if INFO = 0, T(1) returns optimal (or either minimal
71*> or optimal, if query is assumed) TSIZE. See TSIZE for details.
72*> Remaining T contains part of the data structure used to represent Q.
73*> If one wants to apply or construct Q, then one needs to keep T
74*> (in addition to A) and pass it to further subroutines.
75*> \endverbatim
76*>
77*> \param[in] TSIZE
78*> \verbatim
79*> TSIZE is INTEGER
80*> If TSIZE >= 5, the dimension of the array T.
81*> If TSIZE = -1 or -2, then a workspace query is assumed. The routine
82*> only calculates the sizes of the T and WORK arrays, returns these
83*> values as the first entries of the T and WORK arrays, and no error
84*> message related to T or WORK is issued by XERBLA.
85*> If TSIZE = -1, the routine calculates optimal size of T for the
86*> optimum performance and returns this value in T(1).
87*> If TSIZE = -2, the routine calculates minimal size of T and
88*> returns this value in T(1).
89*> \endverbatim
90*>
91*> \param[out] WORK
92*> \verbatim
93*> (workspace) DOUBLE PRECISION array, dimension (MAX(1,LWORK))
94*> On exit, if INFO = 0, WORK(1) contains optimal (or either minimal
95*> or optimal, if query was assumed) LWORK.
96*> See LWORK for details.
97*> \endverbatim
98*>
99*> \param[in] LWORK
100*> \verbatim
101*> LWORK is INTEGER
102*> The dimension of the array WORK.
103*> If LWORK = -1 or -2, then a workspace query is assumed. The routine
104*> only calculates the sizes of the T and WORK arrays, returns these
105*> values as the first entries of the T and WORK arrays, and no error
106*> message related to T or WORK is issued by XERBLA.
107*> If LWORK = -1, the routine calculates optimal size of WORK for the
108*> optimal performance and returns this value in WORK(1).
109*> If LWORK = -2, the routine calculates minimal size of WORK and
110*> returns this value in WORK(1).
111*> \endverbatim
112*>
113*> \param[out] INFO
114*> \verbatim
115*> INFO is INTEGER
116*> = 0: successful exit
117*> < 0: if INFO = -i, the i-th argument had an illegal value
118*> \endverbatim
119*
120* Authors:
121* ========
122*
123*> \author Univ. of Tennessee
124*> \author Univ. of California Berkeley
125*> \author Univ. of Colorado Denver
126*> \author NAG Ltd.
127*
128*> \par Further Details
129* ====================
130*>
131*> \verbatim
132*>
133*> The goal of the interface is to give maximum freedom to the developers for
134*> creating any QR factorization algorithm they wish. The triangular
135*> (trapezoidal) R has to be stored in the upper part of A. The lower part of A
136*> and the array T can be used to store any relevant information for applying or
137*> constructing the Q factor. The WORK array can safely be discarded after exit.
138*>
139*> Caution: One should not expect the sizes of T and WORK to be the same from one
140*> LAPACK implementation to the other, or even from one execution to the other.
141*> A workspace query (for T and WORK) is needed at each execution. However,
142*> for a given execution, the size of T and WORK are fixed and will not change
143*> from one query to the next.
144*>
145*> \endverbatim
146*>
147*> \par Further Details particular to this LAPACK implementation:
148* ==============================================================
149*>
150*> \verbatim
151*>
152*> These details are particular for this LAPACK implementation. Users should not
153*> take them for granted. These details may change in the future, and are not likely
154*> true for another LAPACK implementation. These details are relevant if one wants
155*> to try to understand the code. They are not part of the interface.
156*>
157*> In this version,
158*>
159*> T(2): row block size (MB)
160*> T(3): column block size (NB)
161*> T(6:TSIZE): data structure needed for Q, computed by
162*> DLATSQR or DGEQRT
163*>
164*> Depending on the matrix dimensions M and N, and row and column
165*> block sizes MB and NB returned by ILAENV, DGEQR will use either
166*> DLATSQR (if the matrix is tall-and-skinny) or DGEQRT to compute
167*> the QR factorization.
168*>
169*> \endverbatim
170*>
171* =====================================================================
172 SUBROUTINE dgeqr( M, N, A, LDA, T, TSIZE, WORK, LWORK,
173 $ INFO )
174*
175* -- LAPACK computational routine --
176* -- LAPACK is a software package provided by Univ. of Tennessee, --
177* -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd. --
178*
179* .. Scalar Arguments ..
180 INTEGER INFO, LDA, M, N, TSIZE, LWORK
181* ..
182* .. Array Arguments ..
183 DOUBLE PRECISION A( LDA, * ), T( * ), WORK( * )
184* ..
185*
186* =====================================================================
187*
188* ..
189* .. Local Scalars ..
190 LOGICAL LQUERY, LMINWS, MINT, MINW
191 INTEGER MB, NB, MINTSZ, NBLCKS
192* ..
193* .. External Functions ..
194 LOGICAL LSAME
195 EXTERNAL lsame
196* ..
197* .. External Subroutines ..
198 EXTERNAL dlatsqr, dgeqrt, xerbla
199* ..
200* .. Intrinsic Functions ..
201 INTRINSIC max, min, mod
202* ..
203* .. External Functions ..
204 INTEGER ILAENV
205 EXTERNAL ilaenv
206* ..
207* .. Executable Statements ..
208*
209* Test the input arguments
210*
211 info = 0
212*
213 lquery = ( tsize.EQ.-1 .OR. tsize.EQ.-2 .OR.
214 $ lwork.EQ.-1 .OR. lwork.EQ.-2 )
215*
216 mint = .false.
217 minw = .false.
218 IF( tsize.EQ.-2 .OR. lwork.EQ.-2 ) THEN
219 IF( tsize.NE.-1 ) mint = .true.
220 IF( lwork.NE.-1 ) minw = .true.
221 END IF
222*
223* Determine the block size
224*
225 IF( min( m, n ).GT.0 ) THEN
226 mb = ilaenv( 1, 'DGEQR ', ' ', m, n, 1, -1 )
227 nb = ilaenv( 1, 'DGEQR ', ' ', m, n, 2, -1 )
228 ELSE
229 mb = m
230 nb = 1
231 END IF
232 IF( mb.GT.m .OR. mb.LE.n ) mb = m
233 IF( nb.GT.min( m, n ) .OR. nb.LT.1 ) nb = 1
234 mintsz = n + 5
235 IF( mb.GT.n .AND. m.GT.n ) THEN
236 IF( mod( m - n, mb - n ).EQ.0 ) THEN
237 nblcks = ( m - n ) / ( mb - n )
238 ELSE
239 nblcks = ( m - n ) / ( mb - n ) + 1
240 END IF
241 ELSE
242 nblcks = 1
243 END IF
244*
245* Determine if the workspace size satisfies minimal size
246*
247 lminws = .false.
248 IF( ( tsize.LT.max( 1, nb*n*nblcks + 5 ) .OR. lwork.LT.nb*n )
249 $ .AND. ( lwork.GE.n ) .AND. ( tsize.GE.mintsz )
250 $ .AND. ( .NOT.lquery ) ) THEN
251 IF( tsize.LT.max( 1, nb*n*nblcks + 5 ) ) THEN
252 lminws = .true.
253 nb = 1
254 mb = m
255 END IF
256 IF( lwork.LT.nb*n ) THEN
257 lminws = .true.
258 nb = 1
259 END IF
260 END IF
261*
262 IF( m.LT.0 ) THEN
263 info = -1
264 ELSE IF( n.LT.0 ) THEN
265 info = -2
266 ELSE IF( lda.LT.max( 1, m ) ) THEN
267 info = -4
268 ELSE IF( tsize.LT.max( 1, nb*n*nblcks + 5 )
269 $ .AND. ( .NOT.lquery ) .AND. ( .NOT.lminws ) ) THEN
270 info = -6
271 ELSE IF( ( lwork.LT.max( 1, n*nb ) ) .AND. ( .NOT.lquery )
272 $ .AND. ( .NOT.lminws ) ) THEN
273 info = -8
274 END IF
275*
276 IF( info.EQ.0 ) THEN
277 IF( mint ) THEN
278 t( 1 ) = mintsz
279 ELSE
280 t( 1 ) = nb*n*nblcks + 5
281 END IF
282 t( 2 ) = mb
283 t( 3 ) = nb
284 IF( minw ) THEN
285 work( 1 ) = max( 1, n )
286 ELSE
287 work( 1 ) = max( 1, nb*n )
288 END IF
289 END IF
290 IF( info.NE.0 ) THEN
291 CALL xerbla( 'DGEQR', -info )
292 RETURN
293 ELSE IF( lquery ) THEN
294 RETURN
295 END IF
296*
297* Quick return if possible
298*
299 IF( min( m, n ).EQ.0 ) THEN
300 RETURN
301 END IF
302*
303* The QR Decomposition
304*
305 IF( ( m.LE.n ) .OR. ( mb.LE.n ) .OR. ( mb.GE.m ) ) THEN
306 CALL dgeqrt( m, n, nb, a, lda, t( 6 ), nb, work, info )
307 ELSE
308 CALL dlatsqr( m, n, mb, nb, a, lda, t( 6 ), nb, work,
309 $ lwork, info )
310 END IF
311*
312 work( 1 ) = max( 1, nb*n )
313*
314 RETURN
315*
316* End of DGEQR
317*
318 END
subroutine dgeqr(M, N, A, LDA, T, TSIZE, WORK, LWORK, INFO)
DGEQR
Definition: dgeqr.f:174
subroutine dlatsqr(M, N, MB, NB, A, LDA, T, LDT, WORK, LWORK, INFO)
DLATSQR
Definition: dlatsqr.f:166
subroutine xerbla(SRNAME, INFO)
XERBLA
Definition: xerbla.f:60
subroutine dgeqrt(M, N, NB, A, LDA, T, LDT, WORK, INFO)
DGEQRT
Definition: dgeqrt.f:141