25 #ifndef OPM_SPLINE_TWO_PHASE_MATERIAL_HPP
26 #define OPM_SPLINE_TWO_PHASE_MATERIAL_HPP
30 #include <opm/common/ErrorMacros.hpp>
31 #include <opm/common/Exceptions.hpp>
44 template <
class TraitsT,
class ParamsT = SplineTwoPhaseMaterialParams<TraitsT> >
47 typedef typename ParamsT::SamplePoints SamplePoints;
57 typedef typename Traits::Scalar
Scalar;
61 static_assert(numPhases == 2,
62 "The piecewise linear two-phase capillary pressure law only"
63 "applies to the case of two fluid phases");
67 static const bool implementsTwoPhaseApi =
true;
71 static const bool implementsTwoPhaseSatApi =
true;
75 static const bool isSaturationDependent =
true;
79 static const bool isPressureDependent =
false;
83 static const bool isTemperatureDependent =
false;
87 static const bool isCompositionDependent =
false;
92 template <
class Container,
class Flu
idState>
93 static void capillaryPressures(Container &values,
const Params ¶ms,
const FluidState &fluidState)
95 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
97 values[Traits::wettingPhaseIdx] = 0.0;
98 values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fluidState);
105 template <
class Container,
class Flu
idState>
106 static void saturations(Container &,
const Params &,
const FluidState &)
107 { OPM_THROW(std::logic_error,
"Not implemented: saturations()"); }
112 template <
class Container,
class Flu
idState>
113 static void relativePermeabilities(Container &values,
const Params ¶ms,
const FluidState &fluidState)
115 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
117 values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fluidState);
118 values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fluidState);
124 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
125 static Evaluation pcnw(
const Params ¶ms,
const FluidState &fluidState)
127 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
129 const Evaluation& Sw =
130 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(Traits::wettingPhaseIdx));
132 return twoPhaseSatPcnw(params, Sw);
138 template <
class Evaluation>
139 static Evaluation twoPhaseSatPcnw(
const Params ¶ms,
const Evaluation& Sw)
142 if (Sw <= params.pcnwSpline().xMin())
143 return Evaluation(params.pcnwSpline().yFirst());
144 if (Sw >= params.pcnwSpline().xMax())
145 return Evaluation(params.pcnwSpline().yLast());
147 return params.pcnwSpline().eval(Sw);
150 template <
class Evaluation>
151 static Evaluation twoPhaseSatPcnwInv(
const Params ¶ms,
const Evaluation& pcnw)
153 static const Evaluation nil(0.0);
156 if (pcnw >= params.pcnwSpline().yFirst())
157 return Evaluation(params.pcnwSpline().xMin());
158 if (pcnw <= params.pcnwSpline().yLast())
159 return Evaluation(params.pcnwSpline().xMax());
163 return params.pcnwSpline().intersect(nil, nil, nil, pcnw);
169 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
170 static Evaluation Sw(
const Params &,
const FluidState &)
171 { OPM_THROW(std::logic_error,
"Not implemented: Sw()"); }
173 template <
class Evaluation>
174 static Evaluation twoPhaseSatSw(
const Params &,
const Evaluation& )
175 { OPM_THROW(std::logic_error,
"Not implemented: twoPhaseSatSw()"); }
181 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
182 static Evaluation Sn(
const Params ¶ms,
const FluidState &fluidState)
183 {
return 1 - Sw<FluidState, Evaluation>(params, fluidState); }
185 template <
class Evaluation>
186 static Evaluation twoPhaseSatSn(
const Params ¶ms,
const Evaluation& pC)
187 {
return 1 - twoPhaseSatSw(params, pC); }
193 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
194 static Evaluation krw(
const Params ¶ms,
const FluidState &fluidState)
196 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
198 const Evaluation& Sw =
199 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(Traits::wettingPhaseIdx));
201 return twoPhaseSatKrw(params, Sw);
204 template <
class Evaluation>
205 static Evaluation twoPhaseSatKrw(
const Params ¶ms,
const Evaluation& Sw)
207 if (Sw <= params.krnSpline().xMin())
208 return Evaluation(params.krwSpline().yFirst());
209 if (Sw >= params.krnSpline().xMax())
210 return Evaluation(params.krwSpline().yLast());
212 return params.krwSpline().eval(Sw);
215 template <
class Evaluation>
216 static Evaluation twoPhaseSatKrwInv(
const Params ¶ms,
const Evaluation& krw)
218 static const Evaluation nil(0.0);
220 if (krw <= params.krwSpline().yFirst())
221 return Evaluation(params.krwSpline().xMin());
222 if (krw >= params.krwSpline().yLast())
223 return Evaluation(params.krwSpline().xMax());
225 return params.krwSpline().intersect(nil, nil, nil, krw);
232 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
233 static Evaluation krn(
const Params ¶ms,
const FluidState &fluidState)
235 typedef MathToolbox<typename FluidState::Scalar> FsToolbox;
237 const Evaluation& Sn =
238 FsToolbox::template toLhs<Evaluation>(fluidState.saturation(Traits::nonWettingPhaseIdx));
240 return twoPhaseSatKrn(params, 1.0 - Sn);
243 template <
class Evaluation>
244 static Evaluation twoPhaseSatKrn(
const Params ¶ms,
const Evaluation& Sw)
246 if (Sw <= params.krnSpline().xMin())
247 return Evaluation(params.krnSpline().yFirst());
248 if (Sw >= params.krnSpline().xMax())
249 return Evaluation(params.krnSpline().yLast());
251 return params.krnSpline().eval(Sw);
254 template <
class Evaluation>
255 static Evaluation twoPhaseSatKrnInv(
const Params ¶ms,
const Evaluation& krn)
257 static const Evaluation nil(0.0);
259 if (krn >= params.krnSpline().yFirst())
260 return Evaluation(params.krnSpline().xMin());
261 if (krn <= params.krnSpline().yLast())
262 return Evaluation(params.krnSpline().xMax());
264 return params.krnSpline().intersect(nil, nil, nil, krn);
Definition: Air_Mesitylene.hpp:31
TraitsT Traits
The traits class for this material law.
Definition: SplineTwoPhaseMaterial.hpp:51
static const int numPhases
The number of fluid phases.
Definition: SplineTwoPhaseMaterial.hpp:60
Specification of the material parameters for a two-phase material law which uses a table and spline-b...
Implementation of a tabulated capillary pressure and relperm law which uses spline curves as interpol...
Definition: SplineTwoPhaseMaterial.hpp:45
ParamsT Params
The type of the parameter objects for this law.
Definition: SplineTwoPhaseMaterial.hpp:54
Traits::Scalar Scalar
The type of the scalar values for this law.
Definition: SplineTwoPhaseMaterial.hpp:57