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_LEVEL3_SYMM_TCC
 34 #define CXXBLAS_LEVEL3_SYMM_TCC 1
 35 
 36 namespace cxxblas {
 37 
 38 template <typename IndexType, typename ALPHA, typename MA, typename MB,
 39           typename BETA, typename MC>
 40 void
 41 symm_generic(StorageOrder order, Side sideA, StorageUpLo upLoA,
 42              IndexType m, IndexType n,
 43              const ALPHA &alpha,
 44              const MA *A, IndexType ldA,
 45              const MB *B, IndexType ldB,
 46              const BETA &beta,
 47              MC *C, IndexType ldC)
 48 {
 49     if (order==ColMajor) {
 50         upLoA = (upLoA==Upper) ? Lower : Upper;
 51         sideA = (sideA==Left) ? Right : Left;
 52         symm_generic(RowMajor, sideA, upLoA, n, m,
 53                      alpha, A, ldA, B, ldB,
 54                      beta,
 55                      C, ldC);
 56         return;
 57     }
 58     gescal(order, m, n, beta, C, ldC);
 59     if (sideA==Right) {
 60         for (IndexType i=0; i<m; ++i) {
 61             symv(order, upLoA, n, alpha, A, ldA, B+i*ldB, IndexType(1),
 62                  BETA(1), C+i*ldC, IndexType(1));
 63         }
 64     }
 65     if (sideA==Left) {
 66         for (IndexType j=0; j<n; ++j) {
 67             symv(order, upLoA, m, alpha, A, ldA, B+j, ldB,
 68                  BETA(1), C+j, ldC);
 69         }
 70     }
 71 }
 72 
 73 template <typename IndexType, typename ALPHA, typename MA, typename MB,
 74           typename BETA, typename MC>
 75 void
 76 symm(StorageOrder order, Side side, StorageUpLo upLo,
 77      IndexType m, IndexType n,
 78      const ALPHA &alpha,
 79      const MA *A, IndexType ldA,
 80      const MB *B, IndexType ldB,
 81      const BETA &beta,
 82      MC *C, IndexType ldC)
 83 {
 84     CXXBLAS_DEBUG_OUT("symm_generic");
 85 
 86     symm_generic(order, side, upLo, m, n, alpha, A, ldA, B, ldB, beta, C, ldC);
 87 }
 88 
 89 #ifdef HAVE_CBLAS
 90 
 91 // ssymm
 92 template <typename IndexType>
 93 typename If<IndexType>::isBlasCompatibleInteger
 94 symm(StorageOrder order, Side side, StorageUpLo upLo,
 95      IndexType m, IndexType n,
 96      float alpha,
 97      const float *A, IndexType ldA,
 98      const float *B, IndexType ldB,
 99      float beta,
100      float *C, IndexType ldC)
101 {
102     CXXBLAS_DEBUG_OUT("[" BLAS_IMPL "] cblas_ssymm");
103 
104     cblas_ssymm(CBLAS::getCblasType(order),
105                 CBLAS::getCblasType(side), CBLAS::getCblasType(upLo),
106                 m, n,
107                 alpha,
108                 A, ldA,
109                 B, ldB,
110                 beta,
111                 C, ldC);
112 }
113 
114 // dsymm
115 template <typename IndexType>
116 typename If<IndexType>::isBlasCompatibleInteger
117 symm(StorageOrder order, Side side, StorageUpLo upLo,
118      IndexType m, IndexType n,
119      double alpha,
120      const double *A, IndexType ldA,
121      const double *B, IndexType ldB,
122      double beta,
123      double *C, IndexType ldC)
124 {
125     CXXBLAS_DEBUG_OUT("[" BLAS_IMPL "] cblas_dsymm");
126 
127     cblas_dsymm(CBLAS::getCblasType(order),
128                 CBLAS::getCblasType(side), CBLAS::getCblasType(upLo),
129                 m, n,
130                 alpha,
131                 A, ldA,
132                 B, ldB,
133                 beta,
134                 C, ldC);
135 }
136 
137 // csymm
138 template <typename IndexType>
139 typename If<IndexType>::isBlasCompatibleInteger
140 symm(StorageOrder order, Side side, StorageUpLo upLo,
141      IndexType m, IndexType n,
142      const ComplexFloat &alpha,
143      const ComplexFloat *A, IndexType ldA,
144      const ComplexFloat *B, IndexType ldB,
145      const ComplexFloat &beta,
146      ComplexFloat *C, IndexType ldC)
147 {
148     CXXBLAS_DEBUG_OUT("[" BLAS_IMPL "] cblas_csymm");
149 
150     cblas_csymm(CBLAS::getCblasType(order),
151                 CBLAS::getCblasType(side), CBLAS::getCblasType(upLo),
152                 m, n,
153                 reinterpret_cast<const float *>(&alpha),
154                 reinterpret_cast<const float *>(A), ldA,
155                 reinterpret_cast<const float *>(B), ldB,
156                 reinterpret_cast<const float *>(&beta),
157                 reinterpret_cast<const float *>(C), ldC);
158 }
159 
160 // zsymm
161 template <typename IndexType>
162 typename If<IndexType>::isBlasCompatibleInteger
163 symm(StorageOrder order, Side side, StorageUpLo upLo,
164      IndexType m, IndexType n,
165      const ComplexDouble &alpha,
166      const ComplexDouble *A, IndexType ldA,
167      const ComplexDouble *B, IndexType ldB,
168      const ComplexDouble &beta,
169      ComplexDouble *C, IndexType ldC)
170 {
171     CXXBLAS_DEBUG_OUT("[" BLAS_IMPL "] cblas_zsymm");
172 
173     cblas_zsymm(CBLAS::getCblasType(order),
174                 CBLAS::getCblasType(side), CBLAS::getCblasType(upLo),
175                 m, n,
176                 reinterpret_cast<const double *>(&alpha),
177                 reinterpret_cast<const double *>(A), ldA,
178                 reinterpret_cast<const double *>(B), ldB,
179                 reinterpret_cast<const double *>(&beta),
180                 reinterpret_cast<const double *>(C), ldC);
181 }
182 
183 #endif // HAVE_CBLAS
184 
185 // namespace cxxblas
186 
187 #endif // CXXBLAS_LEVEL3_SYMM_TCC