| blob | 1e952d3477b273d28e254a06fa4bf247d1a42111 |
1 SUBROUTINE DLASR( SIDE, PIVOT, DIRECT, M, N, C, S, A, LDA )
2 !
3 ! -- LAPACK auxiliary routine (version 3.1) --
4 ! Univ. of Tennessee, Univ. of California Berkeley and NAG Ltd..
5 ! November 2006
6 !
7 ! .. Scalar Arguments ..
8 CHARACTER DIRECT, PIVOT, SIDE
9 INTEGER LDA, M, N
10 ! ..
11 ! .. Array Arguments ..
12 DOUBLE PRECISION A( LDA, * ), C( * ), S( * )
13 ! ..
14 !
15 ! Purpose
16 ! =======
17 !
18 ! DLASR applies a sequence of plane rotations to a real matrix A,
19 ! from either the left or the right.
20 !
21 ! When SIDE = 'L', the transformation takes the form
22 !
23 ! A := P*A
24 !
25 ! and when SIDE = 'R', the transformation takes the form
26 !
27 ! A := A*P**T
28 !
29 ! where P is an orthogonal matrix consisting of a sequence of z plane
30 ! rotations, with z = M when SIDE = 'L' and z = N when SIDE = 'R',
31 ! and P**T is the transpose of P.
32 !
33 ! When DIRECT = 'F' (Forward sequence), then
34 !
35 ! P = P(z-1) * ... * P(2) * P(1)
36 !
37 ! and when DIRECT = 'B' (Backward sequence), then
38 !
39 ! P = P(1) * P(2) * ... * P(z-1)
40 !
41 ! where P(k) is a plane rotation matrix defined by the 2-by-2 rotation
42 !
43 ! R(k) = ( c(k) s(k) )
44 ! = ( -s(k) c(k) ).
45 !
46 ! When PIVOT = 'V' (Variable pivot), the rotation is performed
47 ! for the plane (k,k+1), i.e., P(k) has the form
48 !
49 ! P(k) = ( 1 )
50 ! ( ... )
51 ! ( 1 )
52 ! ( c(k) s(k) )
53 ! ( -s(k) c(k) )
54 ! ( 1 )
55 ! ( ... )
56 ! ( 1 )
57 !
58 ! where R(k) appears as a rank-2 modification to the identity matrix in
59 ! rows and columns k and k+1.
60 !
61 ! When PIVOT = 'T' (Top pivot), the rotation is performed for the
62 ! plane (1,k+1), so P(k) has the form
63 !
64 ! P(k) = ( c(k) s(k) )
65 ! ( 1 )
66 ! ( ... )
67 ! ( 1 )
68 ! ( -s(k) c(k) )
69 ! ( 1 )
70 ! ( ... )
71 ! ( 1 )
72 !
73 ! where R(k) appears in rows and columns 1 and k+1.
74 !
75 ! Similarly, when PIVOT = 'B' (Bottom pivot), the rotation is
76 ! performed for the plane (k,z), giving P(k) the form
77 !
78 ! P(k) = ( 1 )
79 ! ( ... )
80 ! ( 1 )
81 ! ( c(k) s(k) )
82 ! ( 1 )
83 ! ( ... )
84 ! ( 1 )
85 ! ( -s(k) c(k) )
86 !
87 ! where R(k) appears in rows and columns k and z. The rotations are
88 ! performed without ever forming P(k) explicitly.
89 !
90 ! Arguments
91 ! =========
92 !
93 ! SIDE (input) CHARACTER*1
94 ! Specifies whether the plane rotation matrix P is applied to
95 ! A on the left or the right.
96 ! = 'L': Left, compute A := P*A
97 ! = 'R': Right, compute A:= A*P**T
98 !
99 ! PIVOT (input) CHARACTER*1
100 ! Specifies the plane for which P(k) is a plane rotation
101 ! matrix.
102 ! = 'V': Variable pivot, the plane (k,k+1)
103 ! = 'T': Top pivot, the plane (1,k+1)
104 ! = 'B': Bottom pivot, the plane (k,z)
105 !
106 ! DIRECT (input) CHARACTER*1
107 ! Specifies whether P is a forward or backward sequence of
108 ! plane rotations.
109 ! = 'F': Forward, P = P(z-1)*...*P(2)*P(1)
110 ! = 'B': Backward, P = P(1)*P(2)*...*P(z-1)
111 !
112 ! M (input) INTEGER
113 ! The number of rows of the matrix A. If m <= 1, an immediate
114 ! return is effected.
115 !
116 ! N (input) INTEGER
117 ! The number of columns of the matrix A. If n <= 1, an
118 ! immediate return is effected.
119 !
120 ! C (input) DOUBLE PRECISION array, dimension
121 ! (M-1) if SIDE = 'L'
122 ! (N-1) if SIDE = 'R'
123 ! The cosines c(k) of the plane rotations.
124 !
125 ! S (input) DOUBLE PRECISION array, dimension
126 ! (M-1) if SIDE = 'L'
127 ! (N-1) if SIDE = 'R'
128 ! The sines s(k) of the plane rotations. The 2-by-2 plane
129 ! rotation part of the matrix P(k), R(k), has the form
130 ! R(k) = ( c(k) s(k) )
131 ! ( -s(k) c(k) ).
132 !
133 ! A (input/output) DOUBLE PRECISION array, dimension (LDA,N)
134 ! The M-by-N matrix A. On exit, A is overwritten by P*A if
135 ! SIDE = 'R' or by A*P**T if SIDE = 'L'.
136 !
137 ! LDA (input) INTEGER
138 ! The leading dimension of the array A. LDA >= max(1,M).
139 !
140 ! =====================================================================
141 !
142 ! .. Parameters ..
143 DOUBLE PRECISION ONE, ZERO
144 PARAMETER ( ONE = 1.0D+0, ZERO = 0.0D+0 )
145 ! ..
146 ! .. Local Scalars ..
147 INTEGER I, INFO, J
148 DOUBLE PRECISION CTEMP, STEMP, TEMP
149 ! ..
150 ! .. External Functions ..
151 ! LOGICAL LSAME
152 ! EXTERNAL LSAME
153 ! ..
154 ! .. External Subroutines ..
155 ! EXTERNAL XERBLA
156 ! ..
157 ! .. Intrinsic Functions ..
158 INTRINSIC MAX
159 ! ..
160 ! .. Executable Statements ..
161 !
162 ! Test the input parameters
163 !
164 INFO = 0
165 IF( .NOT.( LSAME( SIDE, 'L' ) .OR. LSAME( SIDE, 'R' ) ) ) THEN
166 INFO = 1
167 ELSE IF( .NOT.( LSAME( PIVOT, 'V' ) .OR. LSAME( PIVOT, &
168 'T' ) .OR. LSAME( PIVOT, 'B' ) ) ) THEN
169 INFO = 2
170 ELSE IF( .NOT.( LSAME( DIRECT, 'F' ) .OR. LSAME( DIRECT, 'B' ) ) ) &
171 THEN
172 INFO = 3
173 ELSE IF( M.LT.0 ) THEN
174 INFO = 4
175 ELSE IF( N.LT.0 ) THEN
176 INFO = 5
177 ELSE IF( LDA.LT.MAX( 1, M ) ) THEN
178 INFO = 9
179 END IF
180 IF( INFO.NE.0 ) THEN
181 CALL XERBLA( 'DLASR ', INFO )
182 RETURN
183 END IF
184 !
185 ! Quick return if possible
186 !
187 IF( ( M.EQ.0 ) .OR. ( N.EQ.0 ) ) &
188 RETURN
189 IF( LSAME( SIDE, 'L' ) ) THEN
190 !
191 ! Form P * A
192 !
193 IF( LSAME( PIVOT, 'V' ) ) THEN
194 IF( LSAME( DIRECT, 'F' ) ) THEN
195 DO 20 J = 1, M - 1
196 CTEMP = C( J )
197 STEMP = S( J )
198 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
199 DO 10 I = 1, N
200 TEMP = A( J+1, I )
201 A( J+1, I ) = CTEMP*TEMP - STEMP*A( J, I )
202 A( J, I ) = STEMP*TEMP + CTEMP*A( J, I )
203 10 CONTINUE
204 END IF
205 20 CONTINUE
206 ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
207 DO 40 J = M - 1, 1, -1
208 CTEMP = C( J )
209 STEMP = S( J )
210 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
211 DO 30 I = 1, N
212 TEMP = A( J+1, I )
213 A( J+1, I ) = CTEMP*TEMP - STEMP*A( J, I )
214 A( J, I ) = STEMP*TEMP + CTEMP*A( J, I )
215 30 CONTINUE
216 END IF
217 40 CONTINUE
218 END IF
219 ELSE IF( LSAME( PIVOT, 'T' ) ) THEN
220 IF( LSAME( DIRECT, 'F' ) ) THEN
221 DO 60 J = 2, M
222 CTEMP = C( J-1 )
223 STEMP = S( J-1 )
224 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
225 DO 50 I = 1, N
226 TEMP = A( J, I )
227 A( J, I ) = CTEMP*TEMP - STEMP*A( 1, I )
228 A( 1, I ) = STEMP*TEMP + CTEMP*A( 1, I )
229 50 CONTINUE
230 END IF
231 60 CONTINUE
232 ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
233 DO 80 J = M, 2, -1
234 CTEMP = C( J-1 )
235 STEMP = S( J-1 )
236 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
237 DO 70 I = 1, N
238 TEMP = A( J, I )
239 A( J, I ) = CTEMP*TEMP - STEMP*A( 1, I )
240 A( 1, I ) = STEMP*TEMP + CTEMP*A( 1, I )
241 70 CONTINUE
242 END IF
243 80 CONTINUE
244 END IF
245 ELSE IF( LSAME( PIVOT, 'B' ) ) THEN
246 IF( LSAME( DIRECT, 'F' ) ) THEN
247 DO 100 J = 1, M - 1
248 CTEMP = C( J )
249 STEMP = S( J )
250 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
251 DO 90 I = 1, N
252 TEMP = A( J, I )
253 A( J, I ) = STEMP*A( M, I ) + CTEMP*TEMP
254 A( M, I ) = CTEMP*A( M, I ) - STEMP*TEMP
255 90 CONTINUE
256 END IF
257 100 CONTINUE
258 ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
259 DO 120 J = M - 1, 1, -1
260 CTEMP = C( J )
261 STEMP = S( J )
262 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
263 DO 110 I = 1, N
264 TEMP = A( J, I )
265 A( J, I ) = STEMP*A( M, I ) + CTEMP*TEMP
266 A( M, I ) = CTEMP*A( M, I ) - STEMP*TEMP
267 110 CONTINUE
268 END IF
269 120 CONTINUE
270 END IF
271 END IF
272 ELSE IF( LSAME( SIDE, 'R' ) ) THEN
273 !
274 ! Form A * P'
275 !
276 IF( LSAME( PIVOT, 'V' ) ) THEN
277 IF( LSAME( DIRECT, 'F' ) ) THEN
278 DO 140 J = 1, N - 1
279 CTEMP = C( J )
280 STEMP = S( J )
281 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
282 DO 130 I = 1, M
283 TEMP = A( I, J+1 )
284 A( I, J+1 ) = CTEMP*TEMP - STEMP*A( I, J )
285 A( I, J ) = STEMP*TEMP + CTEMP*A( I, J )
286 130 CONTINUE
287 END IF
288 140 CONTINUE
289 ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
290 DO 160 J = N - 1, 1, -1
291 CTEMP = C( J )
292 STEMP = S( J )
293 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
294 DO 150 I = 1, M
295 TEMP = A( I, J+1 )
296 A( I, J+1 ) = CTEMP*TEMP - STEMP*A( I, J )
297 A( I, J ) = STEMP*TEMP + CTEMP*A( I, J )
298 150 CONTINUE
299 END IF
300 160 CONTINUE
301 END IF
302 ELSE IF( LSAME( PIVOT, 'T' ) ) THEN
303 IF( LSAME( DIRECT, 'F' ) ) THEN
304 DO 180 J = 2, N
305 CTEMP = C( J-1 )
306 STEMP = S( J-1 )
307 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
308 DO 170 I = 1, M
309 TEMP = A( I, J )
310 A( I, J ) = CTEMP*TEMP - STEMP*A( I, 1 )
311 A( I, 1 ) = STEMP*TEMP + CTEMP*A( I, 1 )
312 170 CONTINUE
313 END IF
314 180 CONTINUE
315 ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
316 DO 200 J = N, 2, -1
317 CTEMP = C( J-1 )
318 STEMP = S( J-1 )
319 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
320 DO 190 I = 1, M
321 TEMP = A( I, J )
322 A( I, J ) = CTEMP*TEMP - STEMP*A( I, 1 )
323 A( I, 1 ) = STEMP*TEMP + CTEMP*A( I, 1 )
324 190 CONTINUE
325 END IF
326 200 CONTINUE
327 END IF
328 ELSE IF( LSAME( PIVOT, 'B' ) ) THEN
329 IF( LSAME( DIRECT, 'F' ) ) THEN
330 DO 220 J = 1, N - 1
331 CTEMP = C( J )
332 STEMP = S( J )
333 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
334 DO 210 I = 1, M
335 TEMP = A( I, J )
336 A( I, J ) = STEMP*A( I, N ) + CTEMP*TEMP
337 A( I, N ) = CTEMP*A( I, N ) - STEMP*TEMP
338 210 CONTINUE
339 END IF
340 220 CONTINUE
341 ELSE IF( LSAME( DIRECT, 'B' ) ) THEN
342 DO 240 J = N - 1, 1, -1
343 CTEMP = C( J )
344 STEMP = S( J )
345 IF( ( CTEMP.NE.ONE ) .OR. ( STEMP.NE.ZERO ) ) THEN
346 DO 230 I = 1, M
347 TEMP = A( I, J )
348 A( I, J ) = STEMP*A( I, N ) + CTEMP*TEMP
349 A( I, N ) = CTEMP*A( I, N ) - STEMP*TEMP
350 230 CONTINUE
351 END IF
352 240 CONTINUE
353 END IF
354 END IF
355 END IF
356 !
357 RETURN
358 !
359 ! End of DLASR
360 !
361 END SUBROUTINE DLASR