GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/test/porousmediumflow/1pnc/nonequilibrium/problem.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 27 37 73.0%
Functions: 4 14 28.6%
Branches: 27 48 56.2%

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 nonequilibrium model for a thermal nonequilibrium problem problem.
11 *
12 * The simulation domain is a tube where warmer water is injected at the right 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
26 namespace Dumux {
27
28 /*!
29 * \ingroup OnePNCTests
30 * \brief Definition of a problem, for the 1pnc problem.
31 *
32 * The model domain is specified in the input file and
33 * we use homogeneous soil properties.
34 * Initially, the domain is filled with water and a specified nitrogen fraction
35 *
36 * At the right side warmer water is injected via a Neumann boundary and at the left side
37 * Dirichlet values are set to the initial conditions.
38 *
39 * This problem uses the \ref OnePNCModel model.
40 */
41 template <class TypeTag>
42 class OnePTwoCThermalNonequilibriumProblem : public PorousMediumFlowProblem<TypeTag>
43 {
44 using ParentType = PorousMediumFlowProblem<TypeTag>;
45
46 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
47 using Indices = typename GetPropType<TypeTag, Properties::ModelTraits>::Indices;
48 using GridView = typename GetPropType<TypeTag, Properties::GridGeometry>::GridView;
49 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
50 using BoundaryTypes = Dumux::BoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>;
51 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
52 using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
53 using FVElementGeometry = typename GetPropType<TypeTag, Properties::GridGeometry>::LocalView;
54 using NumEqVector = Dumux::NumEqVector<PrimaryVariables>;
55
56 using GridVariables = GetPropType<TypeTag, Properties::GridVariables>;
57 using ElementVolumeVariables = typename GridVariables::GridVolumeVariables::LocalView;
58 using ElementFluxVariablesCache = typename GridVariables::GridFluxVariablesCache::LocalView;
59
60 using SubControlVolumeFace = typename FVElementGeometry::SubControlVolumeFace;
61 using Element = typename GridView::template Codim<0>::Entity;
62 using SolutionVector = GetPropType<TypeTag, Properties::SolutionVector>;
63 using VolumeVariables = GetPropType<TypeTag, Properties::VolumeVariables>;
64 using IapwsH2O = Components::H2O<Scalar>;
65
66 // copy some indices for convenience
67 enum
68 {
69 // indices of the primary variables
70 pressureIdx = Indices::pressureIdx,
71 temperatureIdx = Indices::temperatureIdx,
72 temperatureSolidIdx = Indices::temperatureSolidIdx,
73
74 // component indices
75 H2OIdx = FluidSystem::compIdx(FluidSystem::MultiPhaseFluidSystem::H2OIdx),
76 N2Idx = FluidSystem::compIdx(FluidSystem::MultiPhaseFluidSystem::N2Idx),
77
78 // indices of the equations
79 contiH2OEqIdx = Indices::conti0EqIdx + H2OIdx,
80 contiN2EqIdx = Indices::conti0EqIdx + N2Idx,
81 energyEqIdx = Indices::energyEqIdx
82 };
83
84 //! Property that defines whether mole or mass fractions are used
85 static constexpr bool useMoles = getPropValue<TypeTag, Properties::UseMoles>();
86 static const int dimWorld = GridView::dimensionworld;
87 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
88
89 public:
90 2 OnePTwoCThermalNonequilibriumProblem(std::shared_ptr<const GridGeometry> gridGeometry)
91
7/20
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 2 times.
✗ Branch 5 not taken.
✓ Branch 8 taken 2 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 2 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 2 times.
✓ Branch 15 taken 2 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 2 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.
6 : ParentType(gridGeometry)
92 {
93 //initialize fluid system
94
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 FluidSystem::init();
95
96 // stating in the console whether mole or mass fractions are used
97 if(useMoles)
98
1/2
✓ Branch 2 taken 2 times.
✗ Branch 3 not taken.
4 std::cout<<"problem uses mole fractions"<<std::endl;
99 else
100 std::cout<<"problem uses mass fractions"<<std::endl;
101
102
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 darcyVelocity_ = getParam<Scalar>("Problem.DarcyVelocity");
103 2 }
104
105 void setGridVariables(std::shared_ptr<GridVariables> gridVariables)
106
1/2
✓ Branch 1 taken 2 times.
✗ Branch 2 not taken.
2 { gridVariables_ = gridVariables; }
107
108 const GridVariables& gridVariables() const
109 2783646 { return *gridVariables_; }
110
111 /*!
112 * \name Problem parameters
113 */
114 // \{
115
116 // \}
117
118 /*!
119 * \name Boundary conditions
120 */
121 // \{
122
123 /*!
124 * \brief Specifies which kind of boundary condition should be
125 * used for which equation on a given boundary segment.
126 *
127 * \param globalPos The position for which the bc type should be evaluated
128 */
129 BoundaryTypes boundaryTypesAtPos(const GlobalPosition &globalPos) const
130 {
131 BoundaryTypes values;
132
133 if(globalPos[0] < eps_)
134 values.setAllDirichlet();
135 else
136 values.setAllNeumann();
137
138 return values;
139 }
140
141 /*!
142 * \brief Evaluates the boundary conditions for a Dirichlet boundary segment.
143 *
144 * \param globalPos The position for which the bc type should be evaluated
145 */
146 PrimaryVariables dirichletAtPos(const GlobalPosition &globalPos) const
147 {
148 721 PrimaryVariables values = initial_(globalPos);
149
150 return values;
151 }
152
153 /*!
154 * \brief Evaluates the boundary conditions for a Neumann boundary segment.
155 *
156 * This is the method for the case where the Neumann condition is
157 * potentially solution dependent and requires some quantities that
158 * are specific to the fully-implicit method.
159 *
160 * \param element The finite element
161 * \param fvGeometry The finite-volume geometry
162 * \param elemVolVars All volume variables for the element
163 * \param elemFluxVarsCache Flux variables caches for all faces in stencil
164 * \param scvf The sub-control volume face
165 *
166 * For this method, the \a values parameter stores the flux
167 * in normal direction of each phase. Negative values mean influx.
168 * E.g. for the mass balance that would the mass flux in \f$ [ kg / (m^2 \cdot s)] \f$.
169 */
170 2756641 NumEqVector neumann(const Element& element,
171 const FVElementGeometry& fvGeometry,
172 const ElementVolumeVariables& elemVolVars,
173 const ElementFluxVariablesCache& elemFluxVarsCache,
174 const SubControlVolumeFace& scvf) const
175 {
176 2756641 NumEqVector flux(0.0);
177
2/2
✓ Branch 0 taken 423289 times.
✓ Branch 1 taken 2333352 times.
2756641 const auto& globalPos = scvf.ipGlobal();
178
4/4
✓ Branch 0 taken 423289 times.
✓ Branch 1 taken 2333352 times.
✓ Branch 2 taken 423289 times.
✓ Branch 3 taken 2333352 times.
5513282 const auto& scv = fvGeometry.scv(scvf.insideScvIdx());
179
180
10/10
✓ Branch 0 taken 6889 times.
✓ Branch 1 taken 2749752 times.
✓ Branch 2 taken 6889 times.
✓ Branch 3 taken 2749752 times.
✓ Branch 4 taken 6889 times.
✓ Branch 5 taken 2749752 times.
✓ Branch 6 taken 6889 times.
✓ Branch 7 taken 2749752 times.
✓ Branch 8 taken 6889 times.
✓ Branch 9 taken 2749752 times.
13783205 if (globalPos[0] > this->gridGeometry().bBoxMax()[0] - eps_)
181 {
182 18585 flux[contiH2OEqIdx] = -darcyVelocity_*elemVolVars[scv].molarDensity();
183 30281 flux[contiN2EqIdx] = -darcyVelocity_*elemVolVars[scv].molarDensity()*elemVolVars[scv].moleFraction(0, N2Idx);
184 6889 flux[energyEqIdx] = -darcyVelocity_
185 12737 *elemVolVars[scv].density()
186 18585 *IapwsH2O::liquidEnthalpy(305, elemVolVars[scv].pressure());
187 }
188 2756641 return flux;
189 }
190
191 // \}
192
193 /*!
194 * \name Volume terms
195 */
196 // \{
197
198 /*!
199 * \brief Evaluates the source term for all phases within a given
200 * sub-control volume.
201 *
202 * For this method, the \a priVars parameter stores the rate mass
203 * of a component is generated or annihilated per volume
204 * unit. Positive values mean that mass is created, negative ones
205 * mean that it vanishes.
206 *
207 * The units must be according to either using mole or mass fractions (mole/(m^3*s) or kg/(m^3*s)).
208 */
209 NumEqVector sourceAtPos(const GlobalPosition &globalPos) const
210 2505200 { return NumEqVector(0.0); }
211
212 /*!
213 * \brief Evaluates the initial value for a control volume.
214 *
215 * \param globalPos The position for which the initial condition should be evaluated
216 *
217 * For this method, the \a values parameter stores primary
218 * variables.
219 */
220 PrimaryVariables initialAtPos(const GlobalPosition &globalPos) const
221 602 { return initial_(globalPos); }
222
223 // \}
224 private:
225
226 // the internal method for the initial condition
227 PrimaryVariables initial_(const GlobalPosition &globalPos) const
228 {
229 1323 PrimaryVariables priVars;
230 1323 priVars[pressureIdx] = 1e5; // initial condition for the pressure
231 1323 priVars[N2Idx] = 1e-5; // initial condition for the N2 molefraction
232 1323 priVars[temperatureIdx] = 285.;
233 2646 priVars[temperatureSolidIdx] = 285.;
234 return priVars;
235 }
236 static constexpr Scalar eps_ = 1e-6;
237 std::shared_ptr<GridVariables> gridVariables_;
238 Scalar darcyVelocity_;
239 };
240
241 } // end namespace Dumux
242
243 #endif
244