27 #ifndef OPM_TWO_PHASE_LET_CURVES_PARAMS_HPP
28 #define OPM_TWO_PHASE_LET_CURVES_PARAMS_HPP
43 template <
class TraitsT>
46 using Scalar =
typename TraitsT::Scalar;
49 using Traits = TraitsT;
51 static constexpr
int wIdx = 0;
52 static constexpr
int nwIdx = 1;
56 Valgrind::SetUndefined(*
this);
75 Scalar
Smin(
const unsigned phaseIdx)
const
76 { EnsureFinalized::check();
if (phaseIdx<Traits::numPhases)
return Smin_[phaseIdx];
return 0.0;}
81 Scalar
dS(
const unsigned phaseIdx)
const
82 { EnsureFinalized::check();
if (phaseIdx<Traits::numPhases)
return dS_[phaseIdx];
return 0.0;}
88 { EnsureFinalized::check();
return Sminpc_; }
94 { EnsureFinalized::check();
return dSpc_; }
99 Scalar
L(
const unsigned phaseIdx)
const
100 { EnsureFinalized::check();
if (phaseIdx<Traits::numPhases)
return L_[phaseIdx];
return 0.0;}
105 Scalar
E(
const unsigned phaseIdx)
const
106 { EnsureFinalized::check();
if (phaseIdx<Traits::numPhases)
return E_[phaseIdx];
return 0.0;}
111 Scalar
T(
const unsigned phaseIdx)
const
112 { EnsureFinalized::check();
if (phaseIdx<Traits::numPhases)
return T_[phaseIdx];
return 0.0;}
117 Scalar
Krt(
const unsigned phaseIdx)
const
118 { EnsureFinalized::check();
if (phaseIdx<Traits::numPhases)
return Krt_[phaseIdx];
return 0.0;}
124 { EnsureFinalized::check();
return Lpc_; }
130 { EnsureFinalized::check();
return Epc_; }
136 { EnsureFinalized::check();
return Tpc_; }
142 { EnsureFinalized::check();
return Pcir_; }
148 { EnsureFinalized::check();
return Pct_; }
156 template <
class Container>
159 setLETCoeffs(wIdx, letProp[2], letProp[3], letProp[4], letProp[5]);
162 Smin_[wIdx] = letProp[0];
163 dS_[wIdx] = letProp[1] - letProp[0];
172 template <
class Container>
175 setLETCoeffs(nwIdx, letProp[2], letProp[3], letProp[4], letProp[5]);
178 Smin_[nwIdx] = letProp[0];
179 dS_[nwIdx] = letProp[1] - letProp[0];
189 template <
class Container>
192 setLETPcCoeffs(letProp[2], letProp[3], letProp[4], letProp[5], letProp[6]);
193 Sminpc_ = letProp[0];
194 dSpc_ = 1.0 - letProp[0] - letProp[1];
202 void setLETCoeffs(
unsigned phaseIdx, Scalar
L, Scalar
E, Scalar
T, Scalar
Krt)
204 if (phaseIdx < Traits::numPhases) {
216 void setLETPcCoeffs(Scalar
L, Scalar
E, Scalar
T, Scalar
Pcir, Scalar
Pct)
229 void printLETCoeffs()
257 Scalar Smin_[Traits::numPhases];
258 Scalar dS_[Traits::numPhases];
260 Scalar L_[Traits::numPhases];
261 Scalar E_[Traits::numPhases];
262 Scalar T_[Traits::numPhases];
263 Scalar Krt_[Traits::numPhases];
Default implementation for asserting finalization of parameter objects.
Some templates to wrap the valgrind client request macros.
Default implementation for asserting finalization of parameter objects.
Definition: EnsureFinalized.hpp:47
void finalize()
Mark the object as finalized.
Definition: EnsureFinalized.hpp:75
Specification of the material parameters for the LET constitutive relations.
Definition: TwoPhaseLETCurvesParams.hpp:45
Scalar Sminpc() const
Returns the Epc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:87
void setPcnwSamples(const Container &letProp, const Container &)
Set the LET-related parameters for the capillary pressure curve of the non-wetting phase.
Definition: TwoPhaseLETCurvesParams.hpp:190
Scalar Smin(const unsigned phaseIdx) const
Returns the Smin_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:75
void finalize()
Calculate all dependent quantities once the independent quantities of the parameter object have been ...
Definition: TwoPhaseLETCurvesParams.hpp:65
Scalar T(const unsigned phaseIdx) const
Returns the T_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:111
Scalar dSpc() const
Returns the Epc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:93
Scalar Epc() const
Returns the Epc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:129
Scalar L(const unsigned phaseIdx) const
Returns the L_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:99
Scalar dS(const unsigned phaseIdx) const
Returns the dS_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:81
Scalar Lpc() const
Returns the Lpc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:123
Scalar E(const unsigned phaseIdx) const
Returns the E_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:105
void setKrwSamples(const Container &letProp, const Container &)
Set the LET-related parameters for the relative permeability curve of the wetting phase.
Definition: TwoPhaseLETCurvesParams.hpp:157
Scalar Krt(const unsigned phaseIdx) const
Returns the Krt_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:117
void setKrnSamples(const Container &letProp, const Container &)
Set the LET-related parameters for the relative permeability curve of the non-wetting phase.
Definition: TwoPhaseLETCurvesParams.hpp:173
Scalar Pct() const
Returns the Pct_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:147
Scalar Pcir() const
Returns the Pcir_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:141
Scalar Tpc() const
Returns the Tpc_ parameter.
Definition: TwoPhaseLETCurvesParams.hpp:135