LAPACK: BLAS/SRC/ssbmv.f Source File

Go to the documentation of this file.
00001       SUBROUTINE SSBMV(UPLO,N,K,ALPHA,A,LDA,X,INCX,BETA,Y,INCY)
00002 *     .. Scalar Arguments ..
00003       REAL ALPHA,BETA
00004       INTEGER INCX,INCY,K,LDA,N
00005       CHARACTER UPLO
00006 *     ..
00007 *     .. Array Arguments ..
00008       REAL A(LDA,*),X(*),Y(*)
00009 *     ..
00010 *
00011 *  Purpose
00012 *  =======
00013 *
00014 *  SSBMV  performs the matrix-vector  operation
00015 *
00016 *     y := alpha*A*x + beta*y,
00017 *
00018 *  where alpha and beta are scalars, x and y are n element vectors and
00019 *  A is an n by n symmetric band matrix, with k super-diagonals.
00020 *
00021 *  Arguments
00022 *  ==========
00023 *
00024 *  UPLO   - CHARACTER*1.
00025 *           On entry, UPLO specifies whether the upper or lower
00026 *           triangular part of the band matrix A is being supplied as
00027 *           follows:
00028 *
00029 *              UPLO = 'U' or 'u'   The upper triangular part of A is
00030 *                                  being supplied.
00031 *
00032 *              UPLO = 'L' or 'l'   The lower triangular part of A is
00033 *                                  being supplied.
00034 *
00035 *           Unchanged on exit.
00036 *
00037 *  N      - INTEGER.
00038 *           On entry, N specifies the order of the matrix A.
00039 *           N must be at least zero.
00040 *           Unchanged on exit.
00041 *
00042 *  K      - INTEGER.
00043 *           On entry, K specifies the number of super-diagonals of the
00044 *           matrix A. K must satisfy  0 .le. K.
00045 *           Unchanged on exit.
00046 *
00047 *  ALPHA  - REAL            .
00048 *           On entry, ALPHA specifies the scalar alpha.
00049 *           Unchanged on exit.
00050 *
00051 *  A      - REAL             array of DIMENSION ( LDA, n ).
00052 *           Before entry with UPLO = 'U' or 'u', the leading ( k + 1 )
00053 *           by n part of the array A must contain the upper triangular
00054 *           band part of the symmetric matrix, supplied column by
00055 *           column, with the leading diagonal of the matrix in row
00056 *           ( k + 1 ) of the array, the first super-diagonal starting at
00057 *           position 2 in row k, and so on. The top left k by k triangle
00058 *           of the array A is not referenced.
00059 *           The following program segment will transfer the upper
00060 *           triangular part of a symmetric band matrix from conventional
00061 *           full matrix storage to band storage:
00062 *
00063 *                 DO 20, J = 1, N
00064 *                    M = K + 1 - J
00065 *                    DO 10, I = MAX( 1, J - K ), J
00066 *                       A( M + I, J ) = matrix( I, J )
00067 *              10    CONTINUE
00068 *              20 CONTINUE
00069 *
00070 *           Before entry with UPLO = 'L' or 'l', the leading ( k + 1 )
00071 *           by n part of the array A must contain the lower triangular
00072 *           band part of the symmetric matrix, supplied column by
00073 *           column, with the leading diagonal of the matrix in row 1 of
00074 *           the array, the first sub-diagonal starting at position 1 in
00075 *           row 2, and so on. The bottom right k by k triangle of the
00076 *           array A is not referenced.
00077 *           The following program segment will transfer the lower
00078 *           triangular part of a symmetric band matrix from conventional
00079 *           full matrix storage to band storage:
00080 *
00081 *                 DO 20, J = 1, N
00082 *                    M = 1 - J
00083 *                    DO 10, I = J, MIN( N, J + K )
00084 *                       A( M + I, J ) = matrix( I, J )
00085 *              10    CONTINUE
00086 *              20 CONTINUE
00087 *
00088 *           Unchanged on exit.
00089 *
00090 *  LDA    - INTEGER.
00091 *           On entry, LDA specifies the first dimension of A as declared
00092 *           in the calling (sub) program. LDA must be at least
00093 *           ( k + 1 ).
00094 *           Unchanged on exit.
00095 *
00096 *  X      - REAL             array of DIMENSION at least
00097 *           ( 1 + ( n - 1 )*abs( INCX ) ).
00098 *           Before entry, the incremented array X must contain the
00099 *           vector x.
00100 *           Unchanged on exit.
00101 *
00102 *  INCX   - INTEGER.
00103 *           On entry, INCX specifies the increment for the elements of
00104 *           X. INCX must not be zero.
00105 *           Unchanged on exit.
00106 *
00107 *  BETA   - REAL            .
00108 *           On entry, BETA specifies the scalar beta.
00109 *           Unchanged on exit.
00110 *
00111 *  Y      - REAL             array of DIMENSION at least
00112 *           ( 1 + ( n - 1 )*abs( INCY ) ).
00113 *           Before entry, the incremented array Y must contain the
00114 *           vector y. On exit, Y is overwritten by the updated vector y.
00115 *
00116 *  INCY   - INTEGER.
00117 *           On entry, INCY specifies the increment for the elements of
00118 *           Y. INCY must not be zero.
00119 *           Unchanged on exit.
00120 *
00121 *  Further Details
00122 *  ===============
00123 *
00124 *  Level 2 Blas routine.
00125 *
00126 *  -- Written on 22-October-1986.
00127 *     Jack Dongarra, Argonne National Lab.
00128 *     Jeremy Du Croz, Nag Central Office.
00129 *     Sven Hammarling, Nag Central Office.
00130 *     Richard Hanson, Sandia National Labs.
00131 *
00132 *  =====================================================================
00133 *
00134 *     .. Parameters ..
00135       REAL ONE,ZERO
00136       PARAMETER (ONE=1.0E+0,ZERO=0.0E+0)
00137 *     ..
00138 *     .. Local Scalars ..
00139       REAL TEMP1,TEMP2
00140       INTEGER I,INFO,IX,IY,J,JX,JY,KPLUS1,KX,KY,L
00141 *     ..
00142 *     .. External Functions ..
00143       LOGICAL LSAME
00144       EXTERNAL LSAME
00145 *     ..
00146 *     .. External Subroutines ..
00147       EXTERNAL XERBLA
00148 *     ..
00149 *     .. Intrinsic Functions ..
00150       INTRINSIC MAX,MIN
00151 *     ..
00152 *
00153 *     Test the input parameters.
00154 *
00155       INFO = 0
00156       IF (.NOT.LSAME(UPLO,'U') .AND. .NOT.LSAME(UPLO,'L')) THEN
00157           INFO = 1
00158       ELSE IF (N.LT.0) THEN
00159           INFO = 2
00160       ELSE IF (K.LT.0) THEN
00161           INFO = 3
00162       ELSE IF (LDA.LT. (K+1)) THEN
00163           INFO = 6
00164       ELSE IF (INCX.EQ.0) THEN
00165           INFO = 8
00166       ELSE IF (INCY.EQ.0) THEN
00167           INFO = 11
00168       END IF
00169       IF (INFO.NE.0) THEN
00170           CALL XERBLA('SSBMV ',INFO)
00171           RETURN
00172       END IF
00173 *
00174 *     Quick return if possible.
00175 *
00176       IF ((N.EQ.0) .OR. ((ALPHA.EQ.ZERO).AND. (BETA.EQ.ONE))) RETURN
00177 *
00178 *     Set up the start points in  X  and  Y.
00179 *
00180       IF (INCX.GT.0) THEN
00181           KX = 1
00182       ELSE
00183           KX = 1 - (N-1)*INCX
00184       END IF
00185       IF (INCY.GT.0) THEN
00186           KY = 1
00187       ELSE
00188           KY = 1 - (N-1)*INCY
00189       END IF
00190 *
00191 *     Start the operations. In this version the elements of the array A
00192 *     are accessed sequentially with one pass through A.
00193 *
00194 *     First form  y := beta*y.
00195 *
00196       IF (BETA.NE.ONE) THEN
00197           IF (INCY.EQ.1) THEN
00198               IF (BETA.EQ.ZERO) THEN
00199                   DO 10 I = 1,N
00200                       Y(I) = ZERO
00201    10             CONTINUE
00202               ELSE
00203                   DO 20 I = 1,N
00204                       Y(I) = BETA*Y(I)
00205    20             CONTINUE
00206               END IF
00207           ELSE
00208               IY = KY
00209               IF (BETA.EQ.ZERO) THEN
00210                   DO 30 I = 1,N
00211                       Y(IY) = ZERO
00212                       IY = IY + INCY
00213    30             CONTINUE
00214               ELSE
00215                   DO 40 I = 1,N
00216                       Y(IY) = BETA*Y(IY)
00217                       IY = IY + INCY
00218    40             CONTINUE
00219               END IF
00220           END IF
00221       END IF
00222       IF (ALPHA.EQ.ZERO) RETURN
00223       IF (LSAME(UPLO,'U')) THEN
00224 *
00225 *        Form  y  when upper triangle of A is stored.
00226 *
00227           KPLUS1 = K + 1
00228           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
00229               DO 60 J = 1,N
00230                   TEMP1 = ALPHA*X(J)
00231                   TEMP2 = ZERO
00232                   L = KPLUS1 - J
00233                   DO 50 I = MAX(1,J-K),J - 1
00234                       Y(I) = Y(I) + TEMP1*A(L+I,J)
00235                       TEMP2 = TEMP2 + A(L+I,J)*X(I)
00236    50             CONTINUE
00237                   Y(J) = Y(J) + TEMP1*A(KPLUS1,J) + ALPHA*TEMP2
00238    60         CONTINUE
00239           ELSE
00240               JX = KX
00241               JY = KY
00242               DO 80 J = 1,N
00243                   TEMP1 = ALPHA*X(JX)
00244                   TEMP2 = ZERO
00245                   IX = KX
00246                   IY = KY
00247                   L = KPLUS1 - J
00248                   DO 70 I = MAX(1,J-K),J - 1
00249                       Y(IY) = Y(IY) + TEMP1*A(L+I,J)
00250                       TEMP2 = TEMP2 + A(L+I,J)*X(IX)
00251                       IX = IX + INCX
00252                       IY = IY + INCY
00253    70             CONTINUE
00254                   Y(JY) = Y(JY) + TEMP1*A(KPLUS1,J) + ALPHA*TEMP2
00255                   JX = JX + INCX
00256                   JY = JY + INCY
00257                   IF (J.GT.K) THEN
00258                       KX = KX + INCX
00259                       KY = KY + INCY
00260                   END IF
00261    80         CONTINUE
00262           END IF
00263       ELSE
00264 *
00265 *        Form  y  when lower triangle of A is stored.
00266 *
00267           IF ((INCX.EQ.1) .AND. (INCY.EQ.1)) THEN
00268               DO 100 J = 1,N
00269                   TEMP1 = ALPHA*X(J)
00270                   TEMP2 = ZERO
00271                   Y(J) = Y(J) + TEMP1*A(1,J)
00272                   L = 1 - J
00273                   DO 90 I = J + 1,MIN(N,J+K)
00274                       Y(I) = Y(I) + TEMP1*A(L+I,J)
00275                       TEMP2 = TEMP2 + A(L+I,J)*X(I)
00276    90             CONTINUE
00277                   Y(J) = Y(J) + ALPHA*TEMP2
00278   100         CONTINUE
00279           ELSE
00280               JX = KX
00281               JY = KY
00282               DO 120 J = 1,N
00283                   TEMP1 = ALPHA*X(JX)
00284                   TEMP2 = ZERO
00285                   Y(JY) = Y(JY) + TEMP1*A(1,J)
00286                   L = 1 - J
00287                   IX = JX
00288                   IY = JY
00289                   DO 110 I = J + 1,MIN(N,J+K)
00290                       IX = IX + INCX
00291                       IY = IY + INCY
00292                       Y(IY) = Y(IY) + TEMP1*A(L+I,J)
00293                       TEMP2 = TEMP2 + A(L+I,J)*X(IX)
00294   110             CONTINUE
00295                   Y(JY) = Y(JY) + ALPHA*TEMP2
00296                   JX = JX + INCX
00297                   JY = JY + INCY
00298   120         CONTINUE
00299           END IF
00300       END IF
00301 *
00302       RETURN
00303 *
00304 *     End of SSBMV .
00305 *
00306       END