VectorBase.hpp

#pragma once
 
#include <iterator>
 
#include <dune/common/rangeutilities.hh>
#include <dune/functions/common/multiindex.hh>
 
#include <amdis/algorithm/ForEach.hpp>
#include <amdis/algorithm/InnerProduct.hpp>
#include <amdis/algorithm/Transform.hpp>
#include <amdis/common/FakeContainer.hpp>
#include <amdis/utility/MappedRangeView.hpp>
 
namespace AMDiS
{
  template <class Derived>
  class VectorBase
  {
  private:
    Derived const& asDerived() const
    {
      return static_cast<Derived const&>(*this);
    }
 
    Derived& asDerived()
    {
      return static_cast<Derived&>(*this);
    }
 
  protected:
    VectorBase() = default;
 
  public:
    void finish() { /* do nothing */ }
    void synchronize() { /* do nothing */ }
 
    template <class MultiIndex, class ValueType, class Assign>
    void insert(MultiIndex const& idx, ValueType const& value, Assign assign)
    {
      assign(value, asDerived().at(idx));
    }
 
    template <class IndexRange, class OutputIterator>
    void gather(IndexRange const& localInd, OutputIterator buffer) const
    {
      for (auto const& idx : localInd)
        *buffer++ = asDerived().at(idx);
    }
 
 
    template <class IndexRange, class LocalVec, class Assign>
    void scatter(IndexRange const& localInd, LocalVec const& vec, FakeContainer<bool,true>, Assign assign)
    {
      auto vec_it = std::begin(vec);
      for (auto const& idx : localInd)
        assign(*vec_it++, asDerived().at(idx));
    }
 
    template <class IndexRange, class LocalVec, class MaskRange, class Assign>
    void scatter(IndexRange const& localInd, LocalVec const& vec, MaskRange const& mask, Assign assign)
    {
      auto vec_it = std::begin(vec);
      auto mask_it = std::begin(mask);
      auto ind_it = std::begin(localInd);
      for (; vec_it != std::end(vec); ++vec_it, ++mask_it, ++ind_it) {
        if (*mask_it)
          assign(*vec_it, asDerived().at(*ind_it));
      }
    }
 
 
    template <class IndexRange, class Func>
    void forEach(IndexRange const& localInd, Func&& f) const
    {
      for (auto const& idx : localInd)
        f(idx, asDerived().at(idx));
    }
 
    template <class IndexRange, class Func>
    void forEach(IndexRange const& localInd, Func&& f)
    {
      for (auto const& idx : localInd)
        f(idx, asDerived().at(idx));
    }
  };
 
 
  namespace Recursive
  {
    template <class D>
    struct ForEach<VectorBase<D>>
    {
      template <class UnaryFunction>
      static void impl (D& vec, UnaryFunction f)
      {
        for (std::size_t i = 0; i < std::size_t(vec.size()); ++i) {
          auto idx = Dune::Functions::StaticMultiIndex<std::size_t,1>{i};
          Recursive::forEach(vec.at(idx), f);
        }
      }
    };
 
    template <class D>
    struct Transform<VectorBase<D>>
    {
      template <class Operation, class... Ds>
      static void impl (D& vecOut, Operation op, VectorBase<Ds> const&... vecIn)
      {
        for (std::size_t i = 0; i < std::size_t(vecOut.size()); ++i) {
          auto idx = Dune::Functions::StaticMultiIndex<std::size_t,1>{i};
          Recursive::transform(vecOut.at(idx), op, static_cast<Ds const&>(vecIn).at(idx)...);
        }
      }
    };
 
    template <class D>
    struct InnerProduct<VectorBase<D>>
    {
      template <class D2, class T, class BinOp1, class BinOp2>
      static T impl (D const& in1, VectorBase<D2> const& in2, T init, BinOp1 op1, BinOp2 op2)
      {
        for (std::size_t i = 0; i < std::size_t(in1.size()); ++i) {
          auto idx = Dune::Functions::StaticMultiIndex<std::size_t,1>{i};
          init = Recursive::innerProduct(in1.at(idx), static_cast<D2 const&>(in2).at(idx), std::move(init), op1, op2);
        }
        return init;
      }
    };
 
  } // end namespace Recursive
} // end namespace AMDiS