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