12#ifndef DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_FV_ELEMENT_GEOMETRY_HH
13#define DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_FV_ELEMENT_GEOMETRY_HH
18#include <dune/common/rangeutilities.hh>
19#include <dune/common/reservedvector.hh>
22#include <dune/common/iteratorrange.hh>
31namespace Detail::FCStaggered {
33template<
class FVElementGeometry,
class SubControlVolume>
34typename SubControlVolume::Traits::Geometry scvGeometry(
const FVElementGeometry& fvGeometry,
35 const SubControlVolume& scv)
37 typename SubControlVolume::Traits::CornerStorage corners{};
39 const auto elementGeometry = (scv.elementIndex() != fvGeometry.elementIndex()) ?
40 fvGeometry.element().geometry() :
41 fvGeometry.gridGeometry().element(scv.elementIndex()).geometry();
43 const auto center = elementGeometry.center();
44 const auto dofAxis = scv.dofAxis();
45 for (
int i = 0; i < corners.size(); ++i)
47 auto& corner = corners[i];
50 corner = elementGeometry.corner(i);
54 if ((corner[dofAxis] -
center[dofAxis]) * scv.directionSign() < 0.0)
55 corner[dofAxis] =
center[dofAxis];
58 return {corners.front(), corners.back()};
61template<
class FVElementGeometry,
class SubControlVolumeFace>
62typename SubControlVolumeFace::Traits::Geometry scvfGeometry(
const FVElementGeometry& fvGeometry,
63 const SubControlVolumeFace& scvf)
66 const auto center = scvf.center();
67 const auto shift = scvf.ipGlobal() -
center;
69 const auto insideElementIndex = (fvGeometry.scv(scvf.insideScvIdx())).elementIndex();
70 const auto elementGeometry = (insideElementIndex != fvGeometry.elementIndex()) ?
71 fvGeometry.element().geometry() :
72 fvGeometry.gridGeometry().element(insideElementIndex).geometry();
74 auto corners = std::array{
75 elementGeometry.corner(0),
76 elementGeometry.corner(elementGeometry.corners() - 1)
80 for (
int i = 0; i < corners.size(); ++i)
82 auto& corner = corners[i];
84 if (scvf.isLateral() && (corner -
center)*shift < 0.0)
85 corner[dofAxis] = elementGeometry.center()[dofAxis];
88 auto inPlaneAxes = std::move(std::bitset<SubControlVolumeFace::Traits::dimWorld>{}.set());
91 return {corners[0], corners[1], inPlaneAxes};
97template<
class GG,
bool cachingEnabled>
109 using GridView =
typename GG::GridView;
111 static constexpr auto numScvsPerElement = GG::StaticInformation::numScvsPerElement;
117 using Element =
typename GridView::template Codim<0>::Entity;
121 static constexpr std::size_t maxNumElementScvs = numScvsPerElement;
124 : gridGeometry_(&gridGeometry)
129 {
return gridGeometry().scv(scvIdx); }
133 {
return gridGeometry().scvf(scvfIdx); }
138 assert(scvf.isLateral());
139 const auto otherGlobalIdx = scvfIndices_()[GridGeometry::GeometryHelper::lateralOrthogonalScvfLocalIndex(scvf.localIndex())];
140 return gridGeometry().scvf(otherGlobalIdx);
147 assert(scv.boundary());
150 auto scvfIter = scvfs(*
this, scv).begin();
151 const auto end = scvfs(*
this, scv).end();
152 while (!(scvfIter->isFrontal() && scvfIter->boundary()) && (scvfIter != end))
155 assert(scvfIter->isFrontal());
156 assert(scvfIter->boundary());
167 {
return fvGeometry.gridGeometry().scvs(fvGeometry); }
177 using IndexContainerType = std::decay_t<
decltype(fvGeometry.scvfIndices_())>;
179 return Dune::IteratorRange<ScvfIterator>(
ScvfIterator(fvGeometry.scvfIndices_().begin(), fvGeometry),
180 ScvfIterator(fvGeometry.scvfIndices_().end(), fvGeometry));
191 using IndexContainerType = std::decay_t<
decltype(fvGeometry.scvfIndices_())>;
193 auto begin = ScvfIterator::makeBegin(fvGeometry.scvfIndices_(), fvGeometry, scv.index());
194 auto end = ScvfIterator::makeEnd(fvGeometry.scvfIndices_(), fvGeometry, scv.index());
195 return Dune::IteratorRange<ScvfIterator>(begin, end);
200 {
return numScvsPerElement; }
204 {
return scvfIndices_().size(); }
208 {
return gridGeometry().hasBoundaryScvf(eIdx_); }
218 return std::move(*
this);
223 { this->bindElement(element); }
232 this->bindElement(element);
233 return std::move(*
this);
239 elementPtr_ = &element;
240 eIdx_ = gridGeometry().elementMapper().index(element);
246 assert(gridGeometry_);
247 return *gridGeometry_;
255 {
return *elementPtr_; }
259 {
return GG::GeometryHelper::scvfIntegrationPointInConcaveCorner(*
this, scvf); }
264 const auto& lateralOrthogonalScvf = this->lateralOrthogonalScvf(scvf);
265 assert(!lateralOrthogonalScvf.boundary());
267 const auto otherLocalIdx = GG::GeometryHelper::localIndexOutsideScvfWithSameIntegrationPoint(scvf, scv(scvf.insideScvIdx()));
269 auto outsideFVGeometry =
localView(gridGeometry());
270 const auto outsideElementIdx = scv(lateralOrthogonalScvf.outsideScvIdx()).elementIndex();
271 outsideFVGeometry.bindElement(gridGeometry().element(outsideElementIdx));
273 for (
const auto& otherScvf : scvfs(outsideFVGeometry))
274 if (otherScvf.localIndex() == otherLocalIdx)
277 DUNE_THROW(Dune::InvalidStateException,
"No outside scvf found");
283 return Detail::FCStaggered::scvGeometry(*
this, scv);
289 return Detail::FCStaggered::scvfGeometry(*
this, scvf);
294 const auto& scvfIndices_()
const
296 return gridGeometry().scvfIndicesOfElement(eIdx_);
299 const Element* elementPtr_;
301 const GridGeometry* gridGeometry_;
314 using GridView =
typename GG::GridView;
319 using StaticInfo =
typename GG::StaticInformation;
320 static constexpr auto dim = StaticInfo::dim;
321 static constexpr auto numScvsPerElement = StaticInfo::numScvsPerElement;
322 static constexpr auto numLateralScvfsPerScv = StaticInfo::numLateralScvfsPerScv;
323 static constexpr auto numLateralScvfsPerElement = StaticInfo::numLateralScvfsPerElement;
324 static constexpr auto minNumScvfsPerElement = StaticInfo::minNumScvfsPerElement;
325 static constexpr auto maxNumScvfsPerElement = StaticInfo::maxNumScvfsPerElement;
333 using Element =
typename GridView::template Codim<0>::Entity;
337 static constexpr std::size_t maxNumElementScvs = numScvsPerElement;
340 : gridGeometry_(&gridGeometry)
341 , geometryHelper_(gridGeometry.gridView())
346 {
return scvfs_[findLocalIndex_(scvfIdx, scvfIndices_())]; }
356 using Iter =
typename decltype(fvGeometry.scvfs_)::const_iterator;
357 return Dune::IteratorRange<Iter>(fvGeometry.scvfs_.begin(), fvGeometry.scvfs_.end());
368 using IndexContainerType = std::decay_t<
decltype(fvGeometry.scvfIndices_())>;
370 auto begin = ScvfIterator::makeBegin(fvGeometry.scvfIndices_(), fvGeometry, scv.index());
371 auto end = ScvfIterator::makeEnd(fvGeometry.scvfIndices_(), fvGeometry, scv.index());
372 return Dune::IteratorRange<ScvfIterator>(begin, end);
378 if (scvIdx >= scvIndicesOfElement_.front() && scvIdx <= scvIndicesOfElement_.back())
379 return scvs_[findLocalIndex_(scvIdx, scvIndicesOfElement_)];
381 return neighborScvs_[findLocalIndex_(scvIdx, neighborScvIndices_)];
392 using IteratorType =
typename std::array<SubControlVolume, 1>::const_iterator;
393 return Dune::IteratorRange<IteratorType>(g.scvs_.begin(), g.scvs_.end());
399 assert(scvf.isLateral());
400 const auto otherLocalIdx = GridGeometry::GeometryHelper::lateralOrthogonalScvfLocalIndex(scvf.localIndex());
401 return scvfs_[otherLocalIdx];
407 assert(scv.boundary());
410 auto scvfIter = scvfs(*
this, scv).begin();
411 const auto end = scvfs(*
this, scv).end();
412 while (!(scvfIter->isFrontal() && scvfIter->boundary()) && (scvfIter != end))
415 assert(scvfIter->isFrontal());
416 assert(scvfIter->boundary());
422 {
return gridGeometry().hasBoundaryScvf(eIdx_); }
426 {
return numScvsPerElement; }
430 {
return scvfIndices_().size(); }
439 this->bind_(element);
440 return std::move(*
this);
444 { this->bind_(element); }
453 typename GG::LocalIntersectionMapper localIsMapper;
454 localIsMapper.update(gridGeometry().gridView(), element);
455 this->bindElement_(element, localIsMapper);
456 return std::move(*
this);
461 typename GG::LocalIntersectionMapper localIsMapper;
462 localIsMapper.update(gridGeometry().gridView(), element);
463 this->bindElement_(element, localIsMapper);
471 {
return *elementPtr_; }
476 assert(gridGeometry_);
477 return *gridGeometry_;
482 {
return GG::GeometryHelper::scvfIntegrationPointInConcaveCorner(*
this, scvf); }
489 assert(!lateralOrthogonalScvf.boundary());
491 const auto otherLocalIdx = GG::GeometryHelper::localIndexOutsideScvfWithSameIntegrationPoint(scvf, scv(scvf.insideScvIdx()));
493 auto outsideFVGeometry =
localView(gridGeometry());
494 const auto outsideElementIdx = scv(lateralOrthogonalScvf.outsideScvIdx()).elementIndex();
495 outsideFVGeometry.bindElement(gridGeometry().element(outsideElementIdx));
497 for (
const auto& otherScvf : scvfs(outsideFVGeometry))
498 if (otherScvf.localIndex() == otherLocalIdx)
501 DUNE_THROW(Dune::InvalidStateException,
"No outside scvf found");
507 return Detail::FCStaggered::scvGeometry(*
this, scv);
513 return Detail::FCStaggered::scvfGeometry(*
this, scvf);
519 void bind_(
const Element& element)
521 typename GG::LocalIntersectionMapper localIsMapper;
522 localIsMapper.update(gridGeometry().gridView(), element);
524 bindElement_(element, localIsMapper);
525 neighborScvIndices_.clear();
526 neighborScvs_.clear();
528 for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
530 if (intersection.neighbor())
532 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
533 const auto localOppositeScvIdx = geometryHelper_.localOppositeIdx(localScvIdx);
534 const auto& neighborElement = intersection.outside();
535 const auto neighborElementIdx = gridGeometry().elementMapper().index(neighborElement);
536 const auto& neighborElementGeometry = neighborElement.geometry();
539 std::array<GridIndexType, numScvsPerElement> globalScvIndicesOfNeighborElement;
540 std::iota(globalScvIndicesOfNeighborElement.begin(), globalScvIndicesOfNeighborElement.end(), neighborElementIdx*numScvsPerElement);
542 typename GG::LocalIntersectionMapper localNeighborIsMapper;
543 localNeighborIsMapper.update(gridGeometry().gridView(), neighborElement);
545 for (
const auto& neighborIntersection : intersections(gridGeometry().gridView(), neighborElement))
547 const auto localNeighborScvIdx = localNeighborIsMapper.realToRefIdx(neighborIntersection.indexInInside());
548 if (localNeighborScvIdx != localScvIdx && localNeighborScvIdx != localOppositeScvIdx)
551 const auto dofIndex = gridGeometry().intersectionMapper().globalIntersectionIndex(neighborElement, neighborIntersection.indexInInside());
552 neighborScvs_.push_back(SubControlVolume(
553 neighborElementGeometry,
554 neighborIntersection.geometry(),
555 globalScvIndicesOfNeighborElement[localNeighborScvIdx],
560 onDomainBoundary_(neighborIntersection)
563 neighborScvIndices_.push_back(globalScvIndicesOfNeighborElement[localNeighborScvIdx]);
571 void bindElement_(
const Element& element,
const typename GG::LocalIntersectionMapper& localIsMapper)
574 eIdx_ = gridGeometry().elementMapper().index(element);
575 std::iota(scvIndicesOfElement_.begin(), scvIndicesOfElement_.end(), eIdx_*numScvsPerElement);
577 scvfs_.resize(minNumScvfsPerElement);
579 const auto& elementGeometry =
element.geometry();
581 for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
583 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
584 auto localScvfIdx = localScvIdx*(1 + numLateralScvfsPerScv);
585 const auto& intersectionGeometry = intersection.geometry();
586 const auto dofIndex = gridGeometry().intersectionMapper().globalIntersectionIndex(element, intersection.indexInInside());
588 scvs_[localScvIdx] = SubControlVolume(
590 intersectionGeometry,
591 scvIndicesOfElement_[localScvIdx],
596 onDomainBoundary_(intersection)
600 const auto localOppositeScvIdx = geometryHelper_.localOppositeIdx(localScvIdx);
601 scvfs_[localScvfIdx] = SubControlVolumeFace(
603 intersectionGeometry,
604 std::array{scvIndicesOfElement_[localScvIdx], scvIndicesOfElement_[localOppositeScvIdx]},
606 scvfIndices_()[localScvfIdx],
607 intersection.centerUnitOuterNormal(),
608 SubControlVolumeFace::FaceType::frontal,
609 SubControlVolumeFace::BoundaryType::interior
614 for (
const auto lateralFacetIndex : Dune::transformedRangeView(geometryHelper_.localLaterFaceIndices(localScvIdx),
615 [&](
auto&& idx) { return localIsMapper.refToRealIdx(idx) ;})
618 const auto& lateralIntersection = geometryHelper_.intersection(lateralFacetIndex, element);
619 const auto& lateralFacet = geometryHelper_.facet(lateralFacetIndex, element);
620 const auto& lateralFacetGeometry = lateralFacet.geometry();
623 const auto globalScvIndicesForLateralFace = [&]
625 const auto globalOutsideScvIdx = [&]
627 if (lateralIntersection.neighbor())
629 const auto parallelElemIdx = gridGeometry().elementMapper().index(lateralIntersection.outside());
631 std::array<GridIndexType, numScvsPerElement> globalScvIndicesOfNeighborElement;
632 std::iota(globalScvIndicesOfNeighborElement.begin(), globalScvIndicesOfNeighborElement.end(), parallelElemIdx*numScvsPerElement);
633 return globalScvIndicesOfNeighborElement[localScvIdx];
636 return gridGeometry().outsideVolVarIndex(scvfIndices_()[localScvfIdx]);
639 return std::array{scvIndicesOfElement_[localScvIdx], globalOutsideScvIdx};
642 const auto boundaryType = [&]
644 if (onProcessorBoundary_(lateralIntersection))
645 return SubControlVolumeFace::BoundaryType::processorBoundary;
646 else if (onDomainBoundary_(lateralIntersection))
647 return SubControlVolumeFace::BoundaryType::physicalBoundary;
649 return SubControlVolumeFace::BoundaryType::interior;
652 scvfs_[localScvfIdx] = SubControlVolumeFace(
654 intersectionGeometry,
655 lateralFacetGeometry,
656 globalScvIndicesForLateralFace,
658 scvfIndices_()[localScvfIdx],
659 lateralIntersection.centerUnitOuterNormal(),
660 SubControlVolumeFace::FaceType::lateral,
668 auto localScvfIdx = minNumScvfsPerElement;
669 for (
const auto& intersection : intersections(gridGeometry().gridView(), element))
672 if (onDomainBoundary_(intersection))
674 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
676 scvfs_.push_back(SubControlVolumeFace(
678 intersection.geometry(),
679 std::array{scvIndicesOfElement_[localScvIdx], scvIndicesOfElement_[localScvIdx]},
681 scvfIndices_()[localScvfIdx],
682 intersection.centerUnitOuterNormal(),
683 SubControlVolumeFace::FaceType::frontal,
684 SubControlVolumeFace::BoundaryType::physicalBoundary
690 if constexpr (!ConsistentlyOrientedGrid<typename GridView::Grid>{})
692 static const bool makeConsistentlyOriented = getParam<bool>(
"Grid.MakeConsistentlyOriented",
true);
693 if (makeConsistentlyOriented && scvfs_.size() > minNumScvfsPerElement)
696 std::sort(scvfs_.begin() + minNumScvfsPerElement, scvfs_.end(),
697 [](
const auto& scvfLeft,
const auto& scvfRight)
698 { return scvfLeft.insideScvIdx() < scvfRight.insideScvIdx(); }
704 const auto& scvfIndices_()
const
705 {
return gridGeometry().scvfIndicesOfElement(eIdx_); }
707 template<
class Entry,
class Container>
708 const LocalIndexType findLocalIndex_(
const Entry& entry,
709 const Container& container)
const
711 auto it = std::find(container.begin(), container.end(), entry);
712 assert(it != container.end() &&
"Could not find the entry! Make sure to properly bind this class!");
716 bool onDomainBoundary_(
const typename GridView::Intersection& intersection)
const
718 return !intersection.neighbor() && intersection.boundary();
721 bool onProcessorBoundary_(
const typename GridView::Intersection& intersection)
const
723 return !intersection.neighbor() && !intersection.boundary();
726 Dune::ReservedVector<SubControlVolumeFace, maxNumScvfsPerElement> scvfs_;
728 Dune::ReservedVector<SubControlVolume, numLateralScvfsPerElement> neighborScvs_;
729 Dune::ReservedVector<GridIndexType, numLateralScvfsPerElement> neighborScvIndices_;
731 std::array<SubControlVolume, numScvsPerElement> scvs_;
732 std::array<GridIndexType, numScvsPerElement> scvIndicesOfElement_;
734 const GridGeometry* gridGeometry_;
735 const Element* elementPtr_;
737 typename GridGeometry::GeometryHelper geometryHelper_;
Stencil-local finite volume geometry (scvs and scvfs) for face-centered staggered models Specializati...
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:311
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:429
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:421
friend auto scvfs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry)
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:354
void bindElement(const Element &element) &
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:459
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume face
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:331
FaceCenteredStaggeredFVElementGeometry 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/facecentered/staggered/fvelementgeometry.hh:437
bool scvfIntegrationPointInConcaveCorner(const SubControlVolumeFace &scvf) const
Returns true if the IP of an scvf lies on a concave corner.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:481
FaceCenteredStaggeredFVElementGeometry 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/facecentered/staggered/fvelementgeometry.hh:451
GG GridGeometry
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:334
SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace &scvf) const
Create the geometry of a given sub control volume face.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:511
GridIndexType elementIndex() const
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:466
SubControlVolumeFace outsideScvfWithSameIntegrationPoint(const SubControlVolumeFace &scvf) const
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:486
typename GridView::template Codim< 0 >::Entity Element
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:333
const SubControlVolumeFace & lateralOrthogonalScvf(const SubControlVolumeFace &scvf) const
Return a the lateral sub control volume face which is orthogonal to the given sub control volume face...
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:397
const SubControlVolume & scv(const GridIndexType scvIdx) const
Get a half sub control volume with a global scv index.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:376
const GridGeometry & gridGeometry() const
The grid geometry we are a restriction of.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:474
friend auto scvs(const FaceCenteredStaggeredFVElementGeometry &g)
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:390
SubControlVolume::Traits::Geometry geometry(const SubControlVolume &scv) const
Create the geometry of a given sub control volume.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:505
void bind(const Element &element) &
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:443
const SubControlVolumeFace & scvf(const GridIndexType scvfIdx) const
Get a sub control volume face with a global scvf index.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:345
const Element & element() const
The bound element.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:470
const SubControlVolumeFace & frontalScvfOnBoundary(const SubControlVolume &scv) const
Return the frontal sub control volume face on a the boundary for a given sub control volume.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:405
friend auto scvfs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry, const SubControlVolume &scv)
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:366
FaceCenteredStaggeredFVElementGeometry(const GridGeometry &gridGeometry)
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:339
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:425
typename GG::SubControlVolumeFace SubControlVolumeFace
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:332
Stencil-local finite volume geometry (scvs and scvfs) for face-centered staggered models Specializati...
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:107
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume face
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:115
FaceCenteredStaggeredFVElementGeometry 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/facecentered/staggered/fvelementgeometry.hh:215
FaceCenteredStaggeredFVElementGeometry 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/facecentered/staggered/fvelementgeometry.hh:230
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:199
const Element & element() const
The bound element.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:254
const GridGeometry & gridGeometry() const
The grid geometry we are a restriction of.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:244
friend auto scvfs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry)
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:175
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Get a sub control volume face with a global scvf index.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:132
SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace &scvf) const
Create the geometry of a given sub control volume face.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:287
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:207
void bind(const Element &element) &
Binding of an element, called by the local jacobian to prepare element assembly.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:222
const SubControlVolumeFace & outsideScvfWithSameIntegrationPoint(const SubControlVolumeFace &scvf) const
Returns the the scvf of neighbor element with the same integration point and unit outer normal.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:262
friend auto scvs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry)
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:166
const SubControlVolume & scv(GridIndexType scvIdx) const
Get a sub control volume with a global scv index.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:128
typename GG::SubControlVolumeFace SubControlVolumeFace
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:116
FaceCenteredStaggeredFVElementGeometry(const GridGeometry &gridGeometry)
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:123
GG GridGeometry
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:118
bool scvfIntegrationPointInConcaveCorner(const SubControlVolumeFace &scvf) const
Returns true if the IP of an scvf lies on a concave corner.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:258
SubControlVolume::Traits::Geometry geometry(const SubControlVolume &scv) const
Create the geometry of a given sub control volume.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:281
typename GridView::template Codim< 0 >::Entity Element
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:117
friend auto scvfs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry, const SubControlVolume &scv)
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:189
const SubControlVolumeFace & lateralOrthogonalScvf(const SubControlVolumeFace &scvf) const
Return a the lateral sub control volume face which is orthogonal to the given sub control volume face...
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:136
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:203
const SubControlVolumeFace & frontalScvfOnBoundary(const SubControlVolume &scv) const
Return the frontal sub control volume face on a the boundary for a given sub control volume.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:145
void bindElement(const Element &element) &
Bind only element-local.
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:237
std::size_t elementIndex() const
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:250
Definition: discretization/facecentered/staggered/fvelementgeometry.hh:98
Iterators over sub control volume faces of an fv geometry.
Definition: scvandscvfiterators.hh:70
Iterators over sub control volume faces of an fv geometry and a given sub control volume.
Definition: scvandscvfiterators.hh:110
GridCache::LocalView localView(const GridCache &gridCache)
Free function to get the local view of a grid cache object.
Definition: localview.hh:26
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
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:282
Corners::value_type center(const Corners &corners)
The center of a given list of corners.
Definition: center.hh:24
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.
typename GridView::IndexSet::IndexType GridIndex
Definition: indextraits.hh:27
unsigned int LocalIndex
Definition: indextraits.hh:28