version 3.9.0
Loading...
Searching...
No Matches
discretization/facecentered/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_FACECENTERED_STAGGERED_FV_ELEMENT_GEOMETRY_HH
13#define DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_FV_ELEMENT_GEOMETRY_HH
14
15#include <utility>
16#include <bitset>
17
18#include <dune/common/rangeutilities.hh>
19#include <dune/common/reservedvector.hh>
20#include <dune/common/iteratorrange.hh>
21#include <dune/common/exceptions.hh>
22
28
29namespace Dumux {
30
31#ifndef DOXYGEN
32namespace Detail::FCStaggered {
33
34template<class FVElementGeometry, class SubControlVolume>
35typename SubControlVolume::Traits::Geometry scvGeometry(const FVElementGeometry& fvGeometry,
36 const SubControlVolume& scv)
37{
38 typename SubControlVolume::Traits::CornerStorage corners{};
39 // select the containing element
40 const auto elementGeometry = (scv.elementIndex() != fvGeometry.elementIndex()) ?
41 fvGeometry.element().geometry() :
42 fvGeometry.gridGeometry().element(scv.elementIndex()).geometry();
43
44 const auto center = elementGeometry.center();
45 const auto dofAxis = scv.dofAxis();
46 for (int i = 0; i < corners.size(); ++i)
47 {
48 auto& corner = corners[i];
49
50 // copy the corner of the corresponding element
51 corner = elementGeometry.corner(i);
52
53 // shift the corner such that the scv covers half of the element
54 // (keep the corner positions at the face with the staggered dof)
55 if ((corner[dofAxis] - center[dofAxis]) * scv.directionSign() < 0.0)
56 corner[dofAxis] = center[dofAxis];
57 }
58
59 return {corners.front(), corners.back()};
60}
61
62template<class FVElementGeometry, class SubControlVolumeFace>
63typename SubControlVolumeFace::Traits::Geometry scvfGeometry(const FVElementGeometry& fvGeometry,
64 const SubControlVolumeFace& scvf)
65{
66 const auto normalAxis = scvf.normalAxis();
67 const auto center = scvf.center();
68 const auto shift = scvf.ipGlobal() - center;
69 const auto dofAxis = scvf.isLateral() ? Dumux::normalAxis(shift) : normalAxis;
70 const auto insideElementIndex = (fvGeometry.scv(scvf.insideScvIdx())).elementIndex();
71 const auto elementGeometry = (insideElementIndex != fvGeometry.elementIndex()) ?
72 fvGeometry.element().geometry() :
73 fvGeometry.gridGeometry().element(insideElementIndex).geometry();
74
75 auto corners = std::array{
76 elementGeometry.corner(0),
77 elementGeometry.corner(elementGeometry.corners() - 1)
78 };
79
80 // shift corners to scvf plane and halve lateral faces
81 for (int i = 0; i < corners.size(); ++i)
82 {
83 auto& corner = corners[i];
84 corner[normalAxis] = center[normalAxis];
85 if (scvf.isLateral() && (corner - center)*shift < 0.0)
86 corner[dofAxis] = elementGeometry.center()[dofAxis];
87 }
88
89 auto inPlaneAxes = std::move(std::bitset<SubControlVolumeFace::Traits::dimWorld>{}.set());
90 inPlaneAxes.set(normalAxis, false);
91
92 return {corners[0], corners[1], inPlaneAxes};
93}
94
96template<class FVElementGeometry, class SubControlVolume>
97const SubControlVolume& outsidePeriodicScv(const FVElementGeometry& fvGeometry,
98 const SubControlVolume& selfScv)
99{
100 assert(fvGeometry.gridGeometry().dofOnPeriodicBoundary(selfScv.dofIndex()));
101
102 auto localOppositeIndex = FVElementGeometry::GridGeometry::GeometryHelper::localOppositeIdx(selfScv.localDofIndex());
103 const auto& normalScvf = std::next(scvfs(fvGeometry, selfScv).begin());
104
105 // at a lateral velocity, find the inner and outer normal velocities
106 const auto& orthogonalScvf = fvGeometry.lateralOrthogonalScvf(*normalScvf);
107 const auto orthogonalOutsideScv = fvGeometry.scv(orthogonalScvf.outsideScvIdx());
108
109 auto outsidePeriodicFVGeometry = localView(fvGeometry.gridGeometry());
110 const auto& periodicElement = fvGeometry.gridGeometry().element(orthogonalOutsideScv.elementIndex());
111 outsidePeriodicFVGeometry.bindElement(periodicElement);
112
113 for (const auto& outerPeriodicScv : scvs(outsidePeriodicFVGeometry))
114 if (outerPeriodicScv.localDofIndex() == localOppositeIndex)
115 return outerPeriodicScv;
116
117 DUNE_THROW(Dune::InvalidStateException, "No outside periodic scv found");
118}
119
120} // end namespace Detail::FCStaggered
121#endif // DOXYGEN
122
123template<class GG, bool cachingEnabled>
125
131template<class GG>
132class FaceCenteredStaggeredFVElementGeometry<GG, /*cachingEnabled*/true>
133{
134 using ThisType = FaceCenteredStaggeredFVElementGeometry<GG, /*cachingEnabled*/true>;
135 using GridView = typename GG::GridView;
136 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
137 static constexpr auto numScvsPerElement = GG::StaticInformation::numScvsPerElement;
138
139public:
141 using SubControlVolume = typename GG::SubControlVolume;
142 using SubControlVolumeFace = typename GG::SubControlVolumeFace;
143 using Element = typename GridView::template Codim<0>::Entity;
144 using GridGeometry = GG;
145
147 static constexpr std::size_t maxNumElementScvs = numScvsPerElement;
148
150 : gridGeometry_(&gridGeometry)
151 {}
152
154 const SubControlVolume& scv(GridIndexType scvIdx) const
155 { return gridGeometry().scv(scvIdx); }
156
158 const SubControlVolumeFace& scvf(GridIndexType scvfIdx) const
159 { return gridGeometry().scvf(scvfIdx); }
160
163 {
164 assert(scvf.isLateral());
165 const auto otherGlobalIdx = scvfIndices_()[GridGeometry::GeometryHelper::lateralOrthogonalScvfLocalIndex(scvf.localIndex())];
166 return gridGeometry().scvf(otherGlobalIdx);
167
168 }
169
172 {
173 assert(scv.boundary());
174
175 // frontal boundary faces are always stored after the lateral faces
176 auto scvfIter = scvfs(*this, scv).begin();
177 const auto end = scvfs(*this, scv).end();
178 while (!(scvfIter->isFrontal() && scvfIter->boundary()) && (scvfIter != end))
179 ++scvfIter;
180
181 assert(scvfIter->isFrontal());
182 assert(scvfIter->boundary());
183 return *scvfIter;
184 }
185
191 friend inline auto
193 { return fvGeometry.gridGeometry().scvs(fvGeometry); }
194
200 friend inline auto
202 {
203 using IndexContainerType = std::decay_t<decltype(fvGeometry.scvfIndices_())>;
205 return Dune::IteratorRange<ScvfIterator>(ScvfIterator(fvGeometry.scvfIndices_().begin(), fvGeometry),
206 ScvfIterator(fvGeometry.scvfIndices_().end(), fvGeometry));
207 }
208
214 friend inline auto
216 {
217 using IndexContainerType = std::decay_t<decltype(fvGeometry.scvfIndices_())>;
219 auto begin = ScvfIterator::makeBegin(fvGeometry.scvfIndices_(), fvGeometry, scv.index());
220 auto end = ScvfIterator::makeEnd(fvGeometry.scvfIndices_(), fvGeometry, scv.index());
221 return Dune::IteratorRange<ScvfIterator>(begin, end);
222 }
223
225 std::size_t numScv() const
226 { return numScvsPerElement; }
227
229 std::size_t numScvf() const
230 { return scvfIndices_().size(); }
231
233 bool hasBoundaryScvf() const
234 { return gridGeometry().hasBoundaryScvf(eIdx_); }
235
242 {
243 this->bind(element);
244 return std::move(*this);
245 }
246
248 void bind(const Element& element) &
249 { this->bindElement(element); }
250
257 {
258 this->bindElement(element);
259 return std::move(*this);
260 }
261
263 void bindElement(const Element& element) &
264 {
265 elementPtr_ = &element;
266 eIdx_ = gridGeometry().elementMapper().index(element);
267 }
268
271 {
272 assert(gridGeometry_);
273 return *gridGeometry_;
274 }
275
276 std::size_t elementIndex() const
277 { return eIdx_; }
278
280 const Element& element() const
281 { return *elementPtr_; }
282
285 { return GG::GeometryHelper::scvfIntegrationPointInConcaveCorner(*this, scvf); }
286
289 {
290 const auto& lateralOrthogonalScvf = this->lateralOrthogonalScvf(scvf);
291 assert(!lateralOrthogonalScvf.boundary());
292
293 const auto otherLocalIdx = GG::GeometryHelper::localIndexOutsideScvfWithSameIntegrationPoint(scvf, scv(scvf.insideScvIdx()));
294
295 auto outsideFVGeometry = localView(gridGeometry());
296 const auto outsideElementIdx = scv(lateralOrthogonalScvf.outsideScvIdx()).elementIndex();
297 outsideFVGeometry.bindElement(gridGeometry().element(outsideElementIdx));
298
299 for (const auto& otherScvf : scvfs(outsideFVGeometry))
300 if (otherScvf.localIndex() == otherLocalIdx)
301 return otherScvf;
302
303 DUNE_THROW(Dune::InvalidStateException, "No outside scvf found");
304 }
305
308 {
309 return Detail::FCStaggered::outsidePeriodicScv(*this, selfScv);
310 }
311
313 typename SubControlVolume::Traits::Geometry geometry(const SubControlVolume& scv) const
314 {
315 return Detail::FCStaggered::scvGeometry(*this, scv);
316 }
317
319 typename SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace& scvf) const
320 {
321 return Detail::FCStaggered::scvfGeometry(*this, scvf);
322 }
323
324private:
325
326 const auto& scvfIndices_() const
327 {
328 return gridGeometry().scvfIndicesOfElement(eIdx_);
329 }
330
331 const Element* elementPtr_;
332 GridIndexType eIdx_;
333 const GridGeometry* gridGeometry_;
334};
335
341template<class GG>
342class FaceCenteredStaggeredFVElementGeometry<GG, /*cachingEnabled*/false>
343{
344 using ThisType = FaceCenteredStaggeredFVElementGeometry<GG, /*cachingEnabled*/false>;
345
346 using GridView = typename GG::GridView;
347 using GridIndexType = typename IndexTraits<GridView>::GridIndex;
348
349 //TODO include assert that checks for quad geometry
350
351 using StaticInfo = typename GG::StaticInformation;
352 static constexpr auto dim = StaticInfo::dim;
353 static constexpr auto numScvsPerElement = StaticInfo::numScvsPerElement;
354 static constexpr auto numLateralScvfsPerScv = StaticInfo::numLateralScvfsPerScv;
355 static constexpr auto numLateralScvfsPerElement = StaticInfo::numLateralScvfsPerElement;
356 static constexpr auto minNumScvfsPerElement = StaticInfo::minNumScvfsPerElement;
357 static constexpr auto maxNumScvfsPerElement = StaticInfo::maxNumScvfsPerElement;
358
359 using LocalIndexType = typename IndexTraits<GridView>::LocalIndex;
360
361public:
363 using SubControlVolume = typename GG::SubControlVolume;
364 using SubControlVolumeFace = typename GG::SubControlVolumeFace;
365 using Element = typename GridView::template Codim<0>::Entity;
366 using GridGeometry = GG;
367
369 static constexpr std::size_t maxNumElementScvs = numScvsPerElement;
370
372 : gridGeometry_(&gridGeometry)
373 , geometryHelper_(gridGeometry.gridView())
374 {}
375
377 const SubControlVolumeFace& scvf(const GridIndexType scvfIdx) const
378 { return scvfs_[findLocalIndex_(scvfIdx, scvfIndices_())]; }
379
385 friend inline auto
387 {
388 using Iter = typename decltype(fvGeometry.scvfs_)::const_iterator;
389 return Dune::IteratorRange<Iter>(fvGeometry.scvfs_.begin(), fvGeometry.scvfs_.end());
390 }
391
397 friend inline auto
399 {
400 using IndexContainerType = std::decay_t<decltype(fvGeometry.scvfIndices_())>;
402 auto begin = ScvfIterator::makeBegin(fvGeometry.scvfIndices_(), fvGeometry, scv.index());
403 auto end = ScvfIterator::makeEnd(fvGeometry.scvfIndices_(), fvGeometry, scv.index());
404 return Dune::IteratorRange<ScvfIterator>(begin, end);
405 }
406
408 const SubControlVolume& scv(const GridIndexType scvIdx) const
409 {
410 if (scvIdx >= scvIndicesOfElement_.front() && scvIdx <= scvIndicesOfElement_.back())
411 return scvs_[findLocalIndex_(scvIdx, scvIndicesOfElement_)];
412 else
413 return neighborScvs_[findLocalIndex_(scvIdx, neighborScvIndices_)];
414 }
415
421 friend inline auto
423 {
424 using IteratorType = typename std::array<SubControlVolume, 1>::const_iterator;
425 return Dune::IteratorRange<IteratorType>(g.scvs_.begin(), g.scvs_.end());
426 }
427
430 {
431 assert(scvf.isLateral());
432 const auto otherLocalIdx = GridGeometry::GeometryHelper::lateralOrthogonalScvfLocalIndex(scvf.localIndex());
433 return scvfs_[otherLocalIdx];
434 }
435
438 {
439 assert(scv.boundary());
440
441 // frontal boundary faces are always stored after the lateral faces
442 auto scvfIter = scvfs(*this, scv).begin();
443 const auto end = scvfs(*this, scv).end();
444 while (!(scvfIter->isFrontal() && scvfIter->boundary()) && (scvfIter != end))
445 ++scvfIter;
446
447 assert(scvfIter->isFrontal());
448 assert(scvfIter->boundary());
449 return *scvfIter;
450 }
451
453 bool hasBoundaryScvf() const
454 { return gridGeometry().hasBoundaryScvf(eIdx_); }
455
457 std::size_t numScv() const
458 { return numScvsPerElement; }
459
461 std::size_t numScvf() const
462 { return scvfIndices_().size(); }
463
470 {
471 this->bind_(element);
472 return std::move(*this);
473 }
474
475 void bind(const Element& element) &
476 { this->bind_(element); }
477
484 {
485 typename GG::LocalIntersectionMapper localIsMapper;
486 localIsMapper.update(gridGeometry().gridView(), element);
487 this->bindElement_(element, localIsMapper);
488 return std::move(*this);
489 }
490
491 void bindElement(const Element& element) &
492 {
493 typename GG::LocalIntersectionMapper localIsMapper;
494 localIsMapper.update(gridGeometry().gridView(), element);
495 this->bindElement_(element, localIsMapper);
496 }
497
498 GridIndexType elementIndex() const
499 { return eIdx_; }
500
502 const Element& element() const
503 { return *elementPtr_; }
504
507 {
508 assert(gridGeometry_);
509 return *gridGeometry_;
510 }
511
514 { return GG::GeometryHelper::scvfIntegrationPointInConcaveCorner(*this, scvf); }
515
519 {
520 const SubControlVolumeFace& lateralOrthogonalScvf = this->lateralOrthogonalScvf(scvf);
521 assert(!lateralOrthogonalScvf.boundary());
522
523 const auto otherLocalIdx = GG::GeometryHelper::localIndexOutsideScvfWithSameIntegrationPoint(scvf, scv(scvf.insideScvIdx()));
524
525 auto outsideFVGeometry = localView(gridGeometry());
526 const auto outsideElementIdx = scv(lateralOrthogonalScvf.outsideScvIdx()).elementIndex();
527 outsideFVGeometry.bindElement(gridGeometry().element(outsideElementIdx));
528
529 for (const auto& otherScvf : scvfs(outsideFVGeometry))
530 if (otherScvf.localIndex() == otherLocalIdx)
531 return otherScvf;
532
533 DUNE_THROW(Dune::InvalidStateException, "No outside scvf found");
534 }
535
538 {
539 return Detail::FCStaggered::outsidePeriodicScv(*this, selfScv);
540 }
541
543 typename SubControlVolume::Traits::Geometry geometry(const SubControlVolume& scv) const
544 {
545 return Detail::FCStaggered::scvGeometry(*this, scv);
546 }
547
549 typename SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace& scvf) const
550 {
551 return Detail::FCStaggered::scvfGeometry(*this, scvf);
552 }
553
554private:
557 void bind_(const Element& element)
558 {
559 typename GG::LocalIntersectionMapper localIsMapper;
560 localIsMapper.update(gridGeometry().gridView(), element);
561
562 bindElement_(element, localIsMapper);
563 neighborScvIndices_.clear();
564 neighborScvs_.clear();
565
566 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
567 {
568 if (intersection.neighbor())
569 {
570 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
571 const auto localOppositeScvIdx = geometryHelper_.localOppositeIdx(localScvIdx);
572 const auto& neighborElement = intersection.outside();
573 const auto neighborElementIdx = gridGeometry().elementMapper().index(neighborElement);
574 const auto& neighborElementGeometry = neighborElement.geometry();
575
576 // todo: could be done easier?
577 std::array<GridIndexType, numScvsPerElement> globalScvIndicesOfNeighborElement;
578 std::iota(globalScvIndicesOfNeighborElement.begin(), globalScvIndicesOfNeighborElement.end(), neighborElementIdx*numScvsPerElement);
579
580 typename GG::LocalIntersectionMapper localNeighborIsMapper;
581 localNeighborIsMapper.update(gridGeometry().gridView(), neighborElement);
582
583 for (const auto& neighborIntersection : intersections(gridGeometry().gridView(), neighborElement))
584 {
585 const auto localNeighborScvIdx = localNeighborIsMapper.realToRefIdx(neighborIntersection.indexInInside());
586 if (localNeighborScvIdx != localScvIdx && localNeighborScvIdx != localOppositeScvIdx)
587 {
588
589 const auto dofIndex = gridGeometry().intersectionMapper().globalIntersectionIndex(neighborElement, neighborIntersection.indexInInside());
590 neighborScvs_.push_back(SubControlVolume(
591 neighborElementGeometry,
592 neighborIntersection.geometry(),
593 globalScvIndicesOfNeighborElement[localNeighborScvIdx],
594 localNeighborScvIdx,
595 dofIndex,
596 Dumux::normalAxis(neighborIntersection.centerUnitOuterNormal()),
597 neighborElementIdx,
598 onDomainBoundary_(neighborIntersection)
599 ));
600
601 neighborScvIndices_.push_back(globalScvIndicesOfNeighborElement[localNeighborScvIdx]);
602 }
603 }
604 }
605 }
606 }
607
609 void bindElement_(const Element& element, const typename GG::LocalIntersectionMapper& localIsMapper)
610 {
611 elementPtr_ = &element;
612 eIdx_ = gridGeometry().elementMapper().index(element);
613 std::iota(scvIndicesOfElement_.begin(), scvIndicesOfElement_.end(), eIdx_*numScvsPerElement);
614 scvfs_.clear();
615 scvfs_.resize(minNumScvfsPerElement);
616
617 const auto& elementGeometry = element.geometry();
618
619 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
620 {
621 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
622 auto localScvfIdx = localScvIdx*(1 + numLateralScvfsPerScv);
623 const auto& intersectionGeometry = intersection.geometry();
624 const auto dofIndex = gridGeometry().intersectionMapper().globalIntersectionIndex(element, intersection.indexInInside());
625
626 scvs_[localScvIdx] = SubControlVolume(
627 elementGeometry,
628 intersectionGeometry,
629 scvIndicesOfElement_[localScvIdx],
630 localScvIdx,
631 dofIndex,
632 Dumux::normalAxis(intersection.centerUnitOuterNormal()),
633 eIdx_,
634 onDomainBoundary_(intersection)
635 );
636
637 // the frontal sub control volume face at the element center
638 const auto localOppositeScvIdx = geometryHelper_.localOppositeIdx(localScvIdx);
639 scvfs_[localScvfIdx] = SubControlVolumeFace(
640 elementGeometry,
641 intersectionGeometry,
642 std::array{scvIndicesOfElement_[localScvIdx], scvIndicesOfElement_[localOppositeScvIdx]},
643 localScvfIdx,
644 scvfIndices_()[localScvfIdx],
645 intersection.centerUnitOuterNormal(),
646 SubControlVolumeFace::FaceType::frontal,
647 SubControlVolumeFace::BoundaryType::interior
648 );
649 ++localScvfIdx;
650
651 // the lateral sub control volume faces
652 for (const auto lateralFacetIndex : Dune::transformedRangeView(geometryHelper_.localLaterFaceIndices(localScvIdx),
653 [&](auto&& idx) { return localIsMapper.refToRealIdx(idx) ;})
654 )
655 {
656 const auto& lateralIntersection = geometryHelper_.intersection(lateralFacetIndex, element);
657 const auto& lateralFacet = geometryHelper_.facet(lateralFacetIndex, element);
658 const auto& lateralFacetGeometry = lateralFacet.geometry();
659
660 // helper lambda to get the lateral scvf's global inside and outside scv indices
661 const auto globalScvIndicesForLateralFace = [&]
662 {
663 const auto globalOutsideScvIdx = [&]
664 {
665 if (lateralIntersection.neighbor())
666 {
667 const auto parallelElemIdx = gridGeometry().elementMapper().index(lateralIntersection.outside());
668 // todo: could be done easier?
669 std::array<GridIndexType, numScvsPerElement> globalScvIndicesOfNeighborElement;
670 std::iota(globalScvIndicesOfNeighborElement.begin(), globalScvIndicesOfNeighborElement.end(), parallelElemIdx*numScvsPerElement);
671 return globalScvIndicesOfNeighborElement[localScvIdx];
672 }
673 else
674 return gridGeometry().outsideVolVarIndex(scvfIndices_()[localScvfIdx]);
675 }();
676
677 return std::array{scvIndicesOfElement_[localScvIdx], globalOutsideScvIdx};
678 }();
679
680 const auto boundaryType = [&]
681 {
682 if (onProcessorBoundary_(lateralIntersection))
683 return SubControlVolumeFace::BoundaryType::processorBoundary;
684 else if (onDomainBoundary_(lateralIntersection))
685 return SubControlVolumeFace::BoundaryType::physicalBoundary;
686 else
687 return SubControlVolumeFace::BoundaryType::interior;
688 }();
689
690 scvfs_[localScvfIdx] = SubControlVolumeFace(
691 elementGeometry,
692 intersectionGeometry,
693 lateralFacetGeometry,
694 globalScvIndicesForLateralFace, // TODO higher order
695 localScvfIdx,
696 scvfIndices_()[localScvfIdx],
697 lateralIntersection.centerUnitOuterNormal(),
698 SubControlVolumeFace::FaceType::lateral,
699 boundaryType
700 );
701 ++localScvfIdx;
702 }
703 }
704
705 // do a second loop over all intersections to add frontal boundary faces
706 auto localScvfIdx = minNumScvfsPerElement;
707 for (const auto& intersection : intersections(gridGeometry().gridView(), element))
708 {
709 // the frontal sub control volume face at a domain boundary (coincides with element face)
710 if (onDomainBoundary_(intersection))
711 {
712 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
713 // the frontal sub control volume face at the boundary
714 scvfs_.push_back(SubControlVolumeFace(
715 element.geometry(),
716 intersection.geometry(),
717 std::array{scvIndicesOfElement_[localScvIdx], scvIndicesOfElement_[localScvIdx]},
718 localScvfIdx,
719 scvfIndices_()[localScvfIdx],
720 intersection.centerUnitOuterNormal(),
721 SubControlVolumeFace::FaceType::frontal,
722 SubControlVolumeFace::BoundaryType::physicalBoundary
723 ));
724 ++localScvfIdx;
725 }
726 }
727
728 if constexpr (!ConsistentlyOrientedGrid<typename GridView::Grid>{})
729 {
730 static const bool makeConsistentlyOriented = getParam<bool>("Grid.MakeConsistentlyOriented", true);
731 if (makeConsistentlyOriented && scvfs_.size() > minNumScvfsPerElement)
732 {
733 // make sure frontal boundary scvfs are sorted correctly
734 std::sort(scvfs_.begin() + minNumScvfsPerElement, scvfs_.end(),
735 [](const auto& scvfLeft, const auto& scvfRight)
736 { return scvfLeft.insideScvIdx() < scvfRight.insideScvIdx(); }
737 );
738 }
739 }
740 }
741
742 const auto& scvfIndices_() const
743 { return gridGeometry().scvfIndicesOfElement(eIdx_); }
744
745 template<class Entry, class Container>
746 const LocalIndexType findLocalIndex_(const Entry& entry,
747 const Container& container) const
748 {
749 auto it = std::find(container.begin(), container.end(), entry);
750 assert(it != container.end() && "Could not find the entry! Make sure to properly bind this class!");
751 return std::distance(container.begin(), it);
752 }
753
754 bool onDomainBoundary_(const typename GridView::Intersection& intersection) const
755 {
756 return !intersection.neighbor() && intersection.boundary();
757 }
758
759 bool onProcessorBoundary_(const typename GridView::Intersection& intersection) const
760 {
761 return !intersection.neighbor() && !intersection.boundary();
762 }
763
764 Dune::ReservedVector<SubControlVolumeFace, maxNumScvfsPerElement> scvfs_;
765
766 Dune::ReservedVector<SubControlVolume, numLateralScvfsPerElement> neighborScvs_;
767 Dune::ReservedVector<GridIndexType, numLateralScvfsPerElement> neighborScvIndices_;
768
769 std::array<SubControlVolume, numScvsPerElement> scvs_;
770 std::array<GridIndexType, numScvsPerElement> scvIndicesOfElement_;
771
772 const GridGeometry* gridGeometry_;
773 const Element* elementPtr_;
774 GridIndexType eIdx_;
775 typename GridGeometry::GeometryHelper geometryHelper_;
776};
777
778} // end namespace Dumux
779
780#endif
Stencil-local finite volume geometry (scvs and scvfs) for face-centered staggered models Specializati...
Definition discretization/facecentered/staggered/fvelementgeometry.hh:343
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition discretization/facecentered/staggered/fvelementgeometry.hh:461
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:453
friend auto scvfs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry)
Definition discretization/facecentered/staggered/fvelementgeometry.hh:386
void bindElement(const Element &element) &
Definition discretization/facecentered/staggered/fvelementgeometry.hh:491
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume face
Definition discretization/facecentered/staggered/fvelementgeometry.hh:363
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:469
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:513
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:483
GG GridGeometry
Definition discretization/facecentered/staggered/fvelementgeometry.hh:366
SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace &scvf) const
Create the geometry of a given sub control volume face.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:549
GridIndexType elementIndex() const
Definition discretization/facecentered/staggered/fvelementgeometry.hh:498
SubControlVolumeFace outsideScvfWithSameIntegrationPoint(const SubControlVolumeFace &scvf) const
Definition discretization/facecentered/staggered/fvelementgeometry.hh:518
typename GridView::template Codim< 0 >::Entity Element
Definition discretization/facecentered/staggered/fvelementgeometry.hh:365
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:429
const SubControlVolume & scv(const GridIndexType scvIdx) const
Get a half sub control volume with a global scv index.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:408
const GridGeometry & gridGeometry() const
The grid geometry we are a restriction of.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:506
friend auto scvs(const FaceCenteredStaggeredFVElementGeometry &g)
Definition discretization/facecentered/staggered/fvelementgeometry.hh:422
SubControlVolume::Traits::Geometry geometry(const SubControlVolume &scv) const
Create the geometry of a given sub control volume.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:543
void bind(const Element &element) &
Definition discretization/facecentered/staggered/fvelementgeometry.hh:475
const SubControlVolumeFace & scvf(const GridIndexType scvfIdx) const
Get a sub control volume face with a global scvf index.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:377
const Element & element() const
The bound element.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:502
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:437
friend auto scvfs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry, const SubControlVolume &scv)
Definition discretization/facecentered/staggered/fvelementgeometry.hh:398
FaceCenteredStaggeredFVElementGeometry(const GridGeometry &gridGeometry)
Definition discretization/facecentered/staggered/fvelementgeometry.hh:371
SubControlVolume outsidePeriodicScv(const SubControlVolume &selfScv) const
Get the scv on the outside side of a periodic boundary.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:537
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition discretization/facecentered/staggered/fvelementgeometry.hh:457
typename GG::SubControlVolumeFace SubControlVolumeFace
Definition discretization/facecentered/staggered/fvelementgeometry.hh:364
Stencil-local finite volume geometry (scvs and scvfs) for face-centered staggered models Specializati...
Definition discretization/facecentered/staggered/fvelementgeometry.hh:133
typename GG::SubControlVolume SubControlVolume
export type of subcontrol volume face
Definition discretization/facecentered/staggered/fvelementgeometry.hh:141
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:241
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:256
std::size_t numScv() const
number of sub control volumes in this fv element geometry
Definition discretization/facecentered/staggered/fvelementgeometry.hh:225
const Element & element() const
The bound element.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:280
const GridGeometry & gridGeometry() const
The grid geometry we are a restriction of.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:270
friend auto scvfs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry)
Definition discretization/facecentered/staggered/fvelementgeometry.hh:201
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Get a sub control volume face with a global scvf index.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:158
SubControlVolumeFace::Traits::Geometry geometry(const SubControlVolumeFace &scvf) const
Create the geometry of a given sub control volume face.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:319
bool hasBoundaryScvf() const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:233
void bind(const Element &element) &
Binding of an element, called by the local jacobian to prepare element assembly.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:248
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:288
friend auto scvs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry)
Definition discretization/facecentered/staggered/fvelementgeometry.hh:192
const SubControlVolume & scv(GridIndexType scvIdx) const
Get a sub control volume with a global scv index.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:154
typename GG::SubControlVolumeFace SubControlVolumeFace
Definition discretization/facecentered/staggered/fvelementgeometry.hh:142
FaceCenteredStaggeredFVElementGeometry(const GridGeometry &gridGeometry)
Definition discretization/facecentered/staggered/fvelementgeometry.hh:149
GG GridGeometry
Definition discretization/facecentered/staggered/fvelementgeometry.hh:144
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:284
SubControlVolume::Traits::Geometry geometry(const SubControlVolume &scv) const
Create the geometry of a given sub control volume.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:313
typename GridView::template Codim< 0 >::Entity Element
Definition discretization/facecentered/staggered/fvelementgeometry.hh:143
friend auto scvfs(const FaceCenteredStaggeredFVElementGeometry &fvGeometry, const SubControlVolume &scv)
Definition discretization/facecentered/staggered/fvelementgeometry.hh:215
const SubControlVolume & outsidePeriodicScv(const SubControlVolume &selfScv) const
Get the scv on the outside side of a periodic boundary.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:307
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:162
std::size_t numScvf() const
number of sub control volumes in this fv element geometry
Definition discretization/facecentered/staggered/fvelementgeometry.hh:229
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:171
void bindElement(const Element &element) &
Bind only element-local.
Definition discretization/facecentered/staggered/fvelementgeometry.hh:263
std::size_t elementIndex() const
Definition discretization/facecentered/staggered/fvelementgeometry.hh:276
Definition discretization/facecentered/staggered/fvelementgeometry.hh:124
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
Helper type to determine whether a grid is guaranteed to be oriented consistently....
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
Corners::value_type center(const Corners &corners)
The center of a given list of corners.
Definition center.hh:24
T getParam(Args &&... args)
A free function to get a parameter from the parameter tree singleton.
Definition parameters.hh:139
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...
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