3.5-git
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
discretization/cellcentered/mpfa/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_CCMPFA_FV_ELEMENT_GEOMETRY_HH
27#define DUMUX_DISCRETIZATION_CCMPFA_FV_ELEMENT_GEOMETRY_HH
28
29#include <optional>
30#include <utility>
31
32#include <dune/common/exceptions.hh>
33#include <dune/common/iteratorrange.hh>
34
38
39namespace Dumux {
40
50template<class GG, bool enableGridGeometryCache>
52
59template<class GG>
61{
63 using GridView = typename GG::GridView;
64 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
65
66 static constexpr int dim = GridView::dimension;
67
68public:
70 using Element = typename GridView::template Codim<0>::Entity;
72 using SubControlVolume = typename GG::SubControlVolume;
74 using SubControlVolumeFace = typename GG::SubControlVolumeFace;
76 using GridGeometry = GG;
78 static constexpr std::size_t maxNumElementScvs = 1;
80 static constexpr std::size_t maxNumElementScvfs = dim == 3 ? 24 : 8;
81
84 : gridGeometryPtr_(&gridGeometry) {}
85
87 const SubControlVolume& scv(GridIndexType scvIdx) const
88 {
89 return gridGeometry().scv(scvIdx);
90 }
91
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 CCMpfaFVElementGeometry& 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
151 {
152 this->bindElement(element);
153 return std::move(*this);
154 }
155
156 void bind(const Element& element) &
157 {
158 this->bindElement(element);
159 }
160
167 {
168 this->bindElement(element);
169 return std::move(*this);
170 }
171
173 void bindElement(const Element& element) &
174 {
175 element_ = element;
176 scvIndices_[0] = gridGeometry().elementMapper().index(element);
177 }
178
180 bool isBound() const
181 { return static_cast<bool>(element_); }
182
184 const Element& element() const
185 { return *element_; }
186
189 { return *gridGeometryPtr_; }
190
192 bool hasBoundaryScvf() const
193 { return gridGeometry().hasBoundaryScvf(scvIndices_[0]); }
194
195private:
196
197 std::optional<Element> element_;
198 std::array<GridIndexType, 1> scvIndices_;
199 const GridGeometry* gridGeometryPtr_;
200};
201
207template<class GG>
209{
211 using GridView = typename GG::GridView;
212 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
213 using MpfaHelper = typename GG::MpfaHelper;
214
215 static const int dim = GridView::dimension;
216 static const int dimWorld = GridView::dimensionworld;
217 using CoordScalar = typename GridView::ctype;
218
219public:
221 using Element = typename GridView::template Codim<0>::Entity;
223 using SubControlVolume = typename GG::SubControlVolume;
225 using SubControlVolumeFace = typename GG::SubControlVolumeFace;
227 using GridGeometry = GG;
229 static constexpr std::size_t maxNumElementScvs = 1;
231 static constexpr std::size_t maxNumElementScvfs = dim == 3 ? 24 : 8;
232
235 : gridGeometryPtr_(&gridGeometry) {}
236
239 const SubControlVolume& scv(GridIndexType scvIdx) const
240 {
241 if (scvIdx == scvIndices_[0])
242 return scvs_[0];
243 else
244 return neighborScvs_[findLocalIndex(scvIdx, neighborScvIndices_)];
245 }
246
249 const SubControlVolumeFace& scvf(GridIndexType scvfIdx) const
250 {
251 auto it = std::find(scvfIndices_.begin(), scvfIndices_.end(), scvfIdx);
252 if (it != scvfIndices_.end())
253 return scvfs_[std::distance(scvfIndices_.begin(), it)];
254 else
255 return neighborScvfs_[findLocalIndex(scvfIdx, neighborScvfIndices_)];
256 }
257
260 const SubControlVolumeFace& flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx = 0) const
261 {
262 return scvf( gridGeometry().flipScvfIdx(scvfIdx, outsideScvIdx) );
263 }
264
270 friend inline Dune::IteratorRange<typename std::array<SubControlVolume, 1>::const_iterator>
271 scvs(const ThisType& g)
272 {
273 using IteratorType = typename std::array<SubControlVolume, 1>::const_iterator;
274 return Dune::IteratorRange<IteratorType>(g.scvs_.begin(), g.scvs_.end());
275 }
276
282 friend inline Dune::IteratorRange<typename std::vector<SubControlVolumeFace>::const_iterator>
283 scvfs(const ThisType& g)
284 {
285 using IteratorType = typename std::vector<SubControlVolumeFace>::const_iterator;
286 return Dune::IteratorRange<IteratorType>(g.scvfs_.begin(), g.scvfs_.end());
287 }
288
290 std::size_t numScv() const
291 { return scvs_.size(); }
292
294 std::size_t numScvf() const
295 { return scvfs_.size(); }
296
303 {
304 this->bindElement_(element);
305 return std::move(*this);
306 }
307
308 void bindElement(const Element& element) &
309 {
310 this->bindElement_(element);
311 }
312
319 {
320 this->bind_(element);
321 return std::move(*this);
322 }
323
324 void bind(const Element& element) &
325 {
326 this->bind_(element);
327 }
328
330 bool isBound() const
331 { return static_cast<bool>(element_); }
332
334 const Element& element() const
335 { return *element_; }
336
339 { return *gridGeometryPtr_; }
340
342 bool hasBoundaryScvf() const
343 { return hasBoundaryScvf_; }
344
345private:
346
349 void bind_(const Element& element)
350 {
351 // make inside geometries
352 bindElement_(element);
353
354 // get some references for convenience
355 const auto globalI = gridGeometry().elementMapper().index(element);
356 const auto& assemblyMapI = gridGeometry().connectivityMap()[globalI];
357
358 // reserve memory
359 const auto numNeighbors = assemblyMapI.size();
360 const auto numNeighborScvfs = numNeighborScvfs_(assemblyMapI);
361 neighborScvs_.reserve(numNeighbors);
362 neighborScvIndices_.reserve(numNeighbors);
363 neighborScvfIndices_.reserve(numNeighborScvfs);
364 neighborScvfs_.reserve(numNeighborScvfs);
365
366 // make neighbor geometries
367 // use the assembly map to determine which faces are necessary
368 for (const auto& dataJ : assemblyMapI)
369 makeNeighborGeometries(gridGeometry().element(dataJ.globalJ),
370 dataJ.globalJ,
371 dataJ.scvfsJ,
372 dataJ.additionalScvfs);
373
374 // //! TODO Check if user added additional DOF dependencies, i.e. the residual of DOF globalI depends
375 // //! on additional DOFs not included in the discretization schemes' occupation pattern
376 // const auto& additionalDofDependencies = problem.getAdditionalDofDependencies(globalI);
377 // if (!additionalDofDependencies.empty())
378 // {
379 // const auto newNumNeighbors = neighborScvs_.size() + additionalDofDependencies.size();
380 // neighborScvs_.reserve(newNumNeighbors);
381 // neighborScvIndices_.reserve(newNumNeighbors);
382 // for (auto globalJ : additionalDofDependencies)
383 // {
384 // neighborScvs_.emplace_back(gridGeometry().element(globalJ).geometry(), globalJ);
385 // neighborScvIndices_.emplace_back(globalJ);
386 // }
387 // }
388 }
389
391 void bindElement_(const Element& element)
392 {
393 clear();
394 element_ = element;
395 makeElementGeometries(element);
396 }
397
399 template<class DataJContainer>
400 std::size_t numNeighborScvfs_(const DataJContainer& dataJContainer)
401 {
402 std::size_t numNeighborScvfs = 0;
403 for (const auto& dataJ : dataJContainer)
404 numNeighborScvfs += dataJ.scvfsJ.size() + dataJ.additionalScvfs.size();
405 return numNeighborScvfs;
406 }
407
409 void makeElementGeometries(const Element& element)
410 {
411 // make the scv
412 const auto eIdx = gridGeometry().elementMapper().index(element);
413 scvs_[0] = SubControlVolume(element.geometry(), eIdx);
414 scvIndices_[0] = eIdx;
415
416 // get data on the scv faces
417 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdx);
418 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdx);
419
420 // the quadrature point parameterizaion to be used on scvfs
421 static const auto q = getParam<CoordScalar>("MPFA.Q");
422
423 // reserve memory for the scv faces
424 const auto numLocalScvf = scvFaceIndices.size();
425 scvfIndices_.reserve(numLocalScvf);
426 scvfs_.reserve(numLocalScvf);
427
428 // for network grids we only want to do one scvf per half facet
429 // this approach assumes conforming grids at branching facets
430 std::vector<bool> finishedFacets;
431 if (dim < dimWorld)
432 finishedFacets.resize(element.subEntities(1), false);
433
434 int scvfCounter = 0;
435 for (const auto& is : intersections(gridGeometry().gridView(), element))
436 {
437 // if we are dealing with a lower dimensional network
438 // only make a new scvf if we haven't handled it yet
439 if (dim < dimWorld)
440 {
441 const auto indexInInside = is.indexInInside();
442 if (finishedFacets[indexInInside])
443 continue;
444 else
445 finishedFacets[indexInInside] = true;
446 }
447
448 // if outside level > inside level, use the outside element in the following
449 bool useNeighbor = is.neighbor() && is.outside().level() > element.level();
450 const auto& e = useNeighbor ? is.outside() : element;
451 const auto indexInElement = useNeighbor ? is.indexInOutside() : is.indexInInside();
452 const auto eg = e.geometry();
453 const auto refElement = referenceElement(eg);
454
455 // Set up a container with all relevant positions for scvf corner computation
456 const auto numCorners = is.geometry().corners();
457 const auto isPositions = MpfaHelper::computeScvfCornersOnIntersection(eg,
458 refElement,
459 indexInElement,
460 numCorners);
461
462 // make the scv faces belonging to each corner of the intersection
463 for (int c = 0; c < numCorners; ++c)
464 {
465 // get the global vertex index the scv face is connected to
466 auto vIdxLocal = refElement.subEntity(indexInElement, 1, c, dim);
467 auto vIdxGlobal = gridGeometry().vertexMapper().subIndex(e, vIdxLocal, dim);
468
469 // do not build scvfs connected to a processor boundary
470 if (gridGeometry().isGhostVertex(vIdxGlobal))
471 continue;
472
473 hasBoundaryScvf_ = (hasBoundaryScvf_ || is.boundary());
474
475 scvfs_.emplace_back(MpfaHelper(),
476 MpfaHelper::getScvfCorners(isPositions, numCorners, c),
477 is,
478 vIdxGlobal,
479 vIdxLocal,
480 scvFaceIndices[scvfCounter],
481 eIdx,
482 neighborVolVarIndices[scvfCounter],
483 q,
484 is.boundary());
485
486 scvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
487 scvfCounter++;
488 }
489 }
490 }
491
493 template<typename IndexVector>
494 void makeNeighborGeometries(const Element& element,
495 GridIndexType eIdxGlobal,
496 const IndexVector& scvfIndices,
497 const IndexVector& additionalScvfs)
498 {
499 // create the neighbor scv if it doesn't exist yet
500 neighborScvs_.emplace_back(element.geometry(), eIdxGlobal);
501 neighborScvIndices_.push_back(eIdxGlobal);
502
503 // get data on the scv faces
504 const auto& scvFaceIndices = gridGeometry().scvfIndicesOfScv(eIdxGlobal);
505 const auto& neighborVolVarIndices = gridGeometry().neighborVolVarIndices(eIdxGlobal);
506
507 // the quadrature point parameterizaion to be used on scvfs
508 static const auto q = getParam<CoordScalar>("MPFA.Q");
509
510 // for network grids we only want to do one scvf per half facet
511 // this approach assumes conforming grids at branching facets
512 std::vector<bool> finishedFacets;
513 if (dim < dimWorld)
514 finishedFacets.resize(element.subEntities(1), false);
515
516 int scvfCounter = 0;
517 for (const auto& is : intersections(gridGeometry().gridView(), element))
518 {
519 // if we are dealing with a lower dimensional network
520 // only make a new scvf if we haven't handled it yet
521 if (dim < dimWorld)
522 {
523 auto indexInInside = is.indexInInside();
524 if(finishedFacets[indexInInside])
525 continue;
526 else
527 finishedFacets[indexInInside] = true;
528 }
529
530 // if outside level > inside level, use the outside element in the following
531 bool useNeighbor = is.neighbor() && is.outside().level() > element.level();
532 const auto& e = useNeighbor ? is.outside() : element;
533 const auto indexInElement = useNeighbor ? is.indexInOutside() : is.indexInInside();
534 const auto eg = e.geometry();
535 const auto refElement = referenceElement(eg);
536
537 // Set up a container with all relevant positions for scvf corner computation
538 const auto numCorners = is.geometry().corners();
539 const auto isPositions = MpfaHelper::computeScvfCornersOnIntersection(eg,
540 refElement,
541 indexInElement,
542 numCorners);
543
544 // make the scv faces belonging to each corner of the intersection
545 for (int c = 0; c < numCorners; ++c)
546 {
547 // get the global vertex index the scv face is connected to
548 auto vIdxLocal = refElement.subEntity(indexInElement, 1, c, dim);
549 auto vIdxGlobal = gridGeometry().vertexMapper().subIndex(e, vIdxLocal, dim);
550
551 // do not build scvfs connected to a processor boundary
552 if (gridGeometry().isGhostVertex(vIdxGlobal))
553 continue;
554
555 // only build the scvf if it is in the list of necessary indices
556 if (!MpfaHelper::vectorContainsValue(scvfIndices, scvFaceIndices[scvfCounter])
557 && !MpfaHelper::vectorContainsValue(additionalScvfs, scvFaceIndices[scvfCounter]))
558 {
559 // increment counter either way
560 scvfCounter++;
561 continue;
562 }
563
564 // build scvf
565 neighborScvfs_.emplace_back(MpfaHelper(),
566 MpfaHelper::getScvfCorners(isPositions, numCorners, c),
567 is,
568 vIdxGlobal,
569 vIdxLocal,
570 scvFaceIndices[scvfCounter],
571 eIdxGlobal,
572 neighborVolVarIndices[scvfCounter],
573 q,
574 is.boundary());
575
576 neighborScvfIndices_.emplace_back(scvFaceIndices[scvfCounter]);
577
578 // increment counter
579 scvfCounter++;
580 }
581 }
582 }
583
585 unsigned int findLocalIndex(const GridIndexType idx,
586 const std::vector<GridIndexType>& indices) const
587 {
588 auto it = std::find(indices.begin(), indices.end(), idx);
589 assert(it != indices.end() && "Could not find the scv/scvf! Make sure to properly bind this class!");
590 return std::distance(indices.begin(), it);
591 }
592
594 void clear()
595 {
596 scvfIndices_.clear();
597 scvfs_.clear();
598
599 neighborScvIndices_.clear();
600 neighborScvfIndices_.clear();
601 neighborScvs_.clear();
602 neighborScvfs_.clear();
603
604 hasBoundaryScvf_ = false;
605 }
606
607 const GridGeometry* gridGeometryPtr_;
608 std::optional<Element> element_;
609
610 // local storage after binding an element
611 std::array<GridIndexType, 1> scvIndices_;
612 std::vector<GridIndexType> scvfIndices_;
613 std::array<SubControlVolume, 1> scvs_;
614 std::vector<SubControlVolumeFace> scvfs_;
615
616 std::vector<GridIndexType> neighborScvIndices_;
617 std::vector<GridIndexType> neighborScvfIndices_;
618 std::vector<SubControlVolume> neighborScvs_;
619 std::vector<SubControlVolumeFace> neighborScvfs_;
620
621 bool hasBoundaryScvf_ = false;
622};
623
624} // end namespace
625
626#endif
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.
static ctype distance(const Dune::FieldVector< ctype, dimWorld > &a, const Dune::FieldVector< ctype, dimWorld > &b)
Compute the shortest distance between two points.
Definition: distance.hh:292
Definition: adapt.hh:29
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:51
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered mpfa models Specialization fo...
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:61
const SubControlVolume & scv(GridIndexType scvIdx) const
Get an element sub control volume with a global scv index.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:87
CCMpfaFVElementGeometry 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/cellcentered/mpfa/fvelementgeometry.hh:166
void bindElement(const Element &element) &
Bind only element-local.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:173
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:140
CCMpfaFVElementGeometry(const GridGeometry &gridGeometry)
Constructor.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:83
const Element & element() const
The bound element.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:184
CCMpfaFVElementGeometry 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/cellcentered/mpfa/fvelementgeometry.hh:150
bool isBound() const
Returns true if bind/bindElement has already been called.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:180
friend Dune::IteratorRange< ScvIterator< SubControlVolume, std::array< GridIndexType, 1 >, ThisType > > scvs(const CCMpfaFVElementGeometry &fvGeometry)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:111
void bind(const Element &element) &
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:156
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:192
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:74
friend Dune::IteratorRange< ScvfIterator< SubControlVolumeFace, std::vector< GridIndexType >, ThisType > > scvfs(const CCMpfaFVElementGeometry &fvGeometry)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:124
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Get an element sub control volume face with a global scvf index.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:93
const SubControlVolumeFace & flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx=0) const
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:100
const GridGeometry & gridGeometry() const
The global finite volume geometry we are a restriction of.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:188
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:134
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:70
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:72
GG GridGeometry
export type of finite volume grid geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:76
Stencil-local finite volume geometry (scvs and scvfs) for cell-centered TPFA models Specialization fo...
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:209
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:249
const SubControlVolume & scv(GridIndexType scvIdx) const
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:239
friend Dune::IteratorRange< typename std::vector< SubControlVolumeFace >::const_iterator > scvfs(const ThisType &g)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:283
void bindElement(const Element &element) &
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:308
const GridGeometry & gridGeometry() const
The global finite volume geometry we are a restriction of.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:338
CCMpfaFVElementGeometry 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/cellcentered/mpfa/fvelementgeometry.hh:318
const SubControlVolumeFace & flipScvf(GridIndexType scvfIdx, unsigned int outsideScvIdx=0) const
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:260
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:223
typename GridView::template Codim< 0 >::Entity Element
export type of the element
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:221
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:290
typename GG::SubControlVolumeFace SubControlVolumeFace
export type of subcontrol volume face
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:225
const Element & element() const
The bound element.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:334
bool isBound() const
Returns true if bind/bindElement has already been called.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:330
void bind(const Element &element) &
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:324
CCMpfaFVElementGeometry 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/cellcentered/mpfa/fvelementgeometry.hh:302
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:294
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:342
GG GridGeometry
export type of finite volume grid geometrys
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:227
friend Dune::IteratorRange< typename std::array< SubControlVolume, 1 >::const_iterator > scvs(const ThisType &g)
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:271
CCMpfaFVElementGeometry(const GridGeometry &gridGeometry)
Constructor.
Definition: discretization/cellcentered/mpfa/fvelementgeometry.hh:234
Iterators over sub control volumes.
Definition: scvandscvfiterators.hh:42
Iterators over sub control volume faces of an fv geometry.
Definition: scvandscvfiterators.hh:82