subentitydofs.hh

// -*- tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 2 -*-
// vi: set et ts=4 sw=2 sts=2:
#ifndef DUNE_FUNCTIONS_FUNCTIONSPACEBASES_SUBENTITYDOFS_HH
#define DUNE_FUNCTIONS_FUNCTIONSPACEBASES_SUBENTITYDOFS_HH
 
#include <vector>
 
#include <dune/geometry/referenceelements.hh>
#include <dune/typetree/traversal.hh>
 
 
 
namespace Dune {
namespace Functions {
 
 
 
template<class GridView>
class SubEntityDOFs
{
  static const int dim = GridView::dimension;
 
public:
 
  template<class LocalView>
  SubEntityDOFs& bind(const LocalView& localView, std::size_t subEntityIndex, std::size_t subEntityCodim)
  {
    // fill vector with local indices of all DOFs contained in subentity
    containedDOFs_.clear();
    dofIsContained_.assign(localView.size(), false);
 
    auto re = Dune::referenceElement<double,dim>(localView.element().type());
 
    Dune::TypeTree::forEachLeafNode(localView.tree(), [&](auto&& node, auto&& treePath) {
      const auto& localCoefficients = node.finiteElement().localCoefficients();
      std::size_t localSize = localCoefficients.size();
      for(std::size_t i=0; i<localSize; ++i)
      {
        auto localKey = localCoefficients.localKey(i);
        if (re.subEntities(subEntityIndex, subEntityCodim, localKey.codim()).contains(localKey.subEntity()))
        {
          containedDOFs_.push_back(node.localIndex(i));
          dofIsContained_[node.localIndex(i)] = true;
        }
      }
    });
    return *this;
  }
 
  template<class LocalView, class Intersection>
  SubEntityDOFs& bind(const LocalView& localView, const Intersection& intersection)
  {
    return bind(localView, intersection.indexInInside(), 1);
  }
 
  auto begin() const
  {
    return containedDOFs_.cbegin();
  }
 
  auto end() const
  {
    return containedDOFs_.cend();
  }
 
  auto size() const
  {
    return containedDOFs_.size();
  }
 
  decltype(auto) operator[](std::size_t i) const
  {
    return containedDOFs_[i];
  }
 
  bool contains(std::size_t localIndex) const
  {
    return dofIsContained_[localIndex];
  }
 
private:
 
  std::vector<std::size_t> containedDOFs_;
  std::vector<bool> dofIsContained_;
};
 
 
 
template<class T>
auto subEntityDOFs(const T& t)
{
  return SubEntityDOFs<typename T::GridView>{};
}
 
 
 
template<class LocalView>
auto subEntityDOFs(const LocalView& localView, std::size_t subEntityIndex, std::size_t subEntityCodim)
{
  using GridView = typename LocalView::GridView;
  SubEntityDOFs<GridView> subEntityDOFs;
  subEntityDOFs.bind(localView, subEntityIndex, subEntityCodim);
  return subEntityDOFs;
}
 
 
 
template<class LocalView, class Intersection>
auto subEntityDOFs(const LocalView& localView, const Intersection& intersection)
{
  using GridView = typename LocalView::GridView;
  SubEntityDOFs<GridView> subEntityDOFs;
  subEntityDOFs.bind(localView, intersection);
  return subEntityDOFs;
}
 
 
 
} // namespace Functions
} // namespace Dune
 
#endif // DUNE_FUNCTIONS_FUNCTIONSPACEBASES_SUBENTITYDOFS_HH