1 /*
  2  *   Copyright (c) 2009, Michael Lehn
  3  *
  4  *   All rights reserved.
  5  *
  6  *   Redistribution and use in source and binary forms, with or without
  7  *   modification, are permitted provided that the following conditions
  8  *   are met:
  9  *
 10  *   1) Redistributions of source code must retain the above copyright
 11  *      notice, this list of conditions and the following disclaimer.
 12  *   2) Redistributions in binary form must reproduce the above copyright
 13  *      notice, this list of conditions and the following disclaimer in
 14  *      the documentation and/or other materials provided with the
 15  *      distribution.
 16  *   3) Neither the name of the FLENS development group nor the names of
 17  *      its contributors may be used to endorse or promote products derived
 18  *      from this software without specific prior written permission.
 19  *
 20  *   THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS AND CONTRIBUTORS
 21  *   "AS IS" AND ANY EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT
 22  *   LIMITED TO, THE IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR
 23  *   A PARTICULAR PURPOSE ARE DISCLAIMED. IN NO EVENT SHALL THE COPYRIGHT
 24  *   OWNER OR CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL,
 25  *   SPECIAL, EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT
 26  *   LIMITED TO, PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE,
 27  *   DATA, OR PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY
 28  *   THEORY OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
 29  *   (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
 30  *   OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.
 31  */
 32 
 33 #ifndef CXXBLAS_LEVEL2_SBMV_TCC
 34 #define CXXBLAS_LEVEL2_SBMV_TCC 1
 35 
 36 #include <complex>
 37 #include <cxxblas/level1/level1.h>
 38 
 39 namespace cxxblas {
 40 
 41 template <typename IndexType, typename ALPHA, typename MA, typename VX,
 42           typename BETA, typename VY>
 43 void
 44 sbmv_generic(StorageOrder order, StorageUpLo upLo,
 45              IndexType n, IndexType k,
 46              const ALPHA &alpha,
 47              const MA *A, IndexType ldA,
 48              const VX *x, IndexType incX,
 49              const BETA &beta,
 50              VY *y, IndexType incY)
 51 {
 52     using std::min;
 53     using std::max;
 54 
 55     if (order==ColMajor) {
 56         upLo = (upLo==Upper) ? Lower : Upper;
 57         sbmv_generic(RowMajor, upLo, n, k, alpha, A, ldA,
 58                      x, incX, beta, y, incY);
 59         return;
 60     }
 61 
 62     scal_generic(n, beta, y, incY);
 63     if (upLo==Upper) {
 64         for (IndexType i=0, iX=0, iY=0; i<n; ++i, iX+=incX, iY+=incY) {
 65             IndexType len = min(k+1, n-i);
 66 
 67             VY _y;
 68             dot_generic(len, A+ldA*i, IndexType(1),
 69                              x+iX, IndexType(incX),
 70                              _y);
 71             y[iY] += alpha*_y;
 72 
 73             axpy_generic(len-1, x[iX] * alpha,
 74                                 A+ldA*i+1, IndexType(1),
 75                                 y+iY+incY, incY);
 76         }
 77     } else {
 78         for (IndexType i=0, iY=0; i<n; ++i, iY+=incY) {
 79             IndexType iA = max(k-i, IndexType(0));
 80             IndexType len = min(k, i) + 1;
 81             IndexType _i = max(i-k, IndexType(0));
 82 
 83             VY _y;
 84             dot_generic(len, A+ldA*i+iA, IndexType(1),
 85                              x+_i*incX, IndexType(incX),
 86                              _y);
 87             y[iY] += alpha*_y;
 88 
 89             axpy_generic(len-1, x[i*incX] * alpha,
 90                                 A+ldA*i+iA, IndexType(1),
 91                                 y+_i*incY, incY);
 92         }
 93     }
 94 }
 95 
 96 //------------------------------------------------------------------------------
 97 
 98 template <typename IndexType, typename ALPHA, typename MA, typename VX,
 99           typename BETA, typename VY>
100 void
101 sbmv(StorageOrder order, StorageUpLo upLo,
102      IndexType n, IndexType k,
103      const ALPHA &alpha,
104      const MA *A, IndexType ldA,
105      const VX *x, IndexType incX,
106      const BETA &beta,
107      VY *y, IndexType incY)
108 {
109     CXXBLAS_DEBUG_OUT("sbmv_generic");
110 
111     if (n==0) {
112         return;
113     }
114     if (incX<0) {
115         x -= incX*(n-1);
116     }
117     if (incY<0) {
118         y -= incY*(n-1);
119     }
120     sbmv_generic(order, upLo, n, k, alpha, A, ldA, x, incX, beta, y, incY);
121 }
122 
123 // namespace cxxblas
124 
125 #endif // CXXBLAS_LEVEL2_SBMV_TCC