1 /* ./src_f77/sormbr.f -- translated by f2c (version 20030320).
2 You must link the resulting object file with the libraries:
3 -lf2c -lm (in that order)
4 */
5
6 #include <punc/vf2c.h>
7
8 /* Table of constant values */
9
10 static integer c__1 = 1;
11 static integer c_n1 = -1;
12 static integer c__2 = 2;
13
sormbr_(char * vect,char * side,char * trans,integer * m,integer * n,integer * k,real * a,integer * lda,real * tau,real * c__,integer * ldc,real * work,integer * lwork,integer * info,ftnlen vect_len,ftnlen side_len,ftnlen trans_len)14 /* Subroutine */ int sormbr_(char *vect, char *side, char *trans, integer *m,
15 integer *n, integer *k, real *a, integer *lda, real *tau, real *c__,
16 integer *ldc, real *work, integer *lwork, integer *info, ftnlen
17 vect_len, ftnlen side_len, ftnlen trans_len)
18 {
19 /* System generated locals */
20 address a__1[2];
21 integer a_dim1, a_offset, c_dim1, c_offset, i__1, i__2, i__3[2];
22 char ch__1[2];
23
24 /* Builtin functions */
25 /* Subroutine */ int s_cat(char *, char **, integer *, integer *, ftnlen);
26
27 /* Local variables */
28 static integer i1, i2, nb, mi, ni, nq, nw;
29 static logical left;
30 extern logical lsame_(char *, char *, ftnlen, ftnlen);
31 static integer iinfo;
32 extern /* Subroutine */ int xerbla_(char *, integer *, ftnlen);
33 extern integer ilaenv_(integer *, char *, char *, integer *, integer *,
34 integer *, integer *, ftnlen, ftnlen);
35 static logical notran, applyq;
36 static char transt[1];
37 extern /* Subroutine */ int sormlq_(char *, char *, integer *, integer *,
38 integer *, real *, integer *, real *, real *, integer *, real *,
39 integer *, integer *, ftnlen, ftnlen);
40 static integer lwkopt;
41 static logical lquery;
42 extern /* Subroutine */ int sormqr_(char *, char *, integer *, integer *,
43 integer *, real *, integer *, real *, real *, integer *, real *,
44 integer *, integer *, ftnlen, ftnlen);
45
46
47 /* -- LAPACK routine (version 3.0) -- */
48 /* Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., */
49 /* Courant Institute, Argonne National Lab, and Rice University */
50 /* June 30, 1999 */
51
52 /* .. Scalar Arguments .. */
53 /* .. */
54 /* .. Array Arguments .. */
55 /* .. */
56
57 /* Purpose */
58 /* ======= */
59
60 /* If VECT = 'Q', SORMBR overwrites the general real M-by-N matrix C */
61 /* with */
62 /* SIDE = 'L' SIDE = 'R' */
63 /* TRANS = 'N': Q * C C * Q */
64 /* TRANS = 'T': Q**T * C C * Q**T */
65
66 /* If VECT = 'P', SORMBR overwrites the general real M-by-N matrix C */
67 /* with */
68 /* SIDE = 'L' SIDE = 'R' */
69 /* TRANS = 'N': P * C C * P */
70 /* TRANS = 'T': P**T * C C * P**T */
71
72 /* Here Q and P**T are the orthogonal matrices determined by SGEBRD when */
73 /* reducing a real matrix A to bidiagonal form: A = Q * B * P**T. Q and */
74 /* P**T are defined as products of elementary reflectors H(i) and G(i) */
75 /* respectively. */
76
77 /* Let nq = m if SIDE = 'L' and nq = n if SIDE = 'R'. Thus nq is the */
78 /* order of the orthogonal matrix Q or P**T that is applied. */
79
80 /* If VECT = 'Q', A is assumed to have been an NQ-by-K matrix: */
81 /* if nq >= k, Q = H(1) H(2) . . . H(k); */
82 /* if nq < k, Q = H(1) H(2) . . . H(nq-1). */
83
84 /* If VECT = 'P', A is assumed to have been a K-by-NQ matrix: */
85 /* if k < nq, P = G(1) G(2) . . . G(k); */
86 /* if k >= nq, P = G(1) G(2) . . . G(nq-1). */
87
88 /* Arguments */
89 /* ========= */
90
91 /* VECT (input) CHARACTER*1 */
92 /* = 'Q': apply Q or Q**T; */
93 /* = 'P': apply P or P**T. */
94
95 /* SIDE (input) CHARACTER*1 */
96 /* = 'L': apply Q, Q**T, P or P**T from the Left; */
97 /* = 'R': apply Q, Q**T, P or P**T from the Right. */
98
99 /* TRANS (input) CHARACTER*1 */
100 /* = 'N': No transpose, apply Q or P; */
101 /* = 'T': Transpose, apply Q**T or P**T. */
102
103 /* M (input) INTEGER */
104 /* The number of rows of the matrix C. M >= 0. */
105
106 /* N (input) INTEGER */
107 /* The number of columns of the matrix C. N >= 0. */
108
109 /* K (input) INTEGER */
110 /* If VECT = 'Q', the number of columns in the original */
111 /* matrix reduced by SGEBRD. */
112 /* If VECT = 'P', the number of rows in the original */
113 /* matrix reduced by SGEBRD. */
114 /* K >= 0. */
115
116 /* A (input) REAL array, dimension */
117 /* (LDA,min(nq,K)) if VECT = 'Q' */
118 /* (LDA,nq) if VECT = 'P' */
119 /* The vectors which define the elementary reflectors H(i) and */
120 /* G(i), whose products determine the matrices Q and P, as */
121 /* returned by SGEBRD. */
122
123 /* LDA (input) INTEGER */
124 /* The leading dimension of the array A. */
125 /* If VECT = 'Q', LDA >= max(1,nq); */
126 /* if VECT = 'P', LDA >= max(1,min(nq,K)). */
127
128 /* TAU (input) REAL array, dimension (min(nq,K)) */
129 /* TAU(i) must contain the scalar factor of the elementary */
130 /* reflector H(i) or G(i) which determines Q or P, as returned */
131 /* by SGEBRD in the array argument TAUQ or TAUP. */
132
133 /* C (input/output) REAL array, dimension (LDC,N) */
134 /* On entry, the M-by-N matrix C. */
135 /* On exit, C is overwritten by Q*C or Q**T*C or C*Q**T or C*Q */
136 /* or P*C or P**T*C or C*P or C*P**T. */
137
138 /* LDC (input) INTEGER */
139 /* The leading dimension of the array C. LDC >= max(1,M). */
140
141 /* WORK (workspace/output) REAL array, dimension (LWORK) */
142 /* On exit, if INFO = 0, WORK(1) returns the optimal LWORK. */
143
144 /* LWORK (input) INTEGER */
145 /* The dimension of the array WORK. */
146 /* If SIDE = 'L', LWORK >= max(1,N); */
147 /* if SIDE = 'R', LWORK >= max(1,M). */
148 /* For optimum performance LWORK >= N*NB if SIDE = 'L', and */
149 /* LWORK >= M*NB if SIDE = 'R', where NB is the optimal */
150 /* blocksize. */
151
152 /* If LWORK = -1, then a workspace query is assumed; the routine */
153 /* only calculates the optimal size of the WORK array, returns */
154 /* this value as the first entry of the WORK array, and no error */
155 /* message related to LWORK is issued by XERBLA. */
156
157 /* INFO (output) INTEGER */
158 /* = 0: successful exit */
159 /* < 0: if INFO = -i, the i-th argument had an illegal value */
160
161 /* ===================================================================== */
162
163 /* .. Local Scalars .. */
164 /* .. */
165 /* .. External Functions .. */
166 /* .. */
167 /* .. External Subroutines .. */
168 /* .. */
169 /* .. Intrinsic Functions .. */
170 /* .. */
171 /* .. Executable Statements .. */
172
173 /* Test the input arguments */
174
175 /* Parameter adjustments */
176 a_dim1 = *lda;
177 a_offset = 1 + a_dim1;
178 a -= a_offset;
179 --tau;
180 c_dim1 = *ldc;
181 c_offset = 1 + c_dim1;
182 c__ -= c_offset;
183 --work;
184
185 /* Function Body */
186 *info = 0;
187 applyq = lsame_(vect, "Q", (ftnlen)1, (ftnlen)1);
188 left = lsame_(side, "L", (ftnlen)1, (ftnlen)1);
189 notran = lsame_(trans, "N", (ftnlen)1, (ftnlen)1);
190 lquery = *lwork == -1;
191
192 /* NQ is the order of Q or P and NW is the minimum dimension of WORK */
193
194 if (left) {
195 nq = *m;
196 nw = *n;
197 } else {
198 nq = *n;
199 nw = *m;
200 }
201 if (! applyq && ! lsame_(vect, "P", (ftnlen)1, (ftnlen)1)) {
202 *info = -1;
203 } else if (! left && ! lsame_(side, "R", (ftnlen)1, (ftnlen)1)) {
204 *info = -2;
205 } else if (! notran && ! lsame_(trans, "T", (ftnlen)1, (ftnlen)1)) {
206 *info = -3;
207 } else if (*m < 0) {
208 *info = -4;
209 } else if (*n < 0) {
210 *info = -5;
211 } else if (*k < 0) {
212 *info = -6;
213 } else /* if(complicated condition) */ {
214 /* Computing MAX */
215 i__1 = 1, i__2 = min(nq,*k);
216 if (applyq && *lda < max(1,nq) || ! applyq && *lda < max(i__1,i__2)) {
217 *info = -8;
218 } else if (*ldc < max(1,*m)) {
219 *info = -11;
220 } else if (*lwork < max(1,nw) && ! lquery) {
221 *info = -13;
222 }
223 }
224
225 if (*info == 0) {
226 if (applyq) {
227 if (left) {
228 /* Writing concatenation */
229 i__3[0] = 1, a__1[0] = side;
230 i__3[1] = 1, a__1[1] = trans;
231 s_cat(ch__1, a__1, i__3, &c__2, (ftnlen)2);
232 i__1 = *m - 1;
233 i__2 = *m - 1;
234 nb = ilaenv_(&c__1, "SORMQR", ch__1, &i__1, n, &i__2, &c_n1, (
235 ftnlen)6, (ftnlen)2);
236 } else {
237 /* Writing concatenation */
238 i__3[0] = 1, a__1[0] = side;
239 i__3[1] = 1, a__1[1] = trans;
240 s_cat(ch__1, a__1, i__3, &c__2, (ftnlen)2);
241 i__1 = *n - 1;
242 i__2 = *n - 1;
243 nb = ilaenv_(&c__1, "SORMQR", ch__1, m, &i__1, &i__2, &c_n1, (
244 ftnlen)6, (ftnlen)2);
245 }
246 } else {
247 if (left) {
248 /* Writing concatenation */
249 i__3[0] = 1, a__1[0] = side;
250 i__3[1] = 1, a__1[1] = trans;
251 s_cat(ch__1, a__1, i__3, &c__2, (ftnlen)2);
252 i__1 = *m - 1;
253 i__2 = *m - 1;
254 nb = ilaenv_(&c__1, "SORMLQ", ch__1, &i__1, n, &i__2, &c_n1, (
255 ftnlen)6, (ftnlen)2);
256 } else {
257 /* Writing concatenation */
258 i__3[0] = 1, a__1[0] = side;
259 i__3[1] = 1, a__1[1] = trans;
260 s_cat(ch__1, a__1, i__3, &c__2, (ftnlen)2);
261 i__1 = *n - 1;
262 i__2 = *n - 1;
263 nb = ilaenv_(&c__1, "SORMLQ", ch__1, m, &i__1, &i__2, &c_n1, (
264 ftnlen)6, (ftnlen)2);
265 }
266 }
267 lwkopt = max(1,nw) * nb;
268 work[1] = (real) lwkopt;
269 }
270
271 if (*info != 0) {
272 i__1 = -(*info);
273 xerbla_("SORMBR", &i__1, (ftnlen)6);
274 return 0;
275 } else if (lquery) {
276 return 0;
277 }
278
279 /* Quick return if possible */
280
281 work[1] = 1.f;
282 if (*m == 0 || *n == 0) {
283 return 0;
284 }
285
286 if (applyq) {
287
288 /* Apply Q */
289
290 if (nq >= *k) {
291
292 /* Q was determined by a call to SGEBRD with nq >= k */
293
294 sormqr_(side, trans, m, n, k, &a[a_offset], lda, &tau[1], &c__[
295 c_offset], ldc, &work[1], lwork, &iinfo, (ftnlen)1, (
296 ftnlen)1);
297 } else if (nq > 1) {
298
299 /* Q was determined by a call to SGEBRD with nq < k */
300
301 if (left) {
302 mi = *m - 1;
303 ni = *n;
304 i1 = 2;
305 i2 = 1;
306 } else {
307 mi = *m;
308 ni = *n - 1;
309 i1 = 1;
310 i2 = 2;
311 }
312 i__1 = nq - 1;
313 sormqr_(side, trans, &mi, &ni, &i__1, &a[a_dim1 + 2], lda, &tau[1]
314 , &c__[i1 + i2 * c_dim1], ldc, &work[1], lwork, &iinfo, (
315 ftnlen)1, (ftnlen)1);
316 }
317 } else {
318
319 /* Apply P */
320
321 if (notran) {
322 *(unsigned char *)transt = 'T';
323 } else {
324 *(unsigned char *)transt = 'N';
325 }
326 if (nq > *k) {
327
328 /* P was determined by a call to SGEBRD with nq > k */
329
330 sormlq_(side, transt, m, n, k, &a[a_offset], lda, &tau[1], &c__[
331 c_offset], ldc, &work[1], lwork, &iinfo, (ftnlen)1, (
332 ftnlen)1);
333 } else if (nq > 1) {
334
335 /* P was determined by a call to SGEBRD with nq <= k */
336
337 if (left) {
338 mi = *m - 1;
339 ni = *n;
340 i1 = 2;
341 i2 = 1;
342 } else {
343 mi = *m;
344 ni = *n - 1;
345 i1 = 1;
346 i2 = 2;
347 }
348 i__1 = nq - 1;
349 sormlq_(side, transt, &mi, &ni, &i__1, &a[(a_dim1 << 1) + 1], lda,
350 &tau[1], &c__[i1 + i2 * c_dim1], ldc, &work[1], lwork, &
351 iinfo, (ftnlen)1, (ftnlen)1);
352 }
353 }
354 work[1] = (real) lwkopt;
355 return 0;
356
357 /* End of SORMBR */
358
359 } /* sormbr_ */
360
361