1: *> \brief \b DORGTSQR
2: *
3: * =========== DOCUMENTATION ===========
4: *
5: * Online html documentation available at
6: * http://www.netlib.org/lapack/explore-html/
7: *
8: *> \htmlonly
9: *> Download DORGTSQR + dependencies
10: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.tgz?format=tgz&filename=/lapack/lapack_routine/dorgtsqr.f">
11: *> [TGZ]</a>
12: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.zip?format=zip&filename=/lapack/lapack_routine/dorgtsqr.f">
13: *> [ZIP]</a>
14: *> <a href="http://www.netlib.org/cgi-bin/netlibfiles.txt?format=txt&filename=/lapack/lapack_routine/dorgtsqr.f">
15: *> [TXT]</a>
16: *>
17: * Definition:
18: * ===========
19: *
20: * SUBROUTINE DORGTSQR( M, N, MB, NB, A, LDA, T, LDT, WORK, LWORK,
21: * $ INFO )
22: *
23: * .. Scalar Arguments ..
24: * INTEGER INFO, LDA, LDT, LWORK, M, N, MB, NB
25: * ..
26: * .. Array Arguments ..
27: * DOUBLE PRECISION A( LDA, * ), T( LDT, * ), WORK( * )
28: * ..
29: *
30: *> \par Purpose:
31: * =============
32: *>
33: *> \verbatim
34: *>
35: *> DORGTSQR generates an M-by-N real matrix Q_out with orthonormal columns,
36: *> which are the first N columns of a product of real orthogonal
37: *> matrices of order M which are returned by DLATSQR
38: *>
39: *> Q_out = first_N_columns_of( Q(1)_in * Q(2)_in * ... * Q(k)_in ).
40: *>
41: *> See the documentation for DLATSQR.
42: *> \endverbatim
43: *
44: * Arguments:
45: * ==========
46: *
47: *> \param[in] M
48: *> \verbatim
49: *> M is INTEGER
50: *> The number of rows of the matrix A. M >= 0.
51: *> \endverbatim
52: *>
53: *> \param[in] N
54: *> \verbatim
55: *> N is INTEGER
56: *> The number of columns of the matrix A. M >= N >= 0.
57: *> \endverbatim
58: *>
59: *> \param[in] MB
60: *> \verbatim
61: *> MB is INTEGER
62: *> The row block size used by DLATSQR to return
63: *> arrays A and T. MB > N.
64: *> (Note that if MB > M, then M is used instead of MB
65: *> as the row block size).
66: *> \endverbatim
67: *>
68: *> \param[in] NB
69: *> \verbatim
70: *> NB is INTEGER
71: *> The column block size used by DLATSQR to return
72: *> arrays A and T. NB >= 1.
73: *> (Note that if NB > N, then N is used instead of NB
74: *> as the column block size).
75: *> \endverbatim
76: *>
77: *> \param[in,out] A
78: *> \verbatim
79: *> A is DOUBLE PRECISION array, dimension (LDA,N)
80: *>
81: *> On entry:
82: *>
83: *> The elements on and above the diagonal are not accessed.
84: *> The elements below the diagonal represent the unit
85: *> lower-trapezoidal blocked matrix V computed by DLATSQR
86: *> that defines the input matrices Q_in(k) (ones on the
87: *> diagonal are not stored) (same format as the output A
88: *> below the diagonal in DLATSQR).
89: *>
90: *> On exit:
91: *>
92: *> The array A contains an M-by-N orthonormal matrix Q_out,
93: *> i.e the columns of A are orthogonal unit vectors.
94: *> \endverbatim
95: *>
96: *> \param[in] LDA
97: *> \verbatim
98: *> LDA is INTEGER
99: *> The leading dimension of the array A. LDA >= max(1,M).
100: *> \endverbatim
101: *>
102: *> \param[in] T
103: *> \verbatim
104: *> T is DOUBLE PRECISION array,
105: *> dimension (LDT, N * NIRB)
106: *> where NIRB = Number_of_input_row_blocks
107: *> = MAX( 1, CEIL((M-N)/(MB-N)) )
108: *> Let NICB = Number_of_input_col_blocks
109: *> = CEIL(N/NB)
110: *>
111: *> The upper-triangular block reflectors used to define the
112: *> input matrices Q_in(k), k=(1:NIRB*NICB). The block
113: *> reflectors are stored in compact form in NIRB block
114: *> reflector sequences. Each of NIRB block reflector sequences
115: *> is stored in a larger NB-by-N column block of T and consists
116: *> of NICB smaller NB-by-NB upper-triangular column blocks.
117: *> (same format as the output T in DLATSQR).
118: *> \endverbatim
119: *>
120: *> \param[in] LDT
121: *> \verbatim
122: *> LDT is INTEGER
123: *> The leading dimension of the array T.
124: *> LDT >= max(1,min(NB1,N)).
125: *> \endverbatim
126: *>
127: *> \param[out] WORK
128: *> \verbatim
129: *> (workspace) DOUBLE PRECISION array, dimension (MAX(2,LWORK))
130: *> On exit, if INFO = 0, WORK(1) returns the optimal LWORK.
131: *> \endverbatim
132: *>
133: *> \param[in] LWORK
134: *> \verbatim
135: *> The dimension of the array WORK. LWORK >= (M+NB)*N.
136: *> If LWORK = -1, then a workspace query is assumed.
137: *> The routine only calculates the optimal size of the WORK
138: *> array, returns this value as the first entry of the WORK
139: *> array, and no error message related to LWORK is issued
140: *> by XERBLA.
141: *> \endverbatim
142: *>
143: *> \param[out] INFO
144: *> \verbatim
145: *> INFO is INTEGER
146: *> = 0: successful exit
147: *> < 0: if INFO = -i, the i-th argument had an illegal value
148: *> \endverbatim
149: *>
150: * Authors:
151: * ========
152: *
153: *> \author Univ. of Tennessee
154: *> \author Univ. of California Berkeley
155: *> \author Univ. of Colorado Denver
156: *> \author NAG Ltd.
157: *
158: *> \date November 2019
159: *
160: *> \ingroup doubleOTHERcomputational
161: *
162: *> \par Contributors:
163: * ==================
164: *>
165: *> \verbatim
166: *>
167: *> November 2019, Igor Kozachenko,
168: *> Computer Science Division,
169: *> University of California, Berkeley
170: *>
171: *> \endverbatim
172: *
173: * =====================================================================
174: SUBROUTINE DORGTSQR( M, N, MB, NB, A, LDA, T, LDT, WORK, LWORK,
175: $ INFO )
176: IMPLICIT NONE
177: *
178: * -- LAPACK computational routine (version 3.9.0) --
179: * -- LAPACK is a software package provided by Univ. of Tennessee, --
180: * -- Univ. of California Berkeley, Univ. of Colorado Denver and NAG Ltd..--
181: * November 2019
182: *
183: * .. Scalar Arguments ..
184: INTEGER INFO, LDA, LDT, LWORK, M, N, MB, NB
185: * ..
186: * .. Array Arguments ..
187: DOUBLE PRECISION A( LDA, * ), T( LDT, * ), WORK( * )
188: * ..
189: *
190: * =====================================================================
191: *
192: * .. Parameters ..
193: DOUBLE PRECISION ONE, ZERO
194: PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
195: * ..
196: * .. Local Scalars ..
197: LOGICAL LQUERY
198: INTEGER IINFO, LDC, LWORKOPT, LC, LW, NBLOCAL, J
199: * ..
200: * .. External Subroutines ..
201: EXTERNAL DCOPY, DLAMTSQR, DLASET, XERBLA
202: * ..
203: * .. Intrinsic Functions ..
204: INTRINSIC DBLE, MAX, MIN
205: * ..
206: * .. Executable Statements ..
207: *
208: * Test the input parameters
209: *
210: LQUERY = LWORK.EQ.-1
211: INFO = 0
212: IF( M.LT.0 ) THEN
213: INFO = -1
214: ELSE IF( N.LT.0 .OR. M.LT.N ) THEN
215: INFO = -2
216: ELSE IF( MB.LE.N ) THEN
217: INFO = -3
218: ELSE IF( NB.LT.1 ) THEN
219: INFO = -4
220: ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
221: INFO = -6
222: ELSE IF( LDT.LT.MAX( 1, MIN( NB, N ) ) ) THEN
223: INFO = -8
224: ELSE
225: *
226: * Test the input LWORK for the dimension of the array WORK.
227: * This workspace is used to store array C(LDC, N) and WORK(LWORK)
228: * in the call to DLAMTSQR. See the documentation for DLAMTSQR.
229: *
230: IF( LWORK.LT.2 .AND. (.NOT.LQUERY) ) THEN
231: INFO = -10
232: ELSE
233: *
234: * Set block size for column blocks
235: *
236: NBLOCAL = MIN( NB, N )
237: *
238: * LWORK = -1, then set the size for the array C(LDC,N)
239: * in DLAMTSQR call and set the optimal size of the work array
240: * WORK(LWORK) in DLAMTSQR call.
241: *
242: LDC = M
243: LC = LDC*N
244: LW = N * NBLOCAL
245: *
246: LWORKOPT = LC+LW
247: *
248: IF( ( LWORK.LT.MAX( 1, LWORKOPT ) ).AND.(.NOT.LQUERY) ) THEN
249: INFO = -10
250: END IF
251: END IF
252: *
253: END IF
254: *
255: * Handle error in the input parameters and return workspace query.
256: *
257: IF( INFO.NE.0 ) THEN
258: CALL XERBLA( 'DORGTSQR', -INFO )
259: RETURN
260: ELSE IF ( LQUERY ) THEN
261: WORK( 1 ) = DBLE( LWORKOPT )
262: RETURN
263: END IF
264: *
265: * Quick return if possible
266: *
267: IF( MIN( M, N ).EQ.0 ) THEN
268: WORK( 1 ) = DBLE( LWORKOPT )
269: RETURN
270: END IF
271: *
272: * (1) Form explicitly the tall-skinny M-by-N left submatrix Q1_in
273: * of M-by-M orthogonal matrix Q_in, which is implicitly stored in
274: * the subdiagonal part of input array A and in the input array T.
275: * Perform by the following operation using the routine DLAMTSQR.
276: *
277: * Q1_in = Q_in * ( I ), where I is a N-by-N identity matrix,
278: * ( 0 ) 0 is a (M-N)-by-N zero matrix.
279: *
280: * (1a) Form M-by-N matrix in the array WORK(1:LDC*N) with ones
281: * on the diagonal and zeros elsewhere.
282: *
283: CALL DLASET( 'F', M, N, ZERO, ONE, WORK, LDC )
284: *
285: * (1b) On input, WORK(1:LDC*N) stores ( I );
286: * ( 0 )
287: *
288: * On output, WORK(1:LDC*N) stores Q1_in.
289: *
290: CALL DLAMTSQR( 'L', 'N', M, N, N, MB, NBLOCAL, A, LDA, T, LDT,
291: $ WORK, LDC, WORK( LC+1 ), LW, IINFO )
292: *
293: * (2) Copy the result from the part of the work array (1:M,1:N)
294: * with the leading dimension LDC that starts at WORK(1) into
295: * the output array A(1:M,1:N) column-by-column.
296: *
297: DO J = 1, N
298: CALL DCOPY( M, WORK( (J-1)*LDC + 1 ), 1, A( 1, J ), 1 )
299: END DO
300: *
301: WORK( 1 ) = DBLE( LWORKOPT )
302: RETURN
303: *
304: * End of DORGTSQR
305: *
306: END
CVSweb interface <joel.bertrand@systella.fr>