26 #ifndef REGULARIZED_BROOKS_COREY_HPP
27 #define REGULARIZED_BROOKS_COREY_HPP
63 template <
class TraitsT,
class ParamsT = RegularizedBrooksCoreyParams<TraitsT> >
71 typedef typename Traits::Scalar
Scalar;
75 static_assert(numPhases == 2,
76 "The regularized Brooks-Corey capillary pressure law only "
77 "applies to the case of two fluid phases");
81 static const bool implementsTwoPhaseApi =
true;
85 static const bool implementsTwoPhaseSatApi =
true;
89 static const bool isSaturationDependent =
true;
93 static const bool isPressureDependent =
false;
97 static const bool isTemperatureDependent =
false;
101 static const bool isCompositionDependent =
false;
103 static_assert(Traits::numPhases == 2,
104 "The number of fluid phases must be two if you want to use "
105 "this material law!");
117 template <
class Container,
class Flu
idState>
118 static void capillaryPressures(Container &values,
const Params ¶ms,
const FluidState &fs)
120 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
122 values[Traits::wettingPhaseIdx] = 0.0;
123 values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fs);
130 template <
class Container,
class Flu
idState>
131 static void saturations(Container &values,
const Params ¶ms,
const FluidState &fs)
133 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
135 values[Traits::wettingPhaseIdx] = Sw<FluidState, Evaluation>(params, fs);
136 values[Traits::nonWettingPhaseIdx] = 1 - values[Traits::wettingPhaseIdx];
149 template <
class Container,
class Flu
idState>
150 static void relativePermeabilities(Container &values,
const Params ¶ms,
const FluidState &fs)
152 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
154 values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fs);
155 values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fs);
172 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
173 static Evaluation pcnw(
const Params ¶ms,
const FluidState &fs)
175 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
177 const auto& Sw = FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
178 return twoPhaseSatPcnw(params, Sw);
181 template <
class Evaluation>
182 static Evaluation twoPhaseSatPcnw(
const Params ¶ms,
const Evaluation& Sw)
184 const Scalar Sthres = params.pcnwLowSw();
187 Scalar m = params.pcnwSlopeLow();
188 Scalar pcnw_SwLow = params.pcnwLow();
189 return pcnw_SwLow + m*(Sw - Sthres);
191 else if (Sw >= 1.0) {
192 Scalar m = params.pcnwSlopeHigh();
193 Scalar pcnw_SwHigh = params.pcnwHigh();
194 return pcnw_SwHigh + m*(Sw - 1.0);
199 return BrooksCorey::twoPhaseSatPcnw(params, Sw);
208 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
209 static Evaluation Sw(
const Params ¶ms,
const FluidState &fs)
211 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
213 const Evaluation& pC =
214 FsToolbox::template toLhs<Evaluation>(fs.pressure(Traits::nonWettingPhaseIdx))
215 - FsToolbox::template toLhs<Evaluation>(fs.pressure(Traits::wettingPhaseIdx));
216 return twoPhaseSatSw(params, pC);
219 template <
class Evaluation>
220 static Evaluation twoPhaseSatSw(
const Params ¶ms,
const Evaluation& pcnw)
222 const Scalar Sthres = params.pcnwLowSw();
230 if (pcnw >= params.entryPressure())
231 Sw = BrooksCorey::twoPhaseSatSw(params, pcnw);
241 Scalar m = params.pcnwSlopeLow();
242 Scalar pcnw_SwLow = params.pcnwLow();
243 return Sthres + (pcnw - pcnw_SwLow)/m;
246 Scalar m = params.pcnwSlopeHigh();
247 Scalar pcnw_SwHigh = params.pcnwHigh();
248 return 1.0 + (pcnw - pcnw_SwHigh)/m;;
251 return BrooksCorey::twoPhaseSatSw(params, pcnw);
258 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
259 static Evaluation Sn(
const Params ¶ms,
const FluidState &fs)
260 {
return 1 - Sw<FluidState, Evaluation>(params, fs); }
262 template <
class Evaluation>
263 static Evaluation twoPhaseSatSn(
const Params ¶ms,
const Evaluation& pcnw)
264 {
return 1 - twoPhaseSatSw(params, pcnw); }
280 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
281 static Evaluation krw(
const Params ¶ms,
const FluidState &fs)
283 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
285 const auto& Sw = FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
286 return twoPhaseSatKrw(params, Sw);
289 template <
class Evaluation>
290 static Evaluation twoPhaseSatKrw(
const Params ¶ms,
const Evaluation& Sw)
292 typedef MathToolbox<Evaluation> Toolbox;
295 return Toolbox::createConstant(0.0);
297 return Toolbox::createConstant(1.0);
299 return BrooksCorey::twoPhaseSatKrw(params, Sw);
302 template <
class Evaluation>
303 static Evaluation twoPhaseSatKrwInv(
const Params ¶ms,
const Evaluation& krw)
305 typedef MathToolbox<Evaluation> Toolbox;
308 return Toolbox::createConstant(0.0);
310 return Toolbox::createConstant(1.0);
312 return BrooksCorey::twoPhaseSatKrwInv(params, krw);
329 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
330 static Evaluation krn(
const Params ¶ms,
const FluidState &fs)
332 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
334 const Evaluation& Sw =
335 1.0 - FsToolbox::template toLhs<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
336 return twoPhaseSatKrn(params, Sw);
339 template <
class Evaluation>
340 static Evaluation twoPhaseSatKrn(
const Params ¶ms,
const Evaluation& Sw)
342 typedef MathToolbox<Evaluation> Toolbox;
345 return Toolbox::createConstant(0.0);
347 return Toolbox::createConstant(1.0);
349 return BrooksCorey::twoPhaseSatKrn(params, Sw);
352 template <
class Evaluation>
353 static Evaluation twoPhaseSatKrnInv(
const Params ¶ms,
const Evaluation& krn)
355 typedef MathToolbox<Evaluation> Toolbox;
358 return Toolbox::createConstant(1.0);
360 return Toolbox::createConstant(0.0);
362 return BrooksCorey::twoPhaseSatKrnInv(params, krn);
Implementation of the Brooks-Corey capillary pressure <-> saturation relation.
Definition: Air_Mesitylene.hpp:31
ParamsT Params
Definition: RegularizedBrooksCorey.hpp:70
Class implementing cubic splines.
Implementation of the Brooks-Corey capillary pressure <-> saturation relation.
Definition: BrooksCorey.hpp:54
Parameters that are necessary for the regularization of the Brooks-Corey capillary pressure model...
Implementation of the regularized Brooks-Corey capillary pressure / relative permeability <-> saturat...
Definition: RegularizedBrooksCorey.hpp:64
static const int numPhases
The number of fluid phases.
Definition: RegularizedBrooksCorey.hpp:74
Traits::Scalar Scalar
Definition: RegularizedBrooksCorey.hpp:71
TraitsT Traits
Definition: RegularizedBrooksCorey.hpp:69