version 3.8.0
Loading...
Searching...
No Matches
discretization/staggered/fvelementgeometry.hh
Go to the documentation of this file.
1// -*- mode: C++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*-
2// vi: set et ts=4 sw=4 sts=4:
3//
4// SPDX-FileCopyrightInfo: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5// SPDX-License-Identifier: GPL-3.0-or-later
6//
12#ifndef DUMUX_DISCRETIZATION_STAGGERED_FV_ELEMENT_GEOMETRY_HH
13#define DUMUX_DISCRETIZATION_STAGGERED_FV_ELEMENT_GEOMETRY_HH
14
15#include <optional>
19
20namespace Dumux {
21
31template<class GG, bool enableGridGeometryCache>
33
42template<class GG>
43class StaggeredFVElementGeometry<GG, true> : public CCTpfaFVElementGeometry<GG, true>
44{
46 using GridView = typename GG::GridView;
47 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
48 using LocalIndexType = typename IndexTraits<GridView>::LocalIndex;
49public:
51 using Element = typename GridView::template Codim<0>::Entity;
53 using SubControlVolumeFace = typename GG::SubControlVolumeFace;
54
55 using ParentType::ParentType;
56
59 template<class CellCenterOrFaceFVGridGeometry>
60 StaggeredFVElementGeometry(const CellCenterOrFaceFVGridGeometry& gridGeometry)
61 : ParentType(gridGeometry.actualGridGeometry()) {}
62
64 using ParentType::scvf;
65 const SubControlVolumeFace& scvf(GridIndexType eIdx, LocalIndexType localScvfIdx) const
66 {
67 return this->gridGeometry().scvf(eIdx, localScvfIdx);
68 }
69};
70
79template<class GG>
81{
83 using GridView = typename GG::GridView;
84 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
85 using LocalIndexType = typename IndexTraits<GridView>::LocalIndex;
86public:
88 using Element = typename GridView::template Codim<0>::Entity;
90 using SubControlVolume = typename GG::SubControlVolume;
92 using SubControlVolumeFace = typename GG::SubControlVolumeFace;
94 using GridGeometry = GG;
95
98 template<class CellCenterOrFaceFVGridGeometry>
99 StaggeredFVElementGeometry(const CellCenterOrFaceFVGridGeometry& gridGeometry)
100 : gridGeometryPtr_(&gridGeometry.actualGridGeometry()) {}
101
104 : gridGeometryPtr_(&gridGeometry) {}
105
107 const SubControlVolumeFace& scvf(GridIndexType eIdx, LocalIndexType localScvfIdx) const
108 {
109 return scvf(this->gridGeometry().localToGlobalScvfIndex(eIdx, localScvfIdx));
110 }
111
114 const SubControlVolume& scv(GridIndexType scvIdx) const
115 {
116 if (scvIdx == scvIndices_[0])
117 return scvs_[0];
118 else
119 return neighborScvs_[findLocalIndex_(scvIdx, neighborScvIndices_)];
120 }
121
124 const SubControlVolumeFace& scvf(GridIndexType scvfIdx) const
125 {
126 auto it = std::find(scvfIndices_.begin(), scvfIndices_.end(), scvfIdx);
127 if (it != scvfIndices_.end())
128 return scvfs_[std::distance(scvfIndices_.begin(), it)];
129 else
130 return neighborScvfs_[findLocalIndex_(scvfIdx, neighborScvfIndices_)];
131 }
132
138 friend inline Dune::IteratorRange<typename std::array<SubControlVolume, 1>::const_iterator>
139 scvs(const ThisType& g)
140 {
141 using IteratorType = typename std::array<SubControlVolume, 1>::const_iterator;
142 return Dune::IteratorRange<IteratorType>(g.scvs_.begin(), g.scvs_.end());
143 }
144
150 friend inline Dune::IteratorRange<typename std::vector<SubControlVolumeFace>::const_iterator>
151 scvfs(const ThisType& g)
152 {
153 using IteratorType = typename std::vector<SubControlVolumeFace>::const_iterator;
154 return Dune::IteratorRange<IteratorType>(g.scvfs_.begin(), g.scvfs_.end());
155 }
156
158 std::size_t numScv() const
159 { return scvs_.size(); }
160
162 std::size_t numScvf() const
163 { return scvfs_.size(); }
164
171 {
172 this->bind_(element);
173 return std::move(*this);
174 }
175
177 void bind(const Element& element) &
178 { this->bind_(element); }
179
186 {
187 this->bindElement_(element);
188 return std::move(*this);
189 }
190
192 void bindElement(const Element& element) &
193 { this->bindElement_(element); }
194
196 bool isBound() const
197 { return static_cast<bool>(element_); }
198
200 const Element& element() const
201 { return *element_; }
202
205 { return *gridGeometryPtr_; }
206
208 bool hasBoundaryScvf() const
209 { return hasBoundaryScvf_; }
210
212 typename SubControlVolume::Traits::Geometry geometry(const SubControlVolume& scv) const
213 { return gridGeometryPtr_->element(scv.dofIndex()).geometry(); }
214
216 typename SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace& scvf) const
217 {
218 const auto insideElementIndex = scvf.insideScvIdx();
219 const auto elementGeometry = (insideElementIndex != scvIndices_[0]) ?
220 element_->geometry() :
221 gridGeometryPtr_->element(insideElementIndex).geometry();
222 const auto center = elementGeometry.center();
223 const auto normalAxis = Dumux::normalAxis(scvf.unitOuterNormal());
224
225 auto lowerLeft = elementGeometry.corner(0);
226 auto upperRight = elementGeometry.corner(elementGeometry.corners()-1);
227
228 // shift corners to scvf plane and halve lateral faces
229 lowerLeft[normalAxis] = center[normalAxis];
230 upperRight[normalAxis] = center[normalAxis];
231
232 auto inPlaneAxes = std::move(std::bitset<SubControlVolumeFace::Traits::dimWorld>{}.set());
233 inPlaneAxes.set(normalAxis, false);
234
235 return {lowerLeft, upperRight, inPlaneAxes};
236 }
237
238private:
240 void bindElement_(const Element& element)
241 {
242 clear_();
243 element_ = element;
244 scvfs_.reserve(element.subEntities(1));
245 scvfIndices_.reserve(element.subEntities(1));
246 makeElementGeometries_();
247 }
248
251 void bind_(const Element& element)
252 {
253 bindElement_(element);
254
255 neighborScvs_.reserve(element.subEntities(1));
256 neighborScvfIndices_.reserve(element.subEntities(1));
257 neighborScvfs_.reserve(element.subEntities(1));
258
259 std::vector<GridIndexType> handledNeighbors;
260 handledNeighbors.reserve(element.subEntities(1));
261 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
262 {
263 if (intersection.neighbor())
264 {
265 const auto outside = intersection.outside();
266 const auto outsideIdx = gridGeometry().elementMapper().index(outside);
267
268 // make outside geometries only if not done yet (could happen on non-conforming grids)
269 if ( std::find(handledNeighbors.begin(), handledNeighbors.end(), outsideIdx) == handledNeighbors.end() )
270 {
271 makeNeighborGeometries_(outside, outsideIdx);
272 handledNeighbors.push_back(outsideIdx);
273 }
274 }
275 }
276 }
277
279 void makeElementGeometries_()
280 {
281 const auto& element = *element_;
282 const auto eIdx = gridGeometry().elementMapper().index(element);
283 scvs_[0] = SubControlVolume(element.geometry(), eIdx);
284 scvIndices_[0] = eIdx;
285
286 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
287 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
288
289 typename GridGeometry::GeometryHelper geometryHelper(element, gridGeometry().gridView());
290
291 int scvfCounter = 0;
292 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
293 {
294 const auto& scvfNeighborVolVarIndex = neighborVolVarIndices[scvfCounter];
295
296 if (intersection.neighbor() || intersection.boundary())
297 {
298 geometryHelper.updateLocalFace(gridGeometry().intersectionMapper(), intersection);
299 std::vector<GridIndexType> scvIndices{eIdx, scvfNeighborVolVarIndex};
300 scvfs_.emplace_back(intersection,
301 intersection.geometry(),
302 scvFaceIndices[scvfCounter],
303 scvIndices,
304 geometryHelper);
305 scvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
306 scvfCounter++;
307
308 if (intersection.boundary())
309 hasBoundaryScvf_ = true;
310 }
311 }
312 }
313
315 void makeNeighborGeometries_(const Element& element, const GridIndexType eIdx)
316 {
317 // using ScvfGridIndexStorage = typename SubControlVolumeFace::Traits::GridIndexStorage;
318
319 // create the neighbor scv
320 neighborScvs_.emplace_back(element.geometry(), eIdx);
321 neighborScvIndices_.push_back(eIdx);
322
323 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
324 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
325
326 typename GridGeometry::GeometryHelper geometryHelper(element, gridGeometry().gridView());
327
328 int scvfCounter = 0;
329 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
330 {
331 const auto& scvfNeighborVolVarIndex = neighborVolVarIndices[scvfCounter];
332 geometryHelper.updateLocalFace(gridGeometry().intersectionMapper(), intersection);
333
334 if (intersection.neighbor())
335 {
336 // only create subcontrol faces where the outside element is the bound element
337 if (intersection.outside() == *element_)
338 {
339 std::vector<GridIndexType> scvIndices{eIdx, scvfNeighborVolVarIndex};
340 neighborScvfs_.emplace_back(intersection,
341 intersection.geometry(),
342 scvFaceIndices[scvfCounter],
343 scvIndices,
344 geometryHelper);
345
346 neighborScvfIndices_.push_back(scvFaceIndices[scvfCounter]);
347 }
348 scvfCounter++;
349 }
350 else if (intersection.boundary())
351 scvfCounter++;
352 }
353 }
354
355 const LocalIndexType findLocalIndex_(const GridIndexType idx,
356 const std::vector<GridIndexType>& indices) const
357 {
358 auto it = std::find(indices.begin(), indices.end(), idx);
359 assert(it != indices.end() && "Could not find the scv/scvf! Make sure to properly bind this class!");
360 return std::distance(indices.begin(), it);
361 }
362
364 void clear_()
365 {
366 scvfIndices_.clear();
367 scvfs_.clear();
368
369 neighborScvIndices_.clear();
370 neighborScvfIndices_.clear();
371 neighborScvs_.clear();
372 neighborScvfs_.clear();
373
374 hasBoundaryScvf_ = false;
375 }
376
377 std::optional<Element> element_;
378 const GridGeometry* gridGeometryPtr_;
379
380 // local storage after binding an element
381 std::array<GridIndexType, 1> scvIndices_;
382 std::array<SubControlVolume, 1> scvs_;
383
384 std::vector<GridIndexType> scvfIndices_;
385 std::vector<SubControlVolumeFace> scvfs_;
386
387 std::vector<GridIndexType> neighborScvIndices_;
388 std::vector<SubControlVolume> neighborScvs_;
389
390 std::vector<GridIndexType> neighborScvfIndices_;
391 std::vector<SubControlVolumeFace> neighborScvfs_;
392
393 bool hasBoundaryScvf_ = false;
394};
395
396
397} // end namespace
398
399#endif
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models Specialization fo...
Definition discretization/cellcentered/tpfa/fvelementgeometry.hh:63
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models This builds up th...
Definition discretization/cellcentered/tpfa/fvelementgeometry.hh:53
Base class for the finite volume geometry vector for staggered models This locally builds up the sub ...
Definition discretization/staggered/fvelementgeometry.hh:81
StaggeredFVElementGeometry bind(const Element &element) &&
bind the local view (r-value overload) This overload is called when an instance of this class is a te...
Definition discretization/staggered/fvelementgeometry.hh:170
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition discretization/staggered/fvelementgeometry.hh:88
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition discretization/staggered/fvelementgeometry.hh:208
void bindElement(const Element &element) &
bind this local view to a specific element
Definition discretization/staggered/fvelementgeometry.hh:192
friend Dune::IteratorRange< typename std::vector< SubControlVolumeFace >::const_iterator > scvfs(const ThisType &g)
Definition discretization/staggered/fvelementgeometry.hh:151
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume
Definition discretization/staggered/fvelementgeometry.hh:90
SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace &scvf) const
Create the geometry of a given sub control volume face.
Definition discretization/staggered/fvelementgeometry.hh:216
const SubControlVolume & scv(GridIndexType scvIdx) const
Definition discretization/staggered/fvelementgeometry.hh:114
StaggeredFVElementGeometry(const CellCenterOrFaceFVGridGeometry &gridGeometry)
Definition discretization/staggered/fvelementgeometry.hh:99
SubControlVolume::Traits::Geometry geometry(const SubControlVolume &scv) const
Create the geometry of a given sub control volume.
Definition discretization/staggered/fvelementgeometry.hh:212
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition discretization/staggered/fvelementgeometry.hh:158
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:107
const GridGeometry & gridGeometry() const
The grid finite volume geometry we are a restriction of.
Definition discretization/staggered/fvelementgeometry.hh:204
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition discretization/staggered/fvelementgeometry.hh:162
StaggeredFVElementGeometry bindElement(const Element &element) &&
bind the local view (r-value overload) This overload is called when an instance of this class is a te...
Definition discretization/staggered/fvelementgeometry.hh:185
StaggeredFVElementGeometry(const GridGeometry &gridGeometry)
Constructor.
Definition discretization/staggered/fvelementgeometry.hh:103
const Element & element() const
The bound element.
Definition discretization/staggered/fvelementgeometry.hh:200
void bind(const Element &element) &
bind this local view to a specific element (full stencil)
Definition discretization/staggered/fvelementgeometry.hh:177
bool isBound() const
Returns true if bind/bindElement has already been called.
Definition discretization/staggered/fvelementgeometry.hh:196
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition discretization/staggered/fvelementgeometry.hh:92
friend Dune::IteratorRange< typename std::array< SubControlVolume, 1 >::const_iterator > scvs(const ThisType &g)
Definition discretization/staggered/fvelementgeometry.hh:139
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Definition discretization/staggered/fvelementgeometry.hh:124
GG GridGeometry
export type of finite volume grid geometry
Definition discretization/staggered/fvelementgeometry.hh:94
const SubControlVolumeFace & scvf(GridIndexType eIdx, LocalIndexType localScvfIdx) const
Definition discretization/staggered/fvelementgeometry.hh:65
StaggeredFVElementGeometry(const CellCenterOrFaceFVGridGeometry &gridGeometry)
Definition discretization/staggered/fvelementgeometry.hh:60
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition discretization/staggered/fvelementgeometry.hh:53
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition discretization/staggered/fvelementgeometry.hh:51
Stencil-local finite volume geometry (scvs and scvfs) for staggered models This builds up the sub con...
Definition discretization/staggered/fvelementgeometry.hh:32
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models This builds up th...
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
Corners::value_type center(const Corners &corners)
The center of a given list of corners.
Definition center.hh:24
Defines the index types used for grid and local indices.
Definition adapt.hh:17
Base class for the finite volume geometry vector for face-centered staggered models This builds up th...
typename GridView::IndexSet::IndexType GridIndex
Definition indextraits.hh:27
unsigned int LocalIndex
Definition indextraits.hh:28