3.3.0
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
discretization/cellcentered/tpfa/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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
26#ifndef DUMUX_DISCRETIZATION_CCTPFA_FV_ELEMENT_GEOMETRY_HH
27#define DUMUX_DISCRETIZATION_CCTPFA_FV_ELEMENT_GEOMETRY_HH
28
29#include <algorithm>
30#include <array>
31#include <vector>
32
33#include <dune/common/exceptions.hh>
35#include <dune/common/iteratorrange.hh>
37
38namespace Dumux {
39
49template<class GG, bool enableGridGeometryCache>
51
58template<class GG>
60{
62 using GridView = typename GG::GridView;
63 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
64
65public:
67 using Element = typename GridView::template Codim<0>::Entity;
69 using SubControlVolume = typename GG::SubControlVolume;
71 using SubControlVolumeFace = typename GG::SubControlVolumeFace;
73 using GridGeometry = GG;
74
76 static constexpr std::size_t maxNumElementScvs = 1;
78 static constexpr std::size_t maxNumElementScvfs = 2*GridView::dimension;
79
82 : gridGeometryPtr_(&gridGeometry) {}
83
86 const SubControlVolume& scv(GridIndexType scvIdx) const
87 {
88 return gridGeometry().scv(scvIdx);
89 }
90
93 const SubControlVolumeFace& scvf(GridIndexType scvfIdx) const
94 {
95 return gridGeometry().scvf(scvfIdx);
96 }
97
100 const SubControlVolumeFace& flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx = 0) const
101 {
102 return gridGeometry().flipScvf(scvfIdx, outsideScvIdx);
103 }
104
110 friend inline Dune::IteratorRange< ScvIterator<SubControlVolume, std::array<GridIndexType, 1>, ThisType> >
111 scvs(const CCTpfaFVElementGeometry& fvGeometry)
112 {
114 return Dune::IteratorRange<ScvIterator>(ScvIterator(fvGeometry.scvIndices_.begin(), fvGeometry),
115 ScvIterator(fvGeometry.scvIndices_.end(), fvGeometry));
116 }
117
123 friend inline Dune::IteratorRange< ScvfIterator<SubControlVolumeFace, std::vector<GridIndexType>, ThisType> >
125 {
126 const auto& g = fvGeometry.gridGeometry();
127 const auto scvIdx = fvGeometry.scvIndices_[0];
129 return Dune::IteratorRange<ScvfIterator>(ScvfIterator(g.scvfIndicesOfScv(scvIdx).begin(), fvGeometry),
130 ScvfIterator(g.scvfIndicesOfScv(scvIdx).end(), fvGeometry));
131 }
132
134 std::size_t numScv() const
135 {
136 return scvIndices_.size();
137 }
138
140 std::size_t numScvf() const
141 {
142 return gridGeometry().scvfIndicesOfScv(scvIndices_[0]).size();
143 }
144
146 void bind(const Element& element)
147 {
148 this->bindElement(element);
149 }
150
152 void bindElement(const Element& element)
153 {
154 elementPtr_ = &element;
155 scvIndices_[0] = gridGeometry().elementMapper().index(*elementPtr_);
156 }
157
160 { return *gridGeometryPtr_; }
161
163 bool hasBoundaryScvf() const
164 { return gridGeometry().hasBoundaryScvf(scvIndices_[0]); }
165
166private:
167
168 const Element* elementPtr_;
169 std::array<GridIndexType, 1> scvIndices_;
170 const GridGeometry* gridGeometryPtr_;
171};
172
178template<class GG>
180{
182 using GridView = typename GG::GridView;
183 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
184 using LocalIndexType = typename IndexTraits<GridView>::LocalIndex;
185
186 static const int dim = GridView::dimension;
187 static const int dimWorld = GridView::dimensionworld;
188
189public:
191 using Element = typename GridView::template Codim<0>::Entity;
193 using SubControlVolume = typename GG::SubControlVolume;
195 using SubControlVolumeFace = typename GG::SubControlVolumeFace;
197 using GridGeometry = GG;
199 static constexpr std::size_t maxNumElementScvs = 1;
201 static constexpr std::size_t maxNumElementScvfs = 2*dim;
202
205 : gridGeometryPtr_(&gridGeometry) {}
206
209 const SubControlVolume& scv(GridIndexType scvIdx) const
210 {
211 if (scvIdx == scvIndices_[0])
212 return scvs_[0];
213 else
214 return neighborScvs_[findLocalIndex(scvIdx, neighborScvIndices_)];
215 }
216
219 const SubControlVolumeFace& scvf(GridIndexType scvfIdx) const
220 {
221 auto it = std::find(scvfIndices_.begin(), scvfIndices_.end(), scvfIdx);
222 if (it != scvfIndices_.end())
223 return scvfs_[std::distance(scvfIndices_.begin(), it)];
224 else
225 return neighborScvfs_[findLocalIndex(scvfIdx, neighborScvfIndices_)];
226 }
227
230 const SubControlVolumeFace& flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx = 0) const
231 {
232 auto it = std::find(scvfIndices_.begin(), scvfIndices_.end(), scvfIdx);
233 if (it != scvfIndices_.end())
234 {
235 const auto localScvfIdx = std::distance(scvfIndices_.begin(), it);
236 return neighborScvfs_[flippedScvfIndices_[localScvfIdx][outsideScvIdx]];
237 }
238 else
239 {
240 const auto localScvfIdx = findLocalIndex(scvfIdx, neighborScvfIndices_);
241 const auto localFlippedIndex = flippedNeighborScvfIndices_[localScvfIdx][outsideScvIdx];
242 if (localFlippedIndex < scvfs_.size())
243 return scvfs_[localFlippedIndex];
244 else
245 return neighborScvfs_[localFlippedIndex - scvfs_.size()];
246 }
247 }
248
254 friend inline Dune::IteratorRange<typename std::array<SubControlVolume, 1>::const_iterator>
255 scvs(const ThisType& g)
256 {
257 using IteratorType = typename std::array<SubControlVolume, 1>::const_iterator;
258 return Dune::IteratorRange<IteratorType>(g.scvs_.begin(), g.scvs_.end());
259 }
260
266 friend inline Dune::IteratorRange<typename std::vector<SubControlVolumeFace>::const_iterator>
267 scvfs(const ThisType& g)
268 {
269 using IteratorType = typename std::vector<SubControlVolumeFace>::const_iterator;
270 return Dune::IteratorRange<IteratorType>(g.scvfs_.begin(), g.scvfs_.end());
271 }
272
274 std::size_t numScv() const
275 { return scvs_.size(); }
276
278 std::size_t numScvf() const
279 { return scvfs_.size(); }
280
283 void bind(const Element& element)
284 {
285 bindElement(element);
286
287 neighborScvs_.reserve(element.subEntities(1));
288 neighborScvfIndices_.reserve(element.subEntities(1));
289 neighborScvfs_.reserve(element.subEntities(1));
290
291 std::vector<GridIndexType> handledNeighbors;
292 handledNeighbors.reserve(element.subEntities(1));
293
294 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
295 {
296 // for inner intersections and periodic (according to grid interface) intersections make neighbor geometry
297 if (intersection.neighbor())
298 {
299 const auto outside = intersection.outside();
300 const auto outsideIdx = gridGeometry().elementMapper().index(outside);
301
302 // make outside geometries only if not done yet (could happen on non-conforming grids)
303 if ( std::find(handledNeighbors.begin(), handledNeighbors.end(), outsideIdx) == handledNeighbors.end() )
304 {
305 makeNeighborGeometries(outside, outsideIdx);
306 handledNeighbors.push_back(outsideIdx);
307 }
308 }
309 }
310
311 // build flip index set for network, surface, and periodic grids
312 if (dim < dimWorld || gridGeometry().isPeriodic())
313 {
314 flippedScvfIndices_.resize(scvfs_.size());
315 for (unsigned int localScvfIdx = 0; localScvfIdx < scvfs_.size(); ++localScvfIdx)
316 {
317 const auto& scvf = scvfs_[localScvfIdx];
318 if (scvf.boundary())
319 continue;
320
321 flippedScvfIndices_[localScvfIdx].resize(scvf.numOutsideScvs());
322 for (unsigned int localOutsideScvIdx = 0; localOutsideScvIdx < scvf.numOutsideScvs(); ++localOutsideScvIdx)
323 {
324 const auto globalOutsideScvIdx = scvf.outsideScvIdx(localOutsideScvIdx);
325 for (unsigned int localNeighborScvfIdx = 0; localNeighborScvfIdx < neighborScvfs_.size(); ++localNeighborScvfIdx)
326 {
327 if (neighborScvfs_[localNeighborScvfIdx].insideScvIdx() == globalOutsideScvIdx)
328 {
329 flippedScvfIndices_[localScvfIdx][localOutsideScvIdx] = localNeighborScvfIdx;
330 break;
331 }
332 }
333 }
334 }
335
336 flippedNeighborScvfIndices_.resize(neighborScvfs_.size());
337 for (unsigned int localScvfIdx = 0; localScvfIdx < neighborScvfs_.size(); ++localScvfIdx)
338 {
339 const auto& neighborScvf = neighborScvfs_[localScvfIdx];
340 flippedNeighborScvfIndices_[localScvfIdx].resize(neighborScvf.numOutsideScvs());
341 for (unsigned int localOutsideScvIdx = 0; localOutsideScvIdx < neighborScvf.numOutsideScvs(); ++localOutsideScvIdx)
342 {
343 flippedNeighborScvfIndices_[localScvfIdx][localOutsideScvIdx] = findFlippedScvfIndex_(neighborScvf.insideScvIdx(), neighborScvf.outsideScvIdx(localOutsideScvIdx));
344 }
345 }
346 }
347 }
348
350 void bindElement(const Element& element)
351 {
352 clear();
353 elementPtr_ = &element;
354 scvfs_.reserve(element.subEntities(1));
355 scvfIndices_.reserve(element.subEntities(1));
356 makeElementGeometries(element);
357 }
358
361 { return *gridGeometryPtr_; }
362
364 bool hasBoundaryScvf() const
365 { return hasBoundaryScvf_; }
366
367private:
368
369 GridIndexType findFlippedScvfIndex_(GridIndexType insideScvIdx, GridIndexType globalOutsideScvIdx)
370 {
371 for (unsigned int localNeighborScvfIdx = 0; localNeighborScvfIdx < neighborScvfs_.size(); ++localNeighborScvfIdx)
372 {
373 if (neighborScvfs_[localNeighborScvfIdx].insideScvIdx() == globalOutsideScvIdx)
374 {
375 return scvfs_.size() + localNeighborScvfIdx;
376 }
377 }
378
379 // go over all potential scvfs of the outside scv
380 for (unsigned int localOutsideScvfIdx = 0; localOutsideScvfIdx < scvfs_.size(); ++localOutsideScvfIdx)
381 {
382 const auto& outsideScvf = scvfs_[localOutsideScvfIdx];
383 for (unsigned int j = 0; j < outsideScvf.numOutsideScvs(); ++j)
384 {
385 if (outsideScvf.outsideScvIdx(j) == insideScvIdx)
386 {
387 return localOutsideScvfIdx;
388 }
389 }
390 }
391
392 DUNE_THROW(Dune::InvalidStateException, "No flipped version of this scvf found!");
393 }
394
396 void makeElementGeometries(const Element& element)
397 {
398 using ScvfGridIndexStorage = typename SubControlVolumeFace::Traits::GridIndexStorage;
399
400 const auto eIdx = gridGeometry().elementMapper().index(element);
401 scvs_[0] = SubControlVolume(element.geometry(), eIdx);
402 scvIndices_[0] = eIdx;
403
404 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
405 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
406
407 // for network grids there might be multiple intersection with the same geometryInInside
408 // we indentify those by the indexInInside for now (assumes conforming grids at branching facets)
409 // here we keep track of them
410 std::vector<bool> handledScvf;
411 if (dim < dimWorld)
412 handledScvf.resize(element.subEntities(1), false);
413
414 int scvfCounter = 0;
415 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
416 {
417 if (dim < dimWorld)
418 if (handledScvf[intersection.indexInInside()])
419 continue;
420
421 const auto& scvfNeighborVolVarIndices = neighborVolVarIndices[scvfCounter];
422 if (intersection.neighbor() || intersection.boundary())
423 {
424 ScvfGridIndexStorage scvIndices;
425 scvIndices.resize(scvfNeighborVolVarIndices.size() + 1);
426 scvIndices[0] = eIdx;
427 std::copy(scvfNeighborVolVarIndices.begin(), scvfNeighborVolVarIndices.end(), scvIndices.begin()+1);
428
429 const bool onBoundary = intersection.boundary() && !intersection.neighbor();
430 hasBoundaryScvf_ = (hasBoundaryScvf_ || onBoundary);
431
432 scvfs_.emplace_back(intersection,
433 intersection.geometry(),
434 scvFaceIndices[scvfCounter],
435 scvIndices,
436 onBoundary);
437 scvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
438 scvfCounter++;
439
440 // for surface and network grids mark that we handled this face
441 if (dim < dimWorld)
442 handledScvf[intersection.indexInInside()] = true;
443 }
444 }
445 }
446
448 void makeNeighborGeometries(const Element& element, const GridIndexType eIdx)
449 {
450 using ScvfGridIndexStorage = typename SubControlVolumeFace::Traits::GridIndexStorage;
451
452 // create the neighbor scv
453 neighborScvs_.emplace_back(element.geometry(), eIdx);
454 neighborScvIndices_.push_back(eIdx);
455
456 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
457 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
458
459 // for network grids there might be multiple intersection with the same geometryInInside
460 // we indentify those by the indexInInside for now (assumes conforming grids at branching facets)
461 // here we keep track of them
462 std::vector<bool> handledScvf;
463 if (dim < dimWorld)
464 handledScvf.resize(element.subEntities(1), false);
465
466 int scvfCounter = 0;
467 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
468 {
469 if (dim < dimWorld)
470 if (handledScvf[intersection.indexInInside()])
471 continue;
472
473 if (intersection.neighbor())
474 {
475 // this catches inner and periodic scvfs
476 const auto& scvfNeighborVolVarIndices = neighborVolVarIndices[scvfCounter];
477 if (scvfNeighborVolVarIndices[0] < gridGeometry().gridView().size(0))
478 {
479 // only create subcontrol faces where the outside element is the bound element
480 if (dim == dimWorld)
481 {
482 if (scvfNeighborVolVarIndices[0] == gridGeometry().elementMapper().index(*elementPtr_))
483 {
484 ScvfGridIndexStorage scvIndices({eIdx, scvfNeighborVolVarIndices[0]});
485 neighborScvfs_.emplace_back(intersection,
486 intersection.geometry(),
487 scvFaceIndices[scvfCounter],
488 scvIndices,
489 false);
490
491 neighborScvfIndices_.push_back(scvFaceIndices[scvfCounter]);
492 }
493 }
494 // for network grids we can't use the intersection.outside() index as we can't assure that the
495 // first intersection with this indexInInside is the one that has our bound element as outside
496 // instead we check if the bound element's index is in the outsideScvIndices of the candidate scvf
497 // (will be optimized away for dim == dimWorld)
498 else
499 {
500 for (unsigned outsideScvIdx = 0; outsideScvIdx < scvfNeighborVolVarIndices.size(); ++outsideScvIdx)
501 {
502 if (scvfNeighborVolVarIndices[outsideScvIdx] == gridGeometry().elementMapper().index(*elementPtr_))
503 {
504 ScvfGridIndexStorage scvIndices;
505 scvIndices.resize(scvfNeighborVolVarIndices.size() + 1);
506 scvIndices[0] = eIdx;
507 std::copy(scvfNeighborVolVarIndices.begin(), scvfNeighborVolVarIndices.end(), scvIndices.begin()+1);
508 neighborScvfs_.emplace_back(intersection,
509 intersection.geometry(),
510 scvFaceIndices[scvfCounter],
511 scvIndices,
512 false);
513
514 neighborScvfIndices_.push_back(scvFaceIndices[scvfCounter]);
515 break;
516 }
517 }
518 }
519
520 // for surface and network grids mark that we handled this face
521 if (dim < dimWorld)
522 handledScvf[intersection.indexInInside()] = true;
523 scvfCounter++;
524 }
525 }
526
527 // only increase counter for boundary intersections
528 // (exclude periodic boundaries according to dune grid interface, they have been handled in neighbor==true)
529 else if (intersection.boundary())
530 {
531 if (dim < dimWorld)
532 handledScvf[intersection.indexInInside()] = true;
533 scvfCounter++;
534 }
535 }
536 }
537
538 const LocalIndexType findLocalIndex(const GridIndexType idx,
539 const std::vector<GridIndexType>& indices) const
540 {
541 auto it = std::find(indices.begin(), indices.end(), idx);
542 assert(it != indices.end() && "Could not find the scv/scvf! Make sure to properly bind this class!");
543 return std::distance(indices.begin(), it);
544
545 }
546
548 void clear()
549 {
550 scvfIndices_.clear();
551 scvfs_.clear();
552 flippedScvfIndices_.clear();
553
554 neighborScvIndices_.clear();
555 neighborScvfIndices_.clear();
556 neighborScvs_.clear();
557 neighborScvfs_.clear();
558 flippedNeighborScvfIndices_.clear();
559
560 hasBoundaryScvf_ = false;
561 }
562
563 const Element* elementPtr_;
564 const GridGeometry* gridGeometryPtr_;
565
566 // local storage after binding an element
567 std::array<GridIndexType, 1> scvIndices_;
568 std::array<SubControlVolume, 1> scvs_;
569
570 std::vector<GridIndexType> scvfIndices_;
571 std::vector<SubControlVolumeFace> scvfs_;
572 std::vector<std::vector<GridIndexType>> flippedScvfIndices_;
573
574 std::vector<GridIndexType> neighborScvIndices_;
575 std::vector<SubControlVolume> neighborScvs_;
576
577 std::vector<GridIndexType> neighborScvfIndices_;
578 std::vector<SubControlVolumeFace> neighborScvfs_;
579 std::vector<std::vector<GridIndexType>> flippedNeighborScvfIndices_;
580
581 bool hasBoundaryScvf_ = false;
582};
583
584} // end namespace Dumux
585
586#endif
Defines the index types used for grid and local indices.
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: geometry/distance.hh:138
Definition: adapt.hh:29
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
CCTpfaFVElementGeometry(const GridGeometry &gridGeometry)
Constructor.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:81
friend Dune::IteratorRange< ScvfIterator< SubControlVolumeFace, std::vector< GridIndexType >, ThisType > > scvfs(const CCTpfaFVElementGeometry &fvGeometry)
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:124
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:134
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:163
void bind(const Element &element)
Binding of an element, called by the local jacobian to prepare element assembly.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:146
GG GridGeometry
export type of finite volume grid geometry
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:73
const GridGeometry & gridGeometry() const
The global finite volume geometry we are a restriction of.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:159
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:140
const SubControlVolume & scv(GridIndexType scvIdx) const
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:86
void bindElement(const Element &element)
Bind only element-local.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:152
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:93
const SubControlVolumeFace & flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx=0) const
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:100
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:69
friend Dune::IteratorRange< ScvIterator< SubControlVolume, std::array< GridIndexType, 1 >, ThisType > > scvs(const CCTpfaFVElementGeometry &fvGeometry)
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:111
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models Specialization fo...
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:180
friend Dune::IteratorRange< typename std::vector< SubControlVolumeFace >::const_iterator > scvfs(const ThisType &g)
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:267
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:195
const SubControlVolume & scv(GridIndexType scvIdx) const
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:209
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:274
CCTpfaFVElementGeometry(const GridGeometry &gridGeometry)
Constructor.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:204
void bindElement(const Element &element)
Binding of an element preparing the geometries only inside the element.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:350
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:278
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:191
const GridGeometry & gridGeometry() const
The global finite volume geometry we are a restriction of.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:360
const SubControlVolumeFace & flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx=0) const
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:230
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:219
void bind(const Element &element)
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:283
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:193
friend Dune::IteratorRange< typename std::array< SubControlVolume, 1 >::const_iterator > scvs(const ThisType &g)
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:255
GG GridGeometry
export type of finite volume grid geometry
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:197
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/cellcentered/tpfa/fvelementgeometry.hh:364
Iterators over sub control volumes.
Definition: scvandscvfiterators.hh:42
Iterators over sub control volume faces of an fv geometry.
Definition: scvandscvfiterators.hh:82