• Home
  • Line#
  • Scopes#
  • Navigate#
  • Raw
  • Download
1      SUBROUTINE SSPMV(UPLO,N,ALPHA,AP,X,INCX,BETA,Y,INCY)
2*     .. Scalar Arguments ..
3      REAL ALPHA,BETA
4      INTEGER INCX,INCY,N
5      CHARACTER UPLO
6*     ..
7*     .. Array Arguments ..
8      REAL AP(*),X(*),Y(*)
9*     ..
10*
11*  Purpose
12*  =======
13*
14*  SSPMV  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 symmetric matrix, supplied in packed form.
20*
21*  Arguments
22*  ==========
23*
24*  UPLO   - CHARACTER*1.
25*           On entry, UPLO specifies whether the upper or lower
26*           triangular part of the matrix A is supplied in the packed
27*           array AP as follows:
28*
29*              UPLO = 'U' or 'u'   The upper triangular part of A is
30*                                  supplied in AP.
31*
32*              UPLO = 'L' or 'l'   The lower triangular part of A is
33*                                  supplied in AP.
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  - REAL            .
43*           On entry, ALPHA specifies the scalar alpha.
44*           Unchanged on exit.
45*
46*  AP     - REAL             array of DIMENSION at least
47*           ( ( n*( n + 1 ) )/2 ).
48*           Before entry with UPLO = 'U' or 'u', the array AP must
49*           contain the upper triangular part of the symmetric matrix
50*           packed sequentially, column by column, so that AP( 1 )
51*           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 1, 2 )
52*           and a( 2, 2 ) respectively, and so on.
53*           Before entry with UPLO = 'L' or 'l', the array AP must
54*           contain the lower triangular part of the symmetric matrix
55*           packed sequentially, column by column, so that AP( 1 )
56*           contains a( 1, 1 ), AP( 2 ) and AP( 3 ) contain a( 2, 1 )
57*           and a( 3, 1 ) respectively, and so on.
58*           Unchanged on exit.
59*
60*  X      - REAL             array of dimension at least
61*           ( 1 + ( n - 1 )*abs( INCX ) ).
62*           Before entry, the incremented array X must contain the n
63*           element vector x.
64*           Unchanged on exit.
65*
66*  INCX   - INTEGER.
67*           On entry, INCX specifies the increment for the elements of
68*           X. INCX must not be zero.
69*           Unchanged on exit.
70*
71*  BETA   - REAL            .
72*           On entry, BETA specifies the scalar beta. When BETA is
73*           supplied as zero then Y need not be set on input.
74*           Unchanged on exit.
75*
76*  Y      - REAL             array of dimension at least
77*           ( 1 + ( n - 1 )*abs( INCY ) ).
78*           Before entry, the incremented array Y must contain the n
79*           element vector y. On exit, Y is overwritten by the updated
80*           vector y.
81*
82*  INCY   - INTEGER.
83*           On entry, INCY specifies the increment for the elements of
84*           Y. INCY must not be zero.
85*           Unchanged on exit.
86*
87*  Further Details
88*  ===============
89*
90*  Level 2 Blas routine.
91*
92*  -- Written on 22-October-1986.
93*     Jack Dongarra, Argonne National Lab.
94*     Jeremy Du Croz, Nag Central Office.
95*     Sven Hammarling, Nag Central Office.
96*     Richard Hanson, Sandia National Labs.
97*
98*  =====================================================================
99*
100*     .. Parameters ..
101      REAL ONE,ZERO
102      PARAMETER (ONE=1.0E+0,ZERO=0.0E+0)
103*     ..
104*     .. Local Scalars ..
105      REAL TEMP1,TEMP2
106      INTEGER I,INFO,IX,IY,J,JX,JY,K,KK,KX,KY
107*     ..
108*     .. External Functions ..
109      LOGICAL LSAME
110      EXTERNAL LSAME
111*     ..
112*     .. External Subroutines ..
113      EXTERNAL XERBLA
114*     ..
115*
116*     Test the input parameters.
117*
118      INFO = 0
119      IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
120          INFO = 1
121      ELSE IF (N.LT.0) THEN
122          INFO = 2
123      ELSE IF (INCX.EQ.0) THEN
124          INFO = 6
125      ELSE IF (INCY.EQ.0) THEN
126          INFO = 9
127      END IF
128      IF (INFO.NE.0) THEN
129          CALL XERBLA('SSPMV ',INFO)
130          RETURN
131      END IF
132*
133*     Quick return if possible.
134*
135      IF ((N.EQ.0) .OR. ((ALPHA.EQ.ZERO).AND. (BETA.EQ.ONE))) RETURN
136*
137*     Set up the start points in  X  and  Y.
138*
139      IF (INCX.GT.0) THEN
140          KX = 1
141      ELSE
142          KX = 1 - (N-1)*INCX
143      END IF
144      IF (INCY.GT.0) THEN
145          KY = 1
146      ELSE
147          KY = 1 - (N-1)*INCY
148      END IF
149*
150*     Start the operations. In this version the elements of the array AP
151*     are accessed sequentially with one pass through AP.
152*
153*     First form  y := beta*y.
154*
155      IF (BETA.NE.ONE) THEN
156          IF (INCY.EQ.1) THEN
157              IF (BETA.EQ.ZERO) THEN
158                  DO 10 I = 1,N
159                      Y(I) = ZERO
160   10             CONTINUE
161              ELSE
162                  DO 20 I = 1,N
163                      Y(I) = BETA*Y(I)
164   20             CONTINUE
165              END IF
166          ELSE
167              IY = KY
168              IF (BETA.EQ.ZERO) THEN
169                  DO 30 I = 1,N
170                      Y(IY) = ZERO
171                      IY = IY + INCY
172   30             CONTINUE
173              ELSE
174                  DO 40 I = 1,N
175                      Y(IY) = BETA*Y(IY)
176                      IY = IY + INCY
177   40             CONTINUE
178              END IF
179          END IF
180      END IF
181      IF (ALPHA.EQ.ZERO) 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 SSPMV .
264*
265      END
266