3.4
DUNE for Multi-{Phase, Component, Scale, Physics, ...} flow and transport in porous media
splinemateriallaw.hh
Go to the documentation of this file.
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 * See the file COPYING for full copying permissions. *
5 * *
6 * This program is free software: you can redistribute it and/or modify *
7 * it under the terms of the GNU General Public License as published by *
8 * the Free Software Foundation, either version 3 of the License, or *
9 * (at your option) any later version. *
10 * *
11 * This program is distributed in the hope that it will be useful, *
12 * but WITHOUT ANY WARRANTY; without even the implied warranty of *
13 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the *
14 * GNU General Public License for more details. *
15 * *
16 * You should have received a copy of the GNU General Public License *
17 * along with this program. If not, see <http://www.gnu.org/licenses/>. *
18 *****************************************************************************/
24#ifndef DUMUX_MATERIAL_FLUIDMATRIX_TWOP_SPLINE_MATERIAL_LAW_HH
25#define DUMUX_MATERIAL_FLUIDMATRIX_TWOP_SPLINE_MATERIAL_LAW_HH
26
27#include <memory> // unique_ptr
28
33
34namespace Dumux::FluidMatrix {
35
46template<class TwoPMaterialLaw, bool approximatePcSwInverse = false>
48: public TwoPMaterialLaw
49, public Adapter<SplineTwoPMaterialLaw<TwoPMaterialLaw, approximatePcSwInverse>, PcKrSw>
50{
51public:
53
57
61 static constexpr int numFluidPhases()
62 { return 2; }
63
68 static constexpr bool isRegularized()
69 { return true; }
70
76
81 explicit SplineTwoPMaterialLaw(const std::string& paramGroup)
82 : TwoPMaterialLaw(paramGroup)
83 {
84 const std::array<Scalar, 2> defaultInterval{{ 0.01, 1.0 }};
85 const auto sweInterval = getParamFromGroup<std::array<Scalar, 2>>(paramGroup, "SplineSweInterval", defaultInterval);
86 swInterval_ = {{ TwoPMaterialLaw::EffToAbs::sweToSw(sweInterval[0], this->effToAbsParams()),
87 TwoPMaterialLaw::EffToAbs::sweToSw(sweInterval[1], this->effToAbsParams()) }};
88 swIntervalPc_ = {{ TwoPMaterialLaw::pc(swInterval_[1]),
89 TwoPMaterialLaw::pc(swInterval_[0]) }};
90 numSwSamples_ = getParamFromGroup<Scalar>(paramGroup, "SplineNumSwSamples", 30);
91
92 pcSpline_ = makeSweSpline_(
93 [&](const Scalar s){ return TwoPMaterialLaw::pc(s); },
94 approximatePcSwInverse
95 );
96
97 krwSpline_ = makeSweSpline_(
98 [&](const Scalar s){ return TwoPMaterialLaw::krw(s); }
99 );
100
101 krnSpline_ = makeSweSpline_(
102 [&](const Scalar s){ return TwoPMaterialLaw::krn(s); }
103 );
104 }
105
110 SplineTwoPMaterialLaw(const std::array<Scalar, 2>& sweInterval,
111 std::size_t numSwSamples,
112 TwoPMaterialLaw&& twoP)
113 : TwoPMaterialLaw(std::move(twoP))
114 , numSwSamples_(numSwSamples)
115 {
116 swInterval_ = {{ TwoPMaterialLaw::EffToAbs::sweToSw(sweInterval[0], this->effToAbsParams()),
117 TwoPMaterialLaw::EffToAbs::sweToSw(sweInterval[1], this->effToAbsParams()) }};
118 swIntervalPc_ = {{ TwoPMaterialLaw::pc(swInterval_[1]),
119 TwoPMaterialLaw::pc(swInterval_[0]) }};
120 }
121
125 Scalar pc(const Scalar sw) const
126 {
127 if (sw > swInterval_[0] && sw < swInterval_[1])
128 {
129 if constexpr (approximatePcSwInverse)
130 return pcSpline_->evalInverse(sw);
131 else
132 return pcSpline_->eval(sw);
133 }
134
135 return TwoPMaterialLaw::pc(sw);
136 }
137
141 Scalar dpc_dsw(const Scalar sw) const
142 {
143 if (sw > swInterval_[0] && sw < swInterval_[1])
144 {
145 if constexpr (approximatePcSwInverse)
146 return 1.0/pcSpline_->evalDerivative(pcSpline_->evalInverse(sw));
147 else
148 return pcSpline_->evalDerivative(sw);
149 }
150
152 }
153
157 Scalar sw(const Scalar pc) const
158 {
159 if (pc > swIntervalPc_[0] && pc < swIntervalPc_[1])
160 {
161 if constexpr (approximatePcSwInverse)
162 return pcSpline_->eval(pc);
163 else
164 return pcSpline_->evalInverse(pc);
165 }
166
167 return TwoPMaterialLaw::sw(pc);
168 }
169
173 Scalar dsw_dpc(const Scalar pc) const
174 {
175 if (pc > swIntervalPc_[0] && pc < swIntervalPc_[1])
176 {
177 if constexpr (approximatePcSwInverse)
178 return pcSpline_->evalDerivative(pc);
179 else
180 return 1.0/pcSpline_->evalDerivative(pcSpline_->evalInverse(pc));
181 }
182
184 }
185
189 Scalar krw(const Scalar sw) const
190 {
191 if (sw > swInterval_[0] && sw < swInterval_[1])
192 return krwSpline_->eval(sw);
193
194 return TwoPMaterialLaw::krw(sw);
195 }
196
201 {
202 if (sw > swInterval_[0] && sw < swInterval_[1])
203 return krwSpline_->evalDerivative(sw);
204
206 }
207
211 Scalar krn(const Scalar sw) const
212 {
213 if (sw > swInterval_[0] && sw < swInterval_[1])
214 return krnSpline_->eval(sw);
215
216 return TwoPMaterialLaw::krn(sw);
217 }
218
223 {
224 if (sw > swInterval_[0] && sw < swInterval_[1])
225 return krnSpline_->evalDerivative(sw);
226
228 }
229
230private:
231 template<class Function>
232 std::unique_ptr<MonotoneCubicSpline<Scalar>>
233 makeSweSpline_(const Function& f, bool invert = false) const
234 {
235 const auto sw = linspace(swInterval_[0], swInterval_[1], numSwSamples_);
236
237 auto values = sw;
238 std::transform(sw.begin(), sw.end(), values.begin(), f);
239
240 if (invert)
241 return std::make_unique<MonotoneCubicSpline<Scalar>>(values, sw);
242 else
243 return std::make_unique<MonotoneCubicSpline<Scalar>>(sw, values);
244 }
245
246 std::unique_ptr<MonotoneCubicSpline<Scalar>> pcSpline_;
247 std::unique_ptr<MonotoneCubicSpline<Scalar>> krwSpline_;
248 std::unique_ptr<MonotoneCubicSpline<Scalar>> krnSpline_;
249
250 std::array<Scalar, 2> swInterval_;
251 std::array<Scalar, 2> swIntervalPc_;
252 std::size_t numSwSamples_;
253};
254
255} // end namespace Dumux::FluidMatrix
256
257#endif
A monotone cubic spline.
The infrastructure to retrieve run-time parameters from Dune::ParameterTrees.
Implementation helper for capillary pressure and relative permeability <-> saturation relations for t...
Wrapper type to combine an arbitrary number of different laws for fluid-matrix interaction (e....
std::vector< Scalar > linspace(const Scalar begin, const Scalar end, std::size_t samples, bool endPoint=true)
Generates linearly spaced vectors.
Definition: math.hh:594
Definition: brookscorey.hh:35
Wrapper class to implement regularized material laws (pc-sw, kr-sw) with a conversion policy between ...
Definition: materiallaw.hh:57
Scalar dsw_dpc(const Scalar pc) const
The partial derivative of the saturation to the capillary pressure.
Definition: materiallaw.hh:174
Scalar dkrn_dsw(const Scalar sw) const
The derivative of the relative permeability for the non-wetting phase w.r.t. saturation.
Definition: materiallaw.hh:241
ScalarType Scalar
Definition: materiallaw.hh:60
Scalar krn(const Scalar sw) const
The relative permeability for the non-wetting phase.
Definition: materiallaw.hh:224
Scalar dkrw_dsw(const Scalar sw) const
The derivative of the relative permeability for the wetting phase w.r.t. saturation.
Definition: materiallaw.hh:207
typename Regularization::template Params< Scalar > RegularizationParams
Definition: materiallaw.hh:64
const EffToAbsParams & effToAbsParams() const
Return the parameters of the EffToAbs policy.
Definition: materiallaw.hh:291
Scalar dpc_dsw(const Scalar sw) const
The partial derivative of the capillary pressure w.r.t. the saturation.
Definition: materiallaw.hh:133
Scalar sw(const Scalar pc) const
The saturation-capillary pressure curve.
Definition: materiallaw.hh:158
typename EffToAbsPolicy::template Params< Scalar > EffToAbsParams
Definition: materiallaw.hh:63
typename BaseLaw::template Params< Scalar > BasicParams
Definition: materiallaw.hh:62
Scalar pc(const Scalar sw) const
The capillary pressure-saturation curve.
Definition: materiallaw.hh:116
Scalar krw(const Scalar sw) const
The relative permeability for the wetting phase.
Definition: materiallaw.hh:190
A spline approximation wrapper for 2p material laws.
Definition: splinemateriallaw.hh:50
static constexpr bool isRegularized()
We are always regularized in the sense that we replace the orginal curve by a cubic spline.
Definition: splinemateriallaw.hh:68
Scalar dkrn_dsw(const Scalar sw) const
The derivative of the relative permeability for the non-wetting phase w.r.t. saturation.
Definition: splinemateriallaw.hh:222
Scalar pc(const Scalar sw) const
The capillary pressure-saturation curve.
Definition: splinemateriallaw.hh:125
static constexpr int numFluidPhases()
Return the number of fluid phases.
Definition: splinemateriallaw.hh:61
typename TwoPMaterialLaw::Scalar Scalar
Definition: splinemateriallaw.hh:52
SplineTwoPMaterialLaw(const std::string &paramGroup)
Construct from a subgroup from the global parameter tree.
Definition: splinemateriallaw.hh:81
SplineTwoPMaterialLaw(const std::array< Scalar, 2 > &sweInterval, std::size_t numSwSamples, TwoPMaterialLaw &&twoP)
Construct from parameter structs.
Definition: splinemateriallaw.hh:110
Scalar dsw_dpc(const Scalar pc) const
The partial derivative of the saturation to the capillary pressure.
Definition: splinemateriallaw.hh:173
Scalar krn(const Scalar sw) const
The relative permeability for the non-wetting phase.
Definition: splinemateriallaw.hh:211
Scalar krw(const Scalar sw) const
The relative permeability for the wetting phase.
Definition: splinemateriallaw.hh:189
Scalar dkrw_dsw(const Scalar sw) const
The derivative of the relative permeability for the wetting phase w.r.t. saturation.
Definition: splinemateriallaw.hh:200
SplineTwoPMaterialLaw()=delete
Deleted default constructor (so we are never in an undefined state)
Scalar dpc_dsw(const Scalar sw) const
The partial derivative of the capillary pressure w.r.t. the saturation.
Definition: splinemateriallaw.hh:141
Scalar sw(const Scalar pc) const
The saturation-capillary pressure curve.
Definition: splinemateriallaw.hh:157
Adapter to inherit from, allowing the inheriting class to be wrapped by the makeFluidMatrixInteractio...
Definition: fluidmatrixinteraction.hh:67