1 /*
  2  *   Copyright (c) 2010, 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 FLENS_BLAS_LEVEL3_RK_TCC
 34 #define FLENS_BLAS_LEVEL3_RK_TCC
 35 
 36 #include <flens/typedefs.h>
 37 
 38 namespace flens { namespace blas {
 39 
 40 //-- herk
 41 template <typename ALPHA, typename MA, typename BETA, typename MC>
 42 void
 43 rk(Transpose trans,
 44    const ALPHA &alpha,
 45    const GeMatrix<MA> &A,
 46    const BETA &beta,
 47    HeMatrix<MC> &C)
 48 {
 49     ASSERT(MC::order==MA::order);
 50     typedef typename HeMatrix<MC>::IndexType IndexType;
 51 
 52     IndexType n = (trans==NoTrans) ? A.numRows()
 53                                             : A.numCols();
 54     IndexType k = (trans==NoTrans) ? A.numCols()
 55                                             : A.numRows();
 56 
 57     ASSERT((beta==static_cast<BETA>(0)) || (C.dim()==n));
 58     if (C.dim()!=n) {
 59         C.resize(n, n);
 60     }
 61 
 62     ASSERT(C.dim()==((trans==NoTrans) ? A.numRows() : A.numCols()));
 63 
 64 #   ifdef HAVE_CXXBLAS_HERK
 65     cxxblas::herk(MC::order, C.upLo(),
 66                   trans, n, k, alpha,
 67                   A.data(), A.leadingDimension(),
 68                   beta,
 69                   C.data(), C.leadingDimension());
 70 #   else
 71     ASSERT(0);
 72 #   endif
 73 }
 74 
 75 //-- syrk
 76 template <typename ALPHA, typename MA, typename BETA, typename MC>
 77 void
 78 rk(Transpose trans,
 79    const ALPHA &alpha,
 80    const GeMatrix<MA> &A,
 81    const BETA &beta,
 82    SyMatrix<MC> &C)
 83 {
 84     ASSERT(MC::order==MA::order);
 85 
 86     typedef typename SyMatrix<MC>::IndexType IndexType;
 87 
 88     IndexType n = (trans==NoTrans) ? A.numRows()
 89                                             : A.numCols();
 90     IndexType k = (trans==NoTrans) ? A.numCols()
 91                                             : A.numRows();
 92 
 93     ASSERT((beta==static_cast<BETA>(0)) || (C.dim()==n));
 94     if (C.dim()!=n) {
 95         C.resize(n, n);
 96     }
 97 
 98     ASSERT(C.dim()==((trans==NoTrans) ? A.numRows() : A.numCols()));
 99 
100 #   ifdef HAVE_CXXBLAS_SYRK
101     cxxblas::syrk(MC::order, C.upLo(),
102                   trans, n, k, alpha,
103                   A.data(), A.leadingDimension(),
104                   beta,
105                   C.data(), C.leadingDimension());
106 #   else
107     ASSERT(0);
108 #   endif
109 }
110 
111 } } // namespace blas, flens
112 
113 #endif // FLENS_BLAS_LEVEL3_RK_TCC