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-FileCopyrightText: 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 CO2Tests | ||
10 | * \brief Definition of the spatial parameters for the heterogeneous | ||
11 | * problem which uses the non-isothermal or isothermal CO2 | ||
12 | * fully implicit model. | ||
13 | */ | ||
14 | |||
15 | #ifndef DUMUX_HETEROGENEOUS_SPATIAL_PARAMS_HH | ||
16 | #define DUMUX_HETEROGENEOUS_SPATIAL_PARAMS_HH | ||
17 | |||
18 | #include <dumux/io/grid/griddata.hh> | ||
19 | #include <dumux/porousmediumflow/fvspatialparamsmp.hh> | ||
20 | #include <dumux/material/fluidmatrixinteractions/2p/brookscorey.hh> | ||
21 | |||
22 | #include <dumux/porousmediumflow/co2/model.hh> | ||
23 | |||
24 | namespace Dumux { | ||
25 | |||
26 | /*! | ||
27 | * \ingroup CO2Tests | ||
28 | * \brief Definition of the spatial parameters for the heterogeneous | ||
29 | * problem which uses the non-isothermal or isothermal CO2 | ||
30 | * fully implicit model. | ||
31 | */ | ||
32 | template<class GridGeometry, class Scalar> | ||
33 | class HeterogeneousSpatialParams | ||
34 | : public FVPorousMediumFlowSpatialParamsMP<GridGeometry, | ||
35 | Scalar, | ||
36 | HeterogeneousSpatialParams<GridGeometry, Scalar>> | ||
37 | { | ||
38 | using Grid = typename GridGeometry::Grid; | ||
39 | using GridView = typename GridGeometry::GridView; | ||
40 | using FVElementGeometry = typename GridGeometry::LocalView; | ||
41 | using SubControlVolume = typename FVElementGeometry::SubControlVolume; | ||
42 | using Element = typename GridView::template Codim<0>::Entity; | ||
43 | using ParentType = FVPorousMediumFlowSpatialParamsMP<GridGeometry, Scalar, HeterogeneousSpatialParams<GridGeometry, Scalar>>; | ||
44 | |||
45 | using GlobalPosition = typename SubControlVolume::GlobalPosition; | ||
46 | |||
47 | using PcKrSwCurve = FluidMatrix::BrooksCoreyDefault<Scalar>; | ||
48 | |||
49 | static constexpr int dimWorld = GridView::dimensionworld; | ||
50 | |||
51 | public: | ||
52 | using PermeabilityType = Scalar; | ||
53 | |||
54 | 14 | HeterogeneousSpatialParams(std::shared_ptr<const GridGeometry> gridGeometry, | |
55 | std::shared_ptr<const GridData<Grid>> gridData) | ||
56 | : ParentType(gridGeometry) | ||
57 | 28 | , gridData_(gridData) | |
58 |
3/6✓ Branch 2 taken 14 times.
✗ Branch 3 not taken.
✓ Branch 6 taken 14 times.
✗ Branch 7 not taken.
✓ Branch 9 taken 14 times.
✗ Branch 10 not taken.
|
28 | , pcKrSwCurve_("SpatialParams") |
59 | { | ||
60 | // Set the permeability for the layers | ||
61 | 14 | barrierTopK_ = 1e-17; //sqm | |
62 | 14 | barrierMiddleK_ = 1e-15; //sqm | |
63 | 14 | reservoirK_ = 1e-14; //sqm | |
64 | |||
65 | // Set the effective porosity of the layers | ||
66 | 14 | barrierTopPorosity_ = 0.001; | |
67 | 14 | barrierMiddlePorosity_ = 0.05; | |
68 | 14 | reservoirPorosity_ = 0.2; | |
69 | |||
70 |
1/2✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
|
14 | depthBOR_ = getParam<Scalar>("Problem.DepthBOR"); |
71 | 14 | } | |
72 | |||
73 | /*! | ||
74 | * \brief Reads layer information from the grid | ||
75 | */ | ||
76 | 14 | void getParamsFromGrid() | |
77 | { | ||
78 | 14 | const auto& gridView = this->gridGeometry().gridView(); | |
79 | 14 | paramIdx_.resize(gridView.size(0)); | |
80 | |||
81 |
5/8✓ Branch 1 taken 14 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 14 times.
✗ Branch 5 not taken.
✓ Branch 7 taken 37462 times.
✗ Branch 8 not taken.
✓ Branch 9 taken 37462 times.
✓ Branch 10 taken 14 times.
|
74980 | for (const auto& element : elements(gridView)) |
82 | { | ||
83 |
1/2✓ Branch 1 taken 37462 times.
✗ Branch 2 not taken.
|
37462 | const auto eIdx = this->gridGeometry().elementMapper().index(element); |
84 |
2/4✓ Branch 1 taken 37462 times.
✗ Branch 2 not taken.
✓ Branch 4 taken 37462 times.
✗ Branch 5 not taken.
|
37462 | paramIdx_[eIdx] = gridData_->parameters(element)[0]; |
85 | } | ||
86 | 14 | } | |
87 | |||
88 | /*! | ||
89 | * \brief Function for defining the (intrinsic) permeability \f$[m^2]\f$ | ||
90 | * \note It is possibly solution dependent. | ||
91 | * | ||
92 | * \param element The current element | ||
93 | * \param scv The sub-control volume inside the element. | ||
94 | * \param elemSol The solution at the dofs connected to the element. | ||
95 | * | ||
96 | * \return The intrinsic permeability | ||
97 | */ | ||
98 | template<class ElementSolution> | ||
99 | 71188461 | PermeabilityType permeability(const Element& element, | |
100 | const SubControlVolume& scv, | ||
101 | const ElementSolution& elemSol) const | ||
102 | { | ||
103 | // Get the global index of the element | ||
104 | 71188461 | const auto eIdx = this->gridGeometry().elementMapper().index(element); | |
105 | 142376922 | return permeability(eIdx); | |
106 | } | ||
107 | |||
108 | /*! | ||
109 | * \brief Function for defining the (intrinsic) permeability \f$[m^2]\f$ | ||
110 | * \note It is possibly solution dependent. | ||
111 | * | ||
112 | * \param eIdx The element index | ||
113 | */ | ||
114 | 71225131 | PermeabilityType permeability(std::size_t eIdx) const | |
115 | { | ||
116 |
4/4✓ Branch 0 taken 6622188 times.
✓ Branch 1 taken 64566273 times.
✓ Branch 2 taken 3470 times.
✓ Branch 3 taken 33200 times.
|
71225131 | if (paramIdx_[eIdx] == barrierTop_) |
117 | 6625658 | return barrierTopK_; | |
118 |
4/4✓ Branch 0 taken 13643020 times.
✓ Branch 1 taken 50923253 times.
✓ Branch 2 taken 6570 times.
✓ Branch 3 taken 26630 times.
|
64599473 | else if (paramIdx_[eIdx] == barrierMiddle_) |
119 | 13649590 | return barrierMiddleK_; | |
120 | else | ||
121 | 50949883 | return reservoirK_; | |
122 | } | ||
123 | |||
124 | /*! | ||
125 | * \brief Returns the volume fraction of the inert component with index compIdx \f$[-]\f$ | ||
126 | * | ||
127 | * \param element The current element | ||
128 | * \param scv The sub-control volume inside the element. | ||
129 | * \param elemSol The element solution | ||
130 | * \param compIdx The solid component index | ||
131 | * \return The solid volume fraction | ||
132 | */ | ||
133 | template<class SolidSystem, class ElementSolution> | ||
134 | 71188461 | Scalar inertVolumeFraction(const Element& element, | |
135 | const SubControlVolume& scv, | ||
136 | const ElementSolution& elemSol, | ||
137 | int compIdx) const | ||
138 | { | ||
139 | // Get the global index of the element | ||
140 | 71188461 | const auto eIdx = this->gridGeometry().elementMapper().index(element); | |
141 | 142376922 | return inertVolumeFraction(eIdx); | |
142 | } | ||
143 | |||
144 | /*! | ||
145 | * \brief Returns the porosity \f$[-]\f$ | ||
146 | * | ||
147 | * \param eIdx The element index | ||
148 | */ | ||
149 | 71225131 | Scalar inertVolumeFraction(std::size_t eIdx) const | |
150 | { | ||
151 |
4/4✓ Branch 0 taken 6622188 times.
✓ Branch 1 taken 64566273 times.
✓ Branch 2 taken 3470 times.
✓ Branch 3 taken 33200 times.
|
71225131 | if (paramIdx_[eIdx] == barrierTop_) |
152 | 6625658 | return 1- barrierTopPorosity_; | |
153 |
4/4✓ Branch 0 taken 13643020 times.
✓ Branch 1 taken 50923253 times.
✓ Branch 2 taken 6570 times.
✓ Branch 3 taken 26630 times.
|
64599473 | else if (paramIdx_[eIdx] == barrierMiddle_) |
154 | 13649590 | return 1- barrierMiddlePorosity_; | |
155 | else | ||
156 | 50949883 | return 1- reservoirPorosity_; | |
157 | |||
158 | } | ||
159 | |||
160 | |||
161 | /*! | ||
162 | * \brief Returns the fluid-matrix interaction law at a given location | ||
163 | * | ||
164 | * \param globalPos The global coordinates for the given location | ||
165 | */ | ||
166 | 142376922 | auto fluidMatrixInteractionAtPos(const GlobalPosition& globalPos) const | |
167 |
1/2✗ Branch 0 not taken.
✓ Branch 1 taken 71188461 times.
|
142376922 | { return makeFluidMatrixInteraction(pcKrSwCurve_); } |
168 | |||
169 | /*! | ||
170 | * \brief Function for defining which phase is to be considered as the wetting phase. | ||
171 | * | ||
172 | * \param globalPos The position of the center of the element | ||
173 | * \return The wetting phase index | ||
174 | */ | ||
175 | template<class FluidSystem> | ||
176 | int wettingPhaseAtPos(const GlobalPosition& globalPos) const | ||
177 | { return FluidSystem::BrineIdx; } | ||
178 | |||
179 | /*! | ||
180 | * \brief Returns the temperature within the domain. | ||
181 | * | ||
182 | * \param globalPos The global position | ||
183 | * | ||
184 | * This problem assumes a geothermal gradient with | ||
185 | * a surface temperature of 10 degrees Celsius. | ||
186 | */ | ||
187 | 25149728 | Scalar temperatureAtPos(const GlobalPosition &globalPos) const | |
188 | 25149728 | { return 283.0 + (depthBOR_ - globalPos[dimWorld-1])*0.03; } | |
189 | |||
190 | private: | ||
191 | |||
192 | std::shared_ptr<const GridData<Grid>> gridData_; | ||
193 | |||
194 | int barrierTop_ = 1; | ||
195 | int barrierMiddle_ = 2; | ||
196 | int reservoir_ = 3; | ||
197 | |||
198 | Scalar barrierTopPorosity_; | ||
199 | Scalar barrierMiddlePorosity_; | ||
200 | Scalar reservoirPorosity_; | ||
201 | |||
202 | Scalar barrierTopK_; | ||
203 | Scalar barrierMiddleK_; | ||
204 | Scalar reservoirK_; | ||
205 | |||
206 | Scalar depthBOR_; | ||
207 | |||
208 | std::vector<int> paramIdx_; | ||
209 | const PcKrSwCurve pcKrSwCurve_; | ||
210 | }; | ||
211 | |||
212 | } // end namespace Dumux | ||
213 | |||
214 | #endif | ||
215 |