|
|
|
|
@@ -33,25 +33,81 @@
|
|
|
|
|
|
|
|
|
|
namespace Opm {
|
|
|
|
|
/*!
|
|
|
|
|
* \ingroup fluidmatrixinteractionslaws
|
|
|
|
|
* \ingroup FluidMatrixInteractions
|
|
|
|
|
*
|
|
|
|
|
* \brief Implementation of the van Genuchten capillary pressure <->
|
|
|
|
|
* \brief Implementation of the van Genuchten capillary pressure -
|
|
|
|
|
* saturation relation.
|
|
|
|
|
*
|
|
|
|
|
* This class only implements the "raw" van-Genuchten curves as static
|
|
|
|
|
* members and doesn't concern itself converting absolute to effective
|
|
|
|
|
* saturations and vice versa.
|
|
|
|
|
*
|
|
|
|
|
* For general info: EffToAbsLaw
|
|
|
|
|
* The converion from and to effective saturations can be done using,
|
|
|
|
|
* e.g. EffToAbsLaw.
|
|
|
|
|
*
|
|
|
|
|
* \see VanGenuchtenParams
|
|
|
|
|
*/
|
|
|
|
|
template <class ScalarT, class ParamsT = VanGenuchtenParams<ScalarT> >
|
|
|
|
|
template <class ScalarT,
|
|
|
|
|
int wPhaseIdxV,
|
|
|
|
|
int nPhaseIdxV,
|
|
|
|
|
class ParamsT = VanGenuchtenParams<ScalarT> >
|
|
|
|
|
class VanGenuchten
|
|
|
|
|
{
|
|
|
|
|
public:
|
|
|
|
|
typedef ParamsT Params;
|
|
|
|
|
typedef typename Params::Scalar Scalar;
|
|
|
|
|
//! The type of the parameter objects for this law
|
|
|
|
|
typedef ParamsT Params;
|
|
|
|
|
|
|
|
|
|
//! The type of the scalar values for this law
|
|
|
|
|
typedef typename Params::Scalar Scalar;
|
|
|
|
|
|
|
|
|
|
//! The number of fluid phases to which this capillary pressure law applies
|
|
|
|
|
enum { numPhases = 2 };
|
|
|
|
|
|
|
|
|
|
//! The index of the wetting phase
|
|
|
|
|
enum { wPhaseIdx = wPhaseIdxV };
|
|
|
|
|
|
|
|
|
|
//! The index of the non-wetting phase
|
|
|
|
|
enum { nPhaseIdx = nPhaseIdxV };
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
* \brief The capillary pressure-saturation curves according to van Genuchten.
|
|
|
|
|
*
|
|
|
|
|
* Van Genuchten's empirical capillary pressure <-> saturation
|
|
|
|
|
* function is given by
|
|
|
|
|
* \f[
|
|
|
|
|
* p_{c,wn} = p_n - p_w = ({S_w}^{-1/m} - 1)^{1/n}/\alpha
|
|
|
|
|
* \f]
|
|
|
|
|
*
|
|
|
|
|
* \param values A random access container which stores the
|
|
|
|
|
* relative pressure of each fluid phase.
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state for which the capillary pressure
|
|
|
|
|
* ought to be calculated
|
|
|
|
|
*/
|
|
|
|
|
template <class Container, class FluidState>
|
|
|
|
|
static void capillaryPressures(Container &values, const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
values[wPhaseIdx] = 0.0; // reference phase
|
|
|
|
|
values[nPhaseIdx] = pcwn(params, fs);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
* \brief The relative permeability-saturation curves according to van Genuchten.
|
|
|
|
|
*
|
|
|
|
|
* \param values A random access container which stores the
|
|
|
|
|
* relative permeability of each fluid phase.
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state for which the relative permeabilities
|
|
|
|
|
* ought to be calculated
|
|
|
|
|
*/
|
|
|
|
|
template <class Container, class FluidState>
|
|
|
|
|
static void relativePermeabilities(Container &values, const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
values[wPhaseIdx] = krw(params, fs);
|
|
|
|
|
values[nPhaseIdx] = krn(params, fs);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
* \brief The capillary pressure-saturation curve according to van Genuchten.
|
|
|
|
|
@@ -59,17 +115,20 @@ public:
|
|
|
|
|
* Van Genuchten's empirical capillary pressure <-> saturation
|
|
|
|
|
* function is given by
|
|
|
|
|
* \f[
|
|
|
|
|
p_C = (\overline{S}_w^{-1/m} - 1)^{1/n}/\alpha
|
|
|
|
|
\f]
|
|
|
|
|
* \param Swe Effective saturation of the wetting phase \f$\overline{S}_w\f$
|
|
|
|
|
* \param params A container object that is populated with the appropriate coefficients for the respective law.
|
|
|
|
|
* Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container
|
|
|
|
|
* is constructed accordingly. Afterwards the values are set there, too.
|
|
|
|
|
* p_{c,wn} = p_n - p_w = ({S_w}^{-1/m} - 1)^{1/n}/\alpha
|
|
|
|
|
* \f]
|
|
|
|
|
*
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state for which the capillary pressure
|
|
|
|
|
* ought to be calculated
|
|
|
|
|
*/
|
|
|
|
|
static Scalar pC(const Params ¶ms, Scalar Swe)
|
|
|
|
|
template <class FluidState>
|
|
|
|
|
static Scalar pcwn(const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
assert(0 <= Swe && Swe <= 1);
|
|
|
|
|
return std::pow(std::pow(Swe, -1.0/params.vgM()) - 1, 1.0/params.vgN())/params.vgAlpha();
|
|
|
|
|
Scalar Sw = fs.saturation(wPhaseIdx);
|
|
|
|
|
assert(0 <= Sw && Sw <= 1);
|
|
|
|
|
return std::pow(std::pow(Sw, -1.0/params.vgM()) - 1, 1.0/params.vgN())/params.vgAlpha();
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
@@ -77,119 +136,112 @@ public:
|
|
|
|
|
*
|
|
|
|
|
* This is the inverse of the capillary pressure-saturation curve:
|
|
|
|
|
* \f[
|
|
|
|
|
\overline{S}_w = {p_C}^{-1} = ((\alpha p_C)^n + 1)^{-m}
|
|
|
|
|
\f]
|
|
|
|
|
* S_w = {p_C}^{-1} = ((\alpha p_C)^n + 1)^{-m}
|
|
|
|
|
* \f]
|
|
|
|
|
*
|
|
|
|
|
* \param pC Capillary pressure \f$p_C\f$
|
|
|
|
|
* \param params A container object that is populated with the appropriate coefficients for the respective law.
|
|
|
|
|
* Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container
|
|
|
|
|
* is constructed accordingly. Afterwards the values are set there, too.
|
|
|
|
|
* \return The effective saturation of the wetting phase \f$\overline{S}_w\f$
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state containing valid phase pressures
|
|
|
|
|
*/
|
|
|
|
|
static Scalar Sw(const Params ¶ms, Scalar pC)
|
|
|
|
|
template <class FluidState>
|
|
|
|
|
static Scalar Sw(const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
Scalar pC = fs.pressure(nPhaseIdx) - fs.pressure(wPhaseIdx);
|
|
|
|
|
assert(pC >= 0);
|
|
|
|
|
|
|
|
|
|
return std::pow(std::pow(params.vgAlpha()*pC, params.vgN()) + 1, -params.vgM());
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
* \brief The partial derivative of the capillary
|
|
|
|
|
* pressure w.r.t. the effective saturation according to van Genuchten.
|
|
|
|
|
* \brief The partial derivative of the capillary pressure with
|
|
|
|
|
* regard to the saturation according to van Genuchten.
|
|
|
|
|
*
|
|
|
|
|
* This is equivalent to
|
|
|
|
|
* \f[
|
|
|
|
|
\frac{\partial p_C}{\partial \overline{S}_w} =
|
|
|
|
|
-\frac{1}{\alpha} (\overline{S}_w^{-1/m} - 1)^{1/n - }
|
|
|
|
|
\overline{S}_w^{-1/m} / \overline{S}_w / m
|
|
|
|
|
\f]
|
|
|
|
|
* \frac{\partial p_C}{\partial S_w} =
|
|
|
|
|
* -\frac{1}{\alpha n} ({S_w}^{-1/m} - 1)^{1/n - 1} {S_w}^{-1/m - 1}
|
|
|
|
|
* \f]
|
|
|
|
|
*
|
|
|
|
|
* \param Swe Effective saturation of the wetting phase \f$\overline{S}_w\f$
|
|
|
|
|
* \param params A container object that is populated with the appropriate coefficients for the respective law.
|
|
|
|
|
* Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container
|
|
|
|
|
* is constructed accordingly. Afterwards the values are set there, too.
|
|
|
|
|
*/
|
|
|
|
|
static Scalar dpC_dSw(const Params ¶ms, Scalar Swe)
|
|
|
|
|
{
|
|
|
|
|
assert(0 <= Swe && Swe <= 1);
|
|
|
|
|
|
|
|
|
|
Scalar powSwe = std::pow(Swe, -1/params.vgM());
|
|
|
|
|
return - 1/params.vgAlpha() * std::pow(powSwe - 1, 1/params.vgN() - 1)/params.vgN()
|
|
|
|
|
* powSwe/Swe/params.vgM();
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
* \brief The partial derivative of the effective
|
|
|
|
|
* saturation to the capillary pressure according to van Genuchten.
|
|
|
|
|
*
|
|
|
|
|
* \param pC Capillary pressure \f$p_C\f$
|
|
|
|
|
* \param params A container object that is populated with the appropriate coefficients for the respective law.
|
|
|
|
|
* Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container
|
|
|
|
|
* is constructed accordingly. Afterwards the values are set there, too.
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state containing valid saturations
|
|
|
|
|
*/
|
|
|
|
|
static Scalar dSw_dpC(const Params ¶ms, Scalar pC)
|
|
|
|
|
template <class FluidState>
|
|
|
|
|
static Scalar dpcwn_dSw(const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
assert(pC >= 0);
|
|
|
|
|
Scalar Sw = fs.saturation(wPhaseIdx);
|
|
|
|
|
|
|
|
|
|
Scalar powAlphaPc = std::pow(params.vgAlpha()*pC, params.vgN());
|
|
|
|
|
return -std::pow(powAlphaPc + 1, -params.vgM()-1)*
|
|
|
|
|
params.vgM()*powAlphaPc/pC*params.vgN();
|
|
|
|
|
assert(0 < Sw && Sw < 1);
|
|
|
|
|
|
|
|
|
|
Scalar powSw = std::pow(Sw, -1/params.vgM());
|
|
|
|
|
return
|
|
|
|
|
- 1/(params.vgAlpha() * params.vgN() * Sw)
|
|
|
|
|
* std::pow(powSw - 1, 1/params.vgN() - 1) * powSw;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
* \brief The relative permeability for the wetting phase of
|
|
|
|
|
* the medium implied by van Genuchten's
|
|
|
|
|
* \brief The relative permeability for the wetting phase of the
|
|
|
|
|
* medium according to van Genuchten's curve with Mualem
|
|
|
|
|
* parameterization.
|
|
|
|
|
*
|
|
|
|
|
* \param Swe The mobile saturation of the wetting phase.
|
|
|
|
|
* \param params A container object that is populated with the appropriate coefficients for the respective law.
|
|
|
|
|
* Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container
|
|
|
|
|
* is constructed accordingly. Afterwards the values are set there, too. */
|
|
|
|
|
static Scalar krw(const Params ¶ms, Scalar Swe)
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state for which the relative permeability
|
|
|
|
|
* ought to be calculated
|
|
|
|
|
*/
|
|
|
|
|
template <class FluidState>
|
|
|
|
|
static Scalar krw(const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
assert(0 <= Swe && Swe <= 1);
|
|
|
|
|
Scalar Sw = fs.saturation(wPhaseIdx);
|
|
|
|
|
|
|
|
|
|
Scalar r = 1. - std::pow(1 - std::pow(Swe, 1/params.vgM()), params.vgM());
|
|
|
|
|
return std::sqrt(Swe)*r*r;
|
|
|
|
|
assert(0 <= Sw && Sw <= 1);
|
|
|
|
|
|
|
|
|
|
Scalar r = 1. - std::pow(1 - std::pow(Sw, 1/params.vgM()), params.vgM());
|
|
|
|
|
return std::sqrt(Sw)*r*r;
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
* \brief The derivative of the relative permeability for the
|
|
|
|
|
* \brief The derivative of the relative permeability of the
|
|
|
|
|
* wetting phase in regard to the wetting saturation of the
|
|
|
|
|
* medium implied by the van Genuchten parameterization.
|
|
|
|
|
* medium implied according to the van Genuchten curve with
|
|
|
|
|
* Mualem parameters.
|
|
|
|
|
*
|
|
|
|
|
* \param Swe The mobile saturation of the wetting phase.
|
|
|
|
|
* \param params A container object that is populated with the appropriate coefficients for the respective law.
|
|
|
|
|
* Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container
|
|
|
|
|
* is constructed accordingly. Afterwards the values are set there, too.
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state for which the derivative
|
|
|
|
|
* ought to be calculated
|
|
|
|
|
*/
|
|
|
|
|
static Scalar dkrw_dSw(const Params ¶ms, Scalar Swe)
|
|
|
|
|
template <class FluidState>
|
|
|
|
|
static Scalar dkrw_dSw(const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
assert(0 <= Swe && Swe <= 1);
|
|
|
|
|
Scalar Sw = fs.saturation(wPhaseIdx);
|
|
|
|
|
|
|
|
|
|
const Scalar x = 1 - std::pow(Swe, 1.0/params.vgM());
|
|
|
|
|
assert(0 <= Sw && Sw <= 1);
|
|
|
|
|
|
|
|
|
|
const Scalar x = 1 - std::pow(Sw, 1.0/params.vgM());
|
|
|
|
|
const Scalar xToM = std::pow(x, params.vgM());
|
|
|
|
|
return (1 - xToM)/std::sqrt(Swe) * ( (1 - xToM)/2 + 2*xToM*(1-x)/x );
|
|
|
|
|
return (1 - xToM)/std::sqrt(Sw) * ( (1 - xToM)/2 + 2*xToM*(1-x)/x );
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
* \brief The relative permeability for the non-wetting phase
|
|
|
|
|
* of the medium implied by van Genuchten's
|
|
|
|
|
* parameterization.
|
|
|
|
|
* of the medium according to van Genuchten.
|
|
|
|
|
*
|
|
|
|
|
* \param Swe The mobile saturation of the wetting phase.
|
|
|
|
|
* \param params A container object that is populated with the appropriate coefficients for the respective law.
|
|
|
|
|
* Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container
|
|
|
|
|
* is constructed accordingly. Afterwards the values are set there, too.
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state for which the derivative
|
|
|
|
|
* ought to be calculated
|
|
|
|
|
*/
|
|
|
|
|
static Scalar krn(const Params ¶ms, Scalar Swe)
|
|
|
|
|
static Scalar krn(const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
assert(0 <= Swe && Swe <= 1);
|
|
|
|
|
Scalar Sw = fs.saturation(wPhaseIdx);
|
|
|
|
|
|
|
|
|
|
assert(0 <= Sw && Sw <= 1);
|
|
|
|
|
|
|
|
|
|
return
|
|
|
|
|
std::pow(1 - Swe, 1.0/3) *
|
|
|
|
|
std::pow(1 - std::pow(Swe, 1/params.vgM()), 2*params.vgM());
|
|
|
|
|
std::pow(1 - Sw, 1.0/3) *
|
|
|
|
|
std::pow(1 - std::pow(Sw, 1/params.vgM()), 2*params.vgM());
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
/*!
|
|
|
|
|
@@ -198,22 +250,24 @@ public:
|
|
|
|
|
* the medium as implied by the van Genuchten
|
|
|
|
|
* parameterization.
|
|
|
|
|
*
|
|
|
|
|
* \param Swe The mobile saturation of the wetting phase.
|
|
|
|
|
* \param params A container object that is populated with the appropriate coefficients for the respective law.
|
|
|
|
|
* Therefore, in the (problem specific) spatialParameters first, the material law is chosen, and then the params container
|
|
|
|
|
* is constructed accordingly. Afterwards the values are set there, too.
|
|
|
|
|
* \param params The parameter object expressing the coefficients
|
|
|
|
|
* required by the van Genuchten law.
|
|
|
|
|
* \param fs The fluid state for which the derivative
|
|
|
|
|
* ought to be calculated
|
|
|
|
|
*/
|
|
|
|
|
static Scalar dkrn_dSw(const Params ¶ms, Scalar Swe)
|
|
|
|
|
template <class FluidState>
|
|
|
|
|
static Scalar dkrn_dSw(const Params ¶ms, const FluidState &fs)
|
|
|
|
|
{
|
|
|
|
|
assert(0 <= Swe && Swe <= 1);
|
|
|
|
|
Scalar fs = fs.saturation(wPhaseIdx);
|
|
|
|
|
|
|
|
|
|
const Scalar x = std::pow(Swe, 1.0/params.vgM());
|
|
|
|
|
assert(0 <= Sw && Sw <= 1);
|
|
|
|
|
|
|
|
|
|
const Scalar x = std::pow(Sw, 1.0/params.vgM());
|
|
|
|
|
return
|
|
|
|
|
-std::pow(1 - x, 2*params.vgM())
|
|
|
|
|
*std::pow(1 - Swe, -2/3)
|
|
|
|
|
*(1.0/3 + 2*x/Swe);
|
|
|
|
|
*std::pow(1 - Sw, -2/3)
|
|
|
|
|
*(1.0/3 + 2*x/Sw);
|
|
|
|
|
}
|
|
|
|
|
|
|
|
|
|
};
|
|
|
|
|
} // namespace Opm
|
|
|
|
|
|
|
|
|
|
|
Andreas Lauser