GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/dumux/discretization/pq1bubble/fvgridgeometry.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 82 108 75.9%
Functions: 24 49 49.0%
Branches: 165 360 45.8%

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-FileCopyrightInfo: 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 namespace Dumux {
39
40 namespace Detail {
41 template<class GV, class T>
42 using PQ1BubbleGeometryHelper_t = Dune::Std::detected_or_t<
43 Dumux::PQ1BubbleGeometryHelper<GV, typename T::SubControlVolume, typename T::SubControlVolumeFace>,
44 SpecifiesGeometryHelper,
45 T
46 >;
47 } // end namespace Detail
48
49 template <class GridView>
50 struct PQ1BubbleMapperTraits :public DefaultMapperTraits<GridView>
51 {
52 using DofMapper = Dune::MultipleCodimMultipleGeomTypeMapper<GridView>;
53
54 /**
55 * \brief layout for elements and vertices
56 *
57 */
58 static Dune::MCMGLayout layout()
59 {
60 return [](Dune::GeometryType gt, int dimgrid) {
61
8/8
✓ Branch 0 taken 56 times.
✓ Branch 1 taken 26 times.
✓ Branch 2 taken 56 times.
✓ Branch 3 taken 26 times.
✓ Branch 4 taken 26 times.
✓ Branch 5 taken 30 times.
✓ Branch 6 taken 26 times.
✓ Branch 7 taken 30 times.
164 return (gt.dim() == dimgrid) || (gt.dim() == 0);
62 26 };
63 }
64 };
65
66 /*!
67 * \ingroup PQ1BubbleDiscretization
68 * \brief The default traits for the pq1bubble finite volume grid geometry
69 * Defines the scv and scvf types and the mapper types
70 * \tparam the grid view type
71 */
72 template<class GridView, class MapperTraits = PQ1BubbleMapperTraits<GridView>>
73 struct PQ1BubbleDefaultGridGeometryTraits
74 : public MapperTraits
75 {
76 using SubControlVolume = PQ1BubbleSubControlVolume<GridView>;
77 using SubControlVolumeFace = PQ1BubbleSubControlVolumeFace<GridView>;
78
79 template<class GridGeometry, bool enableCache>
80 using LocalView = PQ1BubbleFVElementGeometry<GridGeometry, enableCache>;
81 };
82
83 /*!
84 * \ingroup PQ1BubbleDiscretization
85 * \brief Base class for the finite volume geometry vector for pq1bubble schemes
86 * This builds up the sub control volumes and sub control volume faces
87 * \note For caching enabled we store the fv geometries for the whole grid view which is memory intensive but faster
88 */
89 template<class Scalar,
90 class GV,
91 bool enableCaching = true,
92 class Traits = PQ1BubbleDefaultGridGeometryTraits<GV>>
93 class PQ1BubbleFVGridGeometry
94 : public BaseGridGeometry<GV, Traits>
95 {
96 using ThisType = PQ1BubbleFVGridGeometry<Scalar, GV, enableCaching, Traits>;
97 using ParentType = BaseGridGeometry<GV, Traits>;
98 using GridIndexType = typename IndexTraits<GV>::GridIndex;
99 using LocalIndexType = typename IndexTraits<GV>::LocalIndex;
100
101 using Element = typename GV::template Codim<0>::Entity;
102 using CoordScalar = typename GV::ctype;
103 static const int dim = GV::dimension;
104 static const int dimWorld = GV::dimensionworld;
105
106 static_assert(dim > 1, "Only implemented for dim > 1");
107
108 public:
109 //! export the discretization method this geometry belongs to
110 using DiscretizationMethod = DiscretizationMethods::PQ1Bubble;
111 static constexpr DiscretizationMethod discMethod{};
112
113 //! export the type of the fv element geometry (the local view type)
114 using LocalView = typename Traits::template LocalView<ThisType, true>;
115 //! export the type of sub control volume
116 using SubControlVolume = typename Traits::SubControlVolume;
117 //! export the type of sub control volume
118 using SubControlVolumeFace = typename Traits::SubControlVolumeFace;
119 //! export the type of extrusion
120 using Extrusion = Extrusion_t<Traits>;
121 //! export dof mapper type
122 using DofMapper = typename Traits::DofMapper;
123 //! export the finite element cache type
124 using FeCache = Dumux::PQ1BubbleFECache<CoordScalar, Scalar, dim>;
125 //! export the grid view type
126 using GridView = GV;
127
128 //! Constructor
129 13 PQ1BubbleFVGridGeometry(const GridView gridView)
130 : ParentType(gridView)
131 , dofMapper_(gridView, Traits::layout())
132
6/16
✓ Branch 2 taken 13 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 13 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 13 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 13 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 13 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 13 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
26 , cache_(*this)
133 {
134 13 update_();
135 13 }
136
137 //! The dofMapper
138 const DofMapper& dofMapper() const
139
7/12
✓ Branch 0 taken 4 times.
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 4 times.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ 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.
4250945 { return dofMapper_; }
140
141 //! The total number of sub control volumes
142 std::size_t numScv() const
143 { return numScv_; }
144
145 //! The total number of sub control volume faces
146 std::size_t numScvf() const
147 { return numScvf_; }
148
149 //! The total number of boundary sub control volume faces
150 std::size_t numBoundaryScvf() const
151 { return numBoundaryScvf_; }
152
153 //! The total number of degrees of freedom
154 std::size_t numDofs() const
155
9/17
✓ Branch 1 taken 3 times.
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 4 taken 5 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 2 times.
✓ Branch 7 taken 3 times.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 3 times.
✓ Branch 11 taken 2 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
62 { return this->dofMapper().size(); }
156
157 //! update all geometries (call this after grid adaption)
158 void update(const GridView& gridView)
159 {
160 ParentType::update(gridView);
161 update_();
162 }
163
164 //! update all geometries (call this after grid adaption)
165 void update(GridView&& gridView)
166 {
167 ParentType::update(std::move(gridView));
168 update_();
169 }
170
171 //! The finite element cache for creating local FE bases
172 const FeCache& feCache() const
173 17939690 { return feCache_; }
174
175 //! If a vertex / d.o.f. is on the boundary
176 bool dofOnBoundary(GridIndexType dofIdx) const
177
4/8
✓ Branch 0 taken 70876 times.
✓ Branch 1 taken 2273824 times.
✓ Branch 2 taken 70876 times.
✓ Branch 3 taken 2273824 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
4689400 { return boundaryDofIndices_[dofIdx]; }
178
179 //! If a vertex / d.o.f. is on a periodic boundary
180 bool dofOnPeriodicBoundary(GridIndexType dofIdx) const
181
1/4
✗ Branch 1 not taken.
✓ Branch 2 taken 68492 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
68492 { return periodicVertexMap_.count(dofIdx); }
182
183 //! The index of the vertex / d.o.f. on the other side of the periodic boundary
184 GridIndexType periodicallyMappedDof(GridIndexType dofIdx) const
185 { return periodicVertexMap_.at(dofIdx); }
186
187 //! Returns the map between dofs across periodic boundaries
188 const std::unordered_map<GridIndexType, GridIndexType>& periodicDofMap() const
189 24 { return periodicVertexMap_; }
190
191 //! Returns the map between dofs across periodic boundaries
192 [[deprecated("Will be removed after release 3.9. Use periodicDofMap() instead.")]]
193 const std::unordered_map<GridIndexType, GridIndexType>& periodicVertexMap() const
194 { return periodicDofMap(); }
195
196 //! local view of this object (constructed with the internal cache)
197 friend inline LocalView localView(const PQ1BubbleFVGridGeometry& gg)
198
22/30
✓ Branch 1 taken 11 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 5 times.
✓ Branch 4 taken 11 times.
✓ Branch 5 taken 1 times.
✓ Branch 6 taken 5 times.
✓ Branch 7 taken 42002 times.
✓ Branch 8 taken 1 times.
✓ Branch 9 taken 6 times.
✓ Branch 10 taken 42002 times.
✓ Branch 11 taken 1 times.
✓ Branch 12 taken 6 times.
✓ Branch 13 taken 903 times.
✓ Branch 14 taken 1 times.
✓ Branch 15 taken 4801 times.
✓ Branch 16 taken 903 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 4801 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✓ Branch 21 taken 21600 times.
✗ Branch 22 not taken.
✓ Branch 23 taken 1 times.
✓ Branch 24 taken 21600 times.
✗ Branch 25 not taken.
✓ Branch 26 taken 1 times.
✓ Branch 27 taken 21604 times.
✗ Branch 28 not taken.
✓ Branch 30 taken 21604 times.
✗ Branch 31 not taken.
1450334 { 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 13 explicit PQ1BubbleGridGeometryCache(const PQ1BubbleFVGridGeometry& gg)
210
5/10
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 13 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 13 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 13 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 13 times.
✗ Branch 14 not taken.
65 : gridGeometry_(&gg)
211 {}
212
213 const PQ1BubbleFVGridGeometry& gridGeometry() const
214 { return *gridGeometry_; }
215
216 //! Get the global sub control volume indices of an element
217 const std::vector<SubControlVolume>& scvs(GridIndexType eIdx) const
218
16/19
✓ Branch 2 taken 7405 times.
✓ Branch 3 taken 834287 times.
✓ Branch 4 taken 7405 times.
✓ Branch 5 taken 834287 times.
✓ Branch 6 taken 64230 times.
✓ Branch 7 taken 256554 times.
✓ Branch 8 taken 56830 times.
✓ Branch 9 taken 216514 times.
✓ Branch 10 taken 32640 times.
✓ Branch 11 taken 43720 times.
✓ Branch 12 taken 14800 times.
✓ Branch 13 taken 14400 times.
✓ Branch 14 taken 58520 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 29200 times.
✓ Branch 17 taken 21600 times.
✗ Branch 18 not taken.
✓ Branch 20 taken 21600 times.
✗ Branch 21 not taken.
285022199 { return scvs_[eIdx]; }
219
220 //! Get the global sub control volume face indices of an element
221 const std::vector<SubControlVolumeFace>& scvfs(GridIndexType eIdx) const
222
2/8
✗ Branch 0 not taken.
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 12000 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
10016671 { return scvfs_[eIdx]; }
223
224 //! Returns whether one of the geometry's scvfs lies on a boundary
225 bool hasBoundaryScvf(GridIndexType eIdx) const
226
4/8
✓ Branch 0 taken 1362 times.
✓ Branch 1 taken 41413 times.
✓ Branch 2 taken 1362 times.
✓ Branch 3 taken 41413 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
85550 { return hasBoundaryScvf_[eIdx]; }
227
228 //! Local mappings necessary to construct geometries of scvfs
229 const std::vector<std::array<LocalIndexType, 2>>& scvfBoundaryGeometryKeys(GridIndexType eIdx) const
230 664 { return scvfBoundaryGeometryKeys_.at(eIdx); }
231
232 private:
233 13 void clear_()
234 {
235 13 scvs_.clear();
236 13 scvfs_.clear();
237 13 hasBoundaryScvf_.clear();
238 13 scvfBoundaryGeometryKeys_.clear();
239 13 }
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 13 void update_()
258 {
259 13 cache_.clear_();
260 26 dofMapper_.update(this->gridView());
261
262 26 auto numElements = this->gridView().size(0);
263 13 cache_.scvs_.resize(numElements);
264 13 cache_.scvfs_.resize(numElements);
265 13 cache_.hasBoundaryScvf_.resize(numElements, false);
266
267 23 boundaryDofIndices_.assign(numDofs(), false);
268
269 13 numScv_ = 0;
270 13 numScvf_ = 0;
271 13 numBoundaryScvf_ = 0;
272
273 // Build the scvs and scv faces
274
10/12
✓ Branch 2 taken 82318 times.
✓ Branch 3 taken 8 times.
✓ Branch 4 taken 61212 times.
✓ Branch 5 taken 5 times.
✓ Branch 6 taken 3 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 12000 times.
✓ Branch 9 taken 3 times.
✓ Branch 10 taken 12000 times.
✓ Branch 11 taken 3 times.
✓ Branch 13 taken 12000 times.
✗ Branch 14 not taken.
188652 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.
188626 auto eIdx = this->elementMapper().index(element);
277
278 // get the element geometry
279
2/4
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 12000 times.
✗ Branch 5 not taken.
106313 auto elementGeometry = element.geometry();
280
1/2
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
94313 const auto refElement = referenceElement(elementGeometry);
281
282 // instantiate the geometry helper
283
1/2
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
94313 GeometryHelper geometryHelper(elementGeometry);
284
285
1/2
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
94313 numScv_ += geometryHelper.numScv();
286 // construct the sub control volumes
287
3/6
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 12000 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 12000 times.
✗ Branch 8 not taken.
188626 cache_.scvs_[eIdx].resize(geometryHelper.numScv());
288
4/4
✓ Branch 1 taken 444021 times.
✓ Branch 2 taken 82313 times.
✓ Branch 3 taken 49600 times.
✓ Branch 4 taken 12000 times.
526334 for (LocalIndexType scvLocalIdx = 0; scvLocalIdx < geometryHelper.numScv(); ++scvLocalIdx)
289 {
290
1/2
✓ Branch 1 taken 49600 times.
✗ Branch 2 not taken.
432021 auto corners = geometryHelper.getScvCorners(scvLocalIdx);
291
4/8
✓ Branch 1 taken 49600 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 49600 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 49600 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 49600 times.
✗ Branch 11 not taken.
2525919 cache_.scvs_[eIdx][scvLocalIdx] = SubControlVolume(
292 geometryHelper.scvVolume(scvLocalIdx, corners),
293 geometryHelper.dofPosition(scvLocalIdx),
294 864042 Dumux::center(corners),
295 scvLocalIdx,
296 eIdx,
297 geometryHelper.dofIndex(this->dofMapper(), element, scvLocalIdx),
298 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.
94313 numScvf_ += geometryHelper.numInteriorScvf();
304
3/6
✓ Branch 1 taken 12000 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 12000 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 12000 times.
✗ Branch 8 not taken.
188626 cache_.scvfs_[eIdx].resize(geometryHelper.numInteriorScvf());
305 94313 LocalIndexType scvfLocalIdx = 0;
306
4/4
✓ Branch 1 taken 696244 times.
✓ Branch 2 taken 82313 times.
✓ Branch 3 taken 75200 times.
✓ Branch 4 taken 12000 times.
778557 for (; scvfLocalIdx < geometryHelper.numInteriorScvf(); ++scvfLocalIdx)
307 {
308
1/2
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
684244 const auto scvPair = geometryHelper.getScvPairForScvf(scvfLocalIdx);
309
1/2
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
684244 const auto corners = geometryHelper.getScvfCorners(scvfLocalIdx);
310
1/2
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
684264 const auto area = Dumux::convexPolytopeVolume<dim-1>(
311 geometryHelper.getInteriorScvfGeometryType(scvfLocalIdx),
312 353120 [&](unsigned int i){ return corners[i]; }
313 );
314
315
1/2
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
2052732 cache_.scvfs_[eIdx][scvfLocalIdx] = SubControlVolumeFace(
316 684244 Dumux::center(corners),
317 area,
318 geometryHelper.normal(corners, scvPair),
319
1/2
✓ Branch 1 taken 75200 times.
✗ Branch 2 not taken.
684244 std::move(scvPair),
320 scvfLocalIdx,
321 geometryHelper.isOverlappingScvf(scvfLocalIdx)
322 );
323 }
324
325 // construct the sub control volume faces on the domain boundary
326
9/15
✓ Branch 1 taken 12000 times.
✓ Branch 2 taken 21101 times.
✓ Branch 3 taken 145618 times.
✓ Branch 4 taken 12000 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 61212 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 326512 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✓ Branch 14 taken 53912 times.
✓ Branch 15 taken 199388 times.
✓ Branch 17 taken 215700 times.
✗ Branch 18 not taken.
947157 for (const auto& intersection : intersections(this->gridView(), element))
327 {
328
5/6
✓ Branch 0 taken 100718 times.
✓ Branch 1 taken 199968 times.
✓ Branch 2 taken 37813 times.
✓ Branch 3 taken 5397 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 5184 times.
348500 if (intersection.boundary() && !intersection.neighbor())
329 {
330
5/6
✓ Branch 0 taken 1368 times.
✓ Branch 1 taken 3816 times.
✓ Branch 2 taken 1368 times.
✓ Branch 3 taken 722 times.
✓ Branch 4 taken 3094 times.
✗ Branch 5 not taken.
13540 cache_.hasBoundaryScvf_[eIdx] = true;
331
332
2/3
✓ Branch 0 taken 1368 times.
✓ Branch 1 taken 3816 times.
✗ Branch 2 not taken.
6770 const auto localFacetIndex = intersection.indexInInside();
333
1/2
✓ Branch 1 taken 5184 times.
✗ Branch 2 not taken.
6770 const auto numBoundaryScvf = geometryHelper.numBoundaryScvf(localFacetIndex);
334 6770 numScvf_ += numBoundaryScvf;
335 6770 numBoundaryScvf_ += numBoundaryScvf;
336
337
2/2
✓ Branch 0 taken 16066 times.
✓ Branch 1 taken 6770 times.
22836 for (unsigned int isScvfLocalIdx = 0; isScvfLocalIdx < numBoundaryScvf; ++isScvfLocalIdx)
338 {
339 // find the scvs this scvf is belonging to
340
1/2
✓ Branch 1 taken 12882 times.
✗ Branch 2 not taken.
16066 const auto scvPair = geometryHelper.getScvPairForBoundaryScvf(localFacetIndex, isScvfLocalIdx);
341
1/2
✓ Branch 1 taken 12882 times.
✗ Branch 2 not taken.
16066 const auto corners = geometryHelper.getBoundaryScvfCorners(localFacetIndex, isScvfLocalIdx);
342 16066 const auto area = Dumux::convexPolytopeVolume<dim-1>(
343 geometryHelper.getBoundaryScvfGeometryType(isScvfLocalIdx),
344 60528 [&](unsigned int i){ return corners[i]; }
345 );
346
4/7
✓ Branch 1 taken 11722 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1160 times.
✓ Branch 4 taken 11722 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 11722 times.
✗ Branch 8 not taken.
48198 cache_.scvfs_[eIdx].emplace_back(
347 16066 Dumux::center(corners),
348 area,
349 intersection.centerUnitOuterNormal(),
350
1/2
✓ Branch 1 taken 12882 times.
✗ Branch 2 not taken.
16066 std::move(scvPair),
351 scvfLocalIdx,
352 typename SubControlVolumeFace::Traits::BoundaryFlag{ intersection },
353 geometryHelper.isOverlappingBoundaryScvf(localFacetIndex)
354 );
355
356 // store look-up map to construct boundary scvf geometries
357
2/4
✓ Branch 1 taken 12882 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 12882 times.
✗ Branch 5 not taken.
16066 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 16066 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 16066 times.
✓ Branch 1 taken 6770 times.
22836 for (int localVIdx = 0; localVIdx < numBoundaryScvf; ++localVIdx)
370 {
371 16066 const auto vIdx = refElement.subEntity(localFacetIndex, 1, localVIdx, dim);
372
1/2
✓ Branch 1 taken 12882 times.
✗ Branch 2 not taken.
16066 const auto vIdxGlobal = this->dofMapper().subIndex(element, vIdx, dim);
373 48198 boundaryDofIndices_[vIdxGlobal] = true;
374 }
375 }
376
377 // inform the grid geometry if we have periodic boundaries
378
2/6
✓ Branch 0 taken 94528 times.
✓ Branch 1 taken 236408 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
330936 else if (intersection.boundary() && intersection.neighbor())
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 (isOutside.boundary() && isOutside.neighbor())
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 periodicVertexMap_[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/10
✗ Branch 0 not taken.
✓ Branch 1 taken 13 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
13 if (this->isPeriodic() && this->gridView().comm().size() > 1)
418 DUNE_THROW(Dune::NotImplemented, "Periodic boundaries for pq1bubble method for parallel simulations!");
419 13 }
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> periodicVertexMap_;
434
435 Cache cache_;
436 };
437
438 } // end namespace Dumux
439
440 #endif
441