24#ifndef DUMUX_DISCRETIZATION_STAGGERED_FV_ELEMENT_GEOMETRY_HH
25#define DUMUX_DISCRETIZATION_STAGGERED_FV_ELEMENT_GEOMETRY_HH
41template<
class GG,
bool enableGr
idGeometryCache>
56 using GridView =
typename GG::GridView;
61 using Element =
typename GridView::template Codim<0>::Entity;
65 using ParentType::ParentType;
69 template<
class CellCenterOrFaceFVGr
idGeometry>
71 :
ParentType(gridGeometry.actualGridGeometry()) {}
74 using ParentType::scvf;
77 return this->gridGeometry().scvf(eIdx, localScvfIdx);
93 using GridView =
typename GG::GridView;
98 using Element =
typename GridView::template Codim<0>::Entity;
108 template<
class CellCenterOrFaceFVGr
idGeometry>
110 : gridGeometryPtr_(&gridGeometry.actualGridGeometry()) {}
114 : gridGeometryPtr_(&gridGeometry) {}
119 return scvf(this->gridGeometry().localToGlobalScvfIndex(eIdx, localScvfIdx));
126 if (scvIdx == scvIndices_[0])
129 return neighborScvs_[findLocalIndex_(scvIdx, neighborScvIndices_)];
136 auto it = std::find(scvfIndices_.begin(), scvfIndices_.end(), scvfIdx);
137 if (it != scvfIndices_.end())
140 return neighborScvfs_[findLocalIndex_(scvfIdx, neighborScvfIndices_)];
148 friend inline Dune::IteratorRange<typename std::array<SubControlVolume, 1>::const_iterator>
151 using IteratorType =
typename std::array<SubControlVolume, 1>::const_iterator;
152 return Dune::IteratorRange<IteratorType>(g.scvs_.begin(), g.scvs_.end());
160 friend inline Dune::IteratorRange<typename std::vector<SubControlVolumeFace>::const_iterator>
163 using IteratorType =
typename std::vector<SubControlVolumeFace>::const_iterator;
164 return Dune::IteratorRange<IteratorType>(g.scvfs_.begin(), g.scvfs_.end());
169 {
return scvs_.size(); }
173 {
return scvfs_.size(); }
179 bindElement(element);
181 neighborScvs_.reserve(element.subEntities(1));
182 neighborScvfIndices_.reserve(element.subEntities(1));
183 neighborScvfs_.reserve(element.subEntities(1));
185 std::vector<GridIndexType> handledNeighbors;
186 handledNeighbors.reserve(element.subEntities(1));
187 for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
189 if (intersection.neighbor())
191 const auto outside = intersection.outside();
192 const auto outsideIdx = gridGeometry().elementMapper().index(outside);
195 if ( std::find(handledNeighbors.begin(), handledNeighbors.end(), outsideIdx) == handledNeighbors.end() )
197 makeNeighborGeometries_(outside, outsideIdx);
198 handledNeighbors.push_back(outsideIdx);
208 elementPtr_ = &element;
209 scvfs_.reserve(element.subEntities(1));
210 scvfIndices_.reserve(element.subEntities(1));
211 makeElementGeometries_(element);
216 {
return *gridGeometryPtr_; }
220 {
return hasBoundaryScvf_; }
225 void makeElementGeometries_(
const Element& element)
227 const auto eIdx = gridGeometry().elementMapper().index(element);
228 scvs_[0] = SubControlVolume(element.geometry(), eIdx);
229 scvIndices_[0] = eIdx;
231 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
232 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
234 typename GridGeometry::GeometryHelper geometryHelper(element, gridGeometry().gridView());
237 for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
239 const auto& scvfNeighborVolVarIndex = neighborVolVarIndices[scvfCounter];
241 if (intersection.neighbor() || intersection.boundary())
243 geometryHelper.updateLocalFace(gridGeometry().intersectionMapper(), intersection);
244 std::vector<GridIndexType> scvIndices{eIdx, scvfNeighborVolVarIndex};
245 scvfs_.emplace_back(intersection,
246 intersection.geometry(),
247 scvFaceIndices[scvfCounter],
250 scvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
253 if (intersection.boundary())
254 hasBoundaryScvf_ =
true;
260 void makeNeighborGeometries_(
const Element& element,
const GridIndexType eIdx)
265 neighborScvs_.emplace_back(element.geometry(), eIdx);
266 neighborScvIndices_.push_back(eIdx);
268 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
269 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
271 typename GridGeometry::GeometryHelper geometryHelper(element, gridGeometry().gridView());
274 for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
276 const auto& scvfNeighborVolVarIndex = neighborVolVarIndices[scvfCounter];
277 geometryHelper.updateLocalFace(gridGeometry().intersectionMapper(), intersection);
279 if (intersection.neighbor())
282 if (intersection.outside() == *elementPtr_)
284 std::vector<GridIndexType> scvIndices{eIdx, scvfNeighborVolVarIndex};
285 neighborScvfs_.emplace_back(intersection,
286 intersection.geometry(),
287 scvFaceIndices[scvfCounter],
291 neighborScvfIndices_.push_back(scvFaceIndices[scvfCounter]);
295 else if (intersection.boundary())
300 const LocalIndexType findLocalIndex_(
const GridIndexType idx,
301 const std::vector<GridIndexType>& indices)
const
303 auto it = std::find(indices.begin(), indices.end(), idx);
304 assert(it != indices.end() &&
"Could not find the scv/scvf! Make sure to properly bind this class!");
311 scvfIndices_.clear();
314 neighborScvIndices_.clear();
315 neighborScvfIndices_.clear();
316 neighborScvs_.clear();
317 neighborScvfs_.clear();
319 hasBoundaryScvf_ =
false;
322 const Element* elementPtr_;
323 const GridGeometry* gridGeometryPtr_;
326 std::array<GridIndexType, 1> scvIndices_;
327 std::array<SubControlVolume, 1> scvs_;
329 std::vector<GridIndexType> scvfIndices_;
330 std::vector<SubControlVolumeFace> scvfs_;
332 std::vector<GridIndexType> neighborScvIndices_;
333 std::vector<SubControlVolume> neighborScvs_;
335 std::vector<GridIndexType> neighborScvfIndices_;
336 std::vector<SubControlVolumeFace> neighborScvfs_;
338 bool hasBoundaryScvf_ =
false;
Defines the index types used for grid and local indices.
ctype distance(const Dune::FieldVector< ctype, dimWorld > &a, const Dune::FieldVector< ctype, dimWorld > &b)
Compute the shortest distance between two points.
Definition: geometry/distance.hh:138
typename GridView::IndexSet::IndexType GridIndex
Definition: indextraits.hh:39
unsigned int LocalIndex
Definition: indextraits.hh:40
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models This builds up th...
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:50
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models Specialization fo...
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:60
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:67
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:71
Stencil-local finite volume geometry (scvs and scvfs) for staggered models This builds up the sub con...
Definition: discretization/staggered/fvelementgeometry.hh:42
const SubControlVolumeFace & scvf(GridIndexType eIdx, LocalIndexType localScvfIdx) const
Definition: discretization/staggered/fvelementgeometry.hh:75
StaggeredFVElementGeometry(const CellCenterOrFaceFVGridGeometry &gridGeometry)
Definition: discretization/staggered/fvelementgeometry.hh:70
Base class for the finite volume geometry vector for staggered models This locally builds up the sub ...
Definition: discretization/staggered/fvelementgeometry.hh:91
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition: discretization/staggered/fvelementgeometry.hh:98
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/staggered/fvelementgeometry.hh:219
friend Dune::IteratorRange< typename std::vector< SubControlVolumeFace >::const_iterator > scvfs(const ThisType &g)
Definition: discretization/staggered/fvelementgeometry.hh:161
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume
Definition: discretization/staggered/fvelementgeometry.hh:100
const SubControlVolume & scv(GridIndexType scvIdx) const
Definition: discretization/staggered/fvelementgeometry.hh:124
StaggeredFVElementGeometry(const CellCenterOrFaceFVGridGeometry &gridGeometry)
Definition: discretization/staggered/fvelementgeometry.hh:109
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/staggered/fvelementgeometry.hh:168
const SubControlVolumeFace & scvf(GridIndexType eIdx, LocalIndexType localScvfIdx) const
Get a sub control volume face with an element index and a local scvf index.
Definition: discretization/staggered/fvelementgeometry.hh:117
const GridGeometry & gridGeometry() const
The grid finite volume geometry we are a restriction of.
Definition: discretization/staggered/fvelementgeometry.hh:215
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/staggered/fvelementgeometry.hh:172
StaggeredFVElementGeometry(const GridGeometry &gridGeometry)
Constructor.
Definition: discretization/staggered/fvelementgeometry.hh:113
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition: discretization/staggered/fvelementgeometry.hh:102
friend Dune::IteratorRange< typename std::array< SubControlVolume, 1 >::const_iterator > scvs(const ThisType &g)
Definition: discretization/staggered/fvelementgeometry.hh:149
void bindElement(const Element &element)
Binding of an element preparing the geometries only inside the element.
Definition: discretization/staggered/fvelementgeometry.hh:205
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Definition: discretization/staggered/fvelementgeometry.hh:134
GG GridGeometry
export type of finite volume grid geometry
Definition: discretization/staggered/fvelementgeometry.hh:104
void bind(const Element &element)
Definition: discretization/staggered/fvelementgeometry.hh:177
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models This builds up th...