27 #ifndef REGULARIZED_BROOKS_COREY_HPP
28 #define REGULARIZED_BROOKS_COREY_HPP
64 template <
class TraitsT,
class ParamsT = RegularizedBrooksCoreyParams<TraitsT> >
67 typedef ::Opm::BrooksCorey<TraitsT, ParamsT>
BrooksCorey;
70 typedef TraitsT Traits;
71 typedef ParamsT Params;
72 typedef typename Traits::Scalar Scalar;
77 "The regularized Brooks-Corey capillary pressure law only "
78 "applies to the case of two fluid phases");
104 static_assert(Traits::numPhases == 2,
105 "The number of fluid phases must be two if you want to use "
106 "this material law!");
118 template <
class Container,
class Flu
idState>
121 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
123 values[Traits::wettingPhaseIdx] = 0.0;
124 values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fs);
131 template <
class Container,
class Flu
idState>
132 static void saturations(Container& values,
const Params& params,
const FluidState& fs)
134 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
136 values[Traits::wettingPhaseIdx] = Sw<FluidState, Evaluation>(params, fs);
137 values[Traits::nonWettingPhaseIdx] = 1 - values[Traits::wettingPhaseIdx];
150 template <
class Container,
class Flu
idState>
153 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
155 values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fs);
156 values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fs);
173 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
174 static Evaluation
pcnw(
const Params& params,
const FluidState& fs)
176 const auto&
Sw = decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
177 return twoPhaseSatPcnw(params,
Sw);
180 template <
class Evaluation>
181 static Evaluation twoPhaseSatPcnw(
const Params& params,
const Evaluation&
Sw)
183 const Scalar Sthres = params.pcnwLowSw();
186 Scalar m = params.pcnwSlopeLow();
187 Scalar pcnw_SwLow = params.pcnwLow();
188 return pcnw_SwLow + m*(
Sw - Sthres);
190 else if (
Sw >= 1.0) {
191 Scalar m = params.pcnwSlopeHigh();
192 Scalar pcnw_SwHigh = params.pcnwHigh();
193 return pcnw_SwHigh + m*(
Sw - 1.0);
198 return BrooksCorey::twoPhaseSatPcnw(params,
Sw);
207 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
208 static Evaluation
Sw(
const Params& params,
const FluidState& fs)
210 const Evaluation& pC =
211 decay<Evaluation>(fs.pressure(Traits::nonWettingPhaseIdx))
212 - decay<Evaluation>(fs.pressure(Traits::wettingPhaseIdx));
213 return twoPhaseSatSw(params, pC);
216 template <
class Evaluation>
217 static Evaluation twoPhaseSatSw(
const Params& params,
const Evaluation&
pcnw)
219 const Scalar Sthres = params.pcnwLowSw();
227 if (
pcnw >= params.entryPressure())
228 Sw = BrooksCorey::twoPhaseSatSw(params,
pcnw);
238 Scalar m = params.pcnwSlopeLow();
239 Scalar pcnw_SwLow = params.pcnwLow();
240 return Sthres + (
pcnw - pcnw_SwLow)/m;
243 Scalar m = params.pcnwSlopeHigh();
244 Scalar pcnw_SwHigh = params.pcnwHigh();
245 return 1.0 + (
pcnw - pcnw_SwHigh)/m;;
248 return BrooksCorey::twoPhaseSatSw(params,
pcnw);
255 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
256 static Evaluation
Sn(
const Params& params,
const FluidState& fs)
257 {
return 1 - Sw<FluidState, Evaluation>(params, fs); }
259 template <
class Evaluation>
260 static Evaluation twoPhaseSatSn(
const Params& params,
const Evaluation&
pcnw)
261 {
return 1 - twoPhaseSatSw(params,
pcnw); }
277 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
278 static Evaluation
krw(
const Params& params,
const FluidState& fs)
280 const auto&
Sw = decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
281 return twoPhaseSatKrw(params,
Sw);
284 template <
class Evaluation>
285 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation&
Sw)
292 return BrooksCorey::twoPhaseSatKrw(params,
Sw);
295 template <
class Evaluation>
296 static Evaluation twoPhaseSatKrwInv(
const Params& params,
const Evaluation&
krw)
303 return BrooksCorey::twoPhaseSatKrwInv(params,
krw);
320 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
321 static Evaluation
krn(
const Params& params,
const FluidState& fs)
323 const Evaluation&
Sw =
324 1.0 - decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
325 return twoPhaseSatKrn(params,
Sw);
328 template <
class Evaluation>
329 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation&
Sw)
336 return BrooksCorey::twoPhaseSatKrn(params,
Sw);
339 template <
class Evaluation>
340 static Evaluation twoPhaseSatKrnInv(
const Params& params,
const Evaluation&
krn)
347 return BrooksCorey::twoPhaseSatKrnInv(params,
krn);
Implementation of the Brooks-Corey capillary pressure <-> saturation relation.
Parameters that are necessary for the regularization of the Brooks-Corey capillary pressure model.
Class implementing cubic splines.
Implementation of the Brooks-Corey capillary pressure <-> saturation relation.
Definition: BrooksCorey.hpp:54
Implementation of the regularized Brooks-Corey capillary pressure / relative permeability <-> saturat...
Definition: RegularizedBrooksCorey.hpp:66
static Evaluation krn(const Params ¶ms, const FluidState &fs)
Regularized version of the relative permeability of the non-wetting phase of the Brooks-Corey curves.
Definition: RegularizedBrooksCorey.hpp:321
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: RegularizedBrooksCorey.hpp:82
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: RegularizedBrooksCorey.hpp:86
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curves depending on absolute saturations.
Definition: RegularizedBrooksCorey.hpp:119
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: RegularizedBrooksCorey.hpp:94
static Evaluation Sn(const Params ¶ms, const FluidState &fs)
Calculate the non-wetting phase saturations depending on the phase pressures.
Definition: RegularizedBrooksCorey.hpp:256
static const int numPhases
The number of fluid phases.
Definition: RegularizedBrooksCorey.hpp:75
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
A regularized Brooks-Corey capillary pressure-saturation curve.
Definition: RegularizedBrooksCorey.hpp:174
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fs)
The relative permeability-saturation curves depending on absolute saturations.
Definition: RegularizedBrooksCorey.hpp:151
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: RegularizedBrooksCorey.hpp:102
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: RegularizedBrooksCorey.hpp:98
static Evaluation krw(const Params ¶ms, const FluidState &fs)
Regularized version of the relative permeability of the wetting phase of the Brooks-Corey curves.
Definition: RegularizedBrooksCorey.hpp:278
static void saturations(Container &values, const Params ¶ms, const FluidState &fs)
Calculate the saturations of the phases starting from their pressure differences.
Definition: RegularizedBrooksCorey.hpp:132
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: RegularizedBrooksCorey.hpp:90
static Evaluation Sw(const Params ¶ms, const FluidState &fs)
A regularized Brooks-Corey saturation-capillary pressure curve.
Definition: RegularizedBrooksCorey.hpp:208