2
3
4
5
6
7
8 INTEGER INFO, KL, KU, LDA, M, N
9
10
11 INTEGER ISEED( 4 )
12 REAL A( LDA, * ), D( * ), WORK( * )
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64 REAL ZERO, ONE
65 parameter( zero = 0.0e+0, one = 1.0e+0 )
66
67
68 INTEGER I, J
69 REAL TAU, WA, WB, WN
70
71
72 EXTERNAL sgemv, sger, slarnv, sscal, xerbla
73
74
76
77
78 REAL SNRM2
79 EXTERNAL snrm2
80
81
82
83
84
85 info = 0
86 IF( m.LT.0 ) THEN
87 info = -1
88 ELSE IF( n.LT.0 ) THEN
89 info = -2
90 ELSE IF( kl.LT.0 .OR. kl.GT.m-1 ) THEN
91 info = -3
92 ELSE IF( ku.LT.0 .OR. ku.GT.n-1 ) THEN
93 info = -4
94 ELSE IF( lda.LT.
max( 1, m ) )
THEN
95 info = -7
96 END IF
97 IF( info.LT.0 ) THEN
98 CALL xerbla( 'SLAGGE', -info )
99 RETURN
100 END IF
101
102
103
104 DO 20 j = 1, n
105 DO 10 i = 1, m
106 a( i, j ) = zero
107 10 CONTINUE
108 20 CONTINUE
109 DO 30 i = 1,
min( m, n )
110 a( i, i ) = d( i )
111 30 CONTINUE
112
113
114
115 DO 40 i =
min( m, n ), 1, -1
116 IF( i.LT.m ) THEN
117
118
119
120 CALL slarnv( 3, iseed, m-i+1, work )
121 wn = snrm2( m-i+1, work, 1 )
122 wa = sign( wn, work( 1 ) )
123 IF( wn.EQ.zero ) THEN
124 tau = zero
125 ELSE
126 wb = work( 1 ) + wa
127 CALL sscal( m-i, one / wb, work( 2 ), 1 )
128 work( 1 ) = one
129 tau = wb / wa
130 END IF
131
132
133
134 CALL sgemv( 'Transpose', m-i+1, n-i+1, one, a( i, i ), lda,
135 $ work, 1, zero, work( m+1 ), 1 )
136 CALL sger( m-i+1, n-i+1, -tau, work, 1, work( m+1 ), 1,
137 $ a( i, i ), lda )
138 END IF
139 IF( i.LT.n ) THEN
140
141
142
143 CALL slarnv( 3, iseed, n-i+1, work )
144 wn = snrm2( n-i+1, work, 1 )
145 wa = sign( wn, work( 1 ) )
146 IF( wn.EQ.zero ) THEN
147 tau = zero
148 ELSE
149 wb = work( 1 ) + wa
150 CALL sscal( n-i, one / wb, work( 2 ), 1 )
151 work( 1 ) = one
152 tau = wb / wa
153 END IF
154
155
156
157 CALL sgemv( 'No transpose', m-i+1, n-i+1, one, a( i, i ),
158 $ lda, work, 1, zero, work( n+1 ), 1 )
159 CALL sger( m-i+1, n-i+1, -tau, work( n+1 ), 1, work, 1,
160 $ a( i, i ), lda )
161 END IF
162 40 CONTINUE
163
164
165
166
167 DO 70 i = 1,
max( m-1-kl, n-1-ku )
168 IF( kl.LE.ku ) THEN
169
170
171
172 IF( i.LE.
min( m-1-kl, n ) )
THEN
173
174
175
176 wn = snrm2( m-kl-i+1, a( kl+i, i ), 1 )
177 wa = sign( wn, a( kl+i, i ) )
178 IF( wn.EQ.zero ) THEN
179 tau = zero
180 ELSE
181 wb = a( kl+i, i ) + wa
182 CALL sscal( m-kl-i, one / wb, a( kl+i+1, i ), 1 )
183 a( kl+i, i ) = one
184 tau = wb / wa
185 END IF
186
187
188
189 CALL sgemv( 'Transpose', m-kl-i+1, n-i, one,
190 $ a( kl+i, i+1 ), lda, a( kl+i, i ), 1, zero,
191 $ work, 1 )
192 CALL sger( m-kl-i+1, n-i, -tau, a( kl+i, i ), 1, work, 1,
193 $ a( kl+i, i+1 ), lda )
194 a( kl+i, i ) = -wa
195 END IF
196
197 IF( i.LE.
min( n-1-ku, m ) )
THEN
198
199
200
201 wn = snrm2( n-ku-i+1, a( i, ku+i ), lda )
202 wa = sign( wn, a( i, ku+i ) )
203 IF( wn.EQ.zero ) THEN
204 tau = zero
205 ELSE
206 wb = a( i, ku+i ) + wa
207 CALL sscal( n-ku-i, one / wb, a( i, ku+i+1 ), lda )
208 a( i, ku+i ) = one
209 tau = wb / wa
210 END IF
211
212
213
214 CALL sgemv( 'No transpose', m-i, n-ku-i+1, one,
215 $ a( i+1, ku+i ), lda, a( i, ku+i ), lda, zero,
216 $ work, 1 )
217 CALL sger( m-i, n-ku-i+1, -tau, work, 1, a( i, ku+i ),
218 $ lda, a( i+1, ku+i ), lda )
219 a( i, ku+i ) = -wa
220 END IF
221 ELSE
222
223
224
225
226 IF( i.LE.
min( n-1-ku, m ) )
THEN
227
228
229
230 wn = snrm2( n-ku-i+1, a( i, ku+i ), lda )
231 wa = sign( wn, a( i, ku+i ) )
232 IF( wn.EQ.zero ) THEN
233 tau = zero
234 ELSE
235 wb = a( i, ku+i ) + wa
236 CALL sscal( n-ku-i, one / wb, a( i, ku+i+1 ), lda )
237 a( i, ku+i ) = one
238 tau = wb / wa
239 END IF
240
241
242
243 CALL sgemv( 'No transpose', m-i, n-ku-i+1, one,
244 $ a( i+1, ku+i ), lda, a( i, ku+i ), lda, zero,
245 $ work, 1 )
246 CALL sger( m-i, n-ku-i+1, -tau, work, 1, a( i, ku+i ),
247 $ lda, a( i+1, ku+i ), lda )
248 a( i, ku+i ) = -wa
249 END IF
250
251 IF( i.LE.
min( m-1-kl, n ) )
THEN
252
253
254
255 wn = snrm2( m-kl-i+1, a( kl+i, i ), 1 )
256 wa = sign( wn, a( kl+i, i ) )
257 IF( wn.EQ.zero ) THEN
258 tau = zero
259 ELSE
260 wb = a( kl+i, i ) + wa
261 CALL sscal( m-kl-i, one / wb, a( kl+i+1, i ), 1 )
262 a( kl+i, i ) = one
263 tau = wb / wa
264 END IF
265
266
267
268 CALL sgemv( 'Transpose', m-kl-i+1, n-i, one,
269 $ a( kl+i, i+1 ), lda, a( kl+i, i ), 1, zero,
270 $ work, 1 )
271 CALL sger( m-kl-i+1, n-i, -tau, a( kl+i, i ), 1, work, 1,
272 $ a( kl+i, i+1 ), lda )
273 a( kl+i, i ) = -wa
274 END IF
275 END IF
276
277 DO 50 j = kl + i + 1, m
278 a( j, i ) = zero
279 50 CONTINUE
280
281 DO 60 j = ku + i + 1, n
282 a( i, j ) = zero
283 60 CONTINUE
284 70 CONTINUE
285 RETURN
286
287
288