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 |