OilfieldToolsNet/opm-common
4
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Add option to store/restore hysteresis parameters

Browse Source
This commit is contained in:
babrodtk 2017-02-24 17:09:40 +01:00
parent fc9885229c
commit 64bdc84a9f
8 changed files with 636 additions and 81 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

62
opm/material/fluidmatrixinteractions/EclDefaultMaterial.hpp
View File

@ -148,6 +148,68 @@ public:
Valgrind::CheckDefined(values[waterPhaseIdx]);
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void oilWaterHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
pcSwMdc = params.oilWaterParams().pcSwMdc();
krnSwMdc = params.oilWaterParams().krnSwMdc();
Valgrind::CheckDefined(pcSwMdc);
Valgrind::CheckDefined(krnSwMdc);
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
const double krwSw = 2.0; //Should not be used
params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void gasOilHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
pcSwMdc = params.gasOilParams().pcSwMdc();
krnSwMdc = params.gasOilParams().krnSwMdc();
Valgrind::CheckDefined(pcSwMdc);
Valgrind::CheckDefined(krnSwMdc);
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
const double krwSw = 2.0; //Should not be used
params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
}
/*!
* \brief Capillary pressure between the gas and the non-wetting
* liquid (i.e., oil) phase.

21
opm/material/fluidmatrixinteractions/EclHysteresisTwoPhaseLawParams.hpp
View File

@ -60,6 +60,9 @@ public:
EclHysteresisTwoPhaseLawParams()
{
// These are initialized to two (even though they represent saturations)
// to signify that the values are larger than physically possible, and force
// using the drainage curve before the first saturation update.
pcSwMdc_ = 2.0;
krnSwMdc_ = 2.0;
// krwSwMdc_ = 2.0;
@ -152,7 +155,7 @@ public:
{ pcSwMdc_ = value; }
/*!
* \brief Set the saturation of the wetting phase where the last switch from the main
* \brief Get the saturation of the wetting phase where the last switch from the main
* drainage curve to imbibition happend on the capillary pressure curve.
*/
Scalar pcSwMdc() const
@ -168,7 +171,7 @@ public:
// { krwSwMdc_ = value; };
/*!
* \brief Set the saturation of the wetting phase where the last switch from the main
* \brief Get the saturation of the wetting phase where the last switch from the main
* drainage curve to imbibition happend on the relperm curve for the
* wetting phase.
*/
@ -185,7 +188,7 @@ public:
{ krnSwMdc_ = value; }
/*!
* \brief Set the saturation of the wetting phase where the last switch from the main
* \brief Get the saturation of the wetting phase where the last switch from the main
* drainage curve to imbibition happend on the relperm curve for the
* non-wetting phase.
*/
@ -286,9 +289,9 @@ private:
Scalar SwKrnMdcImbibition = EffLawT::twoPhaseSatKrnInv(imbibitionParams(), krnMdcDrainage);
deltaSwImbKrn_ = SwKrnMdcImbibition - krnSwMdc_;
Scalar pcMdcDrainage = EffLawT::twoPhaseSatPcnw(drainageParams(), pcSwMdc_);
Scalar SwPcMdcImbibition = EffLawT::twoPhaseSatPcnwInv(imbibitionParams(), pcMdcDrainage);
deltaSwImbPc_ = SwPcMdcImbibition - pcSwMdc_;
// Scalar pcMdcDrainage = EffLawT::twoPhaseSatPcnw(drainageParams(), pcSwMdc_);
// Scalar SwPcMdcImbibition = EffLawT::twoPhaseSatPcnwInv(imbibitionParams(), pcMdcDrainage);
// deltaSwImbPc_ = SwPcMdcImbibition - pcSwMdc_;
// assert(std::abs(EffLawT::twoPhaseSatPcnw(imbibitionParams(), pcSwMdc_ + deltaSwImbPc_)
// - EffLawT::twoPhaseSatPcnw(drainageParams(), pcSwMdc_)) < 1e-8);
@ -311,16 +314,16 @@ private:
// largest wettinging phase saturation which is on the main-drainage curve. These are
// three different values because the sourounding code can choose to use different
// definitions for the saturations for different quantities
// Scalar krwSwMdc_;
//Scalar krwSwMdc_;
Scalar krnSwMdc_;
Scalar pcSwMdc_;
// offsets added to wetting phase saturation uf using the imbibition curves need to
// be used to calculate the wetting phase relperm, the non-wetting phase relperm and
// the capillary pressure
// Scalar deltaSwImbKrw_;
//Scalar deltaSwImbKrw_;
Scalar deltaSwImbKrn_;
Scalar deltaSwImbPc_;
//Scalar deltaSwImbPc_;
// trapped non-wetting phase saturation
//Scalar Sncrt_;

44
opm/material/fluidmatrixinteractions/EclMaterialLawManager.hpp
View File

@ -242,6 +242,50 @@ public:
MaterialLaw::updateHysteresis(*threePhaseParams, fluidState);
}
void oilWaterHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
unsigned elemIdx) const
{
if (!enableHysteresis()) {
OPM_THROW(std::runtime_error, "Cannot get hysteresis parameters if hysteresis not enabled.");
}
const auto& params = materialLawParams(elemIdx);
MaterialLaw::oilWaterHysteresisParams(pcSwMdc, krnSwMdc, params);
}
void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
unsigned elemIdx)
{
if (!enableHysteresis()) {
OPM_THROW(std::runtime_error, "Cannot set hysteresis parameters if hysteresis not enabled.");
}
auto& params = materialLawParams(elemIdx);
MaterialLaw::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc, params);
}
void gasOilHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
unsigned elemIdx) const
{
if (!enableHysteresis()) {
OPM_THROW(std::runtime_error, "Cannot get hysteresis parameters if hysteresis not enabled.");
}
const auto& params = materialLawParams(elemIdx);
MaterialLaw::gasOilHysteresisParams(pcSwMdc, krnSwMdc, params);
}
void setGasOilHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
unsigned elemIdx)
{
if (!enableHysteresis()) {
OPM_THROW(std::runtime_error, "Cannot set hysteresis parameters if hysteresis not enabled.");
}
auto& params = materialLawParams(elemIdx);
MaterialLaw::setGasOilHysteresisParams(pcSwMdc, krnSwMdc, params);
}
EclEpsScalingPoints<Scalar>& oilWaterScaledEpsPointsDrainage(unsigned elemIdx)
{
auto& materialParams = *materialLawParams_[elemIdx];

132
opm/material/fluidmatrixinteractions/EclMultiplexerMaterial.hpp
View File

@ -159,6 +159,138 @@ public:
}
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void oilWaterHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
switch (params.approach()) {
case EclStone1Approach:
Stone1Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclStone1Approach>());
break;
case EclStone2Approach:
Stone2Material::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclStone2Approach>());
break;
case EclDefaultApproach:
DefaultMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclDefaultApproach>());
break;
case EclTwoPhaseApproach:
TwoPhaseMaterial::oilWaterHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclTwoPhaseApproach>());
break;
}
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
switch (params.approach()) {
case EclStone1Approach:
Stone1Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclStone1Approach>());
break;
case EclStone2Approach:
Stone2Material::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclStone2Approach>());
break;
case EclDefaultApproach:
DefaultMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclDefaultApproach>());
break;
case EclTwoPhaseApproach:
TwoPhaseMaterial::setOilWaterHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclTwoPhaseApproach>());
break;
}
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void gasOilHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
switch (params.approach()) {
case EclStone1Approach:
Stone1Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclStone1Approach>());
break;
case EclStone2Approach:
Stone2Material::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclStone2Approach>());
break;
case EclDefaultApproach:
DefaultMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclDefaultApproach>());
break;
case EclTwoPhaseApproach:
TwoPhaseMaterial::gasOilHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclTwoPhaseApproach>());
break;
}
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
switch (params.approach()) {
case EclStone1Approach:
Stone1Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclStone1Approach>());
break;
case EclStone2Approach:
Stone2Material::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclStone2Approach>());
break;
case EclDefaultApproach:
DefaultMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclDefaultApproach>());
break;
case EclTwoPhaseApproach:
TwoPhaseMaterial::setGasOilHysteresisParams(pcSwMdc, krnSwMdc,
params.template getRealParams<EclTwoPhaseApproach>());
break;
}
}
/*!
* \brief Capillary pressure between the gas and the non-wetting
* liquid (i.e., oil) phase.

62
opm/material/fluidmatrixinteractions/EclStone1Material.hpp
View File

@ -146,6 +146,68 @@ public:
Valgrind::CheckDefined(values[waterPhaseIdx]);
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void oilWaterHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
pcSwMdc = params.oilWaterParams().pcSwMdc();
krnSwMdc = params.oilWaterParams().krnSwMdc();
Valgrind::CheckDefined(pcSwMdc);
Valgrind::CheckDefined(krnSwMdc);
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
const double krwSw = 2.0; //Should not be used
params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void gasOilHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
pcSwMdc = params.gasOilParams().pcSwMdc();
krnSwMdc = params.gasOilParams().krnSwMdc();
Valgrind::CheckDefined(pcSwMdc);
Valgrind::CheckDefined(krnSwMdc);
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
const double krwSw = 2.0; //Should not be used
params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
}
/*!
* \brief Capillary pressure between the gas and the non-wetting
* liquid (i.e., oil) phase.

62
opm/material/fluidmatrixinteractions/EclStone2Material.hpp
View File

@ -147,6 +147,68 @@ public:
Valgrind::CheckDefined(values[waterPhaseIdx]);
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void oilWaterHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
pcSwMdc = params.oilWaterParams().pcSwMdc();
krnSwMdc = params.oilWaterParams().krnSwMdc();
Valgrind::CheckDefined(pcSwMdc);
Valgrind::CheckDefined(krnSwMdc);
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
const double krwSw = 2.0; //Should not be used
params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void gasOilHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
pcSwMdc = params.gasOilParams().pcSwMdc();
krnSwMdc = params.gasOilParams().krnSwMdc();
Valgrind::CheckDefined(pcSwMdc);
Valgrind::CheckDefined(krnSwMdc);
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
const double krwSw = 2.0; //Should not be used
params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
}
/*!
* \brief Capillary pressure between the gas and the non-wetting
* liquid (i.e., oil) phase.

62
opm/material/fluidmatrixinteractions/EclTwoPhaseMaterial.hpp
View File

@ -163,6 +163,68 @@ public:
}
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void oilWaterHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
pcSwMdc = params.oilWaterParams().pcSwMdc();
krnSwMdc = params.oilWaterParams().krnSwMdc();
Valgrind::CheckDefined(pcSwMdc);
Valgrind::CheckDefined(krnSwMdc);
}
/*
* Hysteresis parameters for oil-water
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setOilWaterHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
const Scalar krwSw = 2.0; //Should not be used
params.oilWaterParams().update(pcSwMdc, krwSw, krnSwMdc);
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void gasOilHysteresisParams(Scalar& pcSwMdc,
Scalar& krnSwMdc,
const Params& params)
{
pcSwMdc = params.gasOilParams().pcSwMdc();
krnSwMdc = params.gasOilParams().krnSwMdc();
Valgrind::CheckDefined(pcSwMdc);
Valgrind::CheckDefined(krnSwMdc);
}
/*
* Hysteresis parameters for gas-oil
* @see EclHysteresisTwoPhaseLawParams::pcSwMdc(...)
* @see EclHysteresisTwoPhaseLawParams::krnSwMdc(...)
* \param params Parameters
*/
static void setGasOilHysteresisParams(const Scalar& pcSwMdc,
const Scalar& krnSwMdc,
Params& params)
{
const Scalar krwSw = 2.0; //Should not be used
params.gasOilParams().update(pcSwMdc, krwSw, krnSwMdc);
}
/*!
* \brief Capillary pressure between the gas and the non-wetting
* liquid (i.e., oil) phase.

272
tests/test_eclmateriallawmanager.cpp
View File

@ -68,52 +68,52 @@ static const char* fam1DeckString =
"GRID\n"
"\n"
"DX\n"
" 300*1000 /\n"
" 300*1000 /\n"
"DY\n"
" 300*1000 /\n"
" 300*1000 /\n"
"DZ\n"
" 100*20 100*30 100*50 /\n"
" 100*20 100*30 100*50 /\n"
"\n"
"TOPS\n"
" 100*8325 /\n"
" 100*8325 /\n"
"\n"
"\n"
"PROPS\n"
"\n"
"SWOF\n"
"0.12 0 1 0\n"
"0.18 4.64876033057851E-008 1 0\n"
"0.24 0.000000186 0.997 0\n"
"0.3 4.18388429752066E-007 0.98 0\n"
"0.36 7.43801652892562E-007 0.7 0\n"
"0.42 1.16219008264463E-006 0.35 0\n"
"0.48 1.67355371900826E-006 0.2 0\n"
"0.54 2.27789256198347E-006 0.09 0\n"
"0.6 2.97520661157025E-006 0.021 0\n"
"0.66 3.7654958677686E-006 0.01 0\n"
"0.72 4.64876033057851E-006 0.001 0\n"
"0.78 0.000005625 0.0001 0\n"
"0.84 6.69421487603306E-006 0 0\n"
"0.91 8.05914256198347E-006 0 0\n"
"1 0.984 0 0 /\n"
"0.12 0 1 0\n"
"0.18 4.64876033057851E-008 1 0\n"
"0.24 0.000000186 0.997 0\n"
"0.3 4.18388429752066E-007 0.98 0\n"
"0.36 7.43801652892562E-007 0.7 0\n"
"0.42 1.16219008264463E-006 0.35 0\n"
"0.48 1.67355371900826E-006 0.2 0\n"
"0.54 2.27789256198347E-006 0.09 0\n"
"0.6 2.97520661157025E-006 0.021 0\n"
"0.66 3.7654958677686E-006 0.01 0\n"
"0.72 4.64876033057851E-006 0.001 0\n"
"0.78 0.000005625 0.0001 0\n"
"0.84 6.69421487603306E-006 0 0\n"
"0.91 8.05914256198347E-006 0 0\n"
"1 0.984 0 0 /\n"
"\n"
"\n"
"SGOF\n"
"0 0 1 0\n"
"0.001 0 1 0\n"
"0.02 0 0.997 0\n"
"0.05 0.005 0.980 0\n"
"0.12 0.025 0.700 0\n"
"0.2 0.075 0.350 0\n"
"0.25 0.125 0.200 0\n"
"0.3 0.190 0.090 0\n"
"0.4 0.410 0.021 0\n"
"0.45 0.60 0.010 0\n"
"0.5 0.72 0.001 0\n"
"0.6 0.87 0.0001 0\n"
"0.7 0.94 0.000 0\n"
"0.85 0.98 0.000 0\n"
"0.88 0.984 0.000 0 /\n";
"0 0 1 0\n"
"0.001 0 1 0\n"
"0.02 0 0.997 0\n"
"0.05 0.005 0.980 0\n"
"0.12 0.025 0.700 0\n"
"0.2 0.075 0.350 0\n"
"0.25 0.125 0.200 0\n"
"0.3 0.190 0.090 0\n"
"0.4 0.410 0.021 0\n"
"0.45 0.60 0.010 0\n"
"0.5 0.72 0.001 0\n"
"0.6 0.87 0.0001 0\n"
"0.7 0.94 0.000 0\n"
"0.85 0.98 0.000 0\n"
"0.88 0.984 0.000 0 /\n";
static const char* fam2DeckString =
"RUNSPEC\n"
@ -135,62 +135,62 @@ static const char* fam2DeckString =
"GRID\n"
"\n"
"DX\n"
" 300*1000 /\n"
" 300*1000 /\n"
"DY\n"
" 300*1000 /\n"
" 300*1000 /\n"
"DZ\n"
" 100*20 100*30 100*50 /\n"
" 100*20 100*30 100*50 /\n"
"\n"
"TOPS\n"
" 100*8325 /\n"
" 100*8325 /\n"
"\n"
"\n"
"PROPS\n"
"\n"
"PVTW\n"
" 4017.55 1.038 3.22E-6 0.318 0.0 /\n"
" 4017.55 1.038 3.22E-6 0.318 0.0 /\n"
"\n"
"\n"
"SWFN\n"
"0.12 0 0\n"
"0.18 4.64876033057851E-008 0\n"
"0.24 0.000000186 0\n"
"0.3 4.18388429752066E-007 0\n"
"0.36 7.43801652892562E-007 0\n"
"0.42 1.16219008264463E-006 0\n"
"0.48 1.67355371900826E-006 0\n"
"0.54 2.27789256198347E-006 0\n"
"0.6 2.97520661157025E-006 0\n"
"0.66 3.7654958677686E-006 0\n"
"0.72 4.64876033057851E-006 0\n"
"0.78 0.000005625 0\n"
"0.84 6.69421487603306E-006 0\n"
"0.91 8.05914256198347E-006 0\n"
"1 0.984 0 /\n"
"0.12 0 0\n"
"0.18 4.64876033057851E-008 0\n"
"0.24 0.000000186 0\n"
"0.3 4.18388429752066E-007 0\n"
"0.36 7.43801652892562E-007 0\n"
"0.42 1.16219008264463E-006 0\n"
"0.48 1.67355371900826E-006 0\n"
"0.54 2.27789256198347E-006 0\n"
"0.6 2.97520661157025E-006 0\n"
"0.66 3.7654958677686E-006 0\n"
"0.72 4.64876033057851E-006 0\n"
"0.78 0.000005625 0\n"
"0.84 6.69421487603306E-006 0\n"
"0.91 8.05914256198347E-006 0\n"
"1 0.984 0 /\n"
"\n"
"\n"
"SGFN\n"
"0 0 0\n"
"0.001 0 0\n"
"0.02 0 0\n"
"0.05 0.005 0\n"
"0.12 0.025 0\n"
"0.2 0.075 0\n"
"0.25 0.125 0\n"
"0.3 0.190 0\n"
"0.4 0.410 0\n"
"0.45 0.60 0\n"
"0.5 0.72 0\n"
"0.6 0.87 0\n"
"0.7 0.94 0\n"
"0.85 0.98 0\n"
"0.88 0.984 0 /\n"
"0 0 0\n"
"0.001 0 0\n"
"0.02 0 0\n"
"0.05 0.005 0\n"
"0.12 0.025 0\n"
"0.2 0.075 0\n"
"0.25 0.125 0\n"
"0.3 0.190 0\n"
"0.4 0.410 0\n"
"0.45 0.60 0\n"
"0.5 0.72 0\n"
"0.6 0.87 0\n"
"0.7 0.94 0\n"
"0.85 0.98 0\n"
"0.88 0.984 0 /\n"
"\n"
"SOF3\n"
" 0 0 0 \n"
" 0.03 0 0 \n"
" 0.09 0 0 \n"
" 0.16 0 0 \n"
" 0.16 0 0 \n"
" 0.18 1* 0 \n"
" 0.22 0.0001 1* \n"
" 0.28 0.001 0.0001 \n"
@ -213,6 +213,82 @@ static const char* fam2DeckString =
" 0.88 1 1 / \n"
"\n";
//Taken as a mix of the SPE1 cases above, and Norne to enable hysteresis
static const char* hysterDeckString =
"RUNSPEC\n"
"\n"
"DIMENS\n"
" 10 10 3 /\n"
"\n"
"TABDIMS\n"
"/\n"
"\n"
"OIL\n"
"GAS\n"
"WATER\n"
"\n"
"DISGAS\n"
"\n"
"FIELD\n"
"\n"
"GRID\n"
"\n"
"DX\n"
" 300*1000 /\n"
"DY\n"
" 300*1000 /\n"
"DZ\n"
" 100*20 100*30 100*50 /\n"
"\n"
"TOPS\n"
" 100*8325 /\n"
"\n"
"\n"
"EHYSTR\n"
"0.1 0 0.1 1* KR /\n"
"\n"
"SATOPTS\n"
"HYSTER /\n"
"\n"
"PROPS\n"
"\n"
"SWOF\n"
"0.12 0 1 0\n"
"0.18 4.64876033057851E-008 1 0\n"
"0.24 0.000000186 0.997 0\n"
"0.3 4.18388429752066E-007 0.98 0\n"
"0.36 7.43801652892562E-007 0.7 0\n"
"0.42 1.16219008264463E-006 0.35 0\n"
"0.48 1.67355371900826E-006 0.2 0\n"
"0.54 2.27789256198347E-006 0.09 0\n"
"0.6 2.97520661157025E-006 0.021 0\n"
"0.66 3.7654958677686E-006 0.01 0\n"
"0.72 4.64876033057851E-006 0.001 0\n"
"0.78 0.000005625 0.0001 0\n"
"0.84 6.69421487603306E-006 0 0\n"
"0.91 8.05914256198347E-006 0 0\n"
"1 0.984 0 0 /\n"
"\n"
"\n"
"SGOF\n"
"0 0 1 0\n"
"0.001 0 1 0\n"
"0.02 0 0.997 0\n"
"0.05 0.005 0.980 0\n"
"0.12 0.025 0.700 0\n"
"0.2 0.075 0.350 0\n"
"0.25 0.125 0.200 0\n"
"0.3 0.190 0.090 0\n"
"0.4 0.410 0.021 0\n"
"0.45 0.60 0.010 0\n"
"0.5 0.72 0.001 0\n"
"0.6 0.87 0.0001 0\n"
"0.7 0.94 0.000 0\n"
"0.85 0.98 0.000 0\n"
"0.88 0.984 0.000 0 /\n";
template <class Scalar>
inline void testAll()
{
@ -280,8 +356,24 @@ inline void testAll()
OPM_THROW(std::logic_error,
"Discrepancy between the deck and the EclMaterialLawManager");
const auto hysterDeck = parser.parseString(hysterDeckString, parseContext);
const Opm::EclipseState hysterEclState(hysterDeck, parseContext);
Opm::EclMaterialLawManager<MaterialTraits> hysterMaterialLawManager;
hysterMaterialLawManager.initFromDeck(hysterDeck, hysterEclState, compressedToCartesianIdx);
if (hysterMaterialLawManager.enableEndPointScaling())
OPM_THROW(std::logic_error,
"Discrepancy between the deck and the EclMaterialLawManager");
if (hysterMaterialLawManager.enableHysteresis() != true)
OPM_THROW(std::logic_error,
"Discrepancy between the deck and the EclMaterialLawManager");
// make sure that the saturation functions for both keyword families are
// identical
// identical, and that setting and getting the hysteresis parameters works
for (unsigned elemIdx = 0; elemIdx < n; ++ elemIdx) {
for (int i = -10; i < 120; ++ i) {
Scalar Sw = Scalar(i)/100;
@ -318,7 +410,43 @@ inline void testAll()
"Discrepancy between capillary pressure of family 1 and family 2 keywords");
if (std::abs(krFam1[phaseIdx] - krFam2[phaseIdx]) > 1e-1)
OPM_THROW(std::logic_error,
"Discrepancy between capillary pressure of family 1 and family 2 keywords");
"Discrepancy between relative permeabilities of family 1 and family 2 keywords");
}
// This should ideally test each of the materials (stone1, stone2, default, two-phase),
// but currently only tests default
const Scalar pcSwMdc_in[2] = { 1.0/2.0, 1.0/3.0 };
const Scalar krnSwMdc_in[2] = { 1.0/5.0, 1.0/7.0 };
hysterMaterialLawManager.setOilWaterHysteresisParams(
pcSwMdc_in[0],
krnSwMdc_in[0],
elemIdx);
hysterMaterialLawManager.setGasOilHysteresisParams(
pcSwMdc_in[1],
krnSwMdc_in[1],
elemIdx);
Scalar pcSwMdc_out[2];
Scalar krnSwMdc_out[2];
Scalar deltaSwImbKrn_out[2];
hysterMaterialLawManager.oilWaterHysteresisParams(
pcSwMdc_out[0],
krnSwMdc_out[0],
elemIdx);
hysterMaterialLawManager.gasOilHysteresisParams(
pcSwMdc_out[1],
krnSwMdc_out[1],
elemIdx);
for (unsigned phasePairIdx = 0; phasePairIdx < 2; ++ phasePairIdx) {
if ((pcSwMdc_in[phasePairIdx] - pcSwMdc_out[phasePairIdx]) != 0.0)
OPM_THROW(std::logic_error,
"Hysteresis parameters did not propagate correctly");
if ((krnSwMdc_in[phasePairIdx] - krnSwMdc_out[phasePairIdx]) != 0.0)
OPM_THROW(std::logic_error,
"Hysteresis parameters did not propagate correctly");
}
}
}