5#ifndef __DUNE_ACFEM_SUBENTITYMAPPER_HH__
6#define __DUNE_ACFEM_SUBENTITYMAPPER_HH__
8#include <dune/geometry/referenceelements.hh>
9#include <dune/grid/common/intersection.hh>
10#include <dune/common/version.hh>
12#include <unordered_map>
13#include "../mpl/foreach.hh"
31 template<
class DofMapper,
int coDimension>
32 struct MapEntityHierarchyHelper
34 typedef DofMapper DofMapperType;
35 typedef typename DofMapperType::GridPartType GridPartType;
36 typedef typename GridPartType::ctype FieldType;
37 typedef typename DofMapperType::ElementType ElementType;
38 typedef typename DofMapperType::GlobalKeyType GlobalKeyType;
39 typedef ReferenceElements<FieldType, ElementType::dimension> ReferenceElementsType;
40 typedef typename ReferenceElementsType::ReferenceElement ReferenceElementType;
44 struct SubSubEntityFunctor
47 subCoDimension = codim,
48 gridDimension = GridPartType::dimension,
49 subSubEntityDimension = gridDimension-subCoDimension
51 template<
class RefElem,
class Functor>
52 static void apply(
const DofMapperType& mapper,
const ElementType& element,
const RefElem& refElem,
int subIndex, Functor f,
53 const std::unordered_map<GlobalKeyType, unsigned>& globalToLocal)
55 const int numSubSubEntities = refElem.size(subIndex, coDimension, subCoDimension);
57 for (
int i = 0; i < numSubSubEntities; ++i) {
58 const unsigned subSubEntityIndex = refElem.subEntity(subIndex, coDimension, i, subCoDimension);
59 const auto subSubEntity = element.template subEntity<subCoDimension>(subSubEntityIndex);
61 mapper.mapEachEntityDof(subSubEntity, [f, globalToLocal] (
unsigned localIndex,
const GlobalKeyType& globalIndex) {
62 f(globalToLocal.at(globalIndex), globalIndex);
92 template<
class DofMapper,
int coDimension,
class Functor>
94 const typename DofMapper::ElementType& element,
96 const std::integral_constant<int, coDimension>& dummy,
99 typedef DofMapper DofMapperType;
100 typedef MapEntityHierarchyHelper<DofMapperType, coDimension> HelperType;
101 typedef typename HelperType::GridPartType GridPartType;
102 typedef typename HelperType::GlobalKeyType GlobalKeyType;
103 typedef typename HelperType::ReferenceElementsType ReferenceElementsType;
105 enum { dimension = GridPartType::dimension };
108 std::unordered_map<GlobalKeyType, unsigned> globalToLocal;
109 mapper.mapEach(element, [&globalToLocal] (
unsigned local,
const GlobalKeyType& global) {
110 globalToLocal[global] = local;
114 const auto& refElem = ReferenceElementsType::general(element.type());
118 forLoop<dimension-coDimension+1>([&](
auto i){
119 HelperType::template SubSubEntityFunctor<i.value+coDimension>::apply(mapper, element, refElem, subIndex, f, globalToLocal);
142 template<
class DofMapper,
class Intersection,
class Functor>
144 const Intersection& intersection,
148 std::integral_constant<int, 1>(), f);
172 template<
class DofMapper,
class Intersection,
class LocalIndices,
class GlobalIndices>
173 unsigned mapFaceDofs(
const DofMapper& mapper,
const Intersection& intersection,
174 LocalIndices& localIndices, GlobalIndices& globalIndices)
176 typedef typename DofMapper::GlobalKeyType GlobalKeyType;
180 [&localIndices, &globalIndices, &count] (
unsigned localIndex, GlobalKeyType globalIndex) {
181 localIndices[count] = localIndex;
182 globalIndices[count] = globalIndex;
void mapEachFaceDof(const DofMapper &mapper, const Intersection &intersection, Functor f)
Iterate over the codim-1 sub-entity hierarchy linked to the given intersection and call the dof-mappe...
Definition: subentitymapper.hh:143
void mapEntityHierarchy(const DofMapper &mapper, const typename DofMapper::ElementType &element, int subIndex, const std::integral_constant< int, coDimension > &dummy, Functor f)
Iterate over the entire sub-entity hierarchy below one given sub-entity (including the given one) and...
Definition: subentitymapper.hh:93
unsigned mapFaceDofs(const DofMapper &mapper, const Intersection &intersection, LocalIndices &localIndices, GlobalIndices &globalIndices)
Fetch all global DoFs of the codim-1 entity intersection belongs to.
Definition: subentitymapper.hh:173