version 3.9.0
Loading...
Searching...
No Matches
discretization/facecentered/staggered/fvgridgeometry.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_GRID_GEOMETRY
13#define DUMUX_DISCRETIZATION_FACECENTERED_STAGGERED_FV_GRID_GEOMETRY
14
15#include <memory>
16
17#include <dune/common/rangeutilities.hh>
18#include <dune/grid/common/scsgmapper.hh>
19
23#include <dumux/common/math.hh>
24
29
37
38namespace Dumux {
39
46template<class GridView>
74
81template<class GridView,
82 bool cachingEnabled = false,
85
92template<class GV, class Traits>
94: public BaseGridGeometry<GV, Traits>
95{
98 using GridIndexType = typename IndexTraits<GV>::GridIndex;
99 using LocalIndexType = typename IndexTraits<GV>::LocalIndex;
100 using SmallLocalIndexType = typename IndexTraits<GV>::SmallLocalIndex;
101 using Element = typename GV::template Codim<0>::Entity;
102
103 using IntersectionMapper = typename Traits::IntersectionMapper;
104 using ConnectivityMap = typename Traits::template ConnectivityMap<ThisType>;
105
106 using Scalar = typename GV::ctype;
107
108 static constexpr auto dim = Traits::StaticInfo::dim;
109 static constexpr auto numScvsPerElement = Traits::StaticInfo::numScvsPerElement;
110 static constexpr auto numLateralScvfsPerScv = Traits::StaticInfo::numLateralScvfsPerScv;
111 static constexpr auto numLateralScvfsPerElement = Traits::StaticInfo::numLateralScvfsPerElement;
112 static constexpr auto minNumScvfsPerElement = Traits::StaticInfo::minNumScvfsPerElement;
113 static constexpr auto maxNumScvfsPerElement = Traits::StaticInfo::maxNumScvfsPerElement;
114
115 using ScvfCornerStorage = typename Traits::SubControlVolumeFace::Traits::CornerStorage;
116 using ScvCornerStorage = typename Traits::SubControlVolume::Traits::CornerStorage;
117
118public:
121 static constexpr DiscretizationMethod discMethod{};
122
123 static constexpr bool cachingEnabled = true;
124
128 using LocalView = typename Traits::template LocalView<ThisType, true>;
130 using SubControlVolume = typename Traits::SubControlVolume;
132 using SubControlVolumeFace = typename Traits::SubControlVolumeFace;
134 using GridView = GV;
136 using GeometryHelper = typename Traits::GeometryHelper;
138 using LocalIntersectionMapper = typename Traits::LocalIntersectionMapper;
140 using StaticInformation = typename Traits::StaticInfo;
143
145 FaceCenteredStaggeredFVGridGeometry(std::shared_ptr<BasicGridGeometry> gg, const std::string& paramGroup = "")
146 : ParentType(std::move(gg))
147 , intersectionMapper_(this->gridView())
148 {
149 // Check if the overlap size is what we expect
151 DUNE_THROW(Dune::InvalidStateException, "The staggered discretization method needs overlap of exactly 1 for parallel computations. "
152 << " Set the parameter \"Grid.Overlap\" in the input file.");
153
154 update_();
155 }
156
158 FaceCenteredStaggeredFVGridGeometry(const GridView& gridView, const std::string& paramGroup = "")
159 : FaceCenteredStaggeredFVGridGeometry(std::make_shared<BasicGridGeometry>(gridView), paramGroup)
160 {}
161
163 std::size_t numScv() const
164 { return scvs_.size(); }
165
167 std::size_t numScvf() const
168 { return scvfs_.size(); }
169
171 std::size_t numBoundaryScv() const
172 { return numBoundaryScv_; }
173
175 std::size_t numBoundaryScvf() const
176 { return numBoundaryScvf_; }
177
179 std::size_t numIntersections() const
180 { return intersectionMapper_.numIntersections(); }
181
183 std::size_t numDofs() const
184 { return this->gridView().size(1); }
185
187 void update(const GridView& gridView)
188 {
189 ParentType::update(gridView);
190 update_();
191 }
192
194 void update(GridView&& gridView)
195 {
196 ParentType::update(std::move(gridView));
197 update_();
198 }
199
201 const SubControlVolume& scv(GridIndexType scvIdx) const
202 { return scvs_[scvIdx]; }
203
206 auto scvs(const LocalView& fvGeometry) const
207 {
208 auto begin = scvs_.cbegin() + numScvsPerElement*fvGeometry.elementIndex();
209 const auto end = begin + numScvsPerElement;
210 return Dune::IteratorRange<std::decay_t<decltype(begin)>>(begin, end);
211 }
212
214 const SubControlVolumeFace& scvf(GridIndexType scvfIdx) const
215 { return scvfs_[scvfIdx]; }
216
218 const std::vector<GridIndexType>& scvfIndicesOfElement(GridIndexType eIdx) const
219 { return scvfIndicesOfElement_[eIdx]; }
220
225 const ConnectivityMap& connectivityMap() const
226 { return connectivityMap_; }
227
229 bool hasBoundaryScvf(GridIndexType eIdx) const
230 { return hasBoundaryScvf_[eIdx]; }
231
233 const IntersectionMapper& intersectionMapper() const
234 { return intersectionMapper_; }
235
237 bool dofOnPeriodicBoundary(GridIndexType dofIdx) const
238 { return periodicFaceMap_.count(dofIdx); }
239
241 GridIndexType periodicallyMappedDof(GridIndexType dofIdx) const
242 { return periodicFaceMap_.at(dofIdx); }
243
245 const std::unordered_map<GridIndexType, GridIndexType>& periodicDofMap() const
246 { return periodicFaceMap_; }
247
249 [[deprecated("Will be removed after release 3.9. Use periodicDofMap() instead.")]]
250 const std::unordered_map<GridIndexType, GridIndexType>& periodicVertexMap() const
251 { return periodicDofMap(); }
252
253private:
254
255 void update_()
256 {
257 // clear containers (necessary after grid refinement)
258 scvs_.clear();
259 scvfs_.clear();
260 scvfIndicesOfElement_.clear();
261 intersectionMapper_.update(this->gridView());
262
263 // determine size of containers
264 const auto numElements = this->gridView().size(0);
265 scvfIndicesOfElement_.resize(numElements);
266 hasBoundaryScvf_.resize(numElements, false);
267
268 outSideBoundaryVolVarIdx_ = 0;
269 numBoundaryScv_ = 0;
270 numBoundaryScvf_ = 0;
271
272 GeometryHelper geometryHelper(this->gridView());
273
274 // get the global scvf indices first
275 GridIndexType numScvfs = 0;
276 for (const auto& element : elements(this->gridView()))
277 {
278 assert(numScvsPerElement == element.subEntities(1));
279
280 for (const auto& intersection : intersections(this->gridView(), element))
281 {
282 // frontal scvf in element center
283 ++numScvfs;
284
285 // lateral scvfs
286 numScvfs += numLateralScvfsPerScv;
287
288 // handle physical domain boundary
289 if (onDomainBoundary_(intersection))
290 {
291 ++numBoundaryScv_; // frontal face
292 numBoundaryScv_ += numLateralScvfsPerScv; // boundary scvs for lateral faces
293
294 // frontal scvf at boundary
295 ++numScvfs;
296 }
297 }
298 }
299
300 // allocate memory
301 const auto numScvs = numElements*numScvsPerElement;
302 scvs_.resize(numScvs);
303 scvfs_.reserve(numScvfs);
304
305 // Build the scvs and scv faces
306 std::size_t globalScvfIdx = 0;
307 for (const auto& element : elements(this->gridView()))
308 {
309 const auto eIdx = this->elementMapper().index(element);
310 auto& globalScvfIndices = scvfIndicesOfElement_[eIdx];
311 globalScvfIndices.resize(minNumScvfsPerElement);
312 globalScvfIndices.reserve(maxNumScvfsPerElement);
313
314 auto getGlobalScvIdx = [&](const auto elementIdx, const auto localScvIdx)
315 { return numScvsPerElement*elementIdx + localScvIdx; };
316
317 LocalIntersectionMapper localIsMapper;
318 localIsMapper.update(this->gridView(), element);
319
320 for (const auto& intersection : intersections(this->gridView(), element))
321 {
322 const auto& intersectionUnitOuterNormal = intersection.centerUnitOuterNormal();
323 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
324 auto localScvfIdx = localScvIdx*(1 + numLateralScvfsPerScv);
325
326 const auto globalScvIdx = getGlobalScvIdx(eIdx, localScvIdx);
327 const auto dofIndex = intersectionMapper().globalIntersectionIndex(element, intersection.indexInInside());
328 const auto localOppositeScvIdx = geometryHelper.localOppositeIdx(localScvIdx);
329 const auto& intersectionGeometry = intersection.geometry();
330 const auto& elementGeometry = element.geometry();
331
332 assert(localIsMapper.refToRealIdx(localScvIdx) == intersection.indexInInside());
333
334 // handle periodic boundaries
335 if (onPeriodicBoundary_(intersection))
336 {
337 this->setPeriodic();
338
339 const auto& otherElement = intersection.outside();
340
341 SmallLocalIndexType otherIntersectionLocalIdx = 0;
342 bool periodicFaceFound = false;
343
344 for (const auto& otherIntersection : intersections(this->gridView(), otherElement))
345 {
346 if (periodicFaceFound)
347 continue;
348
349 if (Dune::FloatCmp::eq(intersectionUnitOuterNormal*otherIntersection.centerUnitOuterNormal(), -1.0, 1e-7))
350 {
351 const auto periodicDofIdx = intersectionMapper().globalIntersectionIndex(otherElement, otherIntersectionLocalIdx);
352 periodicFaceMap_[dofIndex] = periodicDofIdx;
353 periodicFaceFound = true;
354 }
355
356 ++otherIntersectionLocalIdx;
357 }
358 }
359
360 // the sub control volume
361 scvs_[globalScvIdx] = SubControlVolume(
362 elementGeometry,
363 intersectionGeometry,
364 globalScvIdx,
365 localScvIdx,
366 dofIndex,
367 Dumux::normalAxis(intersectionUnitOuterNormal),
368 this->elementMapper().index(element),
369 onDomainBoundary_(intersection)
370 );
371
372 // the frontal sub control volume face at the element center
373 scvfs_.emplace_back(elementGeometry,
374 intersectionGeometry,
375 std::array{globalScvIdx, getGlobalScvIdx(eIdx, localOppositeScvIdx)},
376 localScvfIdx,
377 globalScvfIdx,
378 intersectionUnitOuterNormal,
379 SubControlVolumeFace::FaceType::frontal,
380 SubControlVolumeFace::BoundaryType::interior
381 );
382
383 globalScvfIndices[localScvfIdx] = globalScvfIdx++;
384 ++localScvfIdx;
385
386 // the lateral sub control volume faces
387 for (const auto lateralFacetIndex : Dune::transformedRangeView(geometryHelper.localLaterFaceIndices(localScvIdx),
388 [&](auto&& idx) { return localIsMapper.refToRealIdx(idx) ;})
389 )
390 {
391 const auto& lateralIntersection = geometryHelper.intersection(lateralFacetIndex, element);
392
393 // helper lambda to get the lateral scvf's global inside and outside scv indices
394 const auto globalScvIndicesForLateralFace = [&]
395 {
396 const auto globalOutsideScvIdx = [&]
397 {
398 if (lateralIntersection.neighbor())
399 {
400 const auto parallelElemIdx = this->elementMapper().index(lateralIntersection.outside());
401 return getGlobalScvIdx(parallelElemIdx, localScvIdx);
402 }
403 else if (onDomainBoundary_(lateralIntersection))
404 return numScvs + outSideBoundaryVolVarIdx_++;
405 else
406 return globalScvIdx; // fallback for parallel, won't be used anyway
407 }();
408
409 return std::array{globalScvIdx, globalOutsideScvIdx};
410 }();
411
412 const auto boundaryType = [&]
413 {
414 if (onProcessorBoundary_(lateralIntersection))
415 return SubControlVolumeFace::BoundaryType::processorBoundary;
416 else if (onDomainBoundary_(lateralIntersection))
417 return SubControlVolumeFace::BoundaryType::physicalBoundary;
418 else
419 return SubControlVolumeFace::BoundaryType::interior;
420 }();
421
422 scvfs_.emplace_back(
423 elementGeometry,
424 intersectionGeometry,
425 geometryHelper.facet(lateralFacetIndex, element).geometry(),
426 globalScvIndicesForLateralFace, // TODO higher order
427 localScvfIdx,
428 globalScvfIdx,
429 lateralIntersection.centerUnitOuterNormal(),
430 SubControlVolumeFace::FaceType::lateral,
431 boundaryType
432 );
433
434 globalScvfIndices[localScvfIdx] = globalScvfIdx++;
435 ++localScvfIdx;
436
437 if (onDomainBoundary_(lateralIntersection))
438 {
439 ++numBoundaryScvf_;
440 hasBoundaryScvf_[eIdx] = true;
441 }
442 } // end loop over lateral facets
443
444 } // end first loop over intersections
445
446 // do a second loop over all intersections to add frontal boundary faces
447 int localScvfIdx = minNumScvfsPerElement;
448 for (const auto& intersection : intersections(this->gridView(), element))
449 {
450 // the frontal sub control volume face at a domain boundary (coincides with element face)
451 if (onDomainBoundary_(intersection))
452 {
453 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
454 const auto globalScvIdx = getGlobalScvIdx(eIdx, localScvIdx);
455 ++numBoundaryScvf_;
456
457 // the frontal sub control volume face at the boundary
458 scvfs_.emplace_back(
459 element.geometry(),
460 intersection.geometry(),
461 std::array{globalScvIdx, globalScvIdx}, // TODO outside boundary, periodic, parallel?
462 localScvfIdx,
463 globalScvfIdx,
464 intersection.centerUnitOuterNormal(),
465 SubControlVolumeFace::FaceType::frontal,
466 SubControlVolumeFace::BoundaryType::physicalBoundary
467 );
468
469 globalScvfIndices.push_back(globalScvfIdx);
470 ++globalScvfIdx;
471 ++localScvfIdx;
472 hasBoundaryScvf_[eIdx] = true;
473 }
474 }
475 }
476
477 connectivityMap_.update(*this);
478 }
479
480 bool onDomainBoundary_(const typename GridView::Intersection& intersection) const
481 {
482 return !intersection.neighbor() && intersection.boundary();
483 }
484
485 bool onProcessorBoundary_(const typename GridView::Intersection& intersection) const
486 {
487 return !intersection.neighbor() && !intersection.boundary();
488 }
489
490 bool onPeriodicBoundary_(const typename GridView::Intersection& intersection) const
491 {
492 return intersection.boundary() && intersection.neighbor();
493 }
494
495 // mappers
496 ConnectivityMap connectivityMap_;
497 IntersectionMapper intersectionMapper_;
498
499 std::vector<SubControlVolume> scvs_;
500 std::vector<SubControlVolumeFace> scvfs_;
501 GridIndexType numBoundaryScv_;
502 GridIndexType numBoundaryScvf_;
503 GridIndexType outSideBoundaryVolVarIdx_;
504 std::vector<bool> hasBoundaryScvf_;
505
506 std::vector<std::vector<GridIndexType>> scvfIndicesOfElement_;
507
508 // a map for periodic boundary vertices
509 std::unordered_map<GridIndexType, GridIndexType> periodicFaceMap_;
510};
511
518template<class GV, class Traits>
520: public BaseGridGeometry<GV, Traits>
521{
524 using GridIndexType = typename IndexTraits<GV>::GridIndex;
525 using LocalIndexType = typename IndexTraits<GV>::LocalIndex;
526 using SmallLocalIndexType = typename IndexTraits<GV>::SmallLocalIndex;
527 using Element = typename GV::template Codim<0>::Entity;
528
529 using IntersectionMapper = typename Traits::IntersectionMapper;
530 using ConnectivityMap = typename Traits::template ConnectivityMap<ThisType>;
531
532 static constexpr auto dim = Traits::StaticInfo::dim;
533 static constexpr auto numScvsPerElement = Traits::StaticInfo::numScvsPerElement;
534 static constexpr auto numLateralScvfsPerScv = Traits::StaticInfo::numLateralScvfsPerScv;
535 static constexpr auto numLateralScvfsPerElement = Traits::StaticInfo::numLateralScvfsPerElement;
536 static constexpr auto minNumScvfsPerElement = Traits::StaticInfo::minNumScvfsPerElement;
537 static constexpr auto maxNumScvfsPerElement = Traits::StaticInfo::maxNumScvfsPerElement;
538
539public:
542 static constexpr DiscretizationMethod discMethod{};
543
544 static constexpr bool cachingEnabled = false;
545
549 using LocalView = typename Traits::template LocalView<ThisType, false>;
551 using SubControlVolume = typename Traits::SubControlVolume;
553 using SubControlVolumeFace = typename Traits::SubControlVolumeFace;
555 using GridView = GV;
557 using GeometryHelper = typename Traits::GeometryHelper;
559 using LocalIntersectionMapper = typename Traits::LocalIntersectionMapper;
561 using StaticInformation = typename Traits::StaticInfo;
564
566 FaceCenteredStaggeredFVGridGeometry(std::shared_ptr<BasicGridGeometry> gg, const std::string& paramGroup = "")
567 : ParentType(std::move(gg))
568 , intersectionMapper_(this->gridView())
569 {
570 // Check if the overlap size is what we expect
572 DUNE_THROW(Dune::InvalidStateException, "The staggered discretization method needs at least an overlap of 1 for parallel computations. "
573 << " Set the parameter \"Grid.Overlap\" in the input file.");
574
575 update_();
576 }
577
579 FaceCenteredStaggeredFVGridGeometry(const GridView& gridView, const std::string& paramGroup = "")
580 : FaceCenteredStaggeredFVGridGeometry(std::make_shared<BasicGridGeometry>(gridView), paramGroup)
581 {}
582
584 std::size_t numScv() const
585 { return numScvs_; }
586
588 std::size_t numScvf() const
589 { return numScvf_; }
590
592 std::size_t numBoundaryScv() const
593 { return numBoundaryScv_; }
594
596 std::size_t numBoundaryScvf() const
597 { return numBoundaryScvf_; }
598
600 std::size_t numIntersections() const
601 { return intersectionMapper_.numIntersections(); }
602
604 std::size_t numDofs() const
605 { return this->gridView().size(1); }
606
611 const ConnectivityMap& connectivityMap() const
612 { return connectivityMap_; }
613
615 bool hasBoundaryScvf(GridIndexType eIdx) const
616 { return hasBoundaryScvf_[eIdx]; }
617
619 const IntersectionMapper& intersectionMapper() const
620 { return intersectionMapper_; }
621
623 const std::vector<GridIndexType>& scvfIndicesOfElement(GridIndexType eIdx) const
624 { return scvfIndicesOfElement_[eIdx]; }
625
627 GridIndexType outsideVolVarIndex(GridIndexType scvfIdx) const
628 { return outsideVolVarIndices_.at(scvfIdx); }
629
631 void update(const GridView& gridView)
632 {
633 ParentType::update(gridView);
634 update_();
635 }
636
638 void update(GridView&& gridView)
639 {
640 ParentType::update(std::move(gridView));
641 update_();
642 }
643
645 bool dofOnPeriodicBoundary(GridIndexType dofIdx) const
646 { return periodicFaceMap_.count(dofIdx); }
647
649 GridIndexType periodicallyMappedDof(GridIndexType dofIdx) const
650 { return periodicFaceMap_.at(dofIdx); }
651
653 const std::unordered_map<GridIndexType, GridIndexType>& periodicDofMap() const
654 { return periodicFaceMap_; }
655
657 [[deprecated("Will be removed after release 3.9. Use periodicDofMap() instead.")]]
658 const std::unordered_map<GridIndexType, GridIndexType>& periodicVertexMap() const
659 { return periodicDofMap(); }
660
661private:
662
663 void update_()
664 {
665 intersectionMapper_.update(this->gridView());
666
667 // clear local data
668 numScvf_ = 0;
669 numBoundaryScv_ = 0;
670 numBoundaryScvf_ = 0;
671 hasBoundaryScvf_.clear();
672 scvfIndicesOfElement_.clear();
673 outsideVolVarIndices_.clear();
674
675 // determine size of containers
676 const auto numElements = this->gridView().size(0);
677 scvfIndicesOfElement_.resize(numElements);
678 hasBoundaryScvf_.resize(numElements, false);
679 numScvs_ = numElements*numScvsPerElement;
680
681 GeometryHelper geometryHelper(this->gridView());
682
683 // get the global scv indices first
684 GridIndexType scvfIdx = 0;
685
686 GridIndexType neighborVolVarIdx = numScvs_;
687
688 for (const auto& element : elements(this->gridView()))
689 {
690 const auto eIdx = this->elementMapper().index(element);
691 assert(numScvsPerElement == element.subEntities(1));
692
693 // the element-wise index sets for finite volume geometry
694 auto& globalScvfIndices = scvfIndicesOfElement_[eIdx];
695 globalScvfIndices.reserve(maxNumScvfsPerElement);
696 globalScvfIndices.resize(minNumScvfsPerElement);
697
698 // keep track of frontal boundary scvfs
699 std::size_t numFrontalBoundaryScvfs = 0;
700
701 using LocalIntersectionIndexMapper = FaceCenteredStaggeredLocalIntersectionIndexMapper<GridView>;
702 LocalIntersectionIndexMapper localIsMapper;
703 localIsMapper.update(this->gridView(), element);
704
705 for (const auto& intersection : intersections(this->gridView(), element))
706 {
707 const auto localScvIdx = localIsMapper.realToRefIdx(intersection.indexInInside());
708 auto localScvfIdx = localScvIdx*(1 + numLateralScvfsPerScv);
709
710 assert(localIsMapper.refToRealIdx(localScvIdx) == intersection.indexInInside());
711 // the frontal sub control volume face at the element center
712 globalScvfIndices[localScvfIdx] = scvfIdx++;
713 ++localScvfIdx;
714
715 if constexpr(dim > 1)
716 {
717 // the lateral sub control volume faces
718 for (const auto lateralFacetIndex : Dune::transformedRangeView(geometryHelper.localLaterFaceIndices(localScvIdx),
719 [&](auto idx) { return localIsMapper.refToRealIdx(idx) ;})
720 )
721 {
722 if (onDomainBoundary_(geometryHelper.intersection(lateralFacetIndex, element)))
723 {
724 outsideVolVarIndices_[scvfIdx] = neighborVolVarIdx++;
725 ++numBoundaryScvf_;
726 hasBoundaryScvf_[eIdx] = true;
727 }
728
729 globalScvfIndices[localScvfIdx] = scvfIdx++;
730 ++localScvfIdx;
731 }
732 }
733
734 // handle physical domain boundary
735 if (onDomainBoundary_(intersection))
736 {
737 ++numBoundaryScv_; // frontal face
738 numBoundaryScv_ += numLateralScvfsPerScv; // boundary scvs for lateral faces
739 ++numFrontalBoundaryScvfs;
740 ++numBoundaryScvf_;
741 hasBoundaryScvf_[eIdx] = true;
742 }
743
744 // handle periodic boundaries
745 if (onPeriodicBoundary_(intersection))
746 {
747 this->setPeriodic();
748
749 const auto& otherElement = intersection.outside();
750
751 SmallLocalIndexType otherIntersectionLocalIdx = 0;
752 bool periodicFaceFound = false;
753
754 for (const auto& otherIntersection : intersections(this->gridView(), otherElement))
755 {
756 if (periodicFaceFound)
757 continue;
758
759 if (Dune::FloatCmp::eq(intersection.centerUnitOuterNormal()*otherIntersection.centerUnitOuterNormal(), -1.0, 1e-7))
760 {
761 const auto periodicDofIdx = intersectionMapper().globalIntersectionIndex(otherElement, otherIntersectionLocalIdx);
762 const auto dofIndex = intersectionMapper().globalIntersectionIndex(element, localScvIdx);
763 periodicFaceMap_[dofIndex] = periodicDofIdx;
764 periodicFaceFound = true;
765 }
766
767 ++otherIntersectionLocalIdx;
768 }
769 }
770 }
771
772 // add global indices of frontal boundary scvfs last
773 for (std::size_t i = 0; i < numFrontalBoundaryScvfs; ++i)
774 globalScvfIndices.push_back(scvfIdx++);
775 }
776
777 // set number of subcontrolvolume faces
778 numScvf_ = scvfIdx;
779
780 connectivityMap_.update(*this);
781 }
782
783 bool onDomainBoundary_(const typename GridView::Intersection& intersection) const
784 {
785 return !intersection.neighbor() && intersection.boundary();
786 }
787
788 bool onProcessorBoundary_(const typename GridView::Intersection& intersection) const
789 {
790 return !intersection.neighbor() && !intersection.boundary();
791 }
792
793 bool onPeriodicBoundary_(const typename GridView::Intersection& intersection) const
794 {
795 return intersection.boundary() && intersection.neighbor();
796 }
797
798 // mappers
799 ConnectivityMap connectivityMap_;
800 IntersectionMapper intersectionMapper_;
801
803 std::size_t numScvs_;
804 std::size_t numScvf_;
805 std::size_t numBoundaryScv_;
806 std::size_t numBoundaryScvf_;
807 std::vector<bool> hasBoundaryScvf_;
808
809 std::vector<std::vector<GridIndexType>> scvfIndicesOfElement_;
810
811 // a map for periodic boundary vertices
812 std::unordered_map<GridIndexType, GridIndexType> periodicFaceMap_;
813 std::unordered_map<GridIndexType, GridIndexType> outsideVolVarIndices_;
814};
815
816} // end namespace Dumux
817
818#endif
Base class for grid geometries.
Check the overlap size for different discretization methods.
Base class for all grid geometries.
Definition basegridgeometry.hh:52
defines a standard intersection mapper for mapping of global DOFs assigned to faces....
Definition intersectionmapper.hh:29
Definition localintersectionindexmapper.hh:29
Stores the dof indices corresponding to the neighboring scvs that contribute to the derivative calcul...
Definition facecentered/staggered/connectivitymap.hh:30
Definition discretization/facecentered/staggered/fvelementgeometry.hh:124
Base class for the finite volume geometry vector for face-centered staggered models This builds up th...
Definition discretization/facecentered/staggered/fvgridgeometry.hh:521
typename Traits::SubControlVolumeFace SubControlVolumeFace
export the type of sub control volume
Definition discretization/facecentered/staggered/fvgridgeometry.hh:553
typename Traits::GeometryHelper GeometryHelper
export the geometry helper type
Definition discretization/facecentered/staggered/fvgridgeometry.hh:557
FaceCenteredStaggeredFVGridGeometry(std::shared_ptr< BasicGridGeometry > gg, const std::string &paramGroup="")
Constructor with basic grid geometry used to share state with another grid geometry on the same grid ...
Definition discretization/facecentered/staggered/fvgridgeometry.hh:566
typename Traits::SubControlVolume SubControlVolume
export the type of sub control volume
Definition discretization/facecentered/staggered/fvgridgeometry.hh:551
std::size_t numDofs() const
the total number of dofs
Definition discretization/facecentered/staggered/fvgridgeometry.hh:604
GridIndexType periodicallyMappedDof(GridIndexType dofIdx) const
The index of the d.o.f. on the other side of the periodic boundary.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:649
const ConnectivityMap & connectivityMap() const
Returns the connectivity map of which dofs have derivatives with respect to a given dof.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:611
std::size_t numBoundaryScv() const
The total number of boundary sub control volumes.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:592
const IntersectionMapper & intersectionMapper() const
Return a reference to the intersection mapper.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:619
void update(GridView &&gridView)
update all fvElementGeometries (call this after grid adaption)
Definition discretization/facecentered/staggered/fvgridgeometry.hh:638
FaceCenteredStaggeredFVGridGeometry(const GridView &gridView, const std::string &paramGroup="")
Constructor from gridView.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:579
bool dofOnPeriodicBoundary(GridIndexType dofIdx) const
If a d.o.f. is on a periodic boundary.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:645
typename Traits::StaticInfo StaticInformation
export static information
Definition discretization/facecentered/staggered/fvgridgeometry.hh:561
const std::vector< GridIndexType > & scvfIndicesOfElement(GridIndexType eIdx) const
Get the global sub control volume face indices of an element.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:623
bool hasBoundaryScvf(GridIndexType eIdx) const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:615
const std::unordered_map< GridIndexType, GridIndexType > & periodicVertexMap() const
Returns the map between dofs across periodic boundaries.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:658
std::size_t numScvf() const
The total number of sub control volume faces.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:588
BasicGridGeometry_t< GV, Traits > BasicGridGeometry
export basic grid geometry type for the alternative constructor
Definition discretization/facecentered/staggered/fvgridgeometry.hh:547
std::size_t numBoundaryScvf() const
The total number of boundary sub control volume faces.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:596
const std::unordered_map< GridIndexType, GridIndexType > & periodicDofMap() const
Returns the map between dofs across periodic boundaries.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:653
std::size_t numScv() const
The total number of sub control volumes.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:584
GV GridView
export the grid view type
Definition discretization/facecentered/staggered/fvgridgeometry.hh:555
typename Traits::LocalIntersectionMapper LocalIntersectionMapper
export the local intersection mapper
Definition discretization/facecentered/staggered/fvgridgeometry.hh:559
GridIndexType outsideVolVarIndex(GridIndexType scvfIdx) const
Get the global sub control volume face indices of an element.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:627
void update(const GridView &gridView)
update all fvElementGeometries (call this after grid adaption)
Definition discretization/facecentered/staggered/fvgridgeometry.hh:631
std::size_t numIntersections() const
The total number of intersections.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:600
Extrusion_t< Traits > Extrusion
export the type of extrusion
Definition discretization/facecentered/staggered/fvgridgeometry.hh:563
typename Traits::template LocalView< ThisType, false > LocalView
export the type of the fv element geometry (the local view type)
Definition discretization/facecentered/staggered/fvgridgeometry.hh:549
Base class for the finite volume geometry vector for staggered models This builds up the sub control ...
Definition discretization/facecentered/staggered/fvgridgeometry.hh:95
const ConnectivityMap & connectivityMap() const
Returns the connectivity map of which dofs have derivatives with respect to a given dof.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:225
std::size_t numScvf() const
The total number of sub control volume faces.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:167
const IntersectionMapper & intersectionMapper() const
Return a reference to the intersection mapper.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:233
bool hasBoundaryScvf(GridIndexType eIdx) const
Returns whether one of the geometry's scvfs lies on a boundary.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:229
void update(const GridView &gridView)
update all fvElementGeometries (call this after grid adaption)
Definition discretization/facecentered/staggered/fvgridgeometry.hh:187
typename Traits::SubControlVolumeFace SubControlVolumeFace
export the type of sub control volume
Definition discretization/facecentered/staggered/fvgridgeometry.hh:132
typename Traits::LocalIntersectionMapper LocalIntersectionMapper
export the local intersection mapper
Definition discretization/facecentered/staggered/fvgridgeometry.hh:138
typename Traits::GeometryHelper GeometryHelper
export the geometry helper type
Definition discretization/facecentered/staggered/fvgridgeometry.hh:136
const std::vector< GridIndexType > & scvfIndicesOfElement(GridIndexType eIdx) const
Get the global sub control volume face indices of an element.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:218
FaceCenteredStaggeredFVGridGeometry(std::shared_ptr< BasicGridGeometry > gg, const std::string &paramGroup="")
Constructor with basic grid geometry used to share state with another grid geometry on the same grid ...
Definition discretization/facecentered/staggered/fvgridgeometry.hh:145
typename Traits::StaticInfo StaticInformation
export static information
Definition discretization/facecentered/staggered/fvgridgeometry.hh:140
std::size_t numDofs() const
the total number of dofs
Definition discretization/facecentered/staggered/fvgridgeometry.hh:183
const std::unordered_map< GridIndexType, GridIndexType > & periodicVertexMap() const
Returns the map between dofs across periodic boundaries.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:250
std::size_t numBoundaryScv() const
The total number of boundary sub control volumes.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:171
FaceCenteredStaggeredFVGridGeometry(const GridView &gridView, const std::string &paramGroup="")
Constructor from gridView.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:158
GridIndexType periodicallyMappedDof(GridIndexType dofIdx) const
The index of the d.o.f. on the other side of the periodic boundary.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:241
GV GridView
export the grid view type
Definition discretization/facecentered/staggered/fvgridgeometry.hh:134
const SubControlVolume & scv(GridIndexType scvIdx) const
Get a sub control volume with a global scv index.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:201
typename Traits::SubControlVolume SubControlVolume
export the type of sub control volume
Definition discretization/facecentered/staggered/fvgridgeometry.hh:130
typename Traits::template LocalView< ThisType, true > LocalView
export the type of the fv element geometry (the local view type)
Definition discretization/facecentered/staggered/fvgridgeometry.hh:128
void update(GridView &&gridView)
update all fvElementGeometries (call this after grid adaption)
Definition discretization/facecentered/staggered/fvgridgeometry.hh:194
BasicGridGeometry_t< GV, Traits > BasicGridGeometry
export basic grid geometry type for the alternative constructor
Definition discretization/facecentered/staggered/fvgridgeometry.hh:126
Extrusion_t< Traits > Extrusion
export the type of extrusion
Definition discretization/facecentered/staggered/fvgridgeometry.hh:142
std::size_t numBoundaryScvf() const
The total number of boundary sub control volume faces.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:175
std::size_t numIntersections() const
The total number of intersections.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:179
const SubControlVolumeFace & scvf(GridIndexType scvfIdx) const
Get a sub control volume face with a global scvf index.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:214
const std::unordered_map< GridIndexType, GridIndexType > & periodicDofMap() const
Returns the map between dofs across periodic boundaries.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:245
bool dofOnPeriodicBoundary(GridIndexType dofIdx) const
If a d.o.f. is on a periodic boundary.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:237
std::size_t numScv() const
The total number of sub control volumes.
Definition discretization/facecentered/staggered/fvgridgeometry.hh:163
auto scvs(const LocalView &fvGeometry) const
Definition discretization/facecentered/staggered/fvgridgeometry.hh:206
Base class for the finite volume geometry vector for face-centered staggered models This builds up th...
Definition discretization/facecentered/staggered/fvgridgeometry.hh:84
Face centered staggered geometry helper.
Definition discretization/facecentered/staggered/geometryhelper.hh:28
Face centered staggered sub control volume face.
Definition discretization/facecentered/staggered/subcontrolvolumeface.hh:54
Face centered staggered sub control volume.
Definition discretization/facecentered/staggered/subcontrolvolume.hh:52
Defines the default element and vertex mapper types.
Geometry helper for face-centered staggered scheme.
Face centered staggered sub control volume.
Face centered staggered sub control volume face.
Helper classes to compute the integration elements.
Stores the dof indices corresponding to the neighboring scvs that contribute to the derivative calcul...
Dune::Std::detected_or_t< Dumux::BasicGridGeometry< GV, typename T::ElementMapper, typename T::VertexMapper >, Detail::SpecifiesBaseGridGeometry, T > BasicGridGeometry_t
Type of the basic grid geometry implementation used as backend.
Definition basegridgeometry.hh:38
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
Defines the index types used for grid and local indices.
defines intersection mappers.
Provides a mapping of local intersection indices (indexInInside) such that the local indices always f...
Define some often used mathematical functions.
The available discretization methods in Dumux.
Definition adapt.hh:17
typename Extrusion< T >::type Extrusion_t
Convenience alias for obtaining the extrusion type.
Definition extrusion.hh:166
Definition common/pdesolver.hh:24
Base class for the finite volume geometry vector for face-centered staggered models This builds up th...
Check if the overlap size is valid for a given discretization method.
Definition checkoverlapsize.hh:28
Definition defaultmappertraits.hh:23
Definition discretization/facecentered/staggered/fvgridgeometry.hh:62
static constexpr auto numLateralScvfsPerElement
Definition discretization/facecentered/staggered/fvgridgeometry.hh:67
static constexpr auto maxNumScvfsPerElement
Definition discretization/facecentered/staggered/fvgridgeometry.hh:70
static constexpr auto numScvsPerElement
Definition discretization/facecentered/staggered/fvgridgeometry.hh:65
static constexpr auto numFacesPerElement
Definition discretization/facecentered/staggered/fvgridgeometry.hh:64
static constexpr auto dim
Definition discretization/facecentered/staggered/fvgridgeometry.hh:63
static constexpr auto numLateralScvfsPerScv
Definition discretization/facecentered/staggered/fvgridgeometry.hh:66
static constexpr auto minNumScvfsPerElement
Definition discretization/facecentered/staggered/fvgridgeometry.hh:68
The default traits for the face-center staggered finite volume grid geometry Defines the scv and scvf...
Definition discretization/facecentered/staggered/fvgridgeometry.hh:48
typename GridView::IndexSet::IndexType GridIndex
Definition indextraits.hh:27
std::uint_least8_t SmallLocalIndex
Definition indextraits.hh:29
unsigned int LocalIndex
Definition indextraits.hh:28