GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/test/freeflow/brinkman/problem.hh
Date: 2024-05-04 19:09:25
Exec Total Coverage
Lines: 49 54 90.7%
Functions: 10 27 37.0%
Branches: 84 154 54.5%

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 NavierStokesTests
10 * \brief Darcy Brinkman model for a single-domain evaluation of coupled freeflow and porous medium flows
11 */
12
13 #ifndef DUMUX_BRINKMAN_TEST_PROBLEM_HH
14 #define DUMUX_BRINKMAN_TEST_PROBLEM_HH
15
16 #include <dumux/common/parameters.hh>
17 #include <dumux/common/properties.hh>
18
19 #include <dumux/freeflow/navierstokes/boundarytypes.hh>
20
21 #include <dumux/freeflow/navierstokes/momentum/fluxhelper.hh>
22 #include <dumux/freeflow/navierstokes/momentum/brinkman.hh>
23 #include <dumux/freeflow/navierstokes/scalarfluxhelper.hh>
24 #include <dumux/freeflow/navierstokes/mass/1p/advectiveflux.hh>
25
26 namespace Dumux {
27
28 /*!
29 * \ingroup NavierStokesTests
30 * \brief Darcy Brinkman model for a single-domain evaluation of coupled freeflow and porous medium flows
31 */
32 template <class TypeTag, class BaseProblem>
33 class BrinkmanProblem : public BaseProblem
34 {
35 using ParentType = BaseProblem;
36 using BoundaryTypes = typename ParentType::BoundaryTypes;
37 using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
38 using FVElementGeometry = typename GridGeometry::LocalView;
39 using SubControlVolume = typename FVElementGeometry::SubControlVolume;
40 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
41 using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>;
42 using InitialValues = typename ParentType::InitialValues;
43 using Sources = typename ParentType::Sources;
44 using DirichletValues = typename ParentType::DirichletValues;
45 using BoundaryFluxes = typename ParentType::BoundaryFluxes;
46 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
47 using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>;
48
49 static constexpr auto dimWorld = GridGeometry::GridView::dimensionworld;
50 static constexpr auto dim = GridGeometry::GridView::dimension;
51 using Element = typename FVElementGeometry::Element;
52 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
53
54 using CouplingManager = GetPropType<TypeTag, Properties::CouplingManager>;
55 using SpatialParams = GetPropType<TypeTag, Properties::SpatialParams>;
56 using PermeabilityType = typename SpatialParams::PermeabilityType;
57
58 public:
59 using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices;
60
61 8 BrinkmanProblem(std::shared_ptr<const GridGeometry> gridGeometry,
62 std::shared_ptr<CouplingManager> couplingManager)
63
10/30
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 4 times.
✗ Branch 5 not taken.
✓ Branch 9 taken 4 times.
✗ Branch 10 not taken.
✓ Branch 11 taken 4 times.
✗ Branch 12 not taken.
✓ Branch 13 taken 4 times.
✗ Branch 14 not taken.
✗ Branch 15 not taken.
✓ Branch 16 taken 4 times.
✓ Branch 18 taken 4 times.
✗ Branch 19 not taken.
✓ Branch 21 taken 4 times.
✗ Branch 22 not taken.
✓ Branch 24 taken 4 times.
✗ Branch 25 not taken.
✓ Branch 27 taken 4 times.
✗ 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 37 not taken.
✗ Branch 38 not taken.
✗ Branch 39 not taken.
✗ Branch 40 not taken.
32 : ParentType(gridGeometry, couplingManager)
64 {
65
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
8 deltaP_ = getParam<Scalar>("Problem.PressureDifference");
66
1/2
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
8 dynamicViscosity_ = getParam<Scalar>("Component.LiquidKinematicViscosity", 1.0)
67
1/4
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
8 * getParam<Scalar>("Component.LiquidDensity", 1.0);
68
2/4
✓ Branch 1 taken 4 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 4 times.
8 problemName_ = getParam<std::string>("Problem.Name");
69
70
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
8 storeBrinkmanParamsForOutput_();
71 8 }
72
73
74 /*!
75 * \brief The problem name.
76 */
77 const std::string& name() const
78
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 { return problemName_ ; }
79
80 /*!
81 * \brief Specifies which kind of boundary condition should be
82 * used for which equation on a given boundary control volume.
83 *
84 * \param globalPos The position of the center of the finite volume
85 */
86 9360 BoundaryTypes boundaryTypesAtPos(const GlobalPosition& globalPos) const
87 {
88
3/4
✓ Branch 0 taken 1548 times.
✓ Branch 1 taken 5200 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 1560 times.
17668 BoundaryTypes values;
89
90 if constexpr (ParentType::isMomentumProblem())
91 {
92 30936 if (onLowerBoundary_(globalPos) || onUpperBoundary_(globalPos))
93 values.setAllDirichlet();
94 else
95 values.setAllNeumann();
96 }
97 else
98
6/8
✓ Branch 0 taken 1548 times.
✓ Branch 1 taken 5200 times.
✓ Branch 2 taken 1548 times.
✓ Branch 3 taken 5200 times.
✗ Branch 4 not taken.
✓ Branch 5 taken 1560 times.
✗ Branch 6 not taken.
✓ Branch 7 taken 1560 times.
16616 values.setNeumann(Indices::conti0EqIdx);
99
100 9360 return values;
101 }
102
103 /*!
104 * \brief Evaluates the boundary conditions for a Dirichlet control volume.
105 *
106 * \param globalPos The center of the finite volume which ought to be set.
107 */
108 DirichletValues dirichletAtPos(const GlobalPosition &globalPos) const
109 {
110 // no-flow/no-slip
111 5064 return DirichletValues(0.0);
112 }
113
114 /*!
115 * \brief Evaluates the boundary conditions for a Neumann control volume.
116 *
117 * \param element The element for which the Neumann boundary condition is set
118 * \param fvGeometry The fvGeometry
119 * \param elemVolVars The element volume variables
120 * \param elemFaceVars The element face variables
121 * \param scvf The boundary sub control volume face
122 */
123 template<class ElementVolumeVariables, class ElementFluxVariablesCache>
124 42640 BoundaryFluxes neumann(const Element& element,
125 const FVElementGeometry& fvGeometry,
126 const ElementVolumeVariables& elemVolVars,
127 const ElementFluxVariablesCache& elemFluxVarsCache,
128 const SubControlVolumeFace& scvf) const
129 {
130
0/2
✗ Branch 1 not taken.
✗ Branch 2 not taken.
47840 BoundaryFluxes values(0.0);
131
2/5
✓ Branch 0 taken 2340 times.
✓ Branch 1 taken 2340 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
53040 const auto& globalPos = scvf.ipGlobal();
132
133 47840 const auto pRef = referencePressure();
134 81640 const auto p = isInlet_(globalPos) ? pRef + deltaP_ : pRef;
135
136 using DM = typename FVElementGeometry::GridGeometry::DiscretizationMethod;
137 if constexpr (ParentType::isMomentumProblem())
138 {
139 if constexpr (DiscretizationMethods::isCVFE<DM>)
140 {
141 9360 const auto mu = this->effectiveViscosity(element, fvGeometry, scvf);
142 9360 const auto& fluxVarCache = elemFluxVarsCache[scvf];
143 9360 Dune::FieldMatrix<Scalar, dimWorld, dimWorld> gradV(0.0);
144
4/4
✓ Branch 0 taken 18720 times.
✓ Branch 1 taken 4680 times.
✓ Branch 2 taken 18720 times.
✓ Branch 3 taken 4680 times.
56160 for (const auto& scv : scvs(fvGeometry))
145 {
146 37440 const auto& volVars = elemVolVars[scv];
147
2/2
✓ Branch 0 taken 37440 times.
✓ Branch 1 taken 18720 times.
112320 for (int dir = 0; dir < dimWorld; ++dir)
148 374400 gradV[dir].axpy(volVars.velocity(dir), fluxVarCache.gradN(scv.indexInElement()));
149 }
150 18720 gradV.mtv(-mu*scvf.unitOuterNormal(), values); // - μ ∇v^T n
151 9360 values.axpy(p, scvf.unitOuterNormal()); // + p n
152 }
153 else
154 values = NavierStokesMomentumBoundaryFluxHelper::fixedPressureMomentumFlux(
155 *this, fvGeometry, scvf, elemVolVars, elemFluxVarsCache, p
156 );
157 }
158 else
159 {
160 if constexpr (DiscretizationMethods::isCVFE<DM>)
161 33280 values = (this->faceVelocity(element, fvGeometry, scvf) * scvf.unitOuterNormal())
162 133120 * elemVolVars[scvf.insideScvIdx()].density();
163 else
164 5200 values = NavierStokesScalarBoundaryFluxHelper<AdvectiveFlux<ModelTraits>>::scalarOutflowFlux(
165 *this, element, fvGeometry, scvf, elemVolVars
166 );
167 }
168
169 42640 return values;
170 }
171
172 template<class ... Args>
173 Scalar referencePressure(Args&&...) const
174 { return 0.0; }
175
176 //! Add the Brinkman term via the source using the helper function addBrinkmanTerm
177 template <class ElementVolumeVariables>
178 Sources source(const Element& element,
179 const FVElementGeometry& fvGeometry,
180 const ElementVolumeVariables& elemVolVars,
181 const SubControlVolume& scv) const
182 {
183
0/4
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✗ Branch 5 not taken.
2705120 Sources source(0.0);
184 if constexpr (ParentType::isMomentumProblem())
185 1178720 addBrinkmanTerm(source, *this, element, fvGeometry, elemVolVars, scv);
186 1178720 return source;
187 }
188
189 const std::vector<PermeabilityType> permeabilityOutput() const
190
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 { return outputPermeability_; }
191
192 const std::vector<Scalar> brinkmanEpsilon() const
193
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 { return outputBrinkmanEpsilon_; }
194
195 private:
196 bool isInlet_(const GlobalPosition& globalPos) const
197
10/22
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2340 times.
✓ Branch 6 taken 2340 times.
✓ Branch 7 taken 2340 times.
✓ Branch 8 taken 2340 times.
✓ Branch 9 taken 2340 times.
✓ Branch 10 taken 2340 times.
✓ Branch 11 taken 2340 times.
✓ Branch 12 taken 2340 times.
✓ Branch 13 taken 2340 times.
✓ Branch 14 taken 2340 times.
✗ Branch 15 not taken.
✗ Branch 17 not taken.
✗ Branch 18 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
132600 { return globalPos[0] < this->gridGeometry().bBoxMin()[0] + eps_; }
198
199 bool isOutlet_(const GlobalPosition& globalPos) const
200 { return globalPos[0] > this->gridGeometry().bBoxMin()[0] + eps_; }
201
202 bool onLowerBoundary_(const GlobalPosition &globalPos) const
203
10/10
✓ Branch 0 taken 3054 times.
✓ Branch 1 taken 1626 times.
✓ Branch 2 taken 3054 times.
✓ Branch 3 taken 1626 times.
✓ Branch 4 taken 3054 times.
✓ Branch 5 taken 1626 times.
✓ Branch 6 taken 3054 times.
✓ Branch 7 taken 1626 times.
✓ Branch 8 taken 3054 times.
✓ Branch 9 taken 1626 times.
23400 { return globalPos[1] < this->gridGeometry().bBoxMin()[1] + eps_; }
204
205 bool onUpperBoundary_(const GlobalPosition &globalPos) const
206
10/10
✓ Branch 0 taken 1428 times.
✓ Branch 1 taken 1626 times.
✓ Branch 2 taken 1428 times.
✓ Branch 3 taken 1626 times.
✓ Branch 4 taken 1428 times.
✓ Branch 5 taken 1626 times.
✓ Branch 6 taken 1428 times.
✓ Branch 7 taken 1626 times.
✓ Branch 8 taken 1428 times.
✓ Branch 9 taken 1626 times.
15270 { return globalPos[1] > this->gridGeometry().bBoxMax()[1] - eps_; }
207
208 4 void storeBrinkmanParamsForOutput_()
209 {
210 if constexpr (!ParentType::isMomentumProblem())
211 {
212 12 outputPermeability_.resize(this->gridGeometry().gridView().size(0));
213 12 outputBrinkmanEpsilon_.resize(this->gridGeometry().gridView().size(0));
214
8/11
✓ Branch 3 taken 3601 times.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 1 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 7200 times.
✓ Branch 10 taken 1 times.
✓ Branch 11 taken 7200 times.
✓ Branch 12 taken 1 times.
✓ Branch 14 taken 7200 times.
✗ Branch 15 not taken.
43216 for (const auto& element : elements(this->gridGeometry().gridView()))
215 {
216
3/6
✓ Branch 1 taken 7200 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 7200 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 7200 times.
✗ Branch 8 not taken.
64800 const auto eIdx = this->gridGeometry().elementMapper().index(element);
217
3/6
✓ Branch 1 taken 7200 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 7200 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 7200 times.
✗ Branch 8 not taken.
50400 outputPermeability_[eIdx] = this->spatialParams().permeabilityAtPos(element.geometry().center());
218
4/8
✓ Branch 1 taken 7200 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 7200 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 7200 times.
✗ Branch 8 not taken.
✓ Branch 13 taken 7200 times.
✗ Branch 14 not taken.
50400 outputBrinkmanEpsilon_[eIdx] = this->spatialParams().brinkmanEpsilonAtPos(element.geometry().center());
219 }
220 }
221 4 }
222
223 static constexpr Scalar eps_ = 1e-6;
224 Scalar dynamicViscosity_;
225 Scalar deltaP_;
226
227 std::vector<PermeabilityType> outputPermeability_;
228 std::vector<Scalar> outputBrinkmanEpsilon_;
229
230 std::string problemName_;
231
232 };
233
234 } // end namespace Dumux
235
236 #endif
237