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