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