GCC Code Coverage Report


Directory: ../../../builds/dumux-repositories/
File: /builds/dumux-repositories/dumux/test/porousmediumflow/2pnc/surfactant/problem.hh
Date: 2024-09-21 20:52:54
Exec Total Coverage
Lines: 42 46 91.3%
Functions: 2 4 50.0%
Branches: 48 78 61.5%

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 TwoPNCTests
10 * \brief Definition of a problem for a surfactant model.
11 */
12 #ifndef DUMUX_TEST_2P3C_SURFACTANT_PROBLEM_HH
13 #define DUMUX_TEST_2P3C_SURFACTANT_PROBLEM_HH
14
15 #include <dumux/common/properties.hh>
16 #include <dumux/common/parameters.hh>
17 #include <dumux/common/boundarytypes.hh>
18 #include <dumux/common/numeqvector.hh>
19
20 #include <dumux/porousmediumflow/problem.hh>
21 #include <dumux/porousmediumflow/2p/formulation.hh>
22
23 namespace Dumux {
24
25 /*!
26 * \ingroup TwoPNCTests
27 * \brief Problem for a surfactant model.
28 *
29 */
30
31 template <class TypeTag>
32 class TestSurfactantProblem
33 : public PorousMediumFlowProblem<TypeTag>
34 {
35 using ParentType = PorousMediumFlowProblem<TypeTag>;
36 using Scalar = GetPropType<TypeTag, Properties::Scalar>;
37 using FluidSystem = GetPropType<TypeTag, Properties::FluidSystem>;
38
39 using GridGeometry = GetPropType<TypeTag, Properties::GridGeometry>;
40 using GridView = typename GridGeometry::GridView;
41 using FVElementGeometry = typename GridGeometry::LocalView;
42 using SubControlVolume = typename GridGeometry::SubControlVolume;
43 using SubControlVolumeFace = typename GridGeometry::SubControlVolumeFace;
44
45 static constexpr int dim = GridView::dimension;
46
47 using PrimaryVariables = GetPropType<TypeTag, Properties::PrimaryVariables>;
48 using NumEqVector = Dumux::NumEqVector<PrimaryVariables>;
49 using BoundaryTypes = Dumux::BoundaryTypes<GetPropType<TypeTag, Properties::ModelTraits>::numEq()>;
50 using Element = typename GridView::Traits::template Codim<0>::Entity;
51 using GlobalPosition = typename Element::Geometry::GlobalCoordinate;
52 using FluidState = GetPropType<TypeTag, Properties::FluidState>;
53
54 public:
55 1 TestSurfactantProblem(std::shared_ptr<const GridGeometry> gridGeometry)
56
8/22
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
✓ Branch 8 taken 1 times.
✗ Branch 9 not taken.
✓ Branch 10 taken 1 times.
✗ Branch 11 not taken.
✗ Branch 12 not taken.
✓ Branch 13 taken 1 times.
✓ Branch 15 taken 1 times.
✗ Branch 16 not taken.
✓ Branch 18 taken 1 times.
✗ Branch 19 not taken.
✓ Branch 21 taken 1 times.
✗ Branch 22 not taken.
✗ Branch 23 not taken.
✗ Branch 24 not taken.
✗ Branch 25 not taken.
✗ Branch 26 not taken.
✗ Branch 29 not taken.
✗ Branch 30 not taken.
3 : ParentType(gridGeometry)
57 {
58
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
1 const auto maxSurfC = getParam<Scalar>("Problem.InjectionSurfactantConcentration");
59 1 injectionFluidState_.setMoleFraction(0, FluidSystem::surfactantCompIdx, maxSurfC);
60
2/4
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 1 times.
✗ Branch 5 not taken.
2 this->spatialParams().setMaxSurfactantConcentration(maxSurfC);
61
62
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
1 initialPressure_ = getParam<Scalar>("Problem.InitialPressure");
63
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
1 initialSw_ = getParam<Scalar>("Problem.InitialSw");
64
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
1 productionWellPressure_ = getParam<Scalar>("Problem.ProductionWellPressure");
65
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
1 injectionWellPressure_ = getParam<Scalar>("Problem.InjectionWellPressure");
66
67
2/4
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 1 times.
1 name_ = getParam<std::string>("Problem.Name");
68 1 }
69
70 /*!
71 * \brief The problem name.
72 *
73 * This is used as a prefix for files generated by the simulation.
74 */
75 const std::string& name() const
76
1/2
✓ Branch 1 taken 1 times.
✗ Branch 2 not taken.
1 { return name_; }
77
78 /*!
79 * \brief Specifies which kind of boundary condition should be
80 * used for which equation on a given boundary segment
81 */
82
83 BoundaryTypes boundaryTypes(const Element& element, const SubControlVolume& scv) const
84 {
85 27048 BoundaryTypes bcTypes;
86
3/6
✓ Branch 0 taken 12084 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 12084 times.
✗ Branch 3 not taken.
✗ Branch 4 not taken.
✓ Branch 5 taken 2880 times.
27048 bcTypes.setAllNeumann();
87 return bcTypes;
88 }
89
90 /*!
91 * \brief Evaluates the boundary conditions for a Neumann boundary segment.
92 */
93
94 template<class ElementVolumeVariables, class ElementFluxVariablesCache>
95 84588 NumEqVector neumann(const Element& element,
96 const FVElementGeometry& fvGeometry,
97 const ElementVolumeVariables& elemVolVars,
98 const ElementFluxVariablesCache& elemFluxVarsCache,
99 const SubControlVolumeFace& scvf) const
100 {
101 84588 NumEqVector values(0.0);
102
103
4/4
✓ Branch 0 taken 42294 times.
✓ Branch 1 taken 42294 times.
✓ Branch 2 taken 42294 times.
✓ Branch 3 taken 42294 times.
169176 const auto& volVars = elemVolVars[scvf.insideScvIdx()];
104
2/2
✓ Branch 0 taken 42294 times.
✓ Branch 1 taken 42294 times.
84588 const Scalar pw = volVars.pressure(0);
105
4/4
✓ Branch 0 taken 42294 times.
✓ Branch 1 taken 42294 times.
✓ Branch 2 taken 42294 times.
✓ Branch 3 taken 42294 times.
169176 const auto K = this->spatialParams().permeability(element);
106
107
8/8
✓ Branch 0 taken 42294 times.
✓ Branch 1 taken 42294 times.
✓ Branch 2 taken 42294 times.
✓ Branch 3 taken 42294 times.
✓ Branch 4 taken 42294 times.
✓ Branch 5 taken 42294 times.
✓ Branch 6 taken 42294 times.
✓ Branch 7 taken 42294 times.
338352 const auto middleX = 0.5*(this->gridGeometry().bBoxMin()[0] + this->gridGeometry().bBoxMax()[0]);
108
4/4
✓ Branch 0 taken 42294 times.
✓ Branch 1 taken 42294 times.
✓ Branch 2 taken 42294 times.
✓ Branch 3 taken 42294 times.
169176 if (scvf.ipGlobal()[0] < middleX)
109 {
110 using std::max;
111
112
2/4
✓ Branch 0 taken 42294 times.
✗ Branch 1 not taken.
✓ Branch 2 taken 42294 times.
✗ Branch 3 not taken.
84588 const Scalar dp0dn = max(0.0, (injectionWellPressure_ - volVars.pressure(0)));
113
114 42294 const Scalar injectionViscosity = FluidSystem::viscosity(injectionFluidState_, 0);
115
116 // m3 / (m2 s) or kg / (m2 s)
117 84588 const auto flux = K * (1 / injectionViscosity) * dp0dn; // volume per second
118 42294 const auto moleflux = flux * FluidSystem::molarDensity(injectionFluidState_, 0);
119
120
2/2
✓ Branch 0 taken 126882 times.
✓ Branch 1 taken 42294 times.
169176 for (int j = 0; j < FluidSystem::numComponents; j++)
121 761292 values[j] = -moleflux * injectionFluidState_.moleFraction(0, j);
122
123 }
124
125 else
126 {
127 using std::min;
128
1/2
✓ Branch 0 taken 42294 times.
✗ Branch 1 not taken.
42294 const auto dp0dn = min(0.0, productionWellPressure_ - pw);
129 42294 const auto dp1dn = dp0dn;
130
131 169176 const Scalar fluxOilphase = K * volVars.mobility(1) * dp1dn; // volume per second
132 169176 const Scalar fluxWaterphase = K * volVars.mobility(0) * dp0dn; // volume per second
133
134 84588 values[FluidSystem::oilCompIdx] = -volVars.molarDensity(1) * fluxOilphase;
135 84588 values[FluidSystem::waterCompIdx] = -volVars.molarDensity(0) * fluxWaterphase;
136
137 42294 const Scalar f = volVars.moleFraction(0, FluidSystem::surfactantCompIdx);
138 253764 values[FluidSystem::surfactantCompIdx] = - volVars.molarDensity(0) * volVars.mobility(0) * K * dp0dn * f;
139 126882 values[FluidSystem::waterCompIdx] -= values[FluidSystem::surfactantCompIdx];
140 }
141
142 84588 return values;
143 }
144
145 /*!
146 * \brief Evaluates the initial values for a control volume.
147 *
148 * \param globalPos The global position
149 */
150
151 PrimaryVariables initialAtPos(const GlobalPosition& globalPos) const
152 {
153 using ModelTraits = GetPropType<TypeTag, Properties::ModelTraits>;
154 static_assert(ModelTraits::priVarFormulation() == TwoPFormulation::p0s1, "must have p0s1 formulation");
155
156 101 PrimaryVariables values(0);
157
158 101 values[0] = initialPressure_;
159 101 values[1] = 1.0 - initialSw_;
160 101 values.setState(ModelTraits::Indices::bothPhases);
161 202 values[FluidSystem::surfactantCompIdx] = 0;
162
163 return values;
164 }
165
166 private:
167 Scalar productionWellPressure_;
168 Scalar injectionWellPressure_;
169 Scalar initialPressure_;
170 Scalar initialSw_;
171
172 FluidState injectionFluidState_;
173 std::string name_;
174 };
175
176 } // end namespace Dumux
177
178 #endif
179