27#ifndef OPM_BROOKS_COREY_HPP
28#define OPM_BROOKS_COREY_HPP
53template <
class TraitsT,
class ParamsT = BrooksCoreyParams<TraitsT> >
57 typedef TraitsT Traits;
58 typedef ParamsT Params;
59 typedef typename Traits::Scalar Scalar;
64 "The Brooks-Corey capillary pressure law only applies "
65 "to the case of two fluid phases");
91 static_assert(Traits::numPhases == 2,
92 "The number of fluid phases must be two if you want to use "
93 "this material law!");
98 template <
class Container,
class Flu
idState>
101 typedef typename std::remove_reference<
decltype(values[0])>::type Evaluation;
103 values[Traits::wettingPhaseIdx] = 0.0;
104 values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fs);
111 template <
class Container,
class Flu
idState>
112 static void saturations(Container& values,
const Params& params,
const FluidState& fs)
114 typedef typename std::remove_reference<
decltype(values[0])>::type Evaluation;
116 values[Traits::wettingPhaseIdx] = Sw<FluidState, Evaluation>(params, fs);
117 values[Traits::nonWettingPhaseIdx] = 1.0 - values[Traits::wettingPhaseIdx];
130 template <
class Container,
class Flu
idState>
133 typedef typename std::remove_reference<
decltype(values[0])>::type Evaluation;
135 values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fs);
136 values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fs);
152 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
153 static Evaluation
pcnw(
const Params& params,
const FluidState& fs)
155 const Evaluation&
Sw =
156 decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
158 assert(0.0 <=
Sw &&
Sw <= 1.0);
160 return twoPhaseSatPcnw(params,
Sw);
163 template <
class Evaluation>
164 static Evaluation twoPhaseSatPcnw(
const Params& params,
const Evaluation&
Sw)
166 assert(0.0 <=
Sw &&
Sw <= 1.0);
168 return params.entryPressure()*pow(
Sw, -1/params.lambda());
171 template <
class Evaluation>
172 static Evaluation twoPhaseSatPcnwInv(
const Params& params,
const Evaluation&
pcnw)
176 return pow(params.entryPressure()/
pcnw, -params.lambda());
191 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
192 static Evaluation
Sw(
const Params& params,
const FluidState& fs)
195 decay<Evaluation>(fs.pressure(Traits::nonWettingPhaseIdx))
196 - decay<Evaluation>(fs.pressure(Traits::wettingPhaseIdx));
197 return twoPhaseSatSw(params, pC);
200 template <
class Evaluation>
201 static Evaluation twoPhaseSatSw(
const Params& params,
const Evaluation& pc)
205 return pow(pc/params.entryPressure(), -params.lambda());
212 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
213 static Evaluation
Sn(
const Params& params,
const FluidState& fs)
214 {
return 1.0 - Sw<FluidState, Evaluation>(params, fs); }
216 template <
class Evaluation>
217 static Evaluation twoPhaseSatSn(
const Params& params,
const Evaluation& pc)
218 {
return 1.0 - twoPhaseSatSw(params, pc); }
227 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
228 static Evaluation
krw(
const Params& params,
const FluidState& fs)
231 decay<Evaluation>(fs.saturation(Traits::wettingPhaseIdx));
233 return twoPhaseSatKrw(params,
Sw);
236 template <
class Evaluation>
237 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation&
Sw)
239 assert(0.0 <=
Sw &&
Sw <= 1.0);
241 return pow(
Sw, 2.0/params.lambda() + 3.0);
244 template <
class Evaluation>
245 static Evaluation twoPhaseSatKrwInv(
const Params& params,
const Evaluation&
krw)
247 return pow(
krw, 1.0/(2.0/params.lambda() + 3.0));
258 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
259 static Evaluation
krn(
const Params& params,
const FluidState& fs)
261 const Evaluation&
Sw =
262 1.0 - decay<Evaluation>(fs.saturation(Traits::nonWettingPhaseIdx));
264 return twoPhaseSatKrn(params,
Sw);
267 template <
class Evaluation>
268 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation&
Sw)
270 assert(0.0 <=
Sw &&
Sw <= 1.0);
272 Scalar exponent = 2.0/params.lambda() + 1.0;
273 const Evaluation
Sn = 1.0 -
Sw;
274 return Sn*
Sn*(1. - pow(
Sw, exponent));
277 template <
class Evaluation>
278 static Evaluation twoPhaseSatKrnInv(
const Params& params,
const Evaluation&
krn)
284 for (
int i = 0; i < 20; ++i) {
285 Evaluation f =
krn - twoPhaseSatKrn(params,
Sw);
286 Evaluation fStar =
krn - twoPhaseSatKrn(params,
Sw + eps);
287 Evaluation fPrime = (fStar - f)/eps;
289 Evaluation delta = f/fPrime;
293 if (abs(delta) < 1e-10)
297 throw NumericalProblem(
"Couldn't invert the Brooks-Corey non-wetting phase"
298 " relperm within 20 iterations");
Specification of the material parameters for the Brooks-Corey constitutive relations.
Provides the OPM specific exception classes.
Implementation of the Brooks-Corey capillary pressure <-> saturation relation.
Definition BrooksCorey.hpp:55
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:112
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure.
Definition BrooksCorey.hpp:81
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition BrooksCorey.hpp:77
static const int numPhases
The number of fluid phases to which this material law applies.
Definition BrooksCorey.hpp:62
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curves.
Definition BrooksCorey.hpp:99
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fs)
The relative permeability-saturation curves.
Definition BrooksCorey.hpp:131
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition BrooksCorey.hpp:69
static Evaluation Sn(const Params ¶ms, const FluidState &fs)
Calculate the non-wetting phase saturations depending on the phase pressures.
Definition BrooksCorey.hpp:213
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition BrooksCorey.hpp:85
static Evaluation Sw(const Params ¶ms, const FluidState &fs)
The saturation-capillary pressure curve according to Brooks & Corey.
Definition BrooksCorey.hpp:192
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:228
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:259
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition.
Definition BrooksCorey.hpp:89
static Evaluation pcnw(const Params ¶ms, const FluidState &fs)
The capillary pressure-saturation curve according to Brooks and Corey.
Definition BrooksCorey.hpp:153
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition BrooksCorey.hpp:73
This class implements a small container which holds the transmissibility mulitpliers for all the face...
Definition Exceptions.hpp:30