GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: dumux/dumux/discretization/pq1bubble/fvgridgeometry.hh
Date: 2025-04-12 19:19:20
Exec Total Coverage
Lines: 106 109 97.2%
Functions: 24 24 100.0%
Branches: 129 238 54.2%

Line Branch Exec Source
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-FileCopyrightText: Copyright © DuMux Project contributors, see AUTHORS.md in root folder
5 // SPDX-License-Identifier: GPL-3.0-or-later
6 //
7 /*!
8 * \file
9 * \ingroup PQ1BubbleDiscretization
10 * \brief Base class for the finite volume geometry vector for the pq1bubble method
11 * This builds up the sub control volumes and sub control volume faces
12 * for each element of the grid partition.
13 */
14 #ifndef DUMUX_DISCRETIZATION_PQ1BUBBLE_GRID_GEOMETRY_HH
15 #define DUMUX_DISCRETIZATION_PQ1BUBBLE_GRID_GEOMETRY_HH
16
17 #include <utility>
18 #include <unordered_map>
19
20 #include <dune/grid/common/mcmgmapper.hh>
21
22 #include <dumux/discretization/method.hh>
23
24 #include <dumux/common/indextraits.hh>
25 #include <dumux/common/defaultmappertraits.hh>
26
27 #include <dumux/geometry/center.hh>
28 #include <dumux/geometry/volume.hh>
29
30 #include <dumux/discretization/basegridgeometry.hh>
31 #include <dumux/discretization/pq1bubble/geometryhelper.hh>
32 #include <dumux/discretization/pq1bubble/fvelementgeometry.hh>
33 #include <dumux/discretization/pq1bubble/subcontrolvolume.hh>
34 #include <dumux/discretization/pq1bubble/subcontrolvolumeface.hh>
35 #include <dumux/discretization/pq1bubble/pq1bubblefecache.hh>
36 #include <dumux/discretization/extrusion.hh>
37
38 #include <dumux/io/grid/periodicgridtraits.hh>
39
40 namespace Dumux {
41
42 namespace Detail {
43 template<class GV, class T>
44 using PQ1BubbleGeometryHelper_t = Dune::Std::detected_or_t<
45 Dumux::PQ1BubbleGeometryHelper<GV, typename T::SubControlVolume, typename T::SubControlVolumeFace>,
46 SpecifiesGeometryHelper,
47 T
48 >;
49 } // end namespace Detail
50
51 template <class GridView>
52 struct PQ1BubbleMapperTraits :public DefaultMapperTraits<GridView>
53 {
54 using DofMapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
55
56 /**
57 * \brief layout for elements and vertices
58 *
59 */
60 14 static Dune::MCMGLayout layout()
61 {
62
2/2
✓ Branch 0 taken 60 times.
✓ Branch 1 taken 28 times.
88 return [](Dune::GeometryType gt, int dimgrid) {
63
4/4
✓ Branch 0 taken 60 times.
✓ Branch 1 taken 28 times.
✓ Branch 2 taken 28 times.
✓ Branch 3 taken 32 times.
88 return (gt.dim() == dimgrid) || (gt.dim() == 0);
64
1/2
✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
14 };
65 }
66 };
67
68 /*!
69 * \ingroup PQ1BubbleDiscretization
70 * \brief The default traits for the pq1bubble finite volume grid geometry
71 * Defines the scv and scvf types and the mapper types
72 * \tparam the grid view type
73 */
74 template<class GridView, class MapperTraits = PQ1BubbleMapperTraits<GridView>>
75 struct PQ1BubbleDefaultGridGeometryTraits
76 : public MapperTraits
77 {
78 using SubControlVolume = PQ1BubbleSubControlVolume<GridView>;
79 using SubControlVolumeFace = PQ1BubbleSubControlVolumeFace<GridView>;
80
81 template<class GridGeometry, bool enableCache>
82 using LocalView = PQ1BubbleFVElementGeometry<GridGeometry, enableCache>;
83 };
84
85 /*!
86 * \ingroup PQ1BubbleDiscretization
87 * \brief Base class for the finite volume geometry vector for pq1bubble schemes
88 * This builds up the sub control volumes and sub control volume faces
89 * \note For caching enabled we store the fv geometries for the whole grid view which is memory intensive but faster
90 */
91 template<class Scalar,
92 class GV,
93 bool enableCaching = true,
94 class Traits = PQ1BubbleDefaultGridGeometryTraits<GV>>
95 class PQ1BubbleFVGridGeometry
96 : public BaseGridGeometry<GV, Traits>
97 {
98 using ThisType = PQ1BubbleFVGridGeometry<Scalar, GV, enableCaching, Traits>;
99 using ParentType = BaseGridGeometry<GV, Traits>;
100 using GridIndexType = typename IndexTraits<GV>::GridIndex;
101 using LocalIndexType = typename IndexTraits<GV>::LocalIndex;
102
103 using Element = typename GV::template Codim<0>::Entity;
104 using CoordScalar = typename GV::ctype;
105 static const int dim = GV::dimension;
106 static const int dimWorld = GV::dimensionworld;
107
108 static_assert(dim > 1, "Only implemented for dim > 1");
109
110 public:
111 //! export the discretization method this geometry belongs to
112 using DiscretizationMethod = DiscretizationMethods::PQ1Bubble;
113 static constexpr DiscretizationMethod discMethod{};
114
115 //! export the type of the fv element geometry (the local view type)
116 using LocalView = typename Traits::template LocalView<ThisType, true>;
117 //! export the type of sub control volume
118 using SubControlVolume = typename Traits::SubControlVolume;
119 //! export the type of sub control volume
120 using SubControlVolumeFace = typename Traits::SubControlVolumeFace;
121 //! export the type of extrusion
122 using Extrusion = Extrusion_t<Traits>;
123 //! export dof mapper type
124 using DofMapper = typename Traits::DofMapper;
125 //! export the finite element cache type
126 using FeCache = Dumux::PQ1BubbleFECache<CoordScalar, Scalar, dim>;
127 //! export the grid view type
128 using GridView = GV;
129 //! export whether the grid(geometry) supports periodicity
130 using SupportsPeriodicity = typename PeriodicGridTraits<typename GV::Grid>::SupportsPeriodicity;
131
132 //! Constructor
133 14 PQ1BubbleFVGridGeometry(const GridView gridView)
134 : ParentType(gridView)
135
1/2
✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
14 , dofMapper_(gridView, Traits::layout())
136 14 , cache_(*this)
137
3/6
✓ Branch 2 taken 14 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 14 times.
✓ Branch 7 taken 14 times.
✗ Branch 8 not taken.
28 , periodicGridTraits_(this->gridView().grid())
138 {
139
1/2
✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
14 update_();
140 14 }
141
142 //! The dofMapper
143 4714103 const DofMapper& dofMapper() const
144
7/12
✓ Branch 0 taken 4 times.
✓ Branch 1 taken 6 times.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 2 times.
✓ Branch 6 taken 6 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 4 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 6 times.
✗ Branch 13 not taken.
✗ Branch 2 not taken.
✗ Branch 5 not taken.
4714113 { return dofMapper_; }
145
146 //! The total number of sub control volumes
147 std::size_t numScv() const
148 { return numScv_; }
149
150 //! The total number of sub control volume faces
151 std::size_t numScvf() const
152 { return numScvf_; }
153
154 //! The total number of boundary sub control volume faces
155 2 std::size_t numBoundaryScvf() const
156
2/4
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
2 { return numBoundaryScvf_; }
157
158 //! The total number of degrees of freedom
159 86 std::size_t numDofs() const
160
6/13
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 11 times.
✗ Branch 5 not taken.
✓ Branch 9 taken 9 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 5 times.
✗ Branch 13 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 15 taken 1 times.
✗ Branch 16 not taken.
✗ Branch 8 not taken.
73 { return this->dofMapper().size(); }
161
162 //! update all geometries (call this after grid adaption)
163 void update(const GridView& gridView)
164 {
165 ParentType::update(gridView);
166 update_();
167 }
168
169 //! update all geometries (call this after grid adaption)
170 void update(GridView&& gridView)
171 {
172 ParentType::update(std::move(gridView));
173 update_();
174 }
175
176 //! The finite element cache for creating local FE bases
177 19787966 const FeCache& feCache() const
178 19787966 { return feCache_; }
179
180 //! If a vertex / d.o.f. is on the boundary
181 2692076 bool dofOnBoundary(GridIndexType dofIdx) const
182
4/4
✓ Branch 0 taken 74970 times.
✓ Branch 1 taken 2501314 times.
✓ Branch 2 taken 2047 times.
✓ Branch 3 taken 113745 times.
2692076 { return boundaryDofIndices_[dofIdx]; }
183
184 //! If a vertex / d.o.f. is on a periodic boundary
185 72336 bool dofOnPeriodicBoundary(GridIndexType dofIdx) const
186
0/4
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
144672 { return periodicDofMap_.count(dofIdx); }
187
188 //! The index of the vertex / d.o.f. on the other side of the periodic boundary
189 GridIndexType periodicallyMappedDof(GridIndexType dofIdx) const
190 { return periodicDofMap_.at(dofIdx); }
191
192 //! Returns the map between dofs across periodic boundaries
193 const std::unordered_map<GridIndexType, GridIndexType>& periodicDofMap() const
194 25 { return periodicDofMap_; }
195
196 //! local view of this object (constructed with the internal cache)
197 1363076 friend inline LocalView localView(const PQ1BubbleFVGridGeometry& gg)
198
13/17
✓ Branch 1 taken 11 times.
✓ Branch 2 taken 6 times.
✓ Branch 4 taken 42003 times.
✓ Branch 5 taken 6 times.
✓ Branch 7 taken 904 times.
✓ Branch 8 taken 4800 times.
✓ Branch 11 taken 21601 times.
✗ Branch 12 not taken.
✓ Branch 14 taken 21604 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✓ Branch 17 taken 1 times.
✓ Branch 3 taken 1 times.
✓ Branch 6 taken 1 times.
✓ Branch 10 taken 4 times.
✓ Branch 13 taken 2 times.
✗ Branch 9 not taken.
840964 { return { gg.cache_ }; }
199
200 private:
201
202 class PQ1BubbleGridGeometryCache
203 {
204 friend class PQ1BubbleFVGridGeometry;
205 public:
206 //! export the geometry helper type
207 using GeometryHelper = Detail::PQ1BubbleGeometryHelper_t<GV, Traits>;
208
209 14 explicit PQ1BubbleGridGeometryCache(const PQ1BubbleFVGridGeometry& gg)
210
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 14 times.
14 : gridGeometry_(&gg)
211 {}
212
213 27970120 const PQ1BubbleFVGridGeometry& gridGeometry() const
214
10/12
✓ Branch 1 taken 66114 times.
✓ Branch 2 taken 126822 times.
✓ Branch 11 taken 4412 times.
✓ Branch 12 taken 10240 times.
✓ Branch 14 taken 4412 times.
✓ Branch 15 taken 1920 times.
✓ Branch 10 taken 7200 times.
✓ Branch 3 taken 2381768 times.
✓ Branch 0 taken 1580 times.
✓ Branch 9 taken 4800 times.
✗ Branch 13 not taken.
✗ Branch 16 not taken.
25109596 { return *gridGeometry_; }
215
216 //! Get the global sub control volume indices of an element
217 134154953 const std::vector<SubControlVolume>& scvs(GridIndexType eIdx) const
218
8/10
✓ Branch 1 taken 7406 times.
✓ Branch 2 taken 921130 times.
✓ Branch 6 taken 43720 times.
✓ Branch 7 taken 14400 times.
✓ Branch 9 taken 21600 times.
✗ Branch 10 not taken.
✓ Branch 3 taken 56830 times.
✓ Branch 4 taken 249154 times.
✓ Branch 5 taken 14800 times.
✗ Branch 8 not taken.
131784151 { return scvs_[eIdx]; }
219
220 //! Get the global sub control volume face indices of an element
221 5645311 const std::vector<SubControlVolumeFace>& scvfs(GridIndexType eIdx) const
222
1/5
✗ Branch 0 not taken.
✓ Branch 1 taken 12000 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 2 not taken.
5500451 { return scvfs_[eIdx]; }
223
224 //! Returns whether one of the geometry's scvfs lies on a boundary
225 42775 bool hasBoundaryScvf(GridIndexType eIdx) const
226
2/4
✓ Branch 0 taken 1362 times.
✓ Branch 1 taken 41413 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
42775 { return hasBoundaryScvf_[eIdx]; }
227
228 //! Local mappings necessary to construct geometries of scvfs
229 664 const std::vector<std::array<LocalIndexType, 2>>& scvfBoundaryGeometryKeys(GridIndexType eIdx) const
230 664 { return scvfBoundaryGeometryKeys_.at(eIdx); }
231
232 private:
233 14 void clear_()
234 {
235
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 14 times.
14 scvs_.clear();
236
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 14 times.
14 scvfs_.clear();
237 14 hasBoundaryScvf_.clear();
238 14 scvfBoundaryGeometryKeys_.clear();
239 14 }
240
241 std::vector<std::vector<SubControlVolume>> scvs_;
242 std::vector<std::vector<SubControlVolumeFace>> scvfs_;
243 std::vector<bool> hasBoundaryScvf_;
244 std::unordered_map<GridIndexType, std::vector<std::array<LocalIndexType, 2>>> scvfBoundaryGeometryKeys_;
245
246 const PQ1BubbleFVGridGeometry* gridGeometry_;
247 };
248
249 public:
250 //! the cache type (only the caching implementation has this)
251 //! this alias should only be used by the local view implementation
252 using Cache = PQ1BubbleGridGeometryCache;
253
254 private:
255 using GeometryHelper = typename Cache::GeometryHelper;
256
257 14 void update_()
258 {
259 14 cache_.clear_();
260 14 dofMapper_.update(this->gridView());
261
262 14 auto numElements = this->gridView().size(0);
263 14 cache_.scvs_.resize(numElements);
264 14 cache_.scvfs_.resize(numElements);
265 14 cache_.hasBoundaryScvf_.resize(numElements, false);
266
267 14 boundaryDofIndices_.assign(numDofs(), false);
268
269 14 numScv_ = 0;
270 14 numScvf_ = 0;
271 14 numBoundaryScvf_ = 0;
272
273 // Build the scvs and scv faces
274
6/9
✓ Branch 2 taken 9 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 3 times.
✗ Branch 6 not taken.
✓ Branch 8 taken 12000 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 12000 times.
✓ Branch 11 taken 3 times.
✓ Branch 1 taken 111266 times.
156385 for (const auto& element : elements(this->gridView()))
275 {
276
2/4
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 12000 times.
✗ Branch 5 not taken.
123261 auto eIdx = this->elementMapper().index(element);
277
278 // get the element geometry
279
1/2
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
123261 auto elementGeometry = element.geometry();
280
1/2
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
123261 const auto refElement = referenceElement(elementGeometry);
281
282 // instantiate the geometry helper
283
1/2
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
123261 GeometryHelper geometryHelper(elementGeometry);
284
285 123261 numScv_ += geometryHelper.numScv();
286 // construct the sub control volumes
287
2/4
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 12000 times.
✗ Branch 5 not taken.
123261 cache_.scvs_[eIdx].resize(geometryHelper.numScv());
288
4/4
✓ Branch 1 taken 559813 times.
✓ Branch 2 taken 111261 times.
✓ Branch 3 taken 49600 times.
✓ Branch 4 taken 12000 times.
671074 for (LocalIndexType scvLocalIdx = 0; scvLocalIdx < geometryHelper.numScv(); ++scvLocalIdx)
289 {
290
1/2
✓ Branch 1 taken 49600 times.
✗ Branch 2 not taken.
547813 auto corners = geometryHelper.getScvCorners(scvLocalIdx);
291 547813 cache_.scvs_[eIdx][scvLocalIdx] = SubControlVolume(
292
1/2
✓ Branch 1 taken 22069 times.
✗ Branch 2 not taken.
547813 geometryHelper.scvVolume(scvLocalIdx, corners),
293 1095626 geometryHelper.dofPosition(scvLocalIdx),
294
1/2
✓ Branch 1 taken 547813 times.
✗ Branch 2 not taken.
1095626 Dumux::center(corners),
295 scvLocalIdx,
296 eIdx,
297 547813 geometryHelper.dofIndex(this->dofMapper(), element, scvLocalIdx),
298
1/2
✓ Branch 1 taken 49600 times.
✗ Branch 2 not taken.
547813 geometryHelper.isOverlappingScv(scvLocalIdx)
299 );
300 }
301
302 // construct the sub control volume faces
303
1/2
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
123261 numScvf_ += geometryHelper.numInteriorScvf();
304
2/4
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 12000 times.
✗ Branch 5 not taken.
123261 cache_.scvfs_[eIdx].resize(geometryHelper.numInteriorScvf());
305 123261 LocalIndexType scvfLocalIdx = 0;
306
4/4
✓ Branch 1 taken 869932 times.
✓ Branch 2 taken 111261 times.
✓ Branch 3 taken 75200 times.
✓ Branch 4 taken 12000 times.
981193 for (; scvfLocalIdx < geometryHelper.numInteriorScvf(); ++scvfLocalIdx)
307 {
308
1/2
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
857932 const auto scvPair = geometryHelper.getScvPairForScvf(scvfLocalIdx);
309
1/2
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
857932 const auto corners = geometryHelper.getScvfCorners(scvfLocalIdx);
310 857932 const auto area = Dumux::convexPolytopeVolume<dim-1>(
311 geometryHelper.getInteriorScvfGeometryType(scvfLocalIdx),
312
1/2
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
1760004 [&](unsigned int i){ return corners[i]; }
313 );
314
315 857932 cache_.scvfs_[eIdx][scvfLocalIdx] = SubControlVolumeFace(
316 857932 Dumux::center(corners),
317 area,
318 857932 geometryHelper.normal(corners, scvPair),
319 std::move(scvPair),
320 scvfLocalIdx,
321
2/4
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 75200 times.
✗ Branch 5 not taken.
1715864 geometryHelper.isOverlappingScvf(scvfLocalIdx)
322 );
323 }
324
325 // construct the sub control volume faces on the domain boundary
326
7/9
✓ Branch 2 taken 174566 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 91746 times.
✓ Branch 6 taken 132420 times.
✓ Branch 7 taken 392704 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 18358 times.
✓ Branch 10 taken 284186 times.
✓ Branch 1 taken 33101 times.
1305696 for (const auto& intersection : intersections(this->gridView(), element))
327 {
328
3/4
✓ Branch 0 taken 7452 times.
✓ Branch 1 taken 119840 times.
✓ Branch 2 taken 580 times.
✗ Branch 3 not taken.
127872 if (intersection.boundary() && !intersection.neighbor())
329 {
330
1/2
✓ Branch 1 taken 5866 times.
✗ Branch 2 not taken.
7452 cache_.hasBoundaryScvf_[eIdx] = true;
331
332 7452 const auto localFacetIndex = intersection.indexInInside();
333
1/2
✓ Branch 1 taken 5866 times.
✗ Branch 2 not taken.
7452 const auto numBoundaryScvf = geometryHelper.numBoundaryScvf(localFacetIndex);
334 7452 numScvf_ += numBoundaryScvf;
335 7452 numBoundaryScvf_ += numBoundaryScvf;
336
337
2/2
✓ Branch 0 taken 17430 times.
✓ Branch 1 taken 7452 times.
24882 for (unsigned int isScvfLocalIdx = 0; isScvfLocalIdx < numBoundaryScvf; ++isScvfLocalIdx)
338 {
339 // find the scvs this scvf is belonging to
340
1/2
✓ Branch 1 taken 14246 times.
✗ Branch 2 not taken.
17430 const auto scvPair = geometryHelper.getScvPairForBoundaryScvf(localFacetIndex, isScvfLocalIdx);
341
1/2
✓ Branch 1 taken 14246 times.
✗ Branch 2 not taken.
17430 const auto corners = geometryHelper.getBoundaryScvfCorners(localFacetIndex, isScvfLocalIdx);
342 17430 const auto area = Dumux::convexPolytopeVolume<dim-1>(
343 geometryHelper.getBoundaryScvfGeometryType(isScvfLocalIdx),
344 42426 [&](unsigned int i){ return corners[i]; }
345 );
346
3/5
✓ Branch 1 taken 14678 times.
✓ Branch 2 taken 2752 times.
✓ Branch 4 taken 5544 times.
✗ Branch 5 not taken.
✗ Branch 3 not taken.
17430 cache_.scvfs_[eIdx].emplace_back(
347 34860 Dumux::center(corners),
348 area,
349
1/2
✓ Branch 1 taken 14246 times.
✗ Branch 2 not taken.
17430 intersection.centerUnitOuterNormal(),
350 std::move(scvPair),
351 scvfLocalIdx,
352 typename SubControlVolumeFace::Traits::BoundaryFlag{ intersection },
353
1/2
✓ Branch 1 taken 13086 times.
✗ Branch 2 not taken.
17430 geometryHelper.isOverlappingBoundaryScvf(localFacetIndex)
354 );
355
356 // store look-up map to construct boundary scvf geometries
357
2/4
✓ Branch 1 taken 14246 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14246 times.
✗ Branch 5 not taken.
17430 cache_.scvfBoundaryGeometryKeys_[eIdx].emplace_back(std::array<LocalIndexType, 2>{{
358 static_cast<LocalIndexType>(localFacetIndex),
359 static_cast<LocalIndexType>(isScvfLocalIdx)
360 }});
361
362 // increment local counter
363 17430 scvfLocalIdx++;
364 }
365
366 // TODO also move this to helper class
367
368 // add all vertices on the intersection to the set of boundary vertices
369
2/2
✓ Branch 0 taken 17430 times.
✓ Branch 1 taken 7452 times.
24882 for (int localVIdx = 0; localVIdx < numBoundaryScvf; ++localVIdx)
370 {
371 17430 const auto vIdx = refElement.subEntity(localFacetIndex, 1, localVIdx, dim);
372
1/2
✓ Branch 1 taken 14246 times.
✗ Branch 2 not taken.
17430 const auto vIdxGlobal = this->dofMapper().subIndex(element, vIdx, dim);
373 17430 boundaryDofIndices_[vIdxGlobal] = true;
374 }
375 }
376
377 // inform the grid geometry if we have periodic boundaries
378
1/2
✓ Branch 1 taken 340144 times.
✗ Branch 2 not taken.
424550 else if (periodicGridTraits_.isPeriodic(intersection))
379 {
380 this->setPeriodic();
381
382 // find the mapped periodic vertex of all vertices on periodic boundaries
383 const auto fIdx = intersection.indexInInside();
384 const auto numFaceVerts = refElement.size(fIdx, 1, dim);
385 const auto eps = 1e-7*(elementGeometry.corner(1) - elementGeometry.corner(0)).two_norm();
386 for (int localVIdx = 0; localVIdx < numFaceVerts; ++localVIdx)
387 {
388 const auto vIdx = refElement.subEntity(fIdx, 1, localVIdx, dim);
389 const auto vIdxGlobal = this->dofMapper().subIndex(element, vIdx, dim);
390 const auto vPos = elementGeometry.corner(vIdx);
391
392 const auto& outside = intersection.outside();
393 const auto outsideGeometry = outside.geometry();
394 for (const auto& isOutside : intersections(this->gridView(), outside))
395 {
396 // only check periodic vertices of the periodic neighbor
397 if (periodicGridTraits_.isPeriodic(isOutside))
398 {
399 const auto fIdxOutside = isOutside.indexInInside();
400 const auto numFaceVertsOutside = refElement.size(fIdxOutside, 1, dim);
401 for (int localVIdxOutside = 0; localVIdxOutside < numFaceVertsOutside; ++localVIdxOutside)
402 {
403 const auto vIdxOutside = refElement.subEntity(fIdxOutside, 1, localVIdxOutside, dim);
404 const auto vPosOutside = outsideGeometry.corner(vIdxOutside);
405 const auto shift = std::abs((this->bBoxMax()-this->bBoxMin())*intersection.centerUnitOuterNormal());
406 if (std::abs((vPosOutside-vPos).two_norm() - shift) < eps)
407 periodicDofMap_[vIdxGlobal] = this->dofMapper().subIndex(outside, vIdxOutside, dim);
408 }
409 }
410 }
411 }
412 }
413 }
414 }
415
416 // error check: periodic boundaries currently don't work for pq1bubble in parallel
417
1/6
✗ Branch 0 not taken.
✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
14 if (this->isPeriodic() && this->gridView().comm().size() > 1)
418 DUNE_THROW(Dune::NotImplemented, "Periodic boundaries for pq1bubble method for parallel simulations!");
419 14 }
420
421 DofMapper dofMapper_;
422
423 const FeCache feCache_;
424
425 std::size_t numScv_;
426 std::size_t numScvf_;
427 std::size_t numBoundaryScvf_;
428
429 // vertices on the boundary
430 std::vector<bool> boundaryDofIndices_;
431
432 // a map for periodic boundary vertices
433 std::unordered_map<GridIndexType, GridIndexType> periodicDofMap_;
434
435 Cache cache_;
436
437 PeriodicGridTraits<typename GridView::Grid> periodicGridTraits_;
438 };
439
440 } // end namespace Dumux
441
442 #endif
443