Merge pull request #860 from daavid00/dyncppres

Support for dynamic scaling of capillary pressure due to salt precipitation
This commit is contained in:
Tor Harald Sandve 2024-01-15 11:12:52 +01:00 committed by GitHub
commit 37612b95db
5 changed files with 60 additions and 6 deletions

View File

@ -35,6 +35,7 @@
#if HAVE_ECL_INPUT
#include <opm/input/eclipse/EclipseState/EclipseState.hpp>
#include <opm/input/eclipse/EclipseState/Tables/PvtwsaltTable.hpp>
#include <opm/input/eclipse/EclipseState/Tables/PcfactTable.hpp>
#include <opm/input/eclipse/EclipseState/Tables/PermfactTable.hpp>
#include <opm/input/eclipse/EclipseState/Tables/SaltSolubilityTable.hpp>
#include <opm/input/eclipse/EclipseState/Tables/TableManager.hpp>
@ -144,6 +145,16 @@ public:
params_.saltdenTable_[pvtRegionIdx] = saltsolTable.getSaltdenColumn().front();
}
}
const TableContainer& pcfactTables = tableManager.getPcfactTables();
if (!pcfactTables.empty()) {
unsigned numSatRegions = tableManager.getTabdims().getNumSatTables();
params_.pcfactTable_.resize(numSatRegions);
for (size_t i = 0; i < pcfactTables.size(); ++i) {
const PcfactTable& pcfactTable = pcfactTables.getTable<PcfactTable>(i);
params_.pcfactTable_[i].setXYContainers(pcfactTable.getPorosityChangeColumn(), pcfactTable.getPcMultiplierColumn());
}
}
}
}
#endif
@ -328,6 +339,11 @@ public:
return params_.bdensityTable_[pvtnumRegionIdx];
}
static const TabulatedFunction& pcfactTable(unsigned satnumRegionIdx)
{
return params_.pcfactTable_[satnumRegionIdx];
}
static const TabulatedFunction& permfactTable(const ElementContext& elemCtx,
unsigned scvIdx,
unsigned timeIdx)
@ -367,6 +383,14 @@ public:
return !params_.saltsolTable_.empty();
}
static bool hasPcfactTables()
{
if constexpr (enableSaltPrecipitation)
return !params_.pcfactTable_.empty();
else
return false;
}
static Scalar saltSol(unsigned regionIdx) {
return params_.saltsolTable_[regionIdx];
}

View File

@ -40,6 +40,7 @@ struct BlackOilBrineParams {
using TabulatedFunction = Tabulated1DFunction<Scalar>;
std::vector<TabulatedFunction> bdensityTable_;
std::vector<TabulatedFunction> pcfactTable_;
std::vector<TabulatedFunction> permfactTable_;
std::vector<Scalar> saltsolTable_;
std::vector<Scalar> saltdenTable_;

View File

@ -130,6 +130,7 @@ class BlackOilIntensiveQuantities
using DispersionIntensiveQuantities = BlackOilDispersionIntensiveQuantities<TypeTag, enableDispersion>;
using DirectionalMobilityPtr = Opm::Utility::CopyablePtr<DirectionalMobility<TypeTag, Evaluation>>;
using BrineModule = BlackOilBrineModule<TypeTag>;
public:
@ -257,6 +258,19 @@ public:
MaterialLaw::capillaryPressures(pC, materialParams, fluidState_);
problem.updateRelperms(mobility_, dirMob_, fluidState_, globalSpaceIdx);
// scaling the capillary pressure due to salt precipitation
if (BrineModule::hasPcfactTables() && priVars.primaryVarsMeaningBrine() == PrimaryVariables::BrineMeaning::Sp) {
unsigned satnumRegionIdx = elemCtx.problem().satnumRegionIndex(elemCtx, dofIdx, timeIdx);
const Evaluation Sp = priVars.makeEvaluation(Indices::saltConcentrationIdx, timeIdx);
const Evaluation porosityFactor = min(1.0 - Sp, 1.0); //phi/phi_0
const auto& pcfactTable = BrineModule::pcfactTable(satnumRegionIdx);
const Evaluation pcFactor = pcfactTable.eval(porosityFactor, /*extrapolation=*/true);
for (unsigned phaseIdx = 0; phaseIdx < numPhases; ++phaseIdx)
if (FluidSystem::phaseIsActive(phaseIdx)) {
pC[phaseIdx] *= pcFactor;
}
}
// oil is the reference phase for pressure
if (priVars.primaryVarsMeaningPressure() == PrimaryVariables::PressureMeaning::Pg) {
const Evaluation& pg = priVars.makeEvaluation(Indices::pressureSwitchIdx, timeIdx);

View File

@ -617,6 +617,17 @@ public:
return changed;
}
if (BrineModule::hasPcfactTables() && primaryVarsMeaningBrine() == BrineMeaning::Sp) {
unsigned satnumRegionIdx = problem.satnumRegionIndex(globalDofIdx);
Scalar Sp = saltConcentration_();
Scalar porosityFactor = min(1.0 - Sp, 1.0); //phi/phi_0
const auto& pcfactTable = BrineModule::pcfactTable(satnumRegionIdx);
pcFactor_ = pcfactTable.eval(porosityFactor, /*extrapolation=*/true);
}
else {
pcFactor_ = 1.0;
}
switch(primaryVarsMeaningWater()) {
case WaterMeaning::Sw:
{
@ -628,7 +639,7 @@ public:
const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
Scalar so = 1.0 - sg - solventSaturation_();
computeCapillaryPressures_(pC, so, sg + solventSaturation_(), /*sw=*/ 0.0, matParams);
p += (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
p += pcFactor_ * (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
}
Scalar rvwSat = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(pvtRegionIdx_,
T,
@ -648,7 +659,7 @@ public:
const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
Scalar so = 1.0 - sw - solventSaturation_();
computeCapillaryPressures_(pC, so, /*sg=*/ 0.0, sw, matParams);
Scalar pw = pg + (pC[waterPhaseIdx] - pC[gasPhaseIdx]);
Scalar pw = pg + pcFactor_ * (pC[waterPhaseIdx] - pC[gasPhaseIdx]);
Scalar rswSat = FluidSystem::waterPvt().saturatedGasDissolutionFactor(pvtRegionIdx_,
T,
pw,
@ -671,7 +682,7 @@ public:
const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
Scalar so = 1.0 - sg - solventSaturation_();
computeCapillaryPressures_(pC, so, sg + solventSaturation_(), /*sw=*/ 0.0, matParams);
p += (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
p += pcFactor_ * (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
}
Scalar rvwSat = FluidSystem::gasPvt().saturatedWaterVaporizationFactor(pvtRegionIdx_,
T,
@ -706,7 +717,7 @@ public:
std::array<Scalar, numPhases> pC = { 0.0 };
const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
computeCapillaryPressures_(pC, /*so=*/ 0.0, /*sg=*/ 0.0, /*sw=*/ 1.0, matParams);
Scalar pg = pw + (pC[gasPhaseIdx] - pC[waterPhaseIdx]);
Scalar pg = pw + pcFactor_ * (pC[gasPhaseIdx] - pC[waterPhaseIdx]);
(*this)[Indices::pressureSwitchIdx] = pg;
changed = true;
}
@ -758,7 +769,7 @@ public:
std::array<Scalar, numPhases> pC = { 0.0 };
const MaterialLawParams& matParams = problem.materialLawParams(globalDofIdx);
computeCapillaryPressures_(pC, /*so=*/0.0, sg + solventSaturation_(), sw, matParams);
Scalar pg = po + (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
Scalar pg = po + pcFactor_ * (pC[gasPhaseIdx] - pC[oilPhaseIdx]);
// we start at the GasMeaning::Rv value that corresponds to that of oil-saturated
// hydrocarbon gas
@ -838,7 +849,7 @@ public:
/*sg=*/sg2 + solventSaturation_(),
sw,
matParams);
Scalar po = pg + (pC[oilPhaseIdx] - pC[gasPhaseIdx]);
Scalar po = pg + pcFactor_ * (pC[oilPhaseIdx] - pC[gasPhaseIdx]);
setPrimaryVarsMeaningGas(GasMeaning::Sg);
setPrimaryVarsMeaningPressure(PressureMeaning::Po);
@ -1084,6 +1095,7 @@ private:
BrineMeaning primaryVarsMeaningBrine_;
SolventMeaning primaryVarsMeaningSolvent_;
unsigned short pvtRegionIdx_;
Scalar pcFactor_;
};

View File

@ -102,6 +102,9 @@ public:
unsigned,
unsigned) const
{ return 0; }
Scalar satnumRegionIndex(unsigned) const
{ return 0; }
/*!
* \brief Returns the index of the relevant region for solvent mixing functions