root / src / blas / cgemv.f @ 11
Historique | Voir | Annoter | Télécharger (7,79 ko)
1 |
SUBROUTINE CGEMV(TRANS,M,N,ALPHA,A,LDA,X,INCX,BETA,Y,INCY) |
---|---|
2 |
* .. Scalar Arguments .. |
3 |
COMPLEX ALPHA,BETA |
4 |
INTEGER INCX,INCY,LDA,M,N |
5 |
CHARACTER TRANS |
6 |
* .. |
7 |
* .. Array Arguments .. |
8 |
COMPLEX A(LDA,*),X(*),Y(*) |
9 |
* .. |
10 |
* |
11 |
* Purpose |
12 |
* ======= |
13 |
* |
14 |
* CGEMV performs one of the matrix-vector operations |
15 |
* |
16 |
* y := alpha*A*x + beta*y, or y := alpha*A'*x + beta*y, or |
17 |
* |
18 |
* y := alpha*conjg( A' )*x + beta*y, |
19 |
* |
20 |
* where alpha and beta are scalars, x and y are vectors and A is an |
21 |
* m by n matrix. |
22 |
* |
23 |
* Arguments |
24 |
* ========== |
25 |
* |
26 |
* TRANS - CHARACTER*1. |
27 |
* On entry, TRANS specifies the operation to be performed as |
28 |
* follows: |
29 |
* |
30 |
* TRANS = 'N' or 'n' y := alpha*A*x + beta*y. |
31 |
* |
32 |
* TRANS = 'T' or 't' y := alpha*A'*x + beta*y. |
33 |
* |
34 |
* TRANS = 'C' or 'c' y := alpha*conjg( A' )*x + beta*y. |
35 |
* |
36 |
* Unchanged on exit. |
37 |
* |
38 |
* M - INTEGER. |
39 |
* On entry, M specifies the number of rows of the matrix A. |
40 |
* M must be at least zero. |
41 |
* Unchanged on exit. |
42 |
* |
43 |
* N - INTEGER. |
44 |
* On entry, N specifies the number of columns of the matrix A. |
45 |
* N must be at least zero. |
46 |
* Unchanged on exit. |
47 |
* |
48 |
* ALPHA - COMPLEX . |
49 |
* On entry, ALPHA specifies the scalar alpha. |
50 |
* Unchanged on exit. |
51 |
* |
52 |
* A - COMPLEX array of DIMENSION ( LDA, n ). |
53 |
* Before entry, the leading m by n part of the array A must |
54 |
* contain the matrix of coefficients. |
55 |
* Unchanged on exit. |
56 |
* |
57 |
* LDA - INTEGER. |
58 |
* On entry, LDA specifies the first dimension of A as declared |
59 |
* in the calling (sub) program. LDA must be at least |
60 |
* max( 1, m ). |
61 |
* Unchanged on exit. |
62 |
* |
63 |
* X - COMPLEX array of DIMENSION at least |
64 |
* ( 1 + ( n - 1 )*abs( INCX ) ) when TRANS = 'N' or 'n' |
65 |
* and at least |
66 |
* ( 1 + ( m - 1 )*abs( INCX ) ) otherwise. |
67 |
* Before entry, the incremented array X must contain the |
68 |
* vector x. |
69 |
* Unchanged on exit. |
70 |
* |
71 |
* INCX - INTEGER. |
72 |
* On entry, INCX specifies the increment for the elements of |
73 |
* X. INCX must not be zero. |
74 |
* Unchanged on exit. |
75 |
* |
76 |
* BETA - COMPLEX . |
77 |
* On entry, BETA specifies the scalar beta. When BETA is |
78 |
* supplied as zero then Y need not be set on input. |
79 |
* Unchanged on exit. |
80 |
* |
81 |
* Y - COMPLEX array of DIMENSION at least |
82 |
* ( 1 + ( m - 1 )*abs( INCY ) ) when TRANS = 'N' or 'n' |
83 |
* and at least |
84 |
* ( 1 + ( n - 1 )*abs( INCY ) ) otherwise. |
85 |
* Before entry with BETA non-zero, the incremented array Y |
86 |
* must contain the vector y. On exit, Y is overwritten by the |
87 |
* updated vector y. |
88 |
* |
89 |
* INCY - INTEGER. |
90 |
* On entry, INCY specifies the increment for the elements of |
91 |
* Y. INCY must not be zero. |
92 |
* Unchanged on exit. |
93 |
* |
94 |
* |
95 |
* Level 2 Blas routine. |
96 |
* |
97 |
* -- Written on 22-October-1986. |
98 |
* Jack Dongarra, Argonne National Lab. |
99 |
* Jeremy Du Croz, Nag Central Office. |
100 |
* Sven Hammarling, Nag Central Office. |
101 |
* Richard Hanson, Sandia National Labs. |
102 |
* |
103 |
* |
104 |
* .. Parameters .. |
105 |
COMPLEX ONE |
106 |
PARAMETER (ONE= (1.0E+0,0.0E+0)) |
107 |
COMPLEX ZERO |
108 |
PARAMETER (ZERO= (0.0E+0,0.0E+0)) |
109 |
* .. |
110 |
* .. Local Scalars .. |
111 |
COMPLEX TEMP |
112 |
INTEGER I,INFO,IX,IY,J,JX,JY,KX,KY,LENX,LENY |
113 |
LOGICAL NOCONJ |
114 |
* .. |
115 |
* .. External Functions .. |
116 |
LOGICAL LSAME |
117 |
EXTERNAL LSAME |
118 |
* .. |
119 |
* .. External Subroutines .. |
120 |
EXTERNAL XERBLA |
121 |
* .. |
122 |
* .. Intrinsic Functions .. |
123 |
INTRINSIC CONJG,MAX |
124 |
* .. |
125 |
* |
126 |
* Test the input parameters. |
127 |
* |
128 |
INFO = 0 |
129 |
IF (.NOT.LSAME(TRANS,'N') .AND. .NOT.LSAME(TRANS,'T') .AND. |
130 |
+ .NOT.LSAME(TRANS,'C')) THEN |
131 |
INFO = 1 |
132 |
ELSE IF (M.LT.0) THEN |
133 |
INFO = 2 |
134 |
ELSE IF (N.LT.0) THEN |
135 |
INFO = 3 |
136 |
ELSE IF (LDA.LT.MAX(1,M)) THEN |
137 |
INFO = 6 |
138 |
ELSE IF (INCX.EQ.0) THEN |
139 |
INFO = 8 |
140 |
ELSE IF (INCY.EQ.0) THEN |
141 |
INFO = 11 |
142 |
END IF |
143 |
IF (INFO.NE.0) THEN |
144 |
CALL XERBLA('CGEMV ',INFO) |
145 |
RETURN |
146 |
END IF |
147 |
* |
148 |
* Quick return if possible. |
149 |
* |
150 |
IF ((M.EQ.0) .OR. (N.EQ.0) .OR. |
151 |
+ ((ALPHA.EQ.ZERO).AND. (BETA.EQ.ONE))) RETURN |
152 |
* |
153 |
NOCONJ = LSAME(TRANS,'T') |
154 |
* |
155 |
* Set LENX and LENY, the lengths of the vectors x and y, and set |
156 |
* up the start points in X and Y. |
157 |
* |
158 |
IF (LSAME(TRANS,'N')) THEN |
159 |
LENX = N |
160 |
LENY = M |
161 |
ELSE |
162 |
LENX = M |
163 |
LENY = N |
164 |
END IF |
165 |
IF (INCX.GT.0) THEN |
166 |
KX = 1 |
167 |
ELSE |
168 |
KX = 1 - (LENX-1)*INCX |
169 |
END IF |
170 |
IF (INCY.GT.0) THEN |
171 |
KY = 1 |
172 |
ELSE |
173 |
KY = 1 - (LENY-1)*INCY |
174 |
END IF |
175 |
* |
176 |
* Start the operations. In this version the elements of A are |
177 |
* accessed sequentially with one pass through A. |
178 |
* |
179 |
* First form y := beta*y. |
180 |
* |
181 |
IF (BETA.NE.ONE) THEN |
182 |
IF (INCY.EQ.1) THEN |
183 |
IF (BETA.EQ.ZERO) THEN |
184 |
DO 10 I = 1,LENY |
185 |
Y(I) = ZERO |
186 |
10 CONTINUE |
187 |
ELSE |
188 |
DO 20 I = 1,LENY |
189 |
Y(I) = BETA*Y(I) |
190 |
20 CONTINUE |
191 |
END IF |
192 |
ELSE |
193 |
IY = KY |
194 |
IF (BETA.EQ.ZERO) THEN |
195 |
DO 30 I = 1,LENY |
196 |
Y(IY) = ZERO |
197 |
IY = IY + INCY |
198 |
30 CONTINUE |
199 |
ELSE |
200 |
DO 40 I = 1,LENY |
201 |
Y(IY) = BETA*Y(IY) |
202 |
IY = IY + INCY |
203 |
40 CONTINUE |
204 |
END IF |
205 |
END IF |
206 |
END IF |
207 |
IF (ALPHA.EQ.ZERO) RETURN |
208 |
IF (LSAME(TRANS,'N')) THEN |
209 |
* |
210 |
* Form y := alpha*A*x + y. |
211 |
* |
212 |
JX = KX |
213 |
IF (INCY.EQ.1) THEN |
214 |
DO 60 J = 1,N |
215 |
IF (X(JX).NE.ZERO) THEN |
216 |
TEMP = ALPHA*X(JX) |
217 |
DO 50 I = 1,M |
218 |
Y(I) = Y(I) + TEMP*A(I,J) |
219 |
50 CONTINUE |
220 |
END IF |
221 |
JX = JX + INCX |
222 |
60 CONTINUE |
223 |
ELSE |
224 |
DO 80 J = 1,N |
225 |
IF (X(JX).NE.ZERO) THEN |
226 |
TEMP = ALPHA*X(JX) |
227 |
IY = KY |
228 |
DO 70 I = 1,M |
229 |
Y(IY) = Y(IY) + TEMP*A(I,J) |
230 |
IY = IY + INCY |
231 |
70 CONTINUE |
232 |
END IF |
233 |
JX = JX + INCX |
234 |
80 CONTINUE |
235 |
END IF |
236 |
ELSE |
237 |
* |
238 |
* Form y := alpha*A'*x + y or y := alpha*conjg( A' )*x + y. |
239 |
* |
240 |
JY = KY |
241 |
IF (INCX.EQ.1) THEN |
242 |
DO 110 J = 1,N |
243 |
TEMP = ZERO |
244 |
IF (NOCONJ) THEN |
245 |
DO 90 I = 1,M |
246 |
TEMP = TEMP + A(I,J)*X(I) |
247 |
90 CONTINUE |
248 |
ELSE |
249 |
DO 100 I = 1,M |
250 |
TEMP = TEMP + CONJG(A(I,J))*X(I) |
251 |
100 CONTINUE |
252 |
END IF |
253 |
Y(JY) = Y(JY) + ALPHA*TEMP |
254 |
JY = JY + INCY |
255 |
110 CONTINUE |
256 |
ELSE |
257 |
DO 140 J = 1,N |
258 |
TEMP = ZERO |
259 |
IX = KX |
260 |
IF (NOCONJ) THEN |
261 |
DO 120 I = 1,M |
262 |
TEMP = TEMP + A(I,J)*X(IX) |
263 |
IX = IX + INCX |
264 |
120 CONTINUE |
265 |
ELSE |
266 |
DO 130 I = 1,M |
267 |
TEMP = TEMP + CONJG(A(I,J))*X(IX) |
268 |
IX = IX + INCX |
269 |
130 CONTINUE |
270 |
END IF |
271 |
Y(JY) = Y(JY) + ALPHA*TEMP |
272 |
JY = JY + INCY |
273 |
140 CONTINUE |
274 |
END IF |
275 |
END IF |
276 |
* |
277 |
RETURN |
278 |
* |
279 |
* End of CGEMV . |
280 |
* |
281 |
END |