27 #ifndef OPM_BROOKS_COREY_HPP
28 #define OPM_BROOKS_COREY_HPP
52 template <
class TraitsT,
class ParamsT = BrooksCoreyParams<TraitsT> >
56 typedef TraitsT Traits;
57 typedef ParamsT Params;
58 typedef typename Traits::Scalar Scalar;
63 "The Brooks-Corey capillary pressure law only applies "
64 "to the case of two fluid phases");
90 static_assert(Traits::numPhases == 2,
91 "The number of fluid phases must be two if you want to use "
92 "this material law!");
97 template <
class Container,
class Flu
idState>
100 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
102 values[Traits::wettingPhaseIdx] = 0.0;
103 values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fs);
110 template <
class Container,
class Flu
idState>
111 static void saturations(Container& values,
const Params& params,
const FluidState& fs)
113 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
115 values[Traits::wettingPhaseIdx] = Sw<FluidState, Evaluation>(params, fs);
116 values[Traits::nonWettingPhaseIdx] = 1.0 - values[Traits::wettingPhaseIdx];
129 template <
class Container,
class Flu
idState>
132 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
134 values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fs);
135 values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fs);
151 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
152 static Evaluation
pcnw(
const Params& params,
const FluidState& fs)
154 const Evaluation&
Sw =
155 decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
157 assert(0.0 <=
Sw &&
Sw <= 1.0);
159 return twoPhaseSatPcnw(params,
Sw);
162 template <
class Evaluation>
163 static Evaluation twoPhaseSatPcnw(
const Params& params,
const Evaluation&
Sw)
165 assert(0.0 <=
Sw &&
Sw <= 1.0);
167 return params.entryPressure()*pow(
Sw, -1/params.lambda());
170 template <
class Evaluation>
171 static Evaluation twoPhaseSatPcnwInv(
const Params& params,
const Evaluation&
pcnw)
175 return pow(params.entryPressure()/
pcnw, -params.lambda());
190 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
191 static Evaluation
Sw(
const Params& params,
const FluidState& fs)
194 decay<Evaluation>(fs.pressure(Traits::nonWettingPhaseIdx))
195 - decay<Evaluation>(fs.pressure(Traits::wettingPhaseIdx));
196 return twoPhaseSatSw(params, pC);
199 template <
class Evaluation>
200 static Evaluation twoPhaseSatSw(
const Params& params,
const Evaluation& pc)
204 return pow(pc/params.entryPressure(), -params.lambda());
211 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
212 static Evaluation
Sn(
const Params& params,
const FluidState& fs)
213 {
return 1.0 - Sw<FluidState, Evaluation>(params, fs); }
215 template <
class Evaluation>
216 static Evaluation twoPhaseSatSn(
const Params& params,
const Evaluation& pc)
217 {
return 1.0 - twoPhaseSatSw(params, pc); }
226 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
227 static Evaluation
krw(
const Params& params,
const FluidState& fs)
230 decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
232 return twoPhaseSatKrw(params,
Sw);
235 template <
class Evaluation>
236 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation&
Sw)
238 assert(0.0 <=
Sw &&
Sw <= 1.0);
240 return pow(
Sw, 2.0/params.lambda() + 3.0);
243 template <
class Evaluation>
244 static Evaluation twoPhaseSatKrwInv(
const Params& params,
const Evaluation&
krw)
246 return pow(
krw, 1.0/(2.0/params.lambda() + 3.0));
257 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
258 static Evaluation
krn(
const Params& params,
const FluidState& fs)
260 const Evaluation&
Sw =
261 1.0 - decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
263 return twoPhaseSatKrn(params,
Sw);
266 template <
class Evaluation>
267 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation&
Sw)
269 assert(0.0 <=
Sw &&
Sw <= 1.0);
271 Scalar exponent = 2.0/params.lambda() + 1.0;
272 const Evaluation
Sn = 1.0 -
Sw;
273 return Sn*
Sn*(1. - pow(
Sw, exponent));
276 template <
class Evaluation>
277 static Evaluation twoPhaseSatKrnInv(
const Params& params,
const Evaluation&
krn)
283 for (
int i = 0; i < 20; ++i) {
284 Evaluation f =
krn - twoPhaseSatKrn(params,
Sw);
285 Evaluation fStar =
krn - twoPhaseSatKrn(params,
Sw + eps);
286 Evaluation fPrime = (fStar - f)/eps;
288 Evaluation delta = f/fPrime;
292 if (abs(delta) < 1e-10)
296 throw NumericalIssue(
"Couldn't invert the Brooks-Corey non-wetting phase"
297 " relperm within 20 iterations");
Specification of the material parameters for the Brooks-Corey constitutive relations.
Provides the opm-material specific exception classes.
Implementation of the Brooks-Corey capillary pressure <-> saturation relation.
Definition: BrooksCorey.hpp:54
static void saturations(Container &values, const Params ¶ms, const FluidState &fs)
Calculate the saturations of the phases starting from their pressure differences.
Definition: BrooksCorey.hpp:111
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition: BrooksCorey.hpp:80
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: BrooksCorey.hpp:76
static const int numPhases
The number of fluid phases to which this material law applies.
Definition: BrooksCorey.hpp:61
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curves.
Definition: BrooksCorey.hpp:98
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fs)
The relative permeability-saturation curves.
Definition: BrooksCorey.hpp:130
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: BrooksCorey.hpp:68
static Evaluation Sn(const Params ¶ms, const FluidState &fs)
Calculate the non-wetting phase saturations depending on the phase pressures.
Definition: BrooksCorey.hpp:212
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: BrooksCorey.hpp:84
static Evaluation Sw(const Params ¶ms, const FluidState &fs)
The saturation-capillary pressure curve according to Brooks & Corey.
Definition: BrooksCorey.hpp:191
static Evaluation krw(const Params ¶ms, const FluidState &fs)
The relative permeability for the wetting phase of the medium implied by the Brooks-Corey parameteriz...
Definition: BrooksCorey.hpp:227
static Evaluation krn(const Params ¶ms, const FluidState &fs)
The relative permeability for the non-wetting phase of the medium as implied by the Brooks-Corey para...
Definition: BrooksCorey.hpp:258
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition: BrooksCorey.hpp:88
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curve according to Brooks and Corey.
Definition: BrooksCorey.hpp:152
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: BrooksCorey.hpp:72