27 #ifndef OPM_SPLINE_TWO_PHASE_MATERIAL_HPP 28 #define OPM_SPLINE_TWO_PHASE_MATERIAL_HPP 32 #include <opm/common/ErrorMacros.hpp> 33 #include <opm/common/Exceptions.hpp> 46 template <
class TraitsT,
class ParamsT = SplineTwoPhaseMaterialParams<TraitsT> >
49 typedef typename ParamsT::SamplePoints SamplePoints;
59 typedef typename Traits::Scalar
Scalar;
64 "The piecewise linear two-phase capillary pressure law only" 65 "applies to the case of two fluid phases");
94 template <
class Container,
class Flu
idState>
97 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
99 values[Traits::wettingPhaseIdx] = 0.0;
100 values[Traits::nonWettingPhaseIdx] = pcnw<FluidState, Evaluation>(params, fluidState);
107 template <
class Container,
class Flu
idState>
109 { OPM_THROW(std::logic_error,
"Not implemented: saturations()"); }
114 template <
class Container,
class Flu
idState>
117 typedef typename std::remove_reference<decltype(values[0])>::type Evaluation;
119 values[Traits::wettingPhaseIdx] = krw<FluidState, Evaluation>(params, fluidState);
120 values[Traits::nonWettingPhaseIdx] = krn<FluidState, Evaluation>(params, fluidState);
126 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
127 static Evaluation
pcnw(
const Params& params,
const FluidState& fluidState)
129 const Evaluation&
Sw =
130 Opm::decay<Evaluation>(fluidState.saturation(Traits::wettingPhaseIdx));
138 template <
class Evaluation>
142 const auto& pcnwSpline = params.pcnwSpline();
143 if (
Sw <= pcnwSpline.xAt(0))
144 return Evaluation(pcnwSpline.valueAt(0));
145 if (
Sw >= pcnwSpline.xAt(pcnwSpline.numSamples() - 1))
146 return Evaluation(pcnwSpline.valueAt(pcnwSpline.numSamples() - 1));
148 return pcnwSpline.eval(
Sw);
151 template <
class Evaluation>
152 static Evaluation twoPhaseSatPcnwInv(
const Params& params,
const Evaluation&
pcnw)
154 static const Evaluation nil(0.0);
157 const auto& pcnwSpline = params.pcnwSpline();
158 if (
pcnw >= pcnwSpline.valueAt(0))
159 return Evaluation(pcnwSpline.xAt(0));
160 if (
pcnw <= pcnwSpline.y(pcnwSpline.numSamples() - 1))
161 return Evaluation(pcnwSpline.xAt(pcnwSpline.numSamples() - 1));
165 return pcnwSpline.intersect(nil, nil, nil,
pcnw);
171 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
172 static Evaluation
Sw(
const Params& ,
const FluidState& )
173 { OPM_THROW(std::logic_error,
"Not implemented: Sw()"); }
175 template <
class Evaluation>
176 static Evaluation twoPhaseSatSw(
const Params& ,
const Evaluation& )
177 { OPM_THROW(std::logic_error,
"Not implemented: twoPhaseSatSw()"); }
183 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
184 static Evaluation
Sn(
const Params& params,
const FluidState& fluidState)
185 {
return 1 - Sw<FluidState, Evaluation>(params, fluidState); }
187 template <
class Evaluation>
188 static Evaluation twoPhaseSatSn(
const Params& params,
const Evaluation& pC)
189 {
return 1 - twoPhaseSatSw(params, pC); }
195 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
196 static Evaluation
krw(
const Params& params,
const FluidState& fluidState)
198 const Evaluation&
Sw =
199 Opm::decay<Evaluation>(fluidState.saturation(Traits::wettingPhaseIdx));
201 return twoPhaseSatKrw(params,
Sw);
204 template <
class Evaluation>
205 static Evaluation twoPhaseSatKrw(
const Params& params,
const Evaluation&
Sw)
207 const auto& krwSpline = params.krwSpline();
208 if (
Sw <= krwSpline.xAt(0))
209 return Evaluation(krwSpline.valueAt(0));
210 if (
Sw >= krwSpline.xAt(krwSpline.numSamples() - 1))
211 return Evaluation(krwSpline.valueAt(krwSpline.numSamples() - 1));
213 return krwSpline.eval(
Sw);
216 template <
class Evaluation>
217 static Evaluation twoPhaseSatKrwInv(
const Params& params,
const Evaluation&
krw)
219 static const Evaluation nil(0.0);
221 const auto& krwSpline = params.krwSpline();
222 if (
krw <= krwSpline.valueAt(0))
223 return Evaluation(krwSpline.xAt(0));
224 if (
krw >= krwSpline.valueAt(krwSpline.numSamples() - 1))
225 return Evaluation(krwSpline.xAt(krwSpline.numSamples() - 1));
227 return krwSpline.intersect(nil, nil, nil,
krw);
234 template <
class Flu
idState,
class Evaluation =
typename Flu
idState::Scalar>
235 static Evaluation
krn(
const Params& params,
const FluidState& fluidState)
237 const Evaluation&
Sn =
238 Opm::decay<Evaluation>(fluidState.saturation(Traits::nonWettingPhaseIdx));
240 return twoPhaseSatKrn(params, 1.0 -
Sn);
243 template <
class Evaluation>
244 static Evaluation twoPhaseSatKrn(
const Params& params,
const Evaluation&
Sw)
246 const auto& krnSpline = params.krnSpline();
247 if (
Sw <= krnSpline.xAt(0))
248 return Evaluation(krnSpline.valueAt(0));
249 if (
Sw >= krnSpline.xAt(krnSpline.numSamples() - 1))
250 return Evaluation(krnSpline.valueAt(krnSpline.numSamples() - 1));
252 return krnSpline.eval(
Sw);
255 template <
class Evaluation>
256 static Evaluation twoPhaseSatKrnInv(
const Params& params,
const Evaluation&
krn)
258 static const Evaluation nil(0.0);
260 const auto& krnSpline = params.krnSpline();
261 if (
krn >= krnSpline.valueAt(0))
262 return Evaluation(krnSpline.xAt(0));
263 if (
krn <= krnSpline.valueAt(krnSpline.numSamples() - 1))
264 return Evaluation(krnSpline.xAt(krnSpline.numSamples() - 1));
266 return krnSpline.intersect(nil, nil, nil,
krn);
static Evaluation Sn(const Params ¶ms, const FluidState &fluidState)
Calculate the non-wetting phase saturations depending on the phase pressures.
Definition: SplineTwoPhaseMaterial.hpp:184
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fluidState)
The relative permeabilities.
Definition: SplineTwoPhaseMaterial.hpp:115
static Evaluation twoPhaseSatPcnw(const Params ¶ms, const Evaluation &Sw)
The saturation-capillary pressure curve.
Definition: SplineTwoPhaseMaterial.hpp:139
static const bool isCompositionDependent
Specify whether the quantities defined by this material law are dependent on the phase composition...
Definition: SplineTwoPhaseMaterial.hpp:89
static const bool isSaturationDependent
Specify whether the quantities defined by this material law are saturation dependent.
Definition: SplineTwoPhaseMaterial.hpp:77
TraitsT Traits
The traits class for this material law.
Definition: SplineTwoPhaseMaterial.hpp:53
Definition: Air_Mesitylene.hpp:33
static Evaluation pcnw(const Params ¶ms, const FluidState &fluidState)
The capillary pressure-saturation curve.
Definition: SplineTwoPhaseMaterial.hpp:127
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fluidState)
The capillary pressure-saturation curve.
Definition: SplineTwoPhaseMaterial.hpp:95
Implementation of a tabulated capillary pressure and relperm law which uses spline curves as interpol...
Definition: SplineTwoPhaseMaterial.hpp:47
static Evaluation krn(const Params ¶ms, const FluidState &fluidState)
The relative permeability for the non-wetting phase of the porous medium.
Definition: SplineTwoPhaseMaterial.hpp:235
static const int numPhases
The number of fluid phases.
Definition: SplineTwoPhaseMaterial.hpp:62
ParamsT Params
The type of the parameter objects for this law.
Definition: SplineTwoPhaseMaterial.hpp:56
Traits::Scalar Scalar
The type of the scalar values for this law.
Definition: SplineTwoPhaseMaterial.hpp:59
Specification of the material parameters for a two-phase material law which uses a table and spline-b...
static const bool isTemperatureDependent
Specify whether the quantities defined by this material law are temperature dependent.
Definition: SplineTwoPhaseMaterial.hpp:85
static const bool isPressureDependent
Specify whether the quantities defined by this material law are dependent on the absolute pressure...
Definition: SplineTwoPhaseMaterial.hpp:81
static void saturations(Container &, const Params &, const FluidState &)
The saturations of the fluid phases starting from their pressure differences.
Definition: SplineTwoPhaseMaterial.hpp:108
static const bool implementsTwoPhaseApi
Specify whether this material law implements the two-phase convenience API.
Definition: SplineTwoPhaseMaterial.hpp:69
static const bool implementsTwoPhaseSatApi
Specify whether this material law implements the two-phase convenience API which only depends on the ...
Definition: SplineTwoPhaseMaterial.hpp:73
static Evaluation krw(const Params ¶ms, const FluidState &fluidState)
The relative permeability for the wetting phase of the porous medium.
Definition: SplineTwoPhaseMaterial.hpp:196
static Evaluation Sw(const Params &, const FluidState &)
The saturation-capillary pressure curve.
Definition: SplineTwoPhaseMaterial.hpp:172