GCC Code Coverage Report

Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/dumux/io/grid/porenetwork/parametersforgeneratedgrid.hh
Date: 2024-05-04 19:09:25
Exec Total Coverage
Lines: 223 255 87.5%
Functions: 52 60 86.7%
Branches: 603 1640 36.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 PoreNetworkModels
10 * \brief Helper class to assign parameters to a generated grid
11 */
12 #ifndef DUMUX_IO_PARAMETERS_FOR_GENERATED_GRID
13 #define DUMUX_IO_PARAMETERS_FOR_GENERATED_GRID
14
15 #include <algorithm>
16 #include <array>
17 #include <numeric>
18 #include <random>
19 #include <string>
20 #include <vector>
21
22 #include <dune/common/exceptions.hh>
23 #include <dune/grid/common/exceptions.hh>
24 #include <dune/geometry/axisalignedcubegeometry.hh>
25
26 #include <dumux/common/parameters.hh>
27 #include <dumux/common/random.hh>
28 #include <dumux/common/stringutilities.hh>
29 #include <dumux/geometry/intersectspointgeometry.hh>
30 #include <dumux/porenetwork/common/throatproperties.hh>
31 #include <dumux/porenetwork/common/poreproperties.hh>
32
33 namespace Dumux::PoreNetwork {
34
35 /*!
36 * \ingroup PoreNetworkModels
37 * \brief Helper class to assign parameters to a generated grid
38 */
39 template <class Grid, class Scalar>
40 ✗ class ParametersForGeneratedGrid
41 {
42 using GridView = typename Grid::LeafGridView;
43 using Element = typename Grid::template Codim<0>::Entity;
44 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
45 using Vertex = typename Grid::template Codim<Grid::dimension>::Entity;
46
47 static constexpr auto dim = Grid::dimension;
48 static constexpr auto dimWorld = Grid::dimensionworld;
49 using BoundaryList = std::array<int, 2*dimWorld>;
50
51 public:
52
53 33 ParametersForGeneratedGrid(const GridView& gridView, const std::string& paramGroup)
54 : gridView_(gridView)
55 , paramGroup_(paramGroup)
56
1/5
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
 99 , priorityList_(getPriorityList_())
57 {
58
1/2
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
 33 computeBoundingBox_();
59
1/2
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
 33 boundaryFaceIndex_ = getBoundaryFacemarkerInput_();
60 33 }
61
62 /*!
63 * \brief Returns the boundary face marker index at given position
64 *
65 * \param pos The current position
66 */
67 int boundaryFaceMarkerAtPos(const GlobalPosition& pos) const
68 {
69 // set the priority which decides the order the vertices on the boundary are indexed
70 // by default, vertices on min/max faces in x direction have the highest priority, followed by y and z
71
6/12
✓ Branch 0 taken 78 times.
✓ Branch 1 taken 36 times.
✓ Branch 2 taken 88616 times.
✓ Branch 3 taken 12659 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✓ Branch 8 taken 568 times.
✓ Branch 9 taken 72 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 118980 for (auto boundaryIdx : priorityList_)
72 {
73
6/12
✓ Branch 0 taken 4 times.
✓ Branch 1 taken 74 times.
✓ Branch 2 taken 4058 times.
✓ Branch 3 taken 84558 times.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✓ Branch 8 taken 122 times.
✓ Branch 9 taken 446 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 89262 if (onBoundary_(pos, boundaryIdx))
74 8368 return boundaryFaceIndex_[boundaryIdx];
75 }
76 return -1;
77 }
78
79 /*!
80 * \brief Computes and returns the label of a given throat
81 *
82 * \param poreLabels The labels of the pores adjacent to the throat.
83 */
84 44289 int throatLabel(const std::array<int, 2>& poreLabels) const
85 {
86
6/6
✓ Branch 0 taken 9795 times.
✓ Branch 1 taken 34040 times.
✓ Branch 2 taken 9795 times.
✓ Branch 3 taken 34040 times.
✓ Branch 4 taken 9795 times.
✓ Branch 5 taken 34040 times.
 132867 if (poreLabels[0] == poreLabels[1]) // both vertices are inside the domain or on the same boundary face
87 return poreLabels[0];
88
4/4
✓ Branch 0 taken 5385 times.
✓ Branch 1 taken 4410 times.
✓ Branch 2 taken 5385 times.
✓ Branch 3 taken 4410 times.
 20042 if (poreLabels[0] == -1) // vertex1 is inside the domain, vertex2 is on a boundary face
89 return poreLabels[1];
90
4/4
✓ Branch 0 taken 898 times.
✓ Branch 1 taken 4487 times.
✓ Branch 2 taken 898 times.
✓ Branch 3 taken 4487 times.
 11054 if (poreLabels[1] == -1) // vertex2 is inside the domain, vertex1 is on a boundary face
91 return poreLabels[0];
92
93 // use the priority list to find out which pore label is favored
94
1/2
✓ Branch 0 taken 1182 times.
✗ Branch 1 not taken.
 3104 for(const auto i : priorityList_)
95 {
96
6/6
✓ Branch 0 taken 700 times.
✓ Branch 1 taken 482 times.
✓ Branch 2 taken 700 times.
✓ Branch 3 taken 482 times.
✓ Branch 4 taken 700 times.
✓ Branch 5 taken 482 times.
 3708 if (poreLabels[0] == boundaryFaceIndex_[i])
97 return poreLabels[0];
98
6/6
✓ Branch 0 taken 284 times.
✓ Branch 1 taken 416 times.
✓ Branch 2 taken 284 times.
✓ Branch 3 taken 416 times.
✓ Branch 4 taken 284 times.
✓ Branch 5 taken 416 times.
 2196 if (poreLabels[1] == boundaryFaceIndex_[i])
99 return poreLabels[1];
100 }
101
102 ✗ DUNE_THROW(Dune::InvalidStateException, "Something went wrong with the throat labels");
103 }
104
105 /*!
106 * \brief Assign parameters for generically created grids
107 */
108 template <class SetParameter, class GetParameter>
109 33 void assignParameters(const SetParameter& setParameter,
110 const GetParameter& getParameter,
111 const std::size_t numSubregions)
112 {
113 using cytpe = typename GlobalPosition::value_type;
114 using InternalBoundingBox = Dune::AxisAlignedCubeGeometry<cytpe, dimWorld, dimWorld>;
115
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 20 times.
 33 std::vector<InternalBoundingBox> internalBoundingBoxes;
116
117 // divide the network into subregions, if specified
118
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 20 times.
 33 if (numSubregions > 0)
119 {
120 // get bounding boxes of subregions
121
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
 2 for (int i = 0; i < numSubregions; ++i)
122 {
123
7/20
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 1 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 1 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 1 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 1 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
 3 auto lowerLeft = getParamFromGroup<GlobalPosition>(paramGroup_, "Grid.Subregion" + std::to_string(i) + ".LowerLeft");
124
7/22
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 1 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 1 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✓ Branch 15 taken 1 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 1 times.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
 3 auto upperRight = getParamFromGroup<GlobalPosition>(paramGroup_, "Grid.Subregion" + std::to_string(i) + ".UpperRight");
125
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 1 internalBoundingBoxes.emplace_back(std::move(lowerLeft), std::move(upperRight));
126 }
127 }
128
129 // get the maximum possible pore body radii such that pore bodies do not intersect
130 // (requires ThroatRegionId, if specified)
131
3/4
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1 times.
✓ Branch 4 taken 20 times.
 66 const std::vector<Scalar> maxPoreRadius = getMaxPoreRadii_(numSubregions, getParameter);
132
3/8
✓ Branch 3 taken 21 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 21 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 21 times.
✗ Branch 9 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
 99 std::vector<bool> poreRadiusLimited(gridView_.size(dim), false);
133
134 // get a helper function for getting the pore radius of a pore body not belonging to a subregion
135
1/4
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 66 auto defaultPoreRadius = poreRadiusGenerator_(-1);
136
137 // get helper functions for pore body radii on subregions
138
1/2
✓ Branch 0 taken 21 times.
✗ Branch 1 not taken.
 66 std::vector<decltype(defaultPoreRadius)> subregionPoreRadius;
139
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 21 times.
 34 for (int i = 0; i < numSubregions; ++i)
140
3/8
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 1 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 1 subregionPoreRadius.emplace_back(poreRadiusGenerator_(i));
141
142 // get helper function for pore volume
143
1/2
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
 33 const auto poreVolume = poreVolumeGenerator_(getParameter);
144
145 // treat the pore body parameters (label, radius and maybe regionId)
146
5/6
✓ Branch 1 taken 16951 times.
✓ Branch 2 taken 21 times.
✓ Branch 3 taken 16951 times.
✓ Branch 4 taken 21 times.
✓ Branch 6 taken 16951 times.
✗ Branch 7 not taken.
 34679 for (const auto& vertex : vertices(gridView_))
147 {
148
2/4
✓ Branch 1 taken 16951 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 16951 times.
 17323 const auto& pos = vertex.geometry().center();
149
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 16951 times.
 17323 const auto vIdxGlobal = gridView_.indexSet().index(vertex);
150 17323 const auto poreLabel = boundaryFaceMarkerAtPos(pos);
151
6/12
✓ Branch 1 taken 16951 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 16951 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 16951 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 16951 times.
✓ Branch 11 taken 15951 times.
✓ Branch 12 taken 1000 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
 34646 setParameter(vertex, "PoreLabel", poreLabel);
152
153 // sets the minimum of the given value and the maximum possible pore body radius
154 // and keeps track of capped radii
155 17820 auto setRadiusAndLogIfCapped = [&](const Scalar value)
156 {
157
8/12
✓ Branch 0 taken 5188 times.
✓ Branch 1 taken 11556 times.
✓ Branch 2 taken 5188 times.
✓ Branch 3 taken 11556 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 138 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 138 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 194 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 194 times.
 34152 if (value > maxPoreRadius[vIdxGlobal])
158 {
159
2/12
✓ Branch 1 taken 5188 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 5188 times.
✗ Branch 5 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
 10376 poreRadiusLimited[vIdxGlobal] = true;
160
4/30
✓ Branch 1 taken 5188 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 5188 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 5188 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 5188 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✗ Branch 34 not taken.
✗ Branch 35 not taken.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
 20752 setParameter(vertex, "PoreInscribedRadius", maxPoreRadius[vIdxGlobal]);
161 }
162 else
163
12/30
✓ Branch 1 taken 11556 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 11556 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 11556 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 11556 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 14 taken 138 times.
✗ Branch 15 not taken.
✓ Branch 17 taken 138 times.
✗ Branch 18 not taken.
✓ Branch 20 taken 138 times.
✗ Branch 21 not taken.
✓ Branch 22 taken 138 times.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✓ Branch 27 taken 194 times.
✗ Branch 28 not taken.
✓ Branch 30 taken 194 times.
✗ Branch 31 not taken.
✓ Branch 33 taken 194 times.
✗ Branch 34 not taken.
✓ Branch 35 taken 194 times.
✗ Branch 36 not taken.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
 47552 setParameter(vertex, "PoreInscribedRadius", value);
164 };
165
166
2/2
✓ Branch 0 taken 15951 times.
✓ Branch 1 taken 1000 times.
 17323 if (numSubregions == 0) // assign radius if no subregions are specified
167
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 15951 times.
✓ Branch 3 taken 15951 times.
✗ Branch 4 not taken.
 32646 setRadiusAndLogIfCapped(defaultPoreRadius(vertex, poreLabel));
168 else // assign region ids and radii to vertices if they are within a subregion
169 {
170 // default value for vertices not belonging to a subregion
171
5/12
✓ Branch 1 taken 1000 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1000 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1000 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 1000 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 1000 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
 2000 setParameter(vertex, "PoreRegionId", -1);
172
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 1000 times.
✓ Branch 3 taken 1000 times.
✗ Branch 4 not taken.
 2000 setRadiusAndLogIfCapped(defaultPoreRadius(vertex, poreLabel));
173
174
2/2
✓ Branch 0 taken 1000 times.
✓ Branch 1 taken 1000 times.
 2000 for (int id = 0; id < numSubregions; ++id)
175 {
176
1/2
✓ Branch 1 taken 1000 times.
✗ Branch 2 not taken.
 1000 const auto& subregion = internalBoundingBoxes[id];
177
3/8
✓ Branch 1 taken 1000 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 125 times.
✓ Branch 6 taken 875 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
 1000 if (intersectsPointGeometry(vertex.geometry().center(), subregion))
178 {
179
5/12
✓ Branch 1 taken 125 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 125 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 125 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 125 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 125 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
 250 setParameter(vertex, "PoreRegionId", id);
180
3/6
✗ Branch 0 not taken.
✓ Branch 1 taken 125 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 125 times.
✓ Branch 5 taken 125 times.
✗ Branch 6 not taken.
 375 setRadiusAndLogIfCapped(subregionPoreRadius[id](vertex, poreLabel));
181 }
182 }
183 }
184
185
5/12
✓ Branch 1 taken 16951 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 16951 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 16951 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 16951 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 16951 times.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
 17323 setParameter(vertex, "PoreVolume", poreVolume(vertex, poreLabel));
186 }
187
188 // treat throat parameters
189
1/2
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
 33 auto defaultThroatInscribedRadius = throatInscribedRadiusGenerator_(-1, getParameter);
190
1/2
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
 33 auto defaultThroatLength = throatLengthGenerator_(-1, getParameter);
191
192 // get helper functions for throat radii and lengths on subregions
193
0/2
✗ Branch 0 not taken.
✗ Branch 1 not taken.
 33 std::vector<decltype(defaultThroatInscribedRadius)> subregionThroatInscribedRadius;
194
4/4
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 20 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 20 times.
 67 std::vector<decltype(defaultThroatLength)> subregionThroatLength;
195
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 21 times.
 34 for (int i = 0; i < numSubregions; ++i)
196 {
197
2/4
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
 1 subregionThroatInscribedRadius.emplace_back(throatInscribedRadiusGenerator_(i, getParameter));
198
2/6
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 1 subregionThroatLength.emplace_back(throatLengthGenerator_(i, getParameter));
199 }
200
201 // set throat parameters
202
6/6
✓ Branch 1 taken 43835 times.
✓ Branch 2 taken 21 times.
✓ Branch 3 taken 43835 times.
✓ Branch 4 taken 21 times.
✓ Branch 5 taken 33359 times.
✓ Branch 6 taken 10476 times.
 88611 for (const auto& element : elements(gridView_))
203 {
204
2/2
✓ Branch 0 taken 33359 times.
✓ Branch 1 taken 10476 times.
 44289 if (numSubregions == 0) // assign values if no subregions are specified
205 {
206
6/14
✓ Branch 1 taken 33359 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 33359 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 33359 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 33359 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 33359 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 33359 times.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 67626 setParameter(element, "ThroatInscribedRadius", defaultThroatInscribedRadius(element));
207
5/12
✓ Branch 1 taken 33359 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 33359 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 33359 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 33359 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 33359 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
 33813 setParameter(element, "ThroatLength", defaultThroatLength(element));
208 }
209 else // assign values to throats if they are within a subregion
210 {
211 // default value for elements not belonging to a subregion
212
5/12
✓ Branch 1 taken 10476 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 10476 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 10476 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 10476 times.
✓ Branch 12 taken 10476 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
 20952 setParameter(element, "ThroatRegionId", -1);
213
6/14
✓ Branch 1 taken 10476 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 10476 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 10476 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 10476 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 10476 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 10476 times.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 20952 setParameter(element, "ThroatInscribedRadius", defaultThroatInscribedRadius(element));
214
5/12
✓ Branch 1 taken 10476 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 10476 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 10476 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 10476 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 10476 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
 10476 setParameter(element, "ThroatLength", defaultThroatLength(element));
215
216
2/2
✓ Branch 0 taken 10476 times.
✓ Branch 1 taken 10476 times.
 20952 for (int id = 0; id < numSubregions; ++id)
217 {
218
1/2
✓ Branch 1 taken 10476 times.
✗ Branch 2 not taken.
 10476 const auto& subregion = internalBoundingBoxes[id];
219
3/8
✓ Branch 1 taken 10476 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 1625 times.
✓ Branch 5 taken 8851 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
 20952 if (intersectsPointGeometry(element.geometry().center(), subregion))
220 {
221
5/12
✓ Branch 1 taken 1625 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1625 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1625 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 1625 times.
✓ Branch 12 taken 1625 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
 3250 setParameter(element, "ThroatRegionId", id);
222
7/16
✓ Branch 1 taken 1625 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1625 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1625 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 1625 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 1625 times.
✗ Branch 14 not taken.
✓ Branch 15 taken 1625 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 1625 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
 4875 setParameter(element, "ThroatInscribedRadius", subregionThroatInscribedRadius[id](element));
223
6/14
✓ Branch 1 taken 1625 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1625 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1625 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 1625 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 1625 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✓ Branch 16 taken 1625 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 3250 setParameter(element, "ThroatLength", subregionThroatLength[id](element));
224 }
225 }
226 }
227
228 // set the throat label
229
1/2
✓ Branch 1 taken 43835 times.
✗ Branch 2 not taken.
 44289 const auto vertex0 = element.template subEntity<dim>(0);
230
1/2
✓ Branch 1 taken 43835 times.
✗ Branch 2 not taken.
 44289 const auto vertex1 = element.template subEntity<dim>(1);
231
5/12
✓ Branch 1 taken 43835 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 43835 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 43835 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 43835 times.
✓ Branch 12 taken 43835 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
 88578 const std::array poreLabels{static_cast<int>(getParameter(vertex0, "PoreLabel")),
232
5/12
✓ Branch 1 taken 43835 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 43835 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 43835 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 43835 times.
✓ Branch 12 taken 43835 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
 88578 static_cast<int>(getParameter(vertex1, "PoreLabel"))};
233
5/12
✓ Branch 1 taken 43835 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 43835 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 43835 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 43835 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 43835 times.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
 44289 setParameter(element, "ThroatLabel", throatLabel(poreLabels));
234 }
235
236 99 const auto numPoreRadiusLimited = std::count(poreRadiusLimited.begin(), poreRadiusLimited.end(), true);
237
2/2
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 19 times.
 33 if (numPoreRadiusLimited > 0)
238
3/6
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
✓ Branch 6 taken 2 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 2 times.
✗ Branch 10 not taken.
 8 std::cout << "*******\nWarning! " << numPoreRadiusLimited << " out of " << poreRadiusLimited.size()
239
1/2
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
 4 << " pore body radii have been capped automatically in order to prevent intersecting pores\n*******" << std::endl;
240 33 }
241
242 private:
243
244 /*!
245 * \brief Returns a list of boundary face priorities from user specified input or default values if no input is given
246 *
247 * This essentially determines the index of a node on an edge or corner. For instance in a 2D case, a list of {0,1,2,3}
248 * will give highest priority to the "x-direction" and lowest to the "diagonal-direction".
249 */
250 BoundaryList getPriorityList_() const
251 {
252 33 const auto list = [&]()
253 {
254 BoundaryList priorityList;
255 99 std::iota(priorityList.begin(), priorityList.end(), 0);
256
257
7/12
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 21 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 21 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 3 times.
✓ Branch 9 taken 18 times.
✓ Branch 10 taken 3 times.
✓ Branch 11 taken 18 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
 165 if (hasParamInGroup(paramGroup_, "Grid.PriorityList"))
258 {
259 try {
260 // priorities can also be set in the input file
261 3 priorityList = getParamFromGroup<BoundaryList>(paramGroup_, "Grid.PriorityList");
262 }
263 // make sure that a priority for each direction is set
264 ✗ catch(Dune::RangeError& e) {
265 ✗ DUNE_THROW(Dumux::ParameterException, "You must specify priorities for all directions (" << dimWorld << ") \n" << e.what());
266 }
267 // make sure each direction is only set once
268
0/2
✗ Branch 0 not taken.
✗ Branch 1 not taken.
 6 if (!isUnique_(priorityList))
269 ✗ DUNE_THROW(Dumux::ParameterException, "You must specify priorities for all directions (duplicate directions)");
270
271 //make sure that the directions are correct (ranging from 0 to dimWorld-1)
272
10/42
✗ Branch 0 not taken.
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 3 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 3 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 3 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 3 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 3 times.
✗ Branch 14 not taken.
✓ Branch 15 taken 3 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 3 times.
✗ Branch 18 not taken.
✓ Branch 19 taken 3 times.
✗ Branch 20 not taken.
✓ Branch 21 taken 3 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
✗ Branch 28 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
✗ Branch 31 not taken.
✗ Branch 32 not taken.
✗ Branch 33 not taken.
✗ Branch 34 not taken.
✗ Branch 35 not taken.
✗ Branch 36 not taken.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
 21 if (std::any_of(priorityList.begin(), priorityList.end(), []( const int i ){ return (i < 0 || i >= 2*dimWorld); }))
273 ✗ DUNE_THROW(Dumux::ParameterException, "You must specify priorities for correct directions (0-" << 2*dimWorld-1 << ")");
274 }
275 33 return priorityList;
276
4/7
✓ Branch 1 taken 23 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 5 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 15 times.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
 57 }();
277
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
 12 return list;
278 }
279
280 33 void computeBoundingBox_()
281 {
282 // calculate the bounding box of the local partition of the grid view
283
4/4
✓ Branch 1 taken 17105 times.
✓ Branch 2 taken 21 times.
✓ Branch 3 taken 17105 times.
✓ Branch 4 taken 21 times.
 17510 for (const auto& vertex : vertices(gridView_))
284 {
285
2/2
✓ Branch 0 taken 51030 times.
✓ Branch 1 taken 17105 times.
 69349 for (int i = 0; i < dimWorld; i++)
286 {
287 using std::min;
288 using std::max;
289
8/8
✓ Branch 1 taken 63 times.
✓ Branch 2 taken 50967 times.
✓ Branch 3 taken 63 times.
✓ Branch 4 taken 50967 times.
✓ Branch 5 taken 63 times.
✓ Branch 6 taken 50967 times.
✓ Branch 7 taken 61 times.
✓ Branch 8 taken 50929 times.
 51963 bBoxMin_[i] = min(bBoxMin_[i], vertex.geometry().corner(0)[i]);
290
8/8
✓ Branch 1 taken 326 times.
✓ Branch 2 taken 50704 times.
✓ Branch 3 taken 326 times.
✓ Branch 4 taken 50704 times.
✓ Branch 5 taken 326 times.
✓ Branch 6 taken 50704 times.
✓ Branch 7 taken 288 times.
✓ Branch 8 taken 50702 times.
 52700 bBoxMax_[i] = max(bBoxMax_[i], vertex.geometry().corner(0)[i]);
291 }
292 }
293 33 }
294
295 /*!
296 * \brief Returns a list of boundary face indices from user specified input or default values if no input is given
297 */
298 33 BoundaryList getBoundaryFacemarkerInput_() const
299 {
300 BoundaryList boundaryFaceMarker;
301 99 std::fill(boundaryFaceMarker.begin(), boundaryFaceMarker.end(), 0);
302
1/2
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
 33 boundaryFaceMarker[0] = 1;
303
1/2
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
 33 boundaryFaceMarker[1] = 1;
304
305
8/14
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 21 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 21 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 21 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 19 times.
✓ Branch 12 taken 2 times.
✓ Branch 13 taken 19 times.
✓ Branch 14 taken 2 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
 165 if (hasParamInGroup(paramGroup_, "Grid.BoundaryPoreLabels"))
306 {
307 62 const auto input = getParamFromGroup<std::vector<std::string>>(paramGroup_, "Grid.BoundaryPoreLabels");
308
4/4
✓ Branch 0 taken 88 times.
✓ Branch 1 taken 19 times.
✓ Branch 2 taken 88 times.
✓ Branch 3 taken 19 times.
 365 for (const auto& entry : input)
309 {
310
2/4
✓ Branch 1 taken 88 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 88 times.
✗ Branch 5 not taken.
 408 const std::string errorMessage = "You must specify BoundaryPoreLabels in the format pos:num, where pos can be xMin, xMax, yMin, yMax, zMin, zMax and num is the corresponding label.\n"
311 "Example (2D, defaults are used for the remaining boundaries): xMin:2 yMax:3\n";
312
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 88 times.
 136 if (entry.find(':') == std::string::npos)
313 ✗ DUNE_THROW(Dumux::ParameterException, errorMessage);
314
315
14/28
✓ Branch 0 taken 19 times.
✓ Branch 1 taken 69 times.
✓ Branch 3 taken 19 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 19 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 19 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 19 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 19 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 19 times.
✗ Branch 19 not taken.
✓ Branch 21 taken 19 times.
✗ Branch 22 not taken.
✗ Branch 24 not taken.
✓ Branch 25 taken 19 times.
✗ Branch 27 not taken.
✓ Branch 28 taken 19 times.
✓ Branch 31 taken 114 times.
✓ Branch 32 taken 19 times.
✗ Branch 33 not taken.
✓ Branch 34 taken 114 times.
✗ Branch 35 not taken.
✗ Branch 36 not taken.
✗ Branch 37 not taken.
✗ Branch 38 not taken.
 322 static const std::map<std::string, int> labels = {{"xMin", 0}, {"xMax", 1}, {"yMin", 2}, {"yMax", 3}, {"zMin", 4}, {"zMax", 5}};
316
3/10
✓ Branch 2 taken 88 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 88 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 88 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
 272 const auto splitEntry = split(entry, ":");
317
5/12
✓ Branch 1 taken 88 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 88 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 88 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 88 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 88 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
 544 const std::string location = std::string(splitEntry[0].begin(), splitEntry[0].end());
318 136 int value = 0;
319 try {
320
7/16
✓ Branch 1 taken 88 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 88 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 88 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 88 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 88 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✓ Branch 16 taken 88 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 88 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
 544 value = std::stoi(std::string(splitEntry[1].begin(), splitEntry[1].end()));
321 }
322 ✗ catch(...) {
323 ✗ DUNE_THROW(Dumux::ParameterException, errorMessage);
324 }
325
326
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 88 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 88 times.
 272 if (splitEntry.size() != 2)
327 ✗ DUNE_THROW(Dumux::ParameterException, errorMessage);
328 136 if (!labels.count(location))
329 ✗ DUNE_THROW(Dumux::ParameterException, errorMessage);
330 else
331
3/6
✓ Branch 1 taken 88 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 88 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 88 times.
 136 boundaryFaceMarker[labels.at(location)] = value;
332 }
333 }
334 33 return boundaryFaceMarker;
335 }
336
337 // returns the maximum possible pore body radii such that pore bodies do not intersect
338 template <class GetParameter>
339 33 std::vector<Scalar> getMaxPoreRadii_(std::size_t numSubregions, const GetParameter& getParameter) const
340 {
341 33 const auto numVertices = gridView_.size(dim);
342
1/2
✓ Branch 3 taken 21 times.
✗ Branch 4 not taken.
 66 std::vector<Scalar> maxPoreRadius(numVertices, std::numeric_limits<Scalar>::max());
343
344
2/4
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 21 times.
✗ Branch 4 not taken.
 33 if (!getParamFromGroup<bool>(paramGroup_, "Grid.CapPoreRadii", true))
345 return maxPoreRadius;
346
347 // check for a user-specified fixed throat length
348
1/4
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 33 const Scalar inputThroatLength = getParamFromGroup<Scalar>(paramGroup_, "Grid.ThroatLength", -1.0);
349
4/4
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 20 times.
✓ Branch 2 taken 1 times.
✓ Branch 3 taken 20 times.
 66 std::vector<Scalar> subregionInputThroatLengths;
350
351
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 20 times.
 33 if (numSubregions > 0)
352 {
353
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 1 times.
 2 for (int i = 0; i < numSubregions; ++i)
354 {
355 // adapt the parameter group if there are subregions
356
5/14
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 1 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 1 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
 3 const std::string paramGroup = paramGroup_ + ".SubRegion" + std::to_string(i);
357
1/4
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 1 const Scalar subregionInputThroatLength = getParamFromGroup<Scalar>(paramGroup, "Grid.ThroatLength", -1.0);
358
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
 1 subregionInputThroatLengths.push_back(subregionInputThroatLength);
359 }
360 }
361
362
6/6
✓ Branch 1 taken 43835 times.
✓ Branch 2 taken 21 times.
✓ Branch 3 taken 43835 times.
✓ Branch 4 taken 21 times.
✓ Branch 5 taken 10476 times.
✓ Branch 6 taken 33359 times.
 44322 for (const auto& element : elements(gridView_))
363 {
364 // do not cap the pore radius if a user-specified throat length if given
365
2/2
✓ Branch 0 taken 10476 times.
✓ Branch 1 taken 33359 times.
 44289 if (numSubregions > 0)
366 {
367 // check subregions
368
5/12
✓ Branch 1 taken 10476 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 10476 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 10476 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✓ Branch 10 taken 10476 times.
✓ Branch 11 taken 10476 times.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
 20952 const auto subregionId = getParameter(element, "ThroatRegionId");
369
1/2
✓ Branch 0 taken 10476 times.
✗ Branch 1 not taken.
 10476 if (subregionId >= 0)
370 {
371 // throat lies within a subregion
372
2/4
✓ Branch 0 taken 10476 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 10476 times.
✗ Branch 3 not taken.
 20952 if (subregionInputThroatLengths[subregionId] > 0.0)
373 246 continue;
374 }
375 ✗ else if (inputThroatLength > 0.0) // throat does not lie within a subregion
376 continue;
377 }
378
2/2
✓ Branch 0 taken 33236 times.
✓ Branch 1 taken 123 times.
 33813 else if (inputThroatLength > 0.0) // no subregions
379 continue;
380
381 // No fixed throat lengths given, check for max. pore radius.
382 // We define this as 1/2 of the length (minus a user specified value) of the shortest pore throat attached to the pore body.
383
3/4
✓ Branch 1 taken 43712 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 18 times.
✓ Branch 4 taken 43694 times.
 44043 const Scalar delta = element.geometry().volume();
384
4/6
✓ Branch 0 taken 18 times.
✓ Branch 1 taken 43694 times.
✓ Branch 3 taken 18 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 18 times.
✗ Branch 7 not taken.
 44043 static const Scalar minThroatLength = getParamFromGroup<Scalar>(paramGroup_, "Grid.MinThroatLength", 1e-6);
385 44043 const Scalar maxRadius = (delta - minThroatLength)/2.0;
386
2/2
✓ Branch 0 taken 87424 times.
✓ Branch 1 taken 43712 times.
 132129 for (int vIdxLocal = 0; vIdxLocal < 2 ; ++vIdxLocal)
387 {
388
2/4
✗ Branch 0 not taken.
✓ Branch 1 taken 87424 times.
✓ Branch 3 taken 87424 times.
✗ Branch 4 not taken.
 88086 const int vIdxGlobal = gridView_.indexSet().subIndex(element, vIdxLocal, dim);
389
4/4
✓ Branch 0 taken 24192 times.
✓ Branch 1 taken 63232 times.
✓ Branch 2 taken 24192 times.
✓ Branch 3 taken 63232 times.
 200633 maxPoreRadius[vIdxGlobal] = std::min(maxPoreRadius[vIdxGlobal], maxRadius);
390 }
391 }
392
393 return maxPoreRadius;
394 }
395
396 // returns a function taking a vertex and a pore label and returning a radius
397 34 std::function<Scalar(const Vertex&, const int)> poreRadiusGenerator_(const int subregionId) const
398 {
399 // adapt the parameter name if there are subregions
400
11/20
✓ Branch 0 taken 21 times.
✓ Branch 1 taken 1 times.
✓ Branch 3 taken 21 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 21 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 1 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 1 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 1 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 1 times.
✓ Branch 18 taken 21 times.
✓ Branch 19 taken 1 times.
✓ Branch 20 taken 21 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
 69 const std::string prefix = subregionId < 0 ? "Grid." : "Grid.Subregion" + std::to_string(subregionId) + ".";
401
402 // prepare random number generation for lognormal parameter distribution
403 34 std::mt19937 generator;
404
405 // check if pores for certain labels should be treated in a special way
406
7/16
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 22 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 22 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✓ Branch 12 taken 22 times.
✓ Branch 13 taken 1 times.
✓ Branch 14 taken 21 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 170 const auto poreLabelsToSetFixedRadius = getParamFromGroup<std::vector<int>>(paramGroup_, prefix + "PoreLabelsToSetFixedRadius", std::vector<int>{});
407
6/16
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 22 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 22 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✓ Branch 12 taken 22 times.
✗ Branch 13 not taken.
✓ Branch 14 taken 22 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 170 const auto poreLabelsToApplyFactorForRadius = getParamFromGroup<std::vector<int>>(paramGroup_, prefix + "PoreLabelsToApplyFactorForRadius", std::vector<int>{});
408
6/16
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 22 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 22 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✓ Branch 12 taken 22 times.
✗ Branch 13 not taken.
✓ Branch 14 taken 22 times.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 170 const auto poreRadiusForLabel = getParamFromGroup<std::vector<Scalar>>(paramGroup_, prefix + "FixedPoreRadiusForLabel", std::vector<Scalar>{});
409
7/18
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 22 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 22 times.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✓ Branch 12 taken 22 times.
✗ Branch 13 not taken.
✓ Branch 14 taken 22 times.
✗ Branch 15 not taken.
✓ Branch 16 taken 22 times.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
 170 const auto poreRadiusFactorForLabel = getParamFromGroup<std::vector<Scalar>>(paramGroup_, prefix + "PoreRadiusFactorForLabel", std::vector<Scalar>{});
410
411 130 const auto generateFunction = [&](auto& poreRadiusDist)
412 {
413 33780 return [=](const auto& vertex, const int poreLabel) mutable
414 {
415 33780 const auto radius = poreRadiusDist(generator);
416
417 // check if pores for certain labels should be treated in a special way
418
4/8
✓ Branch 0 taken 16704 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 16704 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 16704 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 16704 times.
 67560 if (poreLabelsToSetFixedRadius.empty() && poreLabelsToApplyFactorForRadius.empty())
419 return radius; // nothing special to be done
420
421 // set a fixed radius for a given label
422
0/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 22 else if (!poreLabelsToSetFixedRadius.empty() || !poreRadiusForLabel.empty())
423 {
424 ✗ if (poreLabelsToSetFixedRadius.size() != poreRadiusForLabel.size())
425 ✗ DUNE_THROW(Dumux::ParameterException, "PoreLabelsToSetFixedRadius must be of same size as FixedPoreRadiusForLabel");
426
427 ✗ if (const auto it = std::find(poreLabelsToSetFixedRadius.begin(), poreLabelsToSetFixedRadius.end(), poreLabel); it != poreLabelsToSetFixedRadius.end())
428 ✗ return poreRadiusForLabel[std::distance(poreLabelsToSetFixedRadius.begin(), it)];
429 }
430
431 // multiply the pore radius by a given value for a given label
432
0/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 22 else if (!poreLabelsToApplyFactorForRadius.empty() || !poreRadiusFactorForLabel.empty())
433 {
434 ✗ if (poreLabelsToApplyFactorForRadius.size() != poreRadiusFactorForLabel.size())
435 ✗ DUNE_THROW(Dumux::ParameterException, "PoreLabelsToApplyFactorForRadius must be of same size as PoreRadiusFactorForLabel");
436
437 ✗ if (const auto it = std::find(poreLabelsToApplyFactorForRadius.begin(), poreLabelsToApplyFactorForRadius.end(), poreLabel); it != poreLabelsToApplyFactorForRadius.end())
438 ✗ return poreRadiusFactorForLabel[std::distance(poreLabelsToApplyFactorForRadius.begin(), it)] * radius;
439 }
440
441 // default
442 return radius;
443 22 };
444 };
445
446
3/8
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 22 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 68 const Scalar fixedPoreRadius = getParamFromGroup<Scalar>(paramGroup_, prefix + "PoreInscribedRadius", -1.0);
447 // return random radius according to a user-specified distribution
448
2/2
✓ Branch 0 taken 6 times.
✓ Branch 1 taken 16 times.
 34 if (fixedPoreRadius <= 0.0)
449 {
450 // allow to specify a seed to get reproducible results
451
8/20
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 6 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 6 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 6 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 6 times.
✗ Branch 17 not taken.
✓ Branch 18 taken 6 times.
✗ Branch 19 not taken.
✗ Branch 21 not taken.
✓ Branch 22 taken 6 times.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
 36 const auto seed = getParamFromGroup<unsigned int>(paramGroup_, prefix + "ParameterRandomNumberSeed", std::random_device{}());
452 9 generator.seed(seed);
453
454
3/6
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✓ Branch 7 taken 6 times.
 18 const auto type = getParamFromGroup<std::string>(paramGroup_, prefix + "ParameterType", "lognormal");
455
1/2
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
 9 if (type == "lognormal")
456 {
457 // if we use a lognormal distribution, get the mean and standard deviation from input file
458
15/36
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 6 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 6 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 6 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 6 times.
✗ Branch 17 not taken.
✓ Branch 19 taken 6 times.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✓ Branch 22 taken 6 times.
✓ Branch 23 taken 6 times.
✗ Branch 24 not taken.
✓ Branch 25 taken 6 times.
✗ Branch 26 not taken.
✓ Branch 27 taken 6 times.
✗ Branch 28 not taken.
✓ Branch 29 taken 6 times.
✗ Branch 30 not taken.
✓ Branch 32 taken 6 times.
✗ Branch 33 not taken.
✓ Branch 35 taken 6 times.
✗ Branch 36 not taken.
✓ Branch 38 taken 6 times.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
✗ Branch 41 not taken.
✗ Branch 42 not taken.
✗ Branch 43 not taken.
✗ Branch 44 not taken.
✗ Branch 45 not taken.
 63 const auto [meanPoreRadius, stddevPoreRadius] = getDistributionInputParams_("Lognormal", prefix,
459 "MeanPoreInscribedRadius",
460 "StandardDeviationPoreInscribedRadius");
461 9 const Scalar variance = stddevPoreRadius*stddevPoreRadius;
462
463 using std::log;
464 using std::sqrt;
465 9 const Scalar mu = log(meanPoreRadius/sqrt(1.0 + variance/(meanPoreRadius*meanPoreRadius)));
466 9 const Scalar sigma = sqrt(log(1.0 + variance/(meanPoreRadius*meanPoreRadius)));
467
468
1/2
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
 9 Dumux::SimpleLogNormalDistribution<> poreRadiusDist(mu, sigma);
469
2/4
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
 9 return generateFunction(poreRadiusDist);
470 }
471 ✗ else if (type == "uniform")
472 {
473 // if we use a uniform distribution, get the min and max from input file
474 ✗ const auto [minPoreRadius, maxPoreRadius] = getDistributionInputParams_("Uniform", prefix,
475 "MinPoreInscribedRadius",
476 "MaxPoreInscribedRadius");
477 ✗ Dumux::SimpleUniformDistribution<> poreRadiusDist(minPoreRadius, maxPoreRadius);
478 ✗ return generateFunction(poreRadiusDist);
479 }
480 else
481 ✗ DUNE_THROW(Dune::InvalidStateException, "Unknown parameter type " << type);
482 }
483
484 // always return a fixed constant radius
485 else
486 {
487 25 auto poreRadiusDist = [fixedPoreRadius](auto& gen){ return fixedPoreRadius; };
488
2/6
✓ Branch 1 taken 16 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 16 times.
✗ Branch 5 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 25 return generateFunction(poreRadiusDist);
489 }
490 }
491
492 // print helpful error message if params are not properly provided
493 9 std::array<Scalar, 2> getDistributionInputParams_(const std::string& distributionName,
494 const std::string& prefix,
495 const std::string& paramName0,
496 const std::string& paramName1) const
497 {
498 try
499 {
500
2/6
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 9 return std::array{getParamFromGroup<Scalar>(paramGroup_, prefix + paramName0),
501
4/12
✓ Branch 1 taken 6 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 6 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 6 times.
✗ Branch 7 not taken.
✓ Branch 8 taken 6 times.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
 36 getParamFromGroup<Scalar>(paramGroup_, prefix + paramName1)};
502 }
503 ✗ catch (const Dumux::ParameterException& e)
504 {
505 std::cout << "\n" << distributionName << " pore-size distribution needs input parameters "
506 << prefix + paramName0 << " and " << prefix + paramName1 << ".\n"
507 ✗ << "Alternatively, use " << prefix << "PoreInscribedRadius to set a fixed inscribed pore radius." << std::endl;
508 ✗ DUNE_THROW(Dumux::ParameterException, e.what());
509 }
510 }
511
512 template <class GetParameter>
513 33 auto poreVolumeGenerator_(const GetParameter& getParameter) const
514 {
515
2/6
✓ Branch 2 taken 21 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 21 times.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 33 const auto geometry = Pore::shapeFromString(getParamFromGroup<std::string>(paramGroup_, "Grid.PoreGeometry"));
516
2/6
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 21 times.
✗ Branch 5 not taken.
✗ Branch 6 not taken.
✗ Branch 7 not taken.
 66 const auto capPoresOnBoundaries = getParamFromGroup<std::vector<int>>(paramGroup_, "Grid.CapPoresOnBoundaries", std::vector<int>{});
517
4/6
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✓ Branch 5 taken 18 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 21 times.
 37 if (!isUnique_(capPoresOnBoundaries))
518 ✗ DUNE_THROW(Dune::InvalidStateException, "CapPoresOnBoundaries must not contain duplicates");
519
520 // automatically determine the pore volume if not provided by the grid file
521 16951 return [=] (const auto& vertex, const auto vIdx)
522 {
523
4/10
✓ Branch 1 taken 16579 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 16579 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 16579 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 16579 times.
✗ Branch 10 not taken.
✗ Branch 12 not taken.
✗ Branch 13 not taken.
 50853 const Scalar r = getParameter(vertex, "PoreInscribedRadius");
524 16951 const Scalar volume = [&]
525 {
526 16951 if (geometry == Pore::Shape::cylinder)
527 {
528
0/4
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 42 static const Scalar fixedHeight = getParamFromGroup<Scalar>(paramGroup_, "Grid.PoreHeight", -1.0);
529
0/10
✗ Branch 0 not taken.
✗ Branch 1 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.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 40 const Scalar h = fixedHeight > 0.0 ? fixedHeight : getParameter(vertex, "PoreHeight");
530 80 return Pore::volume(geometry, r, h);
531 }
532 else
533 16911 return Pore::volume(geometry, r);
534
1/4
✗ Branch 0 not taken.
✓ Branch 1 taken 16579 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
 18109 }();
535
536
4/4
✓ Branch 0 taken 14778 times.
✓ Branch 1 taken 1801 times.
✓ Branch 2 taken 14778 times.
✓ Branch 3 taken 1801 times.
 33902 if (capPoresOnBoundaries.empty())
537 return volume;
538 else
539 {
540 14847 std::size_t numCaps = 0;
541 14847 Scalar factor = 1.0;
542 14847 const auto& pos = vertex.geometry().center();
543
4/4
✓ Branch 0 taken 88668 times.
✓ Branch 1 taken 14778 times.
✓ Branch 2 taken 88668 times.
✓ Branch 3 taken 14778 times.
 133554 for (auto boundaryIdx : capPoresOnBoundaries)
544 {
545
2/2
✓ Branch 0 taken 3218 times.
✓ Branch 1 taken 85450 times.
 89013 if (onBoundary_(pos, boundaryIdx))
546 {
547 3266 factor *= 0.5;
548 3266 ++numCaps;
549 }
550 }
551
552
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 14778 times.
 14847 if (numCaps > dimWorld)
553 ✗ DUNE_THROW(Dune::InvalidStateException, "Pore " << vIdx << " at " << pos << " capped " << numCaps << " times. Capping should not happen more than " << dimWorld << " times");
554
555 14847 return volume * factor;
556 }
557
10/14
✓ Branch 0 taken 2 times.
✓ Branch 1 taken 21 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1 times.
✓ Branch 4 taken 11 times.
✓ Branch 5 taken 1 times.
✓ Branch 6 taken 1 times.
✓ Branch 7 taken 9 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 2 times.
✓ Branch 10 taken 7 times.
✓ Branch 11 taken 8 times.
✗ Branch 13 not taken.
✗ Branch 14 not taken.
 17012 };
558 }
559
560 // returns a lambda taking a element and returning a radius
561 template <class GetParameter>
562 34 auto throatInscribedRadiusGenerator_(const int subregionId, const GetParameter& getParameter) const
563 {
564 // adapt the parameter name if there are subregions
565
11/20
✓ Branch 0 taken 21 times.
✓ Branch 1 taken 1 times.
✓ Branch 3 taken 21 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 21 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 1 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 1 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 1 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 1 times.
✓ Branch 18 taken 21 times.
✓ Branch 19 taken 1 times.
✓ Branch 20 taken 21 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
 35 const std::string prefix = subregionId < 0 ? "Grid." : "Grid.Subregion" + std::to_string(subregionId) + ".";
566
567 // check for a user-specified fixed throat radius
568
3/8
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 22 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 68 const Scalar inputThroatInscribedRadius = getParamFromGroup<Scalar>(paramGroup_, prefix + "ThroatInscribedRadius", -1.0);
569
570 // shape parameter for calculation of throat radius
571
3/10
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 22 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 68 const Scalar throatN = getParamFromGroup<Scalar>(paramGroup_, prefix + "ThroatInscribedRadiusN", 0.1);
572
573 45460 return [=](const Element& element)
574 {
575
2/2
✓ Branch 1 taken 43592 times.
✓ Branch 2 taken 1414 times.
 45460 const Scalar delta = element.geometry().volume();
576
4/4
✓ Branch 0 taken 43592 times.
✓ Branch 1 taken 1414 times.
✓ Branch 2 taken 43592 times.
✓ Branch 3 taken 1414 times.
 90920 const std::array<Vertex, 2> vertices = {element.template subEntity<dim>(0), element.template subEntity<dim>(1)};
577
578 // the element parameters (throat radius and length)
579
2/2
✓ Branch 0 taken 43592 times.
✓ Branch 1 taken 1414 times.
 45460 if (inputThroatInscribedRadius > 0.0)
580 return inputThroatInscribedRadius;
581 else
582 {
583
6/14
✓ Branch 1 taken 43592 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 43592 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 43592 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 43592 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 43592 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 43592 times.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
 131769 const Scalar poreRadius0 = getParameter(vertices[0], "PoreInscribedRadius");
584
5/12
✓ Branch 1 taken 43592 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 43592 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 43592 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 43592 times.
✗ Branch 11 not taken.
✓ Branch 12 taken 43592 times.
✗ Branch 13 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
 131769 const Scalar poreRadius1 = getParameter(vertices[1], "PoreInscribedRadius");
585 43923 return Throat::averagedRadius(poreRadius0, poreRadius1, delta, throatN);
586 }
587
2/2
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 21 times.
 35 };
588 }
589
590 // returns a lambda taking a element and returning a length
591 template <class GetParameter>
592 34 auto throatLengthGenerator_(int subregionId, const GetParameter& getParameter) const
593 {
594 // adapt the parameter name if there are subregions
595
11/20
✓ Branch 0 taken 21 times.
✓ Branch 1 taken 1 times.
✓ Branch 3 taken 21 times.
✗ Branch 4 not taken.
✓ Branch 6 taken 21 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 1 times.
✗ Branch 10 not taken.
✓ Branch 12 taken 1 times.
✗ Branch 13 not taken.
✓ Branch 15 taken 1 times.
✗ Branch 16 not taken.
✓ Branch 17 taken 1 times.
✓ Branch 18 taken 21 times.
✓ Branch 19 taken 1 times.
✓ Branch 20 taken 21 times.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
 35 const std::string prefix = subregionId < 0 ? "Grid." : "Grid.Subregion" + std::to_string(subregionId) + ".";
596
3/8
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 22 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
 68 const Scalar inputThroatLength = getParamFromGroup<Scalar>(paramGroup_, prefix + "ThroatLength", -1.0);
597 // decide whether to subtract the pore radii from the throat length or not
598
3/10
✓ Branch 1 taken 22 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 22 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 22 times.
✗ Branch 7 not taken.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
 68 const bool subtractRadiiFromThroatLength = getParamFromGroup<bool>(paramGroup_, prefix + "SubtractPoreInscribedRadiiFromThroatLength", true);
599
600 45460 return [=](const Element& element)
601 {
602
1/2
✓ Branch 0 taken 45006 times.
✗ Branch 1 not taken.
 45460 if (inputThroatLength > 0.0)
603 return inputThroatLength;
604
605
1/2
✓ Branch 1 taken 45006 times.
✗ Branch 2 not taken.
 45337 const Scalar delta = element.geometry().volume();
606 45337 if (subtractRadiiFromThroatLength)
607 {
608
2/4
✓ Branch 1 taken 45006 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 45006 times.
✗ Branch 5 not taken.
 90674 const std::array<Vertex, 2> vertices = {element.template subEntity<dim>(0), element.template subEntity<dim>(1)};
609
12/28
✓ Branch 1 taken 45006 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 45006 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 45006 times.
✗ Branch 8 not taken.
✓ Branch 10 taken 45006 times.
✗ Branch 11 not taken.
✓ Branch 13 taken 45006 times.
✗ Branch 14 not taken.
✓ Branch 16 taken 45006 times.
✗ Branch 17 not taken.
✓ Branch 19 taken 45006 times.
✗ Branch 20 not taken.
✓ Branch 22 taken 45006 times.
✗ Branch 23 not taken.
✓ Branch 24 taken 45006 times.
✗ Branch 25 not taken.
✓ Branch 26 taken 45006 times.
✗ Branch 27 not taken.
✓ Branch 28 taken 45006 times.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
✓ Branch 31 taken 45006 times.
✗ Branch 32 not taken.
✗ Branch 33 not taken.
✗ Branch 34 not taken.
✗ Branch 35 not taken.
 226685 const Scalar result = delta - getParameter(vertices[0], "PoreInscribedRadius") - getParameter(vertices[1], "PoreInscribedRadius");
610
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 45006 times.
 45337 if (result <= 0.0)
611 ✗ DUNE_THROW(Dune::GridError, "Pore radii are so large they intersect! Something went wrong at element " << gridView_.indexSet().index(element));
612 else
613 return result;
614 }
615 else
616 return delta;
617
3/6
✓ Branch 0 taken 1 times.
✓ Branch 1 taken 21 times.
✓ Branch 2 taken 45006 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
 45372 };
618 }
619
620 bool onBoundary_(const GlobalPosition& pos, std::size_t boundaryIdx) const
621 {
622 178275 constexpr auto eps = 1e-8; //TODO
623
624 static constexpr std::array boundaryIdxToCoordinateIdx{0, 0, 1, 1, 2, 2};
625
10/18
✓ Branch 0 taken 38 times.
✓ Branch 1 taken 40 times.
✓ Branch 2 taken 43034 times.
✓ Branch 3 taken 45022 times.
✓ Branch 4 taken 44592 times.
✓ Branch 5 taken 44636 times.
✓ Branch 6 taken 138 times.
✓ Branch 7 taken 207 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✓ Branch 12 taken 254 times.
✓ Branch 13 taken 314 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
 178275 const auto coordinateIdx = boundaryIdxToCoordinateIdx[boundaryIdx];
626 178275 auto isMaxBoundary = [] (int n) { return (n % 2 != 0); };
627
628
10/18
✓ Branch 0 taken 38 times.
✓ Branch 1 taken 40 times.
✓ Branch 2 taken 43034 times.
✓ Branch 3 taken 45022 times.
✓ Branch 4 taken 44592 times.
✓ Branch 5 taken 44636 times.
✓ Branch 6 taken 138 times.
✓ Branch 7 taken 207 times.
✗ Branch 8 not taken.
✗ Branch 9 not taken.
✗ Branch 10 not taken.
✗ Branch 11 not taken.
✓ Branch 12 taken 254 times.
✓ Branch 13 taken 314 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✗ Branch 16 not taken.
✗ Branch 17 not taken.
 178275 if (isMaxBoundary(boundaryIdx))
629 264092 return pos[coordinateIdx] > bBoxMax_[coordinateIdx] - eps;
630 else
631 270577 return pos[coordinateIdx] < bBoxMin_[coordinateIdx] + eps;
632 }
633
634 // check if a container is unique
635 // copy the container such that is does not get altered by std::sort
636 template <class T>
637 ✗ bool isUnique_(T t) const
638 {
639 9 std::sort(t.begin(), t.end());
640
2/4
✗ Branch 3 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
✓ Branch 6 taken 3 times.
 9 return (std::unique(t.begin(), t.end()) == t.end());
641 }
642
643 const GridView gridView_;
644
645 const std::string paramGroup_;
646 BoundaryList priorityList_;
647 BoundaryList boundaryFaceIndex_;
648
649 //! the bounding box of the whole domain
650 GlobalPosition bBoxMin_ = GlobalPosition(std::numeric_limits<double>::max());
651 GlobalPosition bBoxMax_ = GlobalPosition(std::numeric_limits<double>::min());
652 };
653
654 } // namespace Dumux::PoreNetwork
655
656 #endif
657