1 /* ./src_f77/slarzt.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 real c_b8 = 0.f;
11 static integer c__1 = 1;
12 
slarzt_(char * direct,char * storev,integer * n,integer * k,real * v,integer * ldv,real * tau,real * t,integer * ldt,ftnlen direct_len,ftnlen storev_len)13 /* Subroutine */ int slarzt_(char *direct, char *storev, integer *n, integer *
14 	k, real *v, integer *ldv, real *tau, real *t, integer *ldt, ftnlen
15 	direct_len, ftnlen storev_len)
16 {
17     /* System generated locals */
18     integer t_dim1, t_offset, v_dim1, v_offset, i__1;
19     real r__1;
20 
21     /* Local variables */
22     static integer i__, j, info;
23     extern logical lsame_(char *, char *, ftnlen, ftnlen);
24     extern /* Subroutine */ int sgemv_(char *, integer *, integer *, real *,
25 	    real *, integer *, real *, integer *, real *, real *, integer *,
26 	    ftnlen), strmv_(char *, char *, char *, integer *, real *,
27 	    integer *, real *, integer *, ftnlen, ftnlen, ftnlen), xerbla_(
28 	    char *, integer *, ftnlen);
29 
30 
31 /*  -- LAPACK routine (version 3.0) -- */
32 /*     Univ. of Tennessee, Univ. of California Berkeley, NAG Ltd., */
33 /*     Courant Institute, Argonne National Lab, and Rice University */
34 /*     June 30, 1999 */
35 
36 /*     .. Scalar Arguments .. */
37 /*     .. */
38 /*     .. Array Arguments .. */
39 /*     .. */
40 
41 /*  Purpose */
42 /*  ======= */
43 
44 /*  SLARZT forms the triangular factor T of a real block reflector */
45 /*  H of order > n, which is defined as a product of k elementary */
46 /*  reflectors. */
47 
48 /*  If DIRECT = 'F', H = H(1) H(2) . . . H(k) and T is upper triangular; */
49 
50 /*  If DIRECT = 'B', H = H(k) . . . H(2) H(1) and T is lower triangular. */
51 
52 /*  If STOREV = 'C', the vector which defines the elementary reflector */
53 /*  H(i) is stored in the i-th column of the array V, and */
54 
55 /*     H  =  I - V * T * V' */
56 
57 /*  If STOREV = 'R', the vector which defines the elementary reflector */
58 /*  H(i) is stored in the i-th row of the array V, and */
59 
60 /*     H  =  I - V' * T * V */
61 
62 /*  Currently, only STOREV = 'R' and DIRECT = 'B' are supported. */
63 
64 /*  Arguments */
65 /*  ========= */
66 
67 /*  DIRECT  (input) CHARACTER*1 */
68 /*          Specifies the order in which the elementary reflectors are */
69 /*          multiplied to form the block reflector: */
70 /*          = 'F': H = H(1) H(2) . . . H(k) (Forward, not supported yet) */
71 /*          = 'B': H = H(k) . . . H(2) H(1) (Backward) */
72 
73 /*  STOREV  (input) CHARACTER*1 */
74 /*          Specifies how the vectors which define the elementary */
75 /*          reflectors are stored (see also Further Details): */
76 /*          = 'C': columnwise                        (not supported yet) */
77 /*          = 'R': rowwise */
78 
79 /*  N       (input) INTEGER */
80 /*          The order of the block reflector H. N >= 0. */
81 
82 /*  K       (input) INTEGER */
83 /*          The order of the triangular factor T (= the number of */
84 /*          elementary reflectors). K >= 1. */
85 
86 /*  V       (input/output) REAL array, dimension */
87 /*                               (LDV,K) if STOREV = 'C' */
88 /*                               (LDV,N) if STOREV = 'R' */
89 /*          The matrix V. See further details. */
90 
91 /*  LDV     (input) INTEGER */
92 /*          The leading dimension of the array V. */
93 /*          If STOREV = 'C', LDV >= max(1,N); if STOREV = 'R', LDV >= K. */
94 
95 /*  TAU     (input) REAL array, dimension (K) */
96 /*          TAU(i) must contain the scalar factor of the elementary */
97 /*          reflector H(i). */
98 
99 /*  T       (output) REAL array, dimension (LDT,K) */
100 /*          The k by k triangular factor T of the block reflector. */
101 /*          If DIRECT = 'F', T is upper triangular; if DIRECT = 'B', T is */
102 /*          lower triangular. The rest of the array is not used. */
103 
104 /*  LDT     (input) INTEGER */
105 /*          The leading dimension of the array T. LDT >= K. */
106 
107 /*  Further Details */
108 /*  =============== */
109 
110 /*  Based on contributions by */
111 /*    A. Petitet, Computer Science Dept., Univ. of Tenn., Knoxville, USA */
112 
113 /*  The shape of the matrix V and the storage of the vectors which define */
114 /*  the H(i) is best illustrated by the following example with n = 5 and */
115 /*  k = 3. The elements equal to 1 are not stored; the corresponding */
116 /*  array elements are modified but restored on exit. The rest of the */
117 /*  array is not used. */
118 
119 /*  DIRECT = 'F' and STOREV = 'C':         DIRECT = 'F' and STOREV = 'R': */
120 
121 /*                                              ______V_____ */
122 /*         ( v1 v2 v3 )                        /            \ */
123 /*         ( v1 v2 v3 )                      ( v1 v1 v1 v1 v1 . . . . 1 ) */
124 /*     V = ( v1 v2 v3 )                      ( v2 v2 v2 v2 v2 . . . 1   ) */
125 /*         ( v1 v2 v3 )                      ( v3 v3 v3 v3 v3 . . 1     ) */
126 /*         ( v1 v2 v3 ) */
127 /*            .  .  . */
128 /*            .  .  . */
129 /*            1  .  . */
130 /*               1  . */
131 /*                  1 */
132 
133 /*  DIRECT = 'B' and STOREV = 'C':         DIRECT = 'B' and STOREV = 'R': */
134 
135 /*                                                        ______V_____ */
136 /*            1                                          /            \ */
137 /*            .  1                           ( 1 . . . . v1 v1 v1 v1 v1 ) */
138 /*            .  .  1                        ( . 1 . . . v2 v2 v2 v2 v2 ) */
139 /*            .  .  .                        ( . . 1 . . v3 v3 v3 v3 v3 ) */
140 /*            .  .  . */
141 /*         ( v1 v2 v3 ) */
142 /*         ( v1 v2 v3 ) */
143 /*     V = ( v1 v2 v3 ) */
144 /*         ( v1 v2 v3 ) */
145 /*         ( v1 v2 v3 ) */
146 
147 /*  ===================================================================== */
148 
149 /*     .. Parameters .. */
150 /*     .. */
151 /*     .. Local Scalars .. */
152 /*     .. */
153 /*     .. External Subroutines .. */
154 /*     .. */
155 /*     .. External Functions .. */
156 /*     .. */
157 /*     .. Executable Statements .. */
158 
159 /*     Check for currently supported options */
160 
161     /* Parameter adjustments */
162     v_dim1 = *ldv;
163     v_offset = 1 + v_dim1;
164     v -= v_offset;
165     --tau;
166     t_dim1 = *ldt;
167     t_offset = 1 + t_dim1;
168     t -= t_offset;
169 
170     /* Function Body */
171     info = 0;
172     if (! lsame_(direct, "B", (ftnlen)1, (ftnlen)1)) {
173 	info = -1;
174     } else if (! lsame_(storev, "R", (ftnlen)1, (ftnlen)1)) {
175 	info = -2;
176     }
177     if (info != 0) {
178 	i__1 = -info;
179 	xerbla_("SLARZT", &i__1, (ftnlen)6);
180 	return 0;
181     }
182 
183     for (i__ = *k; i__ >= 1; --i__) {
184 	if (tau[i__] == 0.f) {
185 
186 /*           H(i)  =  I */
187 
188 	    i__1 = *k;
189 	    for (j = i__; j <= i__1; ++j) {
190 		t[j + i__ * t_dim1] = 0.f;
191 /* L10: */
192 	    }
193 	} else {
194 
195 /*           general case */
196 
197 	    if (i__ < *k) {
198 
199 /*              T(i+1:k,i) = - tau(i) * V(i+1:k,1:n) * V(i,1:n)' */
200 
201 		i__1 = *k - i__;
202 		r__1 = -tau[i__];
203 		sgemv_("No transpose", &i__1, n, &r__1, &v[i__ + 1 + v_dim1],
204 			ldv, &v[i__ + v_dim1], ldv, &c_b8, &t[i__ + 1 + i__ *
205 			t_dim1], &c__1, (ftnlen)12);
206 
207 /*              T(i+1:k,i) = T(i+1:k,i+1:k) * T(i+1:k,i) */
208 
209 		i__1 = *k - i__;
210 		strmv_("Lower", "No transpose", "Non-unit", &i__1, &t[i__ + 1
211 			+ (i__ + 1) * t_dim1], ldt, &t[i__ + 1 + i__ * t_dim1]
212 			, &c__1, (ftnlen)5, (ftnlen)12, (ftnlen)8);
213 	    }
214 	    t[i__ + i__ * t_dim1] = tau[i__];
215 	}
216 /* L20: */
217     }
218     return 0;
219 
220 /*     End of SLARZT */
221 
222 } /* slarzt_ */
223 
224