Add support for miscibilty in the relperms

This commit is contained in:
Tor Harald Sandve 2015-12-02 14:38:49 +01:00
parent bbee43ef7a
commit 238e7c19f3
6 changed files with 88 additions and 8 deletions

View File

@ -63,6 +63,7 @@ namespace Opm {
/// \param[in] has_vapoil turn on vaporized oil feature
/// \param[in] terminal_output request output to cout/cerr
/// \param[in] has_solvent turn on solvent feature
/// \param[in] is_miscible turn on miscible feature
BlackoilSolventModel(const typename Base::ModelParameters& param,
const Grid& grid,
const BlackoilPropsAdInterface& fluid,
@ -75,7 +76,8 @@ namespace Opm {
const bool has_disgas,
const bool has_vapoil,
const bool terminal_output,
const bool has_solvent);
const bool has_solvent,
const bool is_miscible);
/// Apply an update to the primary variables, chopped if appropriate.
/// \param[in] dx updates to apply to primary variables
@ -106,6 +108,8 @@ namespace Opm {
const bool has_solvent_;
const int solvent_pos_;
const SolventPropsAdFromDeck& solvent_props_;
const bool is_miscible_;
// Need to declare Base members we want to use here.
using Base::grid_;

View File

@ -82,12 +82,15 @@ namespace Opm {
const bool has_disgas,
const bool has_vapoil,
const bool terminal_output,
const bool has_solvent)
const bool has_solvent,
const bool is_miscible)
: Base(param, grid, fluid, geo, rock_comp_props, wells_arg, linsolver,
eclState, has_disgas, has_vapoil, terminal_output),
has_solvent_(has_solvent),
solvent_pos_(detail::solventPos(fluid.phaseUsage())),
solvent_props_(solvent_props)
solvent_props_(solvent_props),
is_miscible_(is_miscible)
{
if (has_solvent_) {
@ -523,25 +526,45 @@ namespace Opm {
const ADB& sg = (active_[ Gas ]
? state.saturation[ pu.phase_pos[ Gas ] ]
: zero);
const ADB& so = (active_[ Oil ]
? state.saturation[ pu.phase_pos[ Oil ] ]
: zero);
const ADB sn = ss + so + sg;
Selector<double> zeroSn_selector(sn.value(), Selector<double>::Zero);
const ADB F_totalGas = zeroSn_selector.select( (ss + sg), (ss + sg) / sn);
Selector<double> zero_selector(ss.value() + sg.value(), Selector<double>::Zero);
ADB F_solvent = zero_selector.select(ss, ss / (ss + sg));
V ones = V::Constant(nc, 1.0);
const int canonicalPhaseIdx = canph_[ actph ];
switch (canonicalPhaseIdx) {
case Gas: {
// gas relperm
const ADB krg = solvent_props_.gasRelPermMultiplier( (ones - F_solvent) , cells_) * kr;
ADB krg = solvent_props_.gasRelPermMultiplier( (ones - F_solvent) , cells_) * kr;
if (is_miscible_) {
const ADB misc = solvent_props_.miscibilityFunction(F_solvent, cells_);
const ADB mkrgt = solvent_props_.miscibleSolventGasRelPermMultiplier(F_totalGas, cells_) * solvent_props_.misicibleHydrocarbonWaterRelPerm(sn, cells_);
krg = misc * mkrgt + (ones - misc) * krg;
}
// compute gas mobility and flux
Base::computeMassFlux(actph, transi, krg, phasePressure, state);
// compute solvent mobility and flux
const ADB tr_mult = transMult(state.pressure);
const ADB mu = solvent_props_.muSolvent(phasePressure,cells_);
rq_[solvent_pos_].mob = solvent_props_.solventRelPermMultiplier(F_solvent, cells_) * tr_mult * kr / mu;
ADB krs = solvent_props_.solventRelPermMultiplier(F_solvent, cells_) * kr;
if (is_miscible_) {
const ADB misc = solvent_props_.miscibilityFunction(F_solvent, cells_);
const ADB mkrgt = solvent_props_.miscibleSolventGasRelPermMultiplier(F_totalGas, cells_) * solvent_props_.misicibleHydrocarbonWaterRelPerm(sn, cells_);
krs = misc * mkrgt + (ones - misc) * krs;
}
rq_[solvent_pos_].mob = krs * tr_mult / mu;
const ADB rho_solvent = solvent_props_.solventSurfaceDensity(cells_) * rq_[solvent_pos_].b;
const ADB rhoavg_solvent = ops_.caver * rho_solvent;
@ -554,7 +577,15 @@ namespace Opm {
}
case Oil: {
ADB kro = kr;
if (is_miscible_) {
const ADB misc = solvent_props_.miscibilityFunction(F_solvent, cells_);
const ADB mkro = solvent_props_.miscibleOilRelPermMultiplier( (ones - F_totalGas), cells_) * solvent_props_.misicibleHydrocarbonWaterRelPerm(sn, cells_);
kro = misc * mkro + (ones - misc) * kr;
}
Base::computeMassFlux(actph, transi, kr, phasePressure, state);
break;
}
case Water: {

View File

@ -131,6 +131,7 @@ namespace Opm
bool has_solvent_;
DeckConstPtr deck_;
SolventPropsAdFromDeck solvent_props_;
bool is_miscible_;
};

View File

@ -54,9 +54,10 @@ namespace Opm
, has_solvent_(has_solvent)
, deck_(deck)
, solvent_props_(solvent_props)
, is_miscible_(false)
{
if(deck->hasKeyword("MISCIBLE")) {
std::cerr << "MISICIBLE keyword is present. Mixing is not currently supported" << std::endl;
is_miscible_ = true;
}
}
@ -80,7 +81,8 @@ namespace Opm
BaseType::has_disgas_,
BaseType::has_vapoil_,
BaseType::terminal_output_,
has_solvent_));
has_solvent_,
is_miscible_));
if (!BaseType::threshold_pressures_by_face_.empty()) {
model->setThresholdPressures(BaseType::threshold_pressures_by_face_);

View File

@ -177,6 +177,34 @@ SolventPropsAdFromDeck::SolventPropsAdFromDeck(DeckConstPtr deck,
} else {
OPM_THROW(std::runtime_error, "MSFN must be specified in MISCIBLE (SOLVENT) runs\n");
}
const TableContainer& miscTables = tables->getMiscTables();
if (!miscTables.empty()) {
int numRegions = miscTables.size();
if(numRegions > 1) {
OPM_THROW(std::runtime_error, "Only single table miscibility function supported for MISC");
}
// resize the attributes of the object
misc_.resize(numRegions);
for (int regionIdx = 0; regionIdx < numRegions; ++regionIdx) {
const Opm::MiscTable& miscTable = miscTables.getTable<MiscTable>(regionIdx);
// Copy data
// solventFraction = Ss / (Ss + Sg);
const std::vector<double>& solventFraction = miscTable.getSolventFractionColumn();
const std::vector<double>& misc = miscTable.getMiscibilityColumn();
misc_[regionIdx] = NonuniformTableLinear<double>(solventFraction, misc);
}
} else {
OPM_THROW(std::runtime_error, "MSFN must be specified in MISCIBLE (SOLVENT) runs\n");
}
}
}
@ -247,6 +275,12 @@ ADB SolventPropsAdFromDeck::miscibleOilRelPermMultiplier(const ADB& So,
return SolventPropsAdFromDeck::makeAD(So, cells, mkro_);
}
ADB SolventPropsAdFromDeck::miscibilityFunction(const ADB& solventFraction,
const Cells& cells) const
{
return SolventPropsAdFromDeck::makeAD(solventFraction, cells, misc_);
}
ADB SolventPropsAdFromDeck::makeAD(const ADB& X_AD, const Cells& cells, std::vector<NonuniformTableLinear<double>> table) const {
const int n = cells.size();
assert(Sn.value().size() == n);

View File

@ -101,6 +101,13 @@ public:
ADB miscibleOilRelPermMultiplier (const ADB& So,
const Cells& cells) const;
/// Miscible function
/// \param[in] solventFraction Array of n solvent fraction Ss / (Sg + Ss) values.
/// \param[in] cells Array of n cell indices to be associated with the fraction values.
/// \return Array of n miscibility values
ADB miscibilityFunction (const ADB& solventFraction,
const Cells& cells) const;
/// Solvent surface density
/// \param[in] cells Array of n cell indices to be associated with the fraction values.
/// \return Array of n solvent density values.
@ -121,6 +128,7 @@ private:
std::vector<NonuniformTableLinear<double> > krn_;
std::vector<NonuniformTableLinear<double> > mkro_;
std::vector<NonuniformTableLinear<double> > mkrsg_;
std::vector<NonuniformTableLinear<double> > misc_;
};
} // namespace OPM