GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/dumux/porousmediumflow/3p/volumevariables.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 52 53 98.1%
Functions: 12 12 100.0%
Branches: 19 26 73.1%

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 ThreePModel
10 * \brief Contains the quantities which are constant within a finite volume in the three-phase model.
11 */
12
13 #ifndef DUMUX_3P_VOLUME_VARIABLES_HH
14 #define DUMUX_3P_VOLUME_VARIABLES_HH
15
16 #include <dumux/material/constants.hh>
17 #include <dumux/material/fluidstates/immiscible.hh>
18 #include <dumux/porousmediumflow/volumevariables.hh>
19 #include <dumux/porousmediumflow/nonisothermal/volumevariables.hh>
20 #include <dumux/material/solidstates/updatesolidvolumefractions.hh>
21
22 namespace Dumux {
23
24 /*!
25 * \ingroup ThreePModel
26 * \brief Contains the quantities which are constant within a finite volume in the three-phase model.
27 */
28 template <class Traits>
29 class ThreePVolumeVariables
30 : public PorousMediumFlowVolumeVariables<Traits>
31 , public EnergyVolumeVariables<Traits, ThreePVolumeVariables<Traits> >
32 {
33 using ParentType = PorousMediumFlowVolumeVariables<Traits>;
34 using EnergyVolVars = EnergyVolumeVariables<Traits, ThreePVolumeVariables<Traits> >;
35
36 using Scalar = typename Traits::PrimaryVariables::value_type;
37 using PermeabilityType = typename Traits::PermeabilityType;
38 using Idx = typename Traits::ModelTraits::Indices;
39 using FS = typename Traits::FluidSystem;
40 static constexpr int numFluidComps = ParentType::numFluidComponents();
41
42 enum {
43 wPhaseIdx = FS::wPhaseIdx,
44 gPhaseIdx = FS::gPhaseIdx,
45 nPhaseIdx = FS::nPhaseIdx,
46
47 swIdx = Idx::swIdx,
48 snIdx = Idx::snIdx,
49 pressureIdx = Idx::pressureIdx
50 };
51
52 public:
53 //! Export fluid state type
54 using FluidState = typename Traits::FluidState;
55 //! Export fluid system type
56 using FluidSystem = typename Traits::FluidSystem;
57 //! Export the indices
58 using Indices = Idx;
59 //! Export type of solid state
60 using SolidState = typename Traits::SolidState;
61 //! Export type of solid system
62 using SolidSystem = typename Traits::SolidSystem;
63
64 /*!
65 * \brief Updates all quantities for a given control volume.
66 *
67 * \param elemSol A vector containing all primary variables connected to the element
68 * \param problem The object specifying the problem which ought to
69 * be simulated
70 * \param element An element which contains part of the control volume
71 * \param scv The sub control volume
72 */
73 template<class ElemSol, class Problem, class Element, class Scv>
74 963240 void update(const ElemSol &elemSol,
75 const Problem &problem,
76 const Element &element,
77 const Scv& scv)
78 {
79 963240 ParentType::update(elemSol, problem, element, scv);
80 963240 completeFluidState(elemSol, problem, element, scv, fluidState_, solidState_);
81
82 963240 const auto sw = fluidState_.saturation(wPhaseIdx);
83 963240 const auto sn = fluidState_.saturation(nPhaseIdx);
84
85 // mobilities
86 1926480 const auto fluidMatrixInteraction = problem.spatialParams().fluidMatrixInteraction(element, scv, elemSol);
87
2/2
✓ Branch 0 taken 2889720 times.
✓ Branch 1 taken 963240 times.
3852960 for (int phaseIdx = 0; phaseIdx < ParentType::numFluidPhases(); ++phaseIdx)
88 {
89 2889720 mobility_[phaseIdx] = fluidMatrixInteraction.kr(phaseIdx, sw, sn)
90 2889720 / fluidState_.viscosity(phaseIdx);
91 }
92
93 // porosity
94
2/2
✓ Branch 0 taken 219274 times.
✓ Branch 1 taken 33626 times.
963240 updateSolidVolumeFractions(elemSol, problem, element, scv, solidState_, numFluidComps);
95 963240 EnergyVolVars::updateSolidEnergyParams(elemSol, problem, element, scv, solidState_);
96
4/4
✓ Branch 0 taken 219274 times.
✓ Branch 1 taken 33626 times.
✓ Branch 2 taken 219274 times.
✓ Branch 3 taken 33626 times.
2636820 permeability_ = problem.spatialParams().permeability(element, scv, elemSol);
97 963240 EnergyVolVars::updateEffectiveThermalConductivity();
98 963240 }
99
100 /*!
101 * \brief Sets complete fluid state.
102 *
103 * \param elemSol A vector containing all primary variables connected to the element
104 * \param problem The object specifying the problem which ought to
105 * be simulated
106 * \param element An element which contains part of the control volume
107 * \param scv The sub-control volume
108 * \param fluidState A container with the current (physical) state of the fluid
109 * \param solidState A container with the current (physical) state of the solid
110 *
111 * Set temperature, saturations, capillary pressures, viscosities, densities and enthalpies.
112 */
113 template<class ElemSol, class Problem, class Element, class Scv>
114 963240 void completeFluidState(const ElemSol& elemSol,
115 const Problem& problem,
116 const Element& element,
117 const Scv& scv,
118 FluidState& fluidState,
119 SolidState& solidState)
120 {
121 963240 EnergyVolVars::updateTemperature(elemSol, problem, element, scv, fluidState, solidState);
122
123 963240 const auto& priVars = elemSol[scv.localDofIndex()];
124
125 1926480 const auto fluidMatrixInteraction = problem.spatialParams().fluidMatrixInteraction(element, scv, elemSol);
126
127 963240 const Scalar sw = priVars[swIdx];
128 963240 const Scalar sn = priVars[snIdx];
129 963240 const Scalar sg = 1.0 - sw - sn;
130
131 963240 fluidState.setSaturation(wPhaseIdx, sw);
132 963240 fluidState.setSaturation(gPhaseIdx, sg);
133 963240 fluidState.setSaturation(nPhaseIdx, sn);
134
135 /* now the pressures */
136 963240 const Scalar pg = priVars[pressureIdx];
137
138 // calculate capillary pressures
139 963240 const Scalar pcgw = fluidMatrixInteraction.pcgw(sw, sn);
140 963240 const Scalar pcnw = fluidMatrixInteraction.pcnw(sw, sn);
141 963240 const Scalar pcgn = fluidMatrixInteraction.pcgn(sw, sn);
142
143
2/2
✓ Branch 0 taken 962772 times.
✓ Branch 1 taken 468 times.
963240 const Scalar pcAlpha = fluidMatrixInteraction.pcAlpha(sw, sn);
144 963240 const Scalar pcNW1 = 0.0; // TODO: this should be possible to assign in the problem file
145
146 963240 const Scalar pn = pg- pcAlpha * pcgn - (1.0 - pcAlpha)*(pcgw - pcNW1);
147 963240 const Scalar pw = pn - pcAlpha * pcnw - (1.0 - pcAlpha)*pcNW1;
148
149 963240 fluidState.setPressure(wPhaseIdx, pw);
150 963240 fluidState.setPressure(gPhaseIdx, pg);
151 963240 fluidState.setPressure(nPhaseIdx, pn);
152
153 typename FluidSystem::ParameterCache paramCache;
154 963240 paramCache.updateAll(fluidState);
155
156
2/2
✓ Branch 0 taken 2889720 times.
✓ Branch 1 taken 963240 times.
3852960 for (int phaseIdx = 0; phaseIdx < ParentType::numFluidPhases(); ++phaseIdx)
157 {
158 // compute and set the viscosity
159 2889720 const Scalar mu = FluidSystem::viscosity(fluidState, paramCache, phaseIdx);
160 2889720 fluidState.setViscosity(phaseIdx,mu);
161
162 // compute and set the density
163 2889720 const Scalar rho = FluidSystem::density(fluidState, paramCache, phaseIdx);
164 2889720 fluidState.setDensity(phaseIdx, rho);
165
166 // compute and set the enthalpy
167 2889720 const Scalar h = EnergyVolVars::enthalpy(fluidState, paramCache, phaseIdx);
168 5779440 fluidState.setEnthalpy(phaseIdx, h);
169 }
170 963240 }
171
172 /*!
173 * \brief Returns the phase state for the control-volume.
174 */
175 const FluidState &fluidState() const
176 21746508 { return fluidState_; }
177
178 /*!
179 * \brief Returns the phase state for the control volume.
180 */
181 const SolidState &solidState() const
182
2/4
✓ Branch 2 taken 48 times.
✗ Branch 3 not taken.
✓ Branch 5 taken 48 times.
✗ Branch 6 not taken.
710656 { return solidState_; }
183
184 /*!
185 * \brief Returns the effective saturation of a given phase within
186 * the control volume.
187 *
188 * \param phaseIdx The phase index
189 */
190 Scalar saturation(const int phaseIdx) const
191 36886512 { return fluidState_.saturation(phaseIdx); }
192
193 /*!
194 * \brief Returns the mass density of a given phase within the
195 * control volume.
196 *
197 * \param phaseIdx The phase index
198 */
199 Scalar density(const int phaseIdx) const
200
3/6
✓ Branch 1 taken 48 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 134 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 86 times.
✗ Branch 8 not taken.
100751554 { return fluidState_.density(phaseIdx); }
201
202 /*!
203 * \brief Returns the effective pressure of a given phase within
204 * the control volume.
205 *
206 * \param phaseIdx The phase index
207 */
208 Scalar pressure(const int phaseIdx) const
209 53916404 { return fluidState_.pressure(phaseIdx); }
210
211 /*!
212 * \brief Returns temperature inside the sub-control volume.
213 *
214 * Note that we assume thermodynamic equilibrium, i.e. the
215 * temperatures of the rock matrix and of all fluid phases are
216 * identical.
217 */
218 Scalar temperature() const
219
2/4
✓ Branch 0 taken 428180 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 428180 times.
✗ Branch 3 not taken.
13593680 { return fluidState_.temperature(/*phaseIdx=*/0); }
220
221 /*!
222 * \brief Returns the effective mobility of a given phase within
223 * the control volume.
224 *
225 * \param phaseIdx The phase index
226 */
227 Scalar mobility(const int phaseIdx) const
228 {
229 26979060 return mobility_[phaseIdx];
230 }
231
232 /*!
233 * \brief Returns the effective capillary pressure within the control volume.
234 */
235 Scalar capillaryPressure() const
236 { return fluidState_.capillaryPressure(); }
237
238 /*!
239 * \brief Returns the average porosity within the control volume.
240 */
241 Scalar porosity() const
242 45022056 { return solidState_.porosity(); }
243
244 /*!
245 * \brief Returns the permeability within the control volume in \f$[m^2]\f$.
246 */
247 const PermeabilityType& permeability() const
248 { return permeability_; }
249
250 protected:
251 FluidState fluidState_;
252 SolidState solidState_;
253
254
255 private:
256 PermeabilityType permeability_;
257 Scalar mobility_[ParentType::numFluidPhases()];
258 };
259
260 } // end namespace Dumux
261
262 #endif
263