26#ifndef DUMUX_DISCRETIZATION_CCMPFA_FV_ELEMENT_GEOMETRY_HH
27#define DUMUX_DISCRETIZATION_CCMPFA_FV_ELEMENT_GEOMETRY_HH
30#include <dune/common/exceptions.hh>
31#include <dune/common/iteratorrange.hh>
48template<
class GG,
bool enableGr
idGeometryCache>
61 using GridView =
typename GG::GridView;
64 static constexpr int dim = GridView::dimension;
68 using Element =
typename GridView::template Codim<0>::Entity;
76 static constexpr std::size_t maxNumElementScvs = 1;
78 static constexpr std::size_t maxNumElementScvfs = dim == 3 ? 24 : 8;
82 : gridGeometryPtr_(&gridGeometry) {}
87 return gridGeometry().scv(scvIdx);
93 return gridGeometry().scvf(scvfIdx);
100 return gridGeometry().flipScvf(scvfIdx, outsideScvIdx);
108 friend inline Dune::IteratorRange< ScvIterator<SubControlVolume, std::array<GridIndexType, 1>, ThisType> >
112 return Dune::IteratorRange<ScvIterator>(
ScvIterator(fvGeometry.scvIndices_.begin(), fvGeometry),
113 ScvIterator(fvGeometry.scvIndices_.end(), fvGeometry));
121 friend inline Dune::IteratorRange< ScvfIterator<SubControlVolumeFace, std::vector<GridIndexType>, ThisType> >
124 const auto& g = fvGeometry.gridGeometry();
125 const auto scvIdx = fvGeometry.scvIndices_[0];
127 return Dune::IteratorRange<ScvfIterator>(
ScvfIterator(g.scvfIndicesOfScv(scvIdx).begin(), fvGeometry),
128 ScvfIterator(g.scvfIndicesOfScv(scvIdx).end(), fvGeometry));
134 return scvIndices_.size();
140 return gridGeometry().scvfIndicesOfScv(scvIndices_[0]).size();
146 this->bindElement(element);
153 scvIndices_[0] = gridGeometry().elementMapper().index(element);
158 {
return static_cast<bool>(element_); }
162 {
return *element_; }
166 {
return *gridGeometryPtr_; }
170 {
return gridGeometry().hasBoundaryScvf(scvIndices_[0]); }
174 std::optional<Element> element_;
175 std::array<GridIndexType, 1> scvIndices_;
176 const GridGeometry* gridGeometryPtr_;
188 using GridView =
typename GG::GridView;
190 using MpfaHelper =
typename GG::MpfaHelper;
192 static const int dim = GridView::dimension;
193 static const int dimWorld = GridView::dimensionworld;
194 using CoordScalar =
typename GridView::ctype;
198 using Element =
typename GridView::template Codim<0>::Entity;
206 static constexpr std::size_t maxNumElementScvs = 1;
208 static constexpr std::size_t maxNumElementScvfs = dim == 3 ? 24 : 8;
212 : gridGeometryPtr_(&gridGeometry) {}
218 if (scvIdx == scvIndices_[0])
221 return neighborScvs_[findLocalIndex(scvIdx, neighborScvIndices_)];
228 auto it = std::find(scvfIndices_.begin(), scvfIndices_.end(), scvfIdx);
229 if (it != scvfIndices_.end())
232 return neighborScvfs_[findLocalIndex(scvfIdx, neighborScvfIndices_)];
239 return scvf( gridGeometry().flipScvfIdx(scvfIdx, outsideScvIdx) );
247 friend inline Dune::IteratorRange<typename std::array<SubControlVolume, 1>::const_iterator>
250 using IteratorType =
typename std::array<SubControlVolume, 1>::const_iterator;
251 return Dune::IteratorRange<IteratorType>(g.scvs_.begin(), g.scvs_.end());
259 friend inline Dune::IteratorRange<typename std::vector<SubControlVolumeFace>::const_iterator>
262 using IteratorType =
typename std::vector<SubControlVolumeFace>::const_iterator;
263 return Dune::IteratorRange<IteratorType>(g.scvfs_.begin(), g.scvfs_.end());
268 {
return scvs_.size(); }
272 {
return scvfs_.size(); }
279 bindElement(element);
282 const auto globalI = gridGeometry().elementMapper().index(element);
283 const auto& assemblyMapI = gridGeometry().connectivityMap()[globalI];
286 const auto numNeighbors = assemblyMapI.size();
287 const auto numNeighborScvfs = numNeighborScvfs_(assemblyMapI);
288 neighborScvs_.reserve(numNeighbors);
289 neighborScvIndices_.reserve(numNeighbors);
290 neighborScvfIndices_.reserve(numNeighborScvfs);
291 neighborScvfs_.reserve(numNeighborScvfs);
295 for (
const auto& dataJ : assemblyMapI)
296 makeNeighborGeometries(gridGeometry().element(dataJ.globalJ),
299 dataJ.additionalScvfs);
322 makeElementGeometries(element);
327 {
return static_cast<bool>(element_); }
331 {
return *element_; }
335 {
return *gridGeometryPtr_; }
339 {
return hasBoundaryScvf_; }
344 template<
class DataJContainer>
345 std::size_t numNeighborScvfs_(
const DataJContainer& dataJContainer)
347 std::size_t numNeighborScvfs = 0;
348 for (
const auto& dataJ : dataJContainer)
349 numNeighborScvfs += dataJ.scvfsJ.size() + dataJ.additionalScvfs.size();
350 return numNeighborScvfs;
354 void makeElementGeometries(
const Element& element)
357 const auto eIdx = gridGeometry().elementMapper().index(element);
358 scvs_[0] = SubControlVolume(element.geometry(), eIdx);
359 scvIndices_[0] = eIdx;
362 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
363 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
366 static const auto q = getParam<CoordScalar>(
"MPFA.Q");
369 const auto numLocalScvf = scvFaceIndices.size();
370 scvfIndices_.reserve(numLocalScvf);
371 scvfs_.reserve(numLocalScvf);
375 std::vector<bool> finishedFacets;
377 finishedFacets.resize(element.subEntities(1),
false);
380 for (
const auto& is : intersections(gridGeometry().gridView(), element))
386 const auto indexInInside = is.indexInInside();
387 if (finishedFacets[indexInInside])
390 finishedFacets[indexInInside] =
true;
394 bool useNeighbor = is.neighbor() && is.outside().level() >
element.level();
395 const auto& e = useNeighbor ? is.outside() :
element;
396 const auto indexInElement = useNeighbor ? is.indexInOutside() : is.indexInInside();
397 const auto eg = e.geometry();
398 const auto refElement = referenceElement(eg);
401 const auto numCorners = is.geometry().corners();
402 const auto isPositions = MpfaHelper::computeScvfCornersOnIntersection(eg,
411 auto vIdxLocal = refElement.subEntity(indexInElement, 1, c, dim);
412 auto vIdxGlobal = gridGeometry().vertexMapper().subIndex(e, vIdxLocal, dim);
415 if (gridGeometry().isGhostVertex(vIdxGlobal))
418 hasBoundaryScvf_ = (hasBoundaryScvf_ || is.boundary());
420 scvfs_.emplace_back(MpfaHelper(),
421 MpfaHelper::getScvfCorners(isPositions,
numCorners, c),
425 scvFaceIndices[scvfCounter],
427 neighborVolVarIndices[scvfCounter],
431 scvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
438 template<
typename IndexVector>
439 void makeNeighborGeometries(
const Element& element,
440 GridIndexType eIdxGlobal,
441 const IndexVector& scvfIndices,
442 const IndexVector& additionalScvfs)
445 neighborScvs_.emplace_back(
element.geometry(), eIdxGlobal);
446 neighborScvIndices_.push_back(eIdxGlobal);
449 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdxGlobal);
450 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdxGlobal);
453 static const auto q = getParam<CoordScalar>(
"MPFA.Q");
457 std::vector<bool> finishedFacets;
459 finishedFacets.resize(
element.subEntities(1),
false);
462 for (
const auto& is : intersections(gridGeometry().gridView(), element))
468 auto indexInInside = is.indexInInside();
469 if(finishedFacets[indexInInside])
472 finishedFacets[indexInInside] =
true;
476 bool useNeighbor = is.neighbor() && is.outside().level() >
element.level();
477 const auto& e = useNeighbor ? is.outside() :
element;
478 const auto indexInElement = useNeighbor ? is.indexInOutside() : is.indexInInside();
479 const auto eg = e.geometry();
480 const auto refElement = referenceElement(eg);
483 const auto numCorners = is.geometry().corners();
484 const auto isPositions = MpfaHelper::computeScvfCornersOnIntersection(eg,
493 auto vIdxLocal = refElement.subEntity(indexInElement, 1, c, dim);
494 auto vIdxGlobal = gridGeometry().vertexMapper().subIndex(e, vIdxLocal, dim);
497 if (gridGeometry().isGhostVertex(vIdxGlobal))
501 if (!MpfaHelper::vectorContainsValue(scvfIndices, scvFaceIndices[scvfCounter])
502 && !MpfaHelper::vectorContainsValue(additionalScvfs, scvFaceIndices[scvfCounter]))
510 neighborScvfs_.emplace_back(MpfaHelper(),
511 MpfaHelper::getScvfCorners(isPositions,
numCorners, c),
515 scvFaceIndices[scvfCounter],
517 neighborVolVarIndices[scvfCounter],
521 neighborScvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
530 unsigned int findLocalIndex(
const GridIndexType idx,
531 const std::vector<GridIndexType>& indices)
const
533 auto it = std::find(indices.begin(), indices.end(), idx);
534 assert(it != indices.end() &&
"Could not find the scv/scvf! Make sure to properly bind this class!");
541 scvfIndices_.clear();
544 neighborScvIndices_.clear();
545 neighborScvfIndices_.clear();
546 neighborScvs_.clear();
547 neighborScvfs_.clear();
549 hasBoundaryScvf_ =
false;
552 const GridGeometry* gridGeometryPtr_;
553 std::optional<Element> element_;
556 std::array<GridIndexType, 1> scvIndices_;
557 std::vector<GridIndexType> scvfIndices_;
558 std::array<SubControlVolume, 1> scvs_;
559 std::vector<SubControlVolumeFace> scvfs_;
561 std::vector<GridIndexType> neighborScvIndices_;
562 std::vector<GridIndexType> neighborScvfIndices_;
563 std::vector<SubControlVolume> neighborScvs_;
564 std::vector<SubControlVolumeFace> neighborScvfs_;
566 bool hasBoundaryScvf_ =
false;
Defines the index types used for grid and local indices.
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Class providing iterators over sub control volumes and sub control volume faces of an element.
ctype distance(const Dune::FieldVector< ctype, dimWorld > &a, const Dune::FieldVector< ctype, dimWorld > &b)
Compute the shortest distance between two points.
Definition: distance.hh:138
std::size_t numCorners(Shape shape)
Returns the number of corners of a given geometry.
Definition: throatproperties.hh:215
typename GridView::IndexSet::IndexType GridIndex
Definition: indextraits.hh:39
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered mpfa models This builds up th...
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:49
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered mpfa models Specialization fo...
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:59
const SubControlVolume & scv(GridIndexType scvIdx) const
Get an element sub control volume with a global scv index.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:85
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:138
CCMpfaFVElementGeometry(const GridGeometry &gridGeometry)
Constructor.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:81
const Element & element() const
The bound element.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:161
void bindElement(const Element &element)
Bind only element-local.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:150
bool isBound() const
Returns true if bind/bindElement has already been called.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:157
friend Dune::IteratorRange< ScvIterator< SubControlVolume, std::array< GridIndexType, 1 >, ThisType > > scvs(const CCMpfaFVElementGeometry &fvGeometry)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:109
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:169
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:72
friend Dune::IteratorRange< ScvfIterator< SubControlVolumeFace, std::vector< GridIndexType >, ThisType > > scvfs(const CCMpfaFVElementGeometry &fvGeometry)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:122
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Get an element sub control volume face with a global scvf index.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:91
const SubControlVolumeFace & flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx=0) const
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:98
const GridGeometry & gridGeometry() const
The global finite volume geometry we are a restriction of.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:165
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:132
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:68
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:70
void bind(const Element &element)
Binding of an element, called by the local assembler to prepare element assembly.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:144
GG GridGeometry
export type of finite volume grid geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:74
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models Specialization fo...
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:186
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:226
const SubControlVolume & scv(GridIndexType scvIdx) const
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:216
friend Dune::IteratorRange< typename std::vector< SubControlVolumeFace >::const_iterator > scvfs(const ThisType &g)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:260
void bind(const Element &element)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:276
const GridGeometry & gridGeometry() const
The global finite volume geometry we are a restriction of.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:334
void bindElement(const Element &element)
Binding of an element preparing the geometries only inside the element.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:318
const SubControlVolumeFace & flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx=0) const
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:237
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:200
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:198
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:267
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:202
const Element & element() const
The bound element.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:330
bool isBound() const
Returns true if bind/bindElement has already been called.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:326
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:271
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:338
GG GridGeometry
export type of finite volume grid geometrys
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:204
friend Dune::IteratorRange< typename std::array< SubControlVolume, 1 >::const_iterator > scvs(const ThisType &g)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:248
CCMpfaFVElementGeometry(const GridGeometry &gridGeometry)
Constructor.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:211
Iterators over sub control volumes.
Definition: scvandscvfiterators.hh:42
Iterators over sub control volume faces of an fv geometry.
Definition: scvandscvfiterators.hh:82