yat/utility/BLAS_level3.h

Code
Comments
Other
Rev Date Author Line
3605 27 Jan 17 peter 1 #ifndef _theplu_yat_utility_blas_level3
3605 27 Jan 17 peter 2 #define _theplu_yat_utility_blas_level3
3605 27 Jan 17 peter 3
3605 27 Jan 17 peter 4 // $Id$
3605 27 Jan 17 peter 5
3605 27 Jan 17 peter 6 /*
3904 07 May 20 peter 7   Copyright (C) 2017, 2020 Peter Johansson
3605 27 Jan 17 peter 8
3605 27 Jan 17 peter 9   This file is part of the yat library, http://dev.thep.lu.se/yat
3605 27 Jan 17 peter 10
3605 27 Jan 17 peter 11   The yat library is free software; you can redistribute it and/or
3605 27 Jan 17 peter 12   modify it under the terms of the GNU General Public License as
3605 27 Jan 17 peter 13   published by the Free Software Foundation; either version 3 of the
3605 27 Jan 17 peter 14   License, or (at your option) any later version.
3605 27 Jan 17 peter 15
3605 27 Jan 17 peter 16   The yat library is distributed in the hope that it will be useful,
3605 27 Jan 17 peter 17   but WITHOUT ANY WARRANTY; without even the implied warranty of
3605 27 Jan 17 peter 18   MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
3605 27 Jan 17 peter 19   General Public License for more details.
3605 27 Jan 17 peter 20
3605 27 Jan 17 peter 21   You should have received a copy of the GNU General Public License
3605 27 Jan 17 peter 22   along with yat. If not, see <http://www.gnu.org/licenses/>.
3605 27 Jan 17 peter 23 */
3605 27 Jan 17 peter 24
3651 28 Jun 17 peter 25 #include "BasicMatrix.h"
3605 27 Jan 17 peter 26 #include "BLAS_utility.h"
3605 27 Jan 17 peter 27
3908 13 May 20 peter 28 #include "expression/MatrixBinary.h"
3908 13 May 20 peter 29 #include "expression/MatrixProduct.h"
3908 13 May 20 peter 30 #include "expression/ScaledMatrix.h"
3908 13 May 20 peter 31 #include "expression/TransposedMatrix.h"
3908 13 May 20 peter 32
3904 07 May 20 peter 33 #include <gsl/gsl_blas.h>
3904 07 May 20 peter 34 #include <gsl/gsl_cblas.h>
3605 27 Jan 17 peter 35 #include <gsl/gsl_matrix.h>
3605 27 Jan 17 peter 36
3605 27 Jan 17 peter 37 namespace theplu {
3605 27 Jan 17 peter 38 namespace yat {
3605 27 Jan 17 peter 39 namespace utility {
3605 27 Jan 17 peter 40
3605 27 Jan 17 peter 41   // This file defines operations using both Matrix (but not Vector)
3605 27 Jan 17 peter 42
3605 27 Jan 17 peter 43   /**
3605 27 Jan 17 peter 44      \brief Matrix addition operator
3605 27 Jan 17 peter 45
3605 27 Jan 17 peter 46      Elementwise addition of the elements ,\f$ result_{ij} = lhs_{ij}
3605 27 Jan 17 peter 47      + rhs_{ij} \; \forall i,j \f$.
3605 27 Jan 17 peter 48
3605 27 Jan 17 peter 49      \relates Matrix
3605 27 Jan 17 peter 50
3605 27 Jan 17 peter 51      \since new in yat 0.15
3605 27 Jan 17 peter 52
3605 27 Jan 17 peter 53      \return A Matrix Expression that can be converted to a Matrix
3605 27 Jan 17 peter 54   */
3605 27 Jan 17 peter 55   template<class Derived1, class Derived2>
3605 27 Jan 17 peter 56   expression::MatrixBinary<Derived1, Derived2, expression::Plus>
3605 27 Jan 17 peter 57   operator+(const BasicMatrix<Derived1>& lhs, const BasicMatrix<Derived2>& rhs)
3605 27 Jan 17 peter 58   {
3605 27 Jan 17 peter 59     return expression::MatrixBinary<Derived1, Derived2,
3605 27 Jan 17 peter 60                                     expression::Plus>(lhs, rhs);
3605 27 Jan 17 peter 61   }
3605 27 Jan 17 peter 62
3605 27 Jan 17 peter 63
3605 27 Jan 17 peter 64   /**
3605 27 Jan 17 peter 65      \brief Matrix subtraction operator
3605 27 Jan 17 peter 66
3605 27 Jan 17 peter 67      Elementwise addition of the elements ,\f$ result_{ij} = lhs_{ij}
3605 27 Jan 17 peter 68      - rhs_{ij} \; \forall i,j \f$.
3605 27 Jan 17 peter 69
3605 27 Jan 17 peter 70      \relates Matrix
3605 27 Jan 17 peter 71
3605 27 Jan 17 peter 72      \since new in yat 0.15
3605 27 Jan 17 peter 73
3605 27 Jan 17 peter 74      \return A Matrix Expression that can be converted to a Matrix
3605 27 Jan 17 peter 75   */
3605 27 Jan 17 peter 76   template<class Derived1, class Derived2>
3605 27 Jan 17 peter 77   expression::MatrixBinary<Derived1, Derived2, expression::Minus>
3605 27 Jan 17 peter 78   operator-(const BasicMatrix<Derived1>& lhs, const BasicMatrix<Derived2>& rhs)
3605 27 Jan 17 peter 79   {
3605 27 Jan 17 peter 80     return expression::MatrixBinary<Derived1, Derived2,
3605 27 Jan 17 peter 81                                     expression::Minus>(lhs, rhs);
3605 27 Jan 17 peter 82   }
3605 27 Jan 17 peter 83
3605 27 Jan 17 peter 84
3605 27 Jan 17 peter 85   /**
3605 27 Jan 17 peter 86      \brief Matrix multiplication operator
3605 27 Jan 17 peter 87
3605 27 Jan 17 peter 88      \relates Matrix
3605 27 Jan 17 peter 89
3605 27 Jan 17 peter 90      \since new in yat 0.15
3605 27 Jan 17 peter 91
3605 27 Jan 17 peter 92      \return A Matrix Expression that can be converted to a Matrix
3605 27 Jan 17 peter 93   */
3605 27 Jan 17 peter 94   template<class Derived1, class Derived2>
3908 13 May 20 peter 95   expression::MatrixProduct<BasicMatrix<Derived1>, BasicMatrix<Derived2> >
3605 27 Jan 17 peter 96   operator*(const BasicMatrix<Derived1>& lhs, const BasicMatrix<Derived2>& rhs)
3605 27 Jan 17 peter 97   {
3605 27 Jan 17 peter 98     YAT_ASSERT(lhs.columns() == rhs.rows());
3908 13 May 20 peter 99     return expression::MatrixProduct<BasicMatrix<Derived1>,
3908 13 May 20 peter 100                                      BasicMatrix<Derived2>
3908 13 May 20 peter 101                                      >(lhs, rhs);
3605 27 Jan 17 peter 102   }
3605 27 Jan 17 peter 103
3605 27 Jan 17 peter 104   /**
3605 27 Jan 17 peter 105      \relates Matrix
3605 27 Jan 17 peter 106
3605 27 Jan 17 peter 107      Calculate a scaled Matrix of \a A as \f$ result_{ij} = k * A_{ij}
3605 27 Jan 17 peter 108      \f$
3605 27 Jan 17 peter 109
3605 27 Jan 17 peter 110      \since new in yat 0.15
3605 27 Jan 17 peter 111
3605 27 Jan 17 peter 112      \return A Matrix Expression that can be converted to a Matrix
3605 27 Jan 17 peter 113   */
3605 27 Jan 17 peter 114   template<class T>
3605 27 Jan 17 peter 115   expression::ScaledMatrix<T>
3605 27 Jan 17 peter 116   operator*(const BasicMatrix<T>& A, double k)
3605 27 Jan 17 peter 117   {
3605 27 Jan 17 peter 118     return expression::ScaledMatrix<T>(k, A);
3605 27 Jan 17 peter 119   }
3605 27 Jan 17 peter 120
3605 27 Jan 17 peter 121
3605 27 Jan 17 peter 122   /**
3605 27 Jan 17 peter 123      \relates Matrix
3605 27 Jan 17 peter 124
3605 27 Jan 17 peter 125      Calculate a scaled Matrix of \a A as \f$ result_{ij} = k * A_{ij}
3605 27 Jan 17 peter 126      \f$
3605 27 Jan 17 peter 127
3605 27 Jan 17 peter 128      \since new in yat 0.15
3605 27 Jan 17 peter 129
3605 27 Jan 17 peter 130      \return A Matrix Expression that can be converted to a Matrix
3605 27 Jan 17 peter 131   */
3605 27 Jan 17 peter 132   template<class T>
3605 27 Jan 17 peter 133   expression::ScaledMatrix<T>
3605 27 Jan 17 peter 134   operator*(double k, const BasicMatrix<T>& A)
3605 27 Jan 17 peter 135   {
3605 27 Jan 17 peter 136     return expression::ScaledMatrix<T>(k, A);
3605 27 Jan 17 peter 137   }
3605 27 Jan 17 peter 138
3610 27 Jan 17 peter 139
3610 27 Jan 17 peter 140   /**
3610 27 Jan 17 peter 141      \brief negation operator
3610 27 Jan 17 peter 142
3624 10 Feb 17 peter 143      \relates Matrix
3624 10 Feb 17 peter 144
3610 27 Jan 17 peter 145      \since new in yat 0.15
3610 27 Jan 17 peter 146    */
3610 27 Jan 17 peter 147   template<typename T>
3610 27 Jan 17 peter 148   expression::ScaledMatrix<T>
3610 27 Jan 17 peter 149   operator-(const BasicMatrix<T>& m)
3610 27 Jan 17 peter 150   {
3610 27 Jan 17 peter 151     return expression::ScaledMatrix<T>(-1.0, m);
3610 27 Jan 17 peter 152   }
3610 27 Jan 17 peter 153
3654 10 Jul 17 peter 154
3654 10 Jul 17 peter 155   /**
3904 07 May 20 peter 156      Specialization for ScaledMatrix
3904 07 May 20 peter 157
3904 07 May 20 peter 158      \since New in yat 0.18
3904 07 May 20 peter 159    */
3904 07 May 20 peter 160   template<class T>
3904 07 May 20 peter 161   expression::ScaledMatrix<T>
3904 07 May 20 peter 162   operator*(const BasicMatrix<MatrixExpression<expression::ScaledMatrix<T> > >& A,
3904 07 May 20 peter 163             double k)
3904 07 May 20 peter 164   {
3904 07 May 20 peter 165     const expression::ScaledMatrix<T>& sm =
3904 07 May 20 peter 166       static_cast<const expression::ScaledMatrix<T>&>(A);
3904 07 May 20 peter 167     return expression::ScaledMatrix<T>(k*sm.factor(), sm.base());
3904 07 May 20 peter 168   }
3904 07 May 20 peter 169
3904 07 May 20 peter 170
3904 07 May 20 peter 171   /**
3904 07 May 20 peter 172      Specialization for ScaledMatrix
3904 07 May 20 peter 173
3904 07 May 20 peter 174      \since New in yat 0.18
3904 07 May 20 peter 175    */
3904 07 May 20 peter 176   template<class T>
3904 07 May 20 peter 177   expression::ScaledMatrix<T>
3904 07 May 20 peter 178   operator*(double k,
3904 07 May 20 peter 179             const BasicMatrix<MatrixExpression<expression::ScaledMatrix<T> > >& A)
3904 07 May 20 peter 180   {
3904 07 May 20 peter 181     const expression::ScaledMatrix<T>& sm =
3904 07 May 20 peter 182       static_cast<const expression::ScaledMatrix<T>&>(A);
3904 07 May 20 peter 183     return expression::ScaledMatrix<T>(k*sm.factor(), sm.base());
3904 07 May 20 peter 184   }
3904 07 May 20 peter 185
3904 07 May 20 peter 186
3904 07 May 20 peter 187   /**
3904 07 May 20 peter 188      Specialization for ScaledMatrix
3904 07 May 20 peter 189
3904 07 May 20 peter 190      \since New in yat 0.18
3904 07 May 20 peter 191    */
3904 07 May 20 peter 192   template<class T>
3904 07 May 20 peter 193   expression::ScaledMatrix<T>
3904 07 May 20 peter 194   operator-(const BasicMatrix<MatrixExpression<expression::ScaledMatrix<T> > >& A)
3904 07 May 20 peter 195   {
3904 07 May 20 peter 196     const expression::ScaledMatrix<T>& sm =
3904 07 May 20 peter 197       static_cast<const expression::ScaledMatrix<T>&>(A);
3904 07 May 20 peter 198     return expression::ScaledMatrix<T>(-sm.factor(), sm.base());
3904 07 May 20 peter 199   }
3904 07 May 20 peter 200
3904 07 May 20 peter 201
3904 07 May 20 peter 202   /**
3654 10 Jul 17 peter 203      \brief transpose function
3654 10 Jul 17 peter 204
3654 10 Jul 17 peter 205      \relates Matrix
3654 10 Jul 17 peter 206
3654 10 Jul 17 peter 207      \since new in yat 0.15
3654 10 Jul 17 peter 208    */
3654 10 Jul 17 peter 209   template<typename T>
3654 10 Jul 17 peter 210   expression::TransposedMatrix<T> transpose(const BasicMatrix<T>& A)
3654 10 Jul 17 peter 211   {
3654 10 Jul 17 peter 212     return expression::TransposedMatrix<T>(A);
3654 10 Jul 17 peter 213   }
3654 10 Jul 17 peter 214
3605 27 Jan 17 peter 215 }}} // of namespace utility, yat, and theplu
3605 27 Jan 17 peter 216
3605 27 Jan 17 peter 217 #endif