12#ifndef DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_GEOMETRY_HELPER_HH
13#define DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_GEOMETRY_HELPER_HH
17#include <dune/common/exceptions.hh>
26template<
class Gr
idView>
31 using Element =
typename GridView::template Codim<0>::Entity;
32 using Facet =
typename GridView::template Codim<1>::Entity;
35 static constexpr auto dim = GridView::Grid::dimension;
42 static constexpr SmallLocalIndexType
localOppositeIdx(
const SmallLocalIndexType ownLocalFaceIndex)
44 return isOdd_(ownLocalFaceIndex) ? (ownLocalFaceIndex - 1) : (ownLocalFaceIndex + 1);
50 if constexpr(GridView::Grid::dimension == 1)
52 assert(
false &&
"There are no lateral scvfs in 1D");
56 if constexpr (GridView::Grid::dimension == 2)
58 switch (ownLocalScvfIndex)
70 assert(
false &&
"No lateral orthogonal scvf found");
77 switch (ownLocalScvfIndex)
106 assert(
false &&
"No lateral orthogonal scvf found");
116 constexpr auto table = []
118 using Table = std::array<std::array<SmallLocalIndexType, numLateralFacesPerScv>,
numElementFaces>;
119 if constexpr (
dim == 1)
121 else if constexpr (
dim == 2)
122 return Table {{ {2,3}, {2,3}, {0,1}, {0,1} }};
124 return Table {{ {2,3,4,5}, {2,3,4,5}, {0,1,4,5}, {0,1,4,5}, {0,1,2,3}, {0,1,2,3} }};
127 return table[ownLocalFaceIndex];
131 static Facet
facet(
const SmallLocalIndexType localFacetIdx,
const Element& element)
133 return element.template subEntity <1> (localFacetIdx);
137 auto intersection(
const SmallLocalIndexType localFacetIdx,
const Element& element)
const
144 DUNE_THROW(Dune::InvalidStateException,
"localFacetIdx " << localFacetIdx <<
" out of range");
148 template<
class FVElementGeometry,
class SubControlVolumeFace>
151 assert (scvf.isLateral());
153 using Grid =
typename FVElementGeometry::GridGeometry::Grid;
161 const auto& insideScv = fvGeometry.scv(scvf.insideScvIdx());
162 const auto& outsideScv = fvGeometry.scv(scvf.outsideScvIdx());
163 int onBoundaryCounter = 0;
164 onBoundaryCounter +=
static_cast<int>(insideScv.boundary());
165 onBoundaryCounter +=
static_cast<int>(outsideScv.boundary());
166 return onBoundaryCounter == 1;
170 template<
class SubControlVolumeFace,
class SubControlVolume>
172 const SubControlVolume& scv)
176 const SmallLocalIndexType offset = (
dim == 2) ? 3 : 5;
179 return isOdd_(scv.indexInElement()) ? scvf.localIndex() - offset : scvf.localIndex() + offset;
183 {
return gridView_; }
187 static constexpr bool isOdd_(
int number)
188 {
return number % 2; }
Face centered staggered geometry helper.
Definition: discretization/facecentered/staggered/geometryhelper.hh:28
static Facet facet(const SmallLocalIndexType localFacetIdx, const Element &element)
Returns an element's facet based on the local facet index.
Definition: discretization/facecentered/staggered/geometryhelper.hh:131
static constexpr auto localLaterFaceIndices(const SmallLocalIndexType ownLocalFaceIndex)
Returns the local indices of the faces lateral to the own one.
Definition: discretization/facecentered/staggered/geometryhelper.hh:114
static bool scvfIntegrationPointInConcaveCorner(const FVElementGeometry &fvGeometry, const SubControlVolumeFace &scvf)
Returns true if the IP of an scvf lies on a concave corner.
Definition: discretization/facecentered/staggered/geometryhelper.hh:149
static constexpr SmallLocalIndexType localOppositeIdx(const SmallLocalIndexType ownLocalFaceIndex)
Returns the local index of the opposing face.
Definition: discretization/facecentered/staggered/geometryhelper.hh:42
static constexpr auto numLateralFacesPerScv
Definition: discretization/facecentered/staggered/geometryhelper.hh:37
const GridView & gridView() const
Definition: discretization/facecentered/staggered/geometryhelper.hh:182
static constexpr auto numElementFaces
Definition: discretization/facecentered/staggered/geometryhelper.hh:36
FaceCenteredStaggeredGeometryHelper(const GridView &gridView)
Definition: discretization/facecentered/staggered/geometryhelper.hh:39
static SmallLocalIndexType localIndexOutsideScvfWithSameIntegrationPoint(const SubControlVolumeFace &scvf, const SubControlVolume &scv)
Definition: discretization/facecentered/staggered/geometryhelper.hh:171
auto intersection(const SmallLocalIndexType localFacetIdx, const Element &element) const
Returns an element's intersection based on the local facet index.
Definition: discretization/facecentered/staggered/geometryhelper.hh:137
static constexpr int lateralOrthogonalScvfLocalIndex(const SmallLocalIndexType ownLocalScvfIndex)
Return the local index of a lateral orthogonal scvf.
Definition: discretization/facecentered/staggered/geometryhelper.hh:48
static constexpr auto dim
Definition: discretization/facecentered/staggered/geometryhelper.hh:35
Type trait to determine if a grid supports concave corners (e.g. by cutting out a hole from the domai...
Definition: gridsupportsconcavecorners.hh:31
typename GridView::IndexSet::IndexType GridIndex
Definition: indextraits.hh:27
std::uint_least8_t SmallLocalIndex
Definition: indextraits.hh:29