GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/test/porousmediumflow/1pnc/1p2c/nonisothermal/conduction/problem.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 46 53 86.8%
Functions: 9 24 37.5%
Branches: 88 158 55.7%

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 OnePNCTests
10 * \brief Test for the OnePNCModel in combination with the NI model for a conduction problem.
11 *
12 * The simulation domain is a tube with an elevated temperature on the left hand side.
13 */
14 #ifndef DUMUX_1P2CNI_CONDUCTION_TEST_PROBLEM_HH
15 #define DUMUX_1P2CNI_CONDUCTION_TEST_PROBLEM_HH
16
17 #include <dumux/common/properties.hh>
18 #include <dumux/common/parameters.hh>
19
20 #include <dumux/common/boundarytypes.hh>
21 #include <dumux/common/numeqvector.hh>
22 #include <dumux/porousmediumflow/problem.hh>
23
24 #include <dumux/material/components/h2o.hh>
25 namespace Dumux {
26
27 /*!
28 * \ingroup OnePNCTests
29 * \brief Definition of a problem, for the 1pnc problem.
30 *
31 * Nitrogen is dissolved in the water phase and
32 * is transported with the water flow from the left side to the right.
33 *
34 * The model domain is specified in the input file and
35 * we use homogeneous soil properties.
36 * Initially, the domain is filled with pure water.
37 *
38 * At the left side, a Dirichlet condition defines the nitrogen mole fraction.
39 * The water phase flows from the left side to the right if the applied pressure
40 * gradient is >0. The nitrogen is transported with the water flow
41 * and leaves the domain at the boundary, where again Dirichlet boundary
42 * conditions are applied.
43 *
44 * This problem uses the \ref OnePNCModel model.
45 */
46 template <class TypeTag>
47 class OnePTwoCNIConductionProblem : public PorousMediumFlowProblem<TypeTag>
48 {
49 using ParentType = PorousMediumFlowProblem<TypeTag>;
50
51 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
52 using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices;
53 using GridView = typename GetPropType<TypeTag, Properties::GridGeometry>::GridView;
54 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
55 using BoundaryTypes = Dumux::BoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>;
56 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
57 using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
58 using NumEqVector = Dumux::NumEqVector<PrimaryVariables>;
59 using Element = typename GridView::template Codim<0>::Entity;
60 using ThermalConductivityModel = GetPropType<TypeTag, Properties::ThermalConductivityModel>;
61 using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>;
62 using VolumeVariables = GetPropType<TypeTag, Properties::VolumeVariables>;
63 using IapwsH2O = Components::H2O<Scalar>;
64
65 // copy some indices for convenience
66 enum
67 {
68 // indices of the primary variables
69 pressureIdx = Indices::pressureIdx,
70 temperatureIdx = Indices::temperatureIdx,
71
72 N2Idx = FluidSystem::compIdx(FluidSystem::MultiPhaseFluidSystem::N2Idx)
73 };
74
75 //! Property that defines whether mole or mass fractions are used
76 static constexpr bool useMoles = getPropValue<TypeTag, Properties::UseMoles>();
77 static const int dimWorld = GridView::dimensionworld;
78 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
79
80 public:
81 3 OnePTwoCNIConductionProblem(std::shared_ptr<const GridGeometry> gridGeometry)
82
7/18
✓ 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.
9 : ParentType(gridGeometry), temperatureHigh_(300.0)
83 {
84 //initialize fluid system
85
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
3 FluidSystem::init();
86
87 // stating in the console whether mole or mass fractions are used
88 if(useMoles)
89
1/2
✓ Branch 2 taken 3 times.
✗ Branch 3 not taken.
6 std::cout<<"problem uses mole fractions"<<std::endl;
90 else
91 std::cout<<"problem uses mass fractions"<<std::endl;
92
93
4/10
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 3 taken 2 times.
✓ 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.
✗ Branch 10 not taken.
7 temperatureExact_.resize(gridGeometry->numDofs(), 290.0);
94 3 }
95
96 //! Get the analytical temperature
97 const std::vector<Scalar>& getExactTemperature()
98 {
99
1/2
✓ Branch 1 taken 3 times.
✗ Branch 2 not taken.
3 return temperatureExact_;
100 }
101
102 //! Update the analytical temperature
103 72 void updateExactTemperature(const SolutionVector& curSol, Scalar time)
104 {
105 216 const auto someElement = *(elements(this->gridGeometry().gridView()).begin());
106
107 144 auto someElemSol = elementSolution(someElement, curSol, this->gridGeometry());
108 144 const auto someInitSol = initialAtPos(someElement.geometry().center());
109
110
2/4
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 72 times.
✗ Branch 5 not taken.
288 const auto someFvGeometry = localView(this->gridGeometry()).bindElement(someElement);
111
1/2
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
72 const auto someScv = *(scvs(someFvGeometry).begin());
112
113
1/2
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
72 VolumeVariables volVars;
114
1/2
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
72 volVars.update(someElemSol, *this, someElement, someScv);
115
116
2/4
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 72 times.
✗ Branch 5 not taken.
144 const auto porosity = this->spatialParams().porosity(someElement, someScv, someElemSol);
117
1/2
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
72 const auto densityW = volVars.density();
118
3/6
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 72 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 72 times.
✗ Branch 8 not taken.
216 const auto heatCapacityW = IapwsH2O::liquidHeatCapacity(someInitSol[temperatureIdx], someInitSol[pressureIdx]);
119
1/2
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
72 const auto densityS = volVars.solidDensity();
120
1/2
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
72 const auto heatCapacityS = volVars.solidHeatCapacity();
121 72 const auto storage = densityW*heatCapacityW*porosity + densityS*heatCapacityS*(1 - porosity);
122
1/2
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
72 const auto effectiveThermalConductivity = ThermalConductivityModel::effectiveThermalConductivity(volVars);
123 using std::max;
124
1/2
✗ Branch 0 not taken.
✓ Branch 1 taken 72 times.
72 time = max(time, 1e-10);
125
126
2/4
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 72 times.
✗ Branch 5 not taken.
216 auto fvGeometry = localView(this->gridGeometry());
127
9/14
✓ Branch 1 taken 72 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 72 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 72 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 14400 times.
✓ Branch 10 taken 72 times.
✓ Branch 11 taken 14400 times.
✓ Branch 12 taken 72 times.
✓ Branch 14 taken 14400 times.
✗ Branch 15 not taken.
✓ Branch 17 taken 14400 times.
✗ Branch 18 not taken.
14616 for (const auto& element : elements(this->gridGeometry().gridView()))
128 {
129
1/2
✓ Branch 1 taken 14400 times.
✗ Branch 2 not taken.
14400 fvGeometry.bindElement(element);
130
4/4
✓ Branch 0 taken 28800 times.
✓ Branch 1 taken 14400 times.
✓ Branch 2 taken 19200 times.
✓ Branch 3 taken 4800 times.
67200 for (auto&& scv : scvs(fvGeometry))
131 {
132 28800 auto globalIdx = scv.dofIndex();
133 28800 const auto& globalPos = scv.dofPosition();
134 using std::erf;
135 using std::sqrt;
136 57600 temperatureExact_[globalIdx] = temperatureHigh_ + (someInitSol[temperatureIdx] - temperatureHigh_)
137 86400 *erf(0.5*sqrt(globalPos[0]*globalPos[0]*storage/time/effectiveThermalConductivity));
138
139 }
140 }
141 72 }
142
143 /*!
144 * \name Problem parameters
145 */
146 // \{
147
148 /*!
149 * \brief Returns the temperature within the domain [K].
150 *
151 * This problem assumes a temperature of 20 degrees Celsius.
152 */
153 Scalar temperature() const
154 { return 273.15 + 20; } // in [K]
155
156 // \}
157
158 /*!
159 * \name Boundary conditions
160 */
161 // \{
162
163 /*!
164 * \brief Specifies which kind of boundary condition should be
165 * used for which equation on a given boundary segment.
166 *
167 * \param globalPos The position for which the bc type should be evaluated
168 */
169 925662 BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
170 {
171 925662 BoundaryTypes values;
172
173
14/14
✓ Branch 0 taken 924323 times.
✓ Branch 1 taken 1339 times.
✓ Branch 2 taken 924323 times.
✓ Branch 3 taken 1339 times.
✓ Branch 4 taken 922984 times.
✓ Branch 5 taken 1339 times.
✓ Branch 6 taken 922984 times.
✓ Branch 7 taken 1339 times.
✓ Branch 8 taken 922984 times.
✓ Branch 9 taken 1339 times.
✓ Branch 10 taken 922984 times.
✓ Branch 11 taken 1339 times.
✓ Branch 12 taken 922984 times.
✓ Branch 13 taken 1339 times.
1851324 if(globalPos[0] < eps_ || globalPos[0] > this->gridGeometry().bBoxMax()[0] - eps_)
174 {
175 values.setAllDirichlet();
176 }
177 else
178 {
179 values.setAllNeumann();
180 }
181
182 925662 return values;
183 }
184
185 /*!
186 * \brief Evaluates the boundary conditions for a Dirichlet boundary segment.
187 *
188 * \param globalPos The position for which the bc type should be evaluated
189 */
190 PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
191 {
192 1412 PrimaryVariables values = initial_(globalPos);
193
194 // condition for the temperature at left boundary
195
4/8
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 207 times.
✓ Branch 5 taken 207 times.
✓ Branch 6 taken 146 times.
✓ Branch 7 taken 146 times.
706 if (globalPos[0] < eps_)
196 706 values[temperatureIdx] = temperatureHigh_;
197
198 return values;
199 }
200
201 /*!
202 * \brief Evaluates the boundary conditions for a Neumann boundary segment.
203 */
204 NumEqVector neumannAtPos(const GlobalPosition &globalPos) const
205
2/2
✓ Branch 0 taken 321600 times.
✓ Branch 1 taken 213328 times.
2470888 { return NumEqVector(0.0); }
206
207 // \}
208
209 /*!
210 * \name Volume terms
211 */
212 // \{
213
214 /*!
215 * \brief Evaluates the source term for all phases within a given
216 * sub-control volume.
217 *
218 * For this method, the \a priVars parameter stores the rate mass
219 * of a component is generated or annihilated per volume
220 * unit. Positive values mean that mass is created, negative ones
221 * mean that it vanishes.
222 *
223 * The units must be according to either using mole or mass fractions (mole/(m^3*s) or kg/(m^3*s)).
224 */
225 NumEqVector sourceAtPos(const GlobalPosition &globalPos) const
226 804000 { return NumEqVector(0.0); }
227
228 /*!
229 * \brief Evaluates the initial value for a control volume.
230 *
231 * \param globalPos The position for which the initial condition should be evaluated
232 *
233 * For this method, the \a values parameter stores primary
234 * variables.
235 */
236 PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
237 874 { return initial_(globalPos); }
238
239 // \}
240 private:
241
242 // the internal method for the initial condition
243 PrimaryVariables initial_(const GlobalPosition &globalPos) const
244 {
245 1580 PrimaryVariables priVars;
246
6/12
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 207 times.
✓ Branch 5 taken 207 times.
✓ Branch 6 taken 146 times.
✓ Branch 7 taken 146 times.
✓ Branch 8 taken 24 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 48 times.
✗ Branch 11 not taken.
1580 priVars[pressureIdx] = 1e5; // initial condition for the pressure
247
6/12
✗ Branch 0 not taken.
✗ Branch 1 not taken.
✗ Branch 2 not taken.
✗ Branch 3 not taken.
✓ Branch 4 taken 207 times.
✓ Branch 5 taken 207 times.
✓ Branch 6 taken 146 times.
✓ Branch 7 taken 146 times.
✓ Branch 8 taken 24 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 48 times.
✗ Branch 11 not taken.
1580 priVars[N2Idx] = 1e-5; // initial condition for the N2 molefraction
248
11/22
✗ 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.
✓ Branch 8 taken 207 times.
✓ Branch 9 taken 207 times.
✓ Branch 10 taken 207 times.
✓ Branch 11 taken 207 times.
✓ Branch 12 taken 146 times.
✓ Branch 13 taken 146 times.
✓ Branch 14 taken 146 times.
✓ Branch 15 taken 170 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 24 times.
✓ Branch 19 taken 48 times.
✗ Branch 20 not taken.
✓ Branch 22 taken 48 times.
✗ Branch 23 not taken.
3160 priVars[temperatureIdx] = 290.;
249 return priVars;
250 }
251 static constexpr Scalar eps_ = 1e-6;
252 Scalar temperatureHigh_;
253 std::vector<Scalar> temperatureExact_;
254 };
255
256 } // end namespace Dumux
257
258 #endif
259