transpose.hh

// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_COMMON_TRANSPOSE_HH
#define DUNE_COMMON_TRANSPOSE_HH
 
#include <cstddef>
 
#include <dune/common/fmatrix.hh>
#include <dune/common/promotiontraits.hh>
 
namespace Dune {
 
namespace Impl {
 
  // Wrapper representing the transposed of a matrix.
  // Creating the wrapper does not compute anything
  // but only serves for tagging the wrapped matrix
  // for transposition.
  template<class M>
  class TransposedMatrixWrapper
  {
  public:
 
    enum {
      rows = M::cols,
      cols = M::rows
    };
 
    TransposedMatrixWrapper(const M& matrix) : matrix_(matrix) {}
    TransposedMatrixWrapper(const TransposedMatrixWrapper&) = delete;
    TransposedMatrixWrapper(TransposedMatrixWrapper&&) = delete;
 
    template<class OtherField, int otherRows>
    friend auto operator* (const FieldMatrix<OtherField, otherRows, rows>& matrixA,
                            const TransposedMatrixWrapper& matrixB)
    {
      using ThisField = typename FieldTraits<M>::field_type;
      using Field = typename PromotionTraits<ThisField, OtherField>::PromotedType;
      FieldMatrix<Field, otherRows, cols> result;
      for (std::size_t j=0; j<otherRows; ++j)
        matrixB.matrix_.mv(matrixA[j], result[j]);
      return result;
    }
 
  private:
 
    const M& matrix_;
  };
 
} // namespace Impl
 
template<class Matrix>
auto transpose(const Matrix& matrix) {
  return Impl::TransposedMatrixWrapper<Matrix>(matrix);
}
 
 
} // namespace Dune
 
#endif // DUNE_COMMON_TRANSPOSE_HH