12#ifndef DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_FV_GRID_GEOMETRY
13#define DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_FV_GRID_GEOMETRY
17#include <dune/common/rangeutilities.hh>
18#include <dune/grid/common/scsgmapper.hh>
46template<
class Gr
idView>
55 template<
class Gr
idGeometry>
58 template<
class Gr
idGeometry,
bool enableCache>
63 static constexpr auto dim = GridView::Grid::dimension;
81template<
class GridView,
82 bool cachingEnabled =
false,
92template<
class GV,
class Traits>
101 using Element =
typename GV::template Codim<0>::Entity;
103 using IntersectionMapper =
typename Traits::IntersectionMapper;
104 using ConnectivityMap =
typename Traits::template ConnectivityMap<ThisType>;
106 using Scalar =
typename GV::ctype;
108 static constexpr auto dim = Traits::StaticInfo::dim;
109 static constexpr auto numScvsPerElement = Traits::StaticInfo::numScvsPerElement;
110 static constexpr auto numLateralScvfsPerScv = Traits::StaticInfo::numLateralScvfsPerScv;
111 static constexpr auto numLateralScvfsPerElement = Traits::StaticInfo::numLateralScvfsPerElement;
112 static constexpr auto minNumScvfsPerElement = Traits::StaticInfo::minNumScvfsPerElement;
113 static constexpr auto maxNumScvfsPerElement = Traits::StaticInfo::maxNumScvfsPerElement;
115 using ScvfCornerStorage =
typename Traits::SubControlVolumeFace::Traits::CornerStorage;
116 using ScvCornerStorage =
typename Traits::SubControlVolume::Traits::CornerStorage;
123 static constexpr bool cachingEnabled =
true;
147 , intersectionMapper_(this->gridView())
151 DUNE_THROW(Dune::InvalidStateException,
"The staggered discretization method needs at least an overlap of 1 for parallel computations. "
152 <<
" Set the parameter \"Grid.Overlap\" in the input file.");
164 {
return scvs_.size(); }
168 {
return scvfs_.size(); }
172 {
return numBoundaryScv_; }
176 {
return numBoundaryScvf_; }
180 {
return intersectionMapper_.numIntersections(); }
184 {
return this->gridView().size(1); }
189 ParentType::update(gridView);
196 ParentType::update(std::move(gridView));
202 {
return scvs_[scvIdx]; }
208 auto begin = scvs_.cbegin() + numScvsPerElement*fvGeometry.elementIndex();
209 const auto end = begin + numScvsPerElement;
210 return Dune::IteratorRange<std::decay_t<
decltype(begin)>>(begin, end);
215 {
return scvfs_[scvfIdx]; }
219 {
return scvfIndicesOfElement_[eIdx]; }
226 {
return connectivityMap_; }
230 {
return hasBoundaryScvf_[eIdx]; }
234 {
return intersectionMapper_; }
238 {
return periodicFaceMap_.count(dofIdx); }
242 {
return periodicFaceMap_.at(dofIdx); }
246 {
return periodicFaceMap_; }
255 scvfIndicesOfElement_.clear();
256 intersectionMapper_.update(this->gridView());
259 const auto numElements = this->gridView().size(0);
260 scvfIndicesOfElement_.resize(numElements);
261 hasBoundaryScvf_.resize(numElements,
false);
263 outSideBoundaryVolVarIdx_ = 0;
265 numBoundaryScvf_ = 0;
267 GeometryHelper geometryHelper(this->gridView());
270 GridIndexType numScvfs = 0;
271 for (
const auto& element : elements(this->gridView()))
273 assert(numScvsPerElement == element.subEntities(1));
275 for (
const auto& intersection : intersections(this->gridView(), element))
281 numScvfs += numLateralScvfsPerScv;
284 if (onDomainBoundary_(intersection))
287 numBoundaryScv_ += numLateralScvfsPerScv;
296 const auto numScvs = numElements*numScvsPerElement;
297 scvs_.resize(numScvs);
298 scvfs_.reserve(numScvfs);
301 std::size_t globalScvfIdx = 0;
302 for (
const auto& element : elements(this->gridView()))
304 const auto eIdx = this->elementMapper().index(element);
305 auto& globalScvfIndices = scvfIndicesOfElement_[eIdx];
306 globalScvfIndices.resize(minNumScvfsPerElement);
307 globalScvfIndices.reserve(maxNumScvfsPerElement);
309 auto getGlobalScvIdx = [&](
const auto elementIdx,
const auto localScvIdx)
310 {
return numScvsPerElement*elementIdx + localScvIdx; };
312 LocalIntersectionMapper localIsMapper;
313 localIsMapper.update(this->gridView(),
element);
315 for (
const auto& intersection : intersections(this->gridView(), element))
317 const auto& intersectionUnitOuterNormal = intersection.centerUnitOuterNormal();
318 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
319 auto localScvfIdx = localScvIdx*(1 + numLateralScvfsPerScv);
321 const auto globalScvIdx = getGlobalScvIdx(eIdx, localScvIdx);
322 const auto dofIndex = intersectionMapper().globalIntersectionIndex(element, intersection.indexInInside());
323 const auto localOppositeScvIdx = geometryHelper.localOppositeIdx(localScvIdx);
324 const auto& intersectionGeometry = intersection.geometry();
325 const auto& elementGeometry =
element.geometry();
327 assert(localIsMapper.refToRealIdx(localScvIdx) == intersection.indexInInside());
330 if (onPeriodicBoundary_(intersection))
334 const auto& otherElement = intersection.outside();
336 SmallLocalIndexType otherIntersectionLocalIdx = 0;
337 bool periodicFaceFound =
false;
339 for (
const auto& otherIntersection : intersections(this->gridView(), otherElement))
341 if (periodicFaceFound)
344 if (Dune::FloatCmp::eq(intersectionUnitOuterNormal*otherIntersection.centerUnitOuterNormal(), -1.0, 1e-7))
346 const auto periodicDofIdx = intersectionMapper().globalIntersectionIndex(otherElement, otherIntersectionLocalIdx);
347 periodicFaceMap_[dofIndex] = periodicDofIdx;
348 periodicFaceFound =
true;
351 ++otherIntersectionLocalIdx;
356 scvs_[globalScvIdx] = SubControlVolume(
358 intersectionGeometry,
363 this->elementMapper().index(element),
364 onDomainBoundary_(intersection)
368 scvfs_.emplace_back(elementGeometry,
369 intersectionGeometry,
370 std::array{globalScvIdx, getGlobalScvIdx(eIdx, localOppositeScvIdx)},
373 intersectionUnitOuterNormal,
374 SubControlVolumeFace::FaceType::frontal,
375 SubControlVolumeFace::BoundaryType::interior
378 globalScvfIndices[localScvfIdx] = globalScvfIdx++;
382 for (
const auto lateralFacetIndex : Dune::transformedRangeView(geometryHelper.localLaterFaceIndices(localScvIdx),
383 [&](
auto&& idx) { return localIsMapper.refToRealIdx(idx) ;})
386 const auto& lateralIntersection = geometryHelper.intersection(lateralFacetIndex, element);
389 const auto globalScvIndicesForLateralFace = [&]
391 const auto globalOutsideScvIdx = [&]
393 if (lateralIntersection.neighbor())
395 const auto parallelElemIdx = this->elementMapper().index(lateralIntersection.outside());
396 return getGlobalScvIdx(parallelElemIdx, localScvIdx);
398 else if (onDomainBoundary_(lateralIntersection))
399 return numScvs + outSideBoundaryVolVarIdx_++;
404 return std::array{globalScvIdx, globalOutsideScvIdx};
407 const auto boundaryType = [&]
409 if (onProcessorBoundary_(lateralIntersection))
410 return SubControlVolumeFace::BoundaryType::processorBoundary;
411 else if (onDomainBoundary_(lateralIntersection))
412 return SubControlVolumeFace::BoundaryType::physicalBoundary;
414 return SubControlVolumeFace::BoundaryType::interior;
419 intersectionGeometry,
420 geometryHelper.facet(lateralFacetIndex, element).geometry(),
421 globalScvIndicesForLateralFace,
424 lateralIntersection.centerUnitOuterNormal(),
425 SubControlVolumeFace::FaceType::lateral,
429 globalScvfIndices[localScvfIdx] = globalScvfIdx++;
432 if (onDomainBoundary_(lateralIntersection))
435 hasBoundaryScvf_[eIdx] =
true;
442 int localScvfIdx = minNumScvfsPerElement;
443 for (
const auto& intersection : intersections(this->gridView(),
element))
446 if (onDomainBoundary_(intersection))
448 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
449 const auto globalScvIdx = getGlobalScvIdx(eIdx, localScvIdx);
455 intersection.geometry(),
456 std::array{globalScvIdx, globalScvIdx},
459 intersection.centerUnitOuterNormal(),
460 SubControlVolumeFace::FaceType::frontal,
461 SubControlVolumeFace::BoundaryType::physicalBoundary
464 globalScvfIndices.push_back(globalScvfIdx);
467 hasBoundaryScvf_[eIdx] =
true;
472 connectivityMap_.update(*
this);
475 bool onDomainBoundary_(
const typename GridView::Intersection& intersection)
const
477 return !intersection.neighbor() && intersection.boundary();
480 bool onProcessorBoundary_(
const typename GridView::Intersection& intersection)
const
482 return !intersection.neighbor() && !intersection.boundary();
485 bool onPeriodicBoundary_(
const typename GridView::Intersection& intersection)
const
487 return intersection.boundary() && intersection.neighbor();
491 ConnectivityMap connectivityMap_;
492 IntersectionMapper intersectionMapper_;
494 std::vector<SubControlVolume> scvs_;
495 std::vector<SubControlVolumeFace> scvfs_;
496 GridIndexType numBoundaryScv_;
497 GridIndexType numBoundaryScvf_;
498 GridIndexType outSideBoundaryVolVarIdx_;
499 std::vector<bool> hasBoundaryScvf_;
501 std::vector<std::vector<GridIndexType>> scvfIndicesOfElement_;
504 std::unordered_map<GridIndexType, GridIndexType> periodicFaceMap_;
513template<
class GV,
class Traits>
522 using Element =
typename GV::template Codim<0>::Entity;
524 using IntersectionMapper =
typename Traits::IntersectionMapper;
525 using ConnectivityMap =
typename Traits::template ConnectivityMap<ThisType>;
527 static constexpr auto dim = Traits::StaticInfo::dim;
528 static constexpr auto numScvsPerElement = Traits::StaticInfo::numScvsPerElement;
529 static constexpr auto numLateralScvfsPerScv = Traits::StaticInfo::numLateralScvfsPerScv;
530 static constexpr auto numLateralScvfsPerElement = Traits::StaticInfo::numLateralScvfsPerElement;
531 static constexpr auto minNumScvfsPerElement = Traits::StaticInfo::minNumScvfsPerElement;
532 static constexpr auto maxNumScvfsPerElement = Traits::StaticInfo::maxNumScvfsPerElement;
539 static constexpr bool cachingEnabled =
false;
563 , intersectionMapper_(this->gridView())
567 DUNE_THROW(Dune::InvalidStateException,
"The staggered discretization method needs at least an overlap of 1 for parallel computations. "
568 <<
" Set the parameter \"Grid.Overlap\" in the input file.");
588 {
return numBoundaryScv_; }
592 {
return numBoundaryScvf_; }
596 {
return intersectionMapper_.numIntersections(); }
600 {
return this->gridView().size(1); }
607 {
return connectivityMap_; }
611 {
return hasBoundaryScvf_[eIdx]; }
615 {
return intersectionMapper_; }
619 {
return scvfIndicesOfElement_[eIdx]; }
623 {
return outsideVolVarIndices_.at(scvfIdx); }
628 ParentType::update(gridView);
635 ParentType::update(std::move(gridView));
641 {
return periodicFaceMap_.count(dofIdx); }
645 {
return periodicFaceMap_.at(dofIdx); }
649 {
return periodicFaceMap_; }
655 intersectionMapper_.update(this->gridView());
660 numBoundaryScvf_ = 0;
661 hasBoundaryScvf_.clear();
662 scvfIndicesOfElement_.clear();
663 outsideVolVarIndices_.clear();
666 const auto numElements = this->gridView().size(0);
667 scvfIndicesOfElement_.resize(numElements);
668 hasBoundaryScvf_.resize(numElements,
false);
669 numScvs_ = numElements*numScvsPerElement;
671 GeometryHelper geometryHelper(this->gridView());
674 GridIndexType scvfIdx = 0;
676 GridIndexType neighborVolVarIdx = numScvs_;
678 for (
const auto& element : elements(this->gridView()))
680 const auto eIdx = this->elementMapper().index(element);
681 assert(numScvsPerElement == element.subEntities(1));
684 auto& globalScvfIndices = scvfIndicesOfElement_[eIdx];
685 globalScvfIndices.reserve(maxNumScvfsPerElement);
686 globalScvfIndices.resize(minNumScvfsPerElement);
689 std::size_t numFrontalBoundaryScvfs = 0;
692 LocalIntersectionIndexMapper localIsMapper;
693 localIsMapper.update(this->gridView(), element);
695 for (
const auto& intersection : intersections(this->gridView(), element))
697 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
698 auto localScvfIdx = localScvIdx*(1 + numLateralScvfsPerScv);
700 assert(localIsMapper.refToRealIdx(localScvIdx) == intersection.indexInInside());
702 globalScvfIndices[localScvfIdx] = scvfIdx++;
705 if constexpr(dim > 1)
708 for (
const auto lateralFacetIndex : Dune::transformedRangeView(geometryHelper.localLaterFaceIndices(localScvIdx),
709 [&](
auto idx) { return localIsMapper.refToRealIdx(idx) ;})
712 if (onDomainBoundary_(geometryHelper.intersection(lateralFacetIndex, element)))
714 outsideVolVarIndices_[scvfIdx] = neighborVolVarIdx++;
716 hasBoundaryScvf_[eIdx] =
true;
719 globalScvfIndices[localScvfIdx] = scvfIdx++;
725 if (onDomainBoundary_(intersection))
728 numBoundaryScv_ += numLateralScvfsPerScv;
729 ++numFrontalBoundaryScvfs;
731 hasBoundaryScvf_[eIdx] =
true;
735 if (onPeriodicBoundary_(intersection))
739 const auto& otherElement = intersection.outside();
741 SmallLocalIndexType otherIntersectionLocalIdx = 0;
742 bool periodicFaceFound =
false;
744 for (
const auto& otherIntersection : intersections(this->gridView(), otherElement))
746 if (periodicFaceFound)
749 if (Dune::FloatCmp::eq(intersection.centerUnitOuterNormal()*otherIntersection.centerUnitOuterNormal(), -1.0, 1e-7))
751 const auto periodicDofIdx = intersectionMapper().globalIntersectionIndex(otherElement, otherIntersectionLocalIdx);
752 const auto dofIndex = intersectionMapper().globalIntersectionIndex(element, localScvIdx);
753 periodicFaceMap_[dofIndex] = periodicDofIdx;
754 periodicFaceFound =
true;
757 ++otherIntersectionLocalIdx;
763 for (std::size_t i = 0; i < numFrontalBoundaryScvfs; ++i)
764 globalScvfIndices.push_back(scvfIdx++);
770 connectivityMap_.update(*
this);
773 bool onDomainBoundary_(
const typename GridView::Intersection& intersection)
const
775 return !intersection.neighbor() && intersection.boundary();
778 bool onProcessorBoundary_(
const typename GridView::Intersection& intersection)
const
780 return !intersection.neighbor() && !intersection.boundary();
783 bool onPeriodicBoundary_(
const typename GridView::Intersection& intersection)
const
785 return intersection.boundary() && intersection.neighbor();
789 ConnectivityMap connectivityMap_;
790 IntersectionMapper intersectionMapper_;
793 std::size_t numScvs_;
794 std::size_t numScvf_;
795 std::size_t numBoundaryScv_;
796 std::size_t numBoundaryScvf_;
797 std::vector<bool> hasBoundaryScvf_;
799 std::vector<std::vector<GridIndexType>> scvfIndicesOfElement_;
802 std::unordered_map<GridIndexType, GridIndexType> periodicFaceMap_;
803 std::unordered_map<GridIndexType, GridIndexType> outsideVolVarIndices_;
Base class for grid geometries.
Check the overlap size for different discretization methods.
Base class for all grid geometries.
Definition: basegridgeometry.hh:52
typename BaseImplementation::GridView GridView
export the grid view type
Definition: basegridgeometry.hh:60
Definition: localintersectionindexmapper.hh:29
Stores the dof indices corresponding to the neighboring scvs that contribute to the derivative calcul...
Definition: facecentered/staggered/connectivitymap.hh:30
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:98
Base class for the finite volume geometry vector for face-centered staggered models This builds up th...
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:516
typename Traits::SubControlVolumeFace SubControlVolumeFace
export the type of sub control volume
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:548
typename Traits::GeometryHelper GeometryHelper
export the geometry helper type
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:552
FaceCenteredStaggeredFVGridGeometry(std::shared_ptr< BasicGridGeometry > gg, const std::string ¶mGroup="")
Constructor with basic grid geometry used to share state with another grid geometry on the same grid ...
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:561
typename Traits::SubControlVolume SubControlVolume
export the type of sub control volume
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:546
std::size_t numDofs() const
the total number of dofs
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:599
GridIndexType periodicallyMappedDof(GridIndexType dofIdx) const
The index of the d.o.f. on the other side of the periodic boundary.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:644
const ConnectivityMap & connectivityMap() const
Returns the connectivity map of which dofs have derivatives with respect to a given dof.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:606
std::size_t numBoundaryScv() const
The total number of boundary sub control volumes.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:587
const IntersectionMapper & intersectionMapper() const
Return a reference to the intersection mapper.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:614
void update(GridView &&gridView)
update all fvElementGeometries (call this after grid adaption)
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:633
FaceCenteredStaggeredFVGridGeometry(const GridView &gridView, const std::string ¶mGroup="")
Constructor from gridView.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:574
bool dofOnPeriodicBoundary(GridIndexType dofIdx) const
If a d.o.f. is on a periodic boundary.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:640
typename Traits::StaticInfo StaticInformation
export static information
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:556
const std::vector< GridIndexType > & scvfIndicesOfElement(GridIndexType eIdx) const
Get the global sub control volume face indices of an element.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:618
bool hasBoundaryScvf(GridIndexType eIdx) const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:610
const std::unordered_map< GridIndexType, GridIndexType > & periodicVertexMap() const
Returns the map between dofs across periodic boundaries // TODO rename to periodic dof map in fvassem...
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:648
std::size_t numScvf() const
The total number of sub control volume faces.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:583
BasicGridGeometry_t< GV, Traits > BasicGridGeometry
export basic grid geometry type for the alternative constructor
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:542
std::size_t numBoundaryScvf() const
The total number of boundary sub control volume faces.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:591
std::size_t numScv() const
The total number of sub control volumes.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:579
typename Traits::LocalIntersectionMapper LocalIntersectionMapper
export the local intersection mapper
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:554
GridIndexType outsideVolVarIndex(GridIndexType scvfIdx) const
Get the global sub control volume face indices of an element.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:622
void update(const GridView &gridView)
update all fvElementGeometries (call this after grid adaption)
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:626
std::size_t numIntersections() const
The total number of intersections.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:595
Extrusion_t< Traits > Extrusion
export the type of extrusion
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:558
typename Traits::template LocalView< ThisType, false > LocalView
export the type of the fv element geometry (the local view type)
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:544
Base class for the finite volume geometry vector for staggered models This builds up the sub control ...
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:95
const ConnectivityMap & connectivityMap() const
Returns the connectivity map of which dofs have derivatives with respect to a given dof.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:225
std::size_t numScvf() const
The total number of sub control volume faces.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:167
const IntersectionMapper & intersectionMapper() const
Return a reference to the intersection mapper.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:233
bool hasBoundaryScvf(GridIndexType eIdx) const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:229
void update(const GridView &gridView)
update all fvElementGeometries (call this after grid adaption)
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:187
typename Traits::SubControlVolumeFace SubControlVolumeFace
export the type of sub control volume
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:132
typename Traits::LocalIntersectionMapper LocalIntersectionMapper
export the local intersection mapper
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:138
typename Traits::GeometryHelper GeometryHelper
export the geometry helper type
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:136
const std::vector< GridIndexType > & scvfIndicesOfElement(GridIndexType eIdx) const
Get the global sub control volume face indices of an element.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:218
FaceCenteredStaggeredFVGridGeometry(std::shared_ptr< BasicGridGeometry > gg, const std::string ¶mGroup="")
Constructor with basic grid geometry used to share state with another grid geometry on the same grid ...
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:145
typename Traits::StaticInfo StaticInformation
export static information
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:140
std::size_t numDofs() const
the total number of dofs
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:183
const std::unordered_map< GridIndexType, GridIndexType > & periodicVertexMap() const
Returns the map between dofs across periodic boundaries // TODO rename to periodic dof map in fvassem...
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:245
std::size_t numBoundaryScv() const
The total number of boundary sub control volumes.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:171
FaceCenteredStaggeredFVGridGeometry(const GridView &gridView, const std::string ¶mGroup="")
Constructor from gridView.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:158
GridIndexType periodicallyMappedDof(GridIndexType dofIdx) const
The index of the d.o.f. on the other side of the periodic boundary.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:241
const SubControlVolume & scv(GridIndexType scvIdx) const
Get a sub control volume with a global scv index.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:201
typename Traits::SubControlVolume SubControlVolume
export the type of sub control volume
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:130
typename Traits::template LocalView< ThisType, true > LocalView
export the type of the fv element geometry (the local view type)
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:128
void update(GridView &&gridView)
update all fvElementGeometries (call this after grid adaption)
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:194
BasicGridGeometry_t< GV, Traits > BasicGridGeometry
export basic grid geometry type for the alternative constructor
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:126
Extrusion_t< Traits > Extrusion
export the type of extrusion
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:142
std::size_t numBoundaryScvf() const
The total number of boundary sub control volume faces.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:175
std::size_t numIntersections() const
The total number of intersections.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:179
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Get a sub control volume face with a global scvf index.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:214
bool dofOnPeriodicBoundary(GridIndexType dofIdx) const
If a d.o.f. is on a periodic boundary.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:237
std::size_t numScv() const
The total number of sub control volumes.
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:163
auto scvs(const LocalView &fvGeometry) const
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:206
Base class for the finite volume geometry vector for face-centered staggered models This builds up th...
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:84
Face centered staggered geometry helper.
Definition: discretization/facecentered/staggered/geometryhelper.hh:28
Face centered staggered sub control volume face.
Definition: discretization/facecentered/staggered/subcontrolvolumeface.hh:54
Face centered staggered sub control volume.
Definition: discretization/facecentered/staggered/subcontrolvolume.hh:52
Defines the default element and vertex mapper types.
Geometry helper for face-centered staggered scheme.
Helper classes to compute the integration elements.
Dune::Std::detected_or_t< Dumux::BasicGridGeometry< GV, typename T::ElementMapper, typename T::VertexMapper >, Detail::SpecifiesBaseGridGeometry, T > BasicGridGeometry_t
Type of the basic grid geometry implementation used as backend.
Definition: basegridgeometry.hh:42
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
defines intersection mappers.
Define some often used mathematical functions.
The available discretization methods in Dumux.
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition: extrusion.hh:166
Check if the overlap size is valid for a given discretization method.
Definition: checkoverlapsize.hh:28
Definition: defaultmappertraits.hh:23
Definition: method.hh:130
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:62
static constexpr auto numLateralScvfsPerElement
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:67
static constexpr auto maxNumScvfsPerElement
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:70
static constexpr auto numScvsPerElement
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:65
static constexpr auto numFacesPerElement
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:64
static constexpr auto dim
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:63
static constexpr auto numLateralScvfsPerScv
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:66
static constexpr auto minNumScvfsPerElement
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:68
The default traits for the face-center staggered finite volume grid geometry Defines the scv and scvf...
Definition: discretization/facecentered/staggered/fvgridgeometry.hh:48
Structure to define the index types used for grid and local indices.
Definition: indextraits.hh:26