GCC Code Coverage Report

Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/test/porousmediumflow/2p2c/waterair/problem.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 46 50 92.0%
Functions: 10 16 62.5%
Branches: 96 130 73.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 TwoPTwoCTests
10 * \brief Non-isothermal gas injection problem where a gas (e.g. air)
11 * is injected into a fully water saturated medium.
12 */
13
14 #ifndef DUMUX_WATER_AIR_PROBLEM_HH
15 #define DUMUX_WATER_AIR_PROBLEM_HH
16
17 #include <dumux/common/properties.hh>
18 #include <dumux/common/parameters.hh>
19 #include <dumux/common/boundarytypes.hh>
20 #include <dumux/common/numeqvector.hh>
21
22 #include <dumux/material/components/n2.hh>
23 #include <dumux/porousmediumflow/problem.hh>
24
25 namespace Dumux {
26
27 /*!
28 * \ingroup TwoPTwoCModel
29 * \brief Non-isothermal gas injection problem where a gas (e.g. air)
30 * is injected into a fully water saturated medium.
31 *
32 * */
33 template <class TypeTag >
34 class WaterAirProblem : public PorousMediumFlowProblem<TypeTag>
35 {
36 using ParentType = PorousMediumFlowProblem<TypeTag>;
37
38 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
39 using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
40 using GridView = typename GridGeometry::GridView;
41 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
42 using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>;
43 using Indices = typename ModelTraits::Indices;
44
45 // primary variable indices
46 enum
47 {
48 pressureIdx = Indices::pressureIdx,
49 switchIdx = Indices::switchIdx,
50 temperatureIdx = Indices::temperatureIdx,
51 energyEqIdx = Indices::energyEqIdx
52 };
53
54 // equation indices
55 enum
56 {
57 contiH2OEqIdx = Indices::conti0EqIdx + FluidSystem::H2OIdx,
58 contiN2EqIdx = Indices::conti0EqIdx + FluidSystem::N2Idx
59 };
60
61 // phase presence
62 enum { wPhaseOnly = Indices::firstPhaseOnly };
63 // component index
64 enum { N2Idx = FluidSystem::N2Idx };
65
66 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
67 using NumEqVector = Dumux::NumEqVector<PrimaryVariables>;
68 using BoundaryTypes = Dumux::BoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>;
69 using Element = typename GridView::template Codim<0>::Entity;
70 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
71
72 //! Property that defines whether mole or mass fractions are used
73 static constexpr bool useMoles = ModelTraits::useMoles();
74
75 public:
76 3 WaterAirProblem(std::shared_ptr<const GridGeometry> gridGeometry)
77
7/20
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
✓ Branch 8 taken 3 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 3 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 3 times.
✓ Branch 15 taken 3 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 3 times.
✗ Branch 19 not taken.
✗ Branch 20 not taken.
✗ Branch 21 not taken.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 26 not taken.
✗ Branch 27 not taken.
 9 : ParentType(gridGeometry)
78 {
79 3 maxDepth_ = 1000.0; // [m]
80
81
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 FluidSystem::init();
82
83
2/4
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 3 times.
 3 name_ = getParam<std::string>("Problem.Name");
84 3 useDirichlet_ = name_.find("buoyancy") != std::string::npos;
85
86 //stating in the console whether mole or mass fractions are used
87 if(useMoles)
88
1/2
✓ Branch 2 taken 3 times.
✗ Branch 3 not taken.
 6 std::cout << "The problem uses mole-fractions" << std::endl;
89 else
90 std::cout << "The problem uses mass-fractions" << std::endl;
91
92
2/4
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 3 times.
✗ Branch 5 not taken.
 6 this->spatialParams().plotMaterialLaw();
93 3 }
94
95 /*!
96 * \brief The problem name.
97 *
98 * This is used as a prefix for files generated by the simulation.
99 */
100 const std::string& name() const
101
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
 3 { return name_; }
102
103
104 /*!
105 * \brief Specifies which kind of boundary condition should be
106 * used for which equation on a given boundary segment.
107 *
108 * \param globalPos The position for which the bc type should be evaluated
109 */
110 248639 BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
111 {
112 248639 BoundaryTypes bcTypes;
113
8/8
✓ Branch 0 taken 181887 times.
✓ Branch 1 taken 66752 times.
✓ Branch 2 taken 181887 times.
✓ Branch 3 taken 66752 times.
✓ Branch 4 taken 115135 times.
✓ Branch 5 taken 66752 times.
✓ Branch 6 taken 115135 times.
✓ Branch 7 taken 66752 times.
 497278 if(globalPos[0] > 40 - eps_ || globalPos[0] < eps_)
114 bcTypes.setAllDirichlet();
115 else
116 bcTypes.setAllNeumann();
117
118
2/2
✓ Branch 0 taken 75519 times.
✓ Branch 1 taken 173120 times.
 248639 if (useDirichlet_)
119 {
120 75519 if (isInjectionArea_(globalPos))
121 bcTypes.setAllDirichlet();
122 }
123
124 248639 return bcTypes;
125 }
126
127 /*!
128 * \brief Evaluates the boundary conditions for a Dirichlet boundary segment.
129 *
130 * \param globalPos The position for which the Dirichlet condition should be evaluated
131 */
132 64321 PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
133 {
134 128642 PrimaryVariables priVars = initial_(globalPos);
135
136 64321 if (useDirichlet_)
137 {
138 28802 if (isInjectionArea_(globalPos))
139 {
140 3009 priVars.setState(Indices::bothPhases);
141 3009 priVars[switchIdx] = 0.2;
142 3009 return priVars;
143 }
144 }
145
146 return priVars;
147 }
148
149 /*!
150 * \brief Evaluates the boundary conditions for a Neumann boundary segment.
151 *
152 * \param globalPos The position for which the Neumann condition should be evaluated
153 * \return the mole/mass flux of each component. Negative values mean influx.
154 *
155 * The units must be according to using either mole or mass fractions (mole/(m^2*s) or kg/(m^2*s)).
156 */
157 479520 NumEqVector neumannAtPos(const GlobalPosition &globalPos) const
158 {
159 479520 NumEqVector values(0.0);
160
161 // we inject pure gasious nitrogen at the initial condition temperature and pressure from the bottom (negative values mean injection)
162 514812 if (isInjectionArea_(globalPos))
163 {
164 35292 values[contiN2EqIdx] = useMoles ? -1e-3/FluidSystem::molarMass(N2Idx) : -1e-3; // kg/(m^2*s) or mole/(m^2*s)
165
166 70584 const auto initialValues = initial_(globalPos);
167 35292 const auto& fluidMatrixInteraction = this->spatialParams().fluidMatrixInteractionAtPos(globalPos);
168 70584 const auto pn = initialValues[pressureIdx] + fluidMatrixInteraction.endPointPc();
169 35292 const auto t = initialValues[temperatureIdx];
170
171 // note: energy equation is always formulated in terms of mass specific quantities, not per mole
172 35292 values[energyEqIdx] = -1e-3*Components::N2<Scalar>::gasEnthalpy(t, pn);
173 }
174
175 479520 return values;
176 }
177
178
179 /*!
180 * \brief Evaluates the initial value for a control volume.
181 *
182 * \param globalPos The position for which the initial condition should be evaluated
183 *
184 * For this method, the \a values parameter stores primary
185 * variables.
186 */
187 12546 PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
188 {
189 25092 auto priVars = initial_(globalPos);
190 // initially there is a heated lens in the domain
191
10/10
✓ Branch 0 taken 6338 times.
✓ Branch 1 taken 6208 times.
✓ Branch 2 taken 3234 times.
✓ Branch 3 taken 3104 times.
✓ Branch 4 taken 3234 times.
✓ Branch 5 taken 3104 times.
✓ Branch 6 taken 2434 times.
✓ Branch 7 taken 800 times.
✓ Branch 8 taken 2434 times.
✓ Branch 9 taken 800 times.
 12546 if (globalPos[0] > 20 - eps_ && globalPos[0] < 30 + eps_ && globalPos[1] < 30 + eps_)
192 4868 priVars[temperatureIdx] = 380.0;
193
194 12546 return priVars;
195 }
196
197 private:
198 // internal method for the initial condition (reused for the
199 // dirichlet conditions!)
200 ✗ PrimaryVariables initial_(const GlobalPosition &globalPos) const
201 {
202 112159 PrimaryVariables priVars(0.0);
203
4/4
✓ Branch 1 taken 25793 times.
✓ Branch 2 taken 38528 times.
✓ Branch 3 taken 6338 times.
✓ Branch 4 taken 6208 times.
 112159 priVars.setState(wPhaseOnly);
204 112159 Scalar densityW = 1000.0;
205
8/8
✓ Branch 2 taken 25793 times.
✓ Branch 3 taken 38528 times.
✓ Branch 4 taken 25793 times.
✓ Branch 5 taken 38528 times.
✓ Branch 6 taken 6338 times.
✓ Branch 7 taken 6208 times.
✓ Branch 8 taken 6338 times.
✓ Branch 9 taken 6208 times.
 224318 priVars[pressureIdx] = 1e5 + (maxDepth_ - globalPos[1])*densityW*9.81;
206
4/4
✓ Branch 1 taken 25793 times.
✓ Branch 2 taken 38528 times.
✓ Branch 3 taken 6338 times.
✓ Branch 4 taken 6208 times.
 112159 priVars[switchIdx] = 0.0;
207
4/4
✓ Branch 1 taken 25793 times.
✓ Branch 2 taken 38528 times.
✓ Branch 3 taken 6338 times.
✓ Branch 4 taken 6208 times.
 112159 priVars[temperatureIdx] = initialTemperatureProfile_(globalPos);
208 ✗ return priVars;
209 }
210
211 ✗ Scalar initialTemperatureProfile_(const GlobalPosition &globalPos) const
212
8/8
✓ Branch 2 taken 25793 times.
✓ Branch 3 taken 38528 times.
✓ Branch 4 taken 25793 times.
✓ Branch 5 taken 38528 times.
✓ Branch 6 taken 6338 times.
✓ Branch 7 taken 6208 times.
✓ Branch 8 taken 6338 times.
✓ Branch 9 taken 6208 times.
 224318 { return 283.0 + (maxDepth_ - globalPos[1])*0.03; }
213
214 ✗ bool isInjectionArea_(const GlobalPosition& globalPos) const
215 {
216
34/48
✓ Branch 0 taken 291848 times.
✓ Branch 1 taken 187672 times.
✓ Branch 2 taken 291848 times.
✓ Branch 3 taken 187672 times.
✓ Branch 4 taken 104176 times.
✓ Branch 5 taken 187672 times.
✓ Branch 6 taken 104176 times.
✓ Branch 7 taken 187672 times.
✓ Branch 8 taken 35292 times.
✓ Branch 9 taken 68884 times.
✓ Branch 10 taken 35292 times.
✓ Branch 11 taken 68884 times.
✓ Branch 12 taken 14401 times.
✓ Branch 13 taken 11392 times.
✓ Branch 14 taken 14401 times.
✓ Branch 15 taken 11392 times.
✓ Branch 16 taken 3009 times.
✓ Branch 17 taken 11392 times.
✓ Branch 18 taken 3009 times.
✓ Branch 19 taken 11392 times.
✓ Branch 20 taken 3009 times.
✗ Branch 21 not taken.
✓ Branch 22 taken 3009 times.
✗ 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 taken 42715 times.
✓ Branch 37 taken 32804 times.
✓ Branch 38 taken 42715 times.
✓ Branch 39 taken 32804 times.
✓ Branch 40 taken 9911 times.
✓ Branch 41 taken 32804 times.
✓ Branch 42 taken 9911 times.
✓ Branch 43 taken 32804 times.
✓ Branch 44 taken 5593 times.
✓ Branch 45 taken 4318 times.
✓ Branch 46 taken 5593 times.
✓ Branch 47 taken 4318 times.
 1161664 return globalPos[0] > 14.8 - eps_ && globalPos[0] < 25.2 + eps_ && globalPos[1] < eps_;
217 }
218
219 Scalar maxDepth_;
220 static constexpr Scalar eps_ = 1e-2;
221 std::string name_;
222 bool useDirichlet_;
223 };
224
225 } // end namespace Dumux
226
227 #endif
228