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