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