12#ifndef DUMUX_DISCRETIZATION_LOCAL_INTERSECTION_INDEX_MAPPER_HH 
   13#define DUMUX_DISCRETIZATION_LOCAL_INTERSECTION_INDEX_MAPPER_HH 
   17#include <dune/common/float_cmp.hh> 
   28template<
class Gr
idView, 
bool consistentlyOrientedGr
id>
 
   37template<
class Gr
idView>
 
   41    using Element = 
typename GridView::template Codim<0>::Entity;
 
   42    static constexpr auto numElementFaces = GridView::Grid::dimension * 2;
 
   44    void update(
const GridView& gv, 
const Element& element)
 
   46        static const bool makeConsistentlyOriented = 
getParam<bool>(
"Grid.MakeConsistentlyOriented", 
true);
 
   47        if (!makeConsistentlyOriented)
 
   49            std::iota(realToRefMap_.begin(), realToRefMap_.end(), 0);
 
   50            refToRealMap_ = realToRefMap_;
 
   54        for (
const auto& is : intersections(gv, element))
 
   56            const auto& otherOuterNormal = is.centerUnitOuterNormal();
 
   58            const int positveOrientation = !std::signbit(otherOuterNormal[idx]);
 
   59            const auto refIdx = idx * 2 + positveOrientation;
 
   60            const auto realIdx = is.indexInInside();
 
   61            realToRefMap_[realIdx] = refIdx;
 
   62            refToRealMap_[refIdx] = realIdx;
 
 
   67    SmallLocalIndexType 
realToRefIdx(
const SmallLocalIndexType localIsIdx)
 const 
   68    { 
return realToRefMap_[localIsIdx]; }
 
 
   71    SmallLocalIndexType 
refToRealIdx(
const SmallLocalIndexType localIsIdx)
 const 
   72    { 
return refToRealMap_[localIsIdx]; }
 
 
   75    std::array<SmallLocalIndexType, numElementFaces> realToRefMap_ = {};
 
   76    std::array<SmallLocalIndexType, numElementFaces> refToRealMap_ = {};
 
 
   87template<
class Gr
idView>
 
   91    using Element = 
typename GridView::template Codim<0>::Entity;
 
   95    void update(
const GridView&, 
const Element&) {}
 
   98    SmallLocalIndexType 
realToRefIdx(
const SmallLocalIndexType localIsIdx)
 const 
   99    { 
return localIsIdx; }
 
 
  102    SmallLocalIndexType 
refToRealIdx(
const SmallLocalIndexType localIsIdx)
 const 
  103    { 
return localIsIdx; }
 
 
 
  115template<
class Gr
idView>
 
SmallLocalIndexType refToRealIdx(const SmallLocalIndexType localIsIdx) const
Return the intersection's local reference indexInElement given an actual local index.
Definition localintersectionindexmapper.hh:71
void update(const GridView &gv, const Element &element)
Definition localintersectionindexmapper.hh:44
SmallLocalIndexType realToRefIdx(const SmallLocalIndexType localIsIdx) const
Return the intersection's actual local indexInElement given a local reference index.
Definition localintersectionindexmapper.hh:67
SmallLocalIndexType realToRefIdx(const SmallLocalIndexType localIsIdx) const
Return the intersection's actual local indexInElement given a local reference index.
Definition localintersectionindexmapper.hh:98
SmallLocalIndexType refToRealIdx(const SmallLocalIndexType localIsIdx) const
Return the intersection's local reference indexInElement given an actual local index.
Definition localintersectionindexmapper.hh:102
void update(const GridView &, const Element &)
Update the map for getting the corresponding local face indices in another element.
Definition localintersectionindexmapper.hh:95
Definition localintersectionindexmapper.hh:29
Helper type to determine whether a grid is guaranteed to be oriented consistently....
Detail::FaceCenteredStaggeredLocalIntersectionIndexMapper< GridView, ConsistentlyOrientedGrid< typename GridView::Grid >{}> FaceCenteredStaggeredLocalIntersectionIndexMapper
Provides a mapping of local intersection indices (indexInInside) such that the local indices always f...
Definition localintersectionindexmapper.hh:116
static std::size_t normalAxis(const Vector &v)
Returns the normal axis index of a unit vector (0 = x, 1 = y, 2 = z)
Definition normalaxis.hh:26
T getParam(Args &&... args)
A free function to get a parameter from the parameter tree singleton.
Definition parameters.hh:139
Define some often used mathematical functions.
Distance implementation details.
Definition cvfelocalresidual.hh:25
Base class for the finite volume geometry vector for face-centered staggered models This builds up th...
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
std::uint_least8_t SmallLocalIndex
Definition indextraits.hh:29