changed: unify StandardWell::getMobility(Eval|Scalar)

This commit is contained in:
Arne Morten Kvarving
2023-05-07 21:56:09 +02:00
parent 6021388045
commit 35c56e4ce4
2 changed files with 37 additions and 102 deletions

View File

@@ -333,17 +333,11 @@ namespace Opm
virtual double getRefDensity() const override;
// get the mobility for specific perforation
void getMobilityEval(const Simulator& ebosSimulator,
const int perf,
std::vector<EvalWell>& mob,
DeferredLogger& deferred_logger) const;
// get the mobility for specific perforation
void getMobilityScalar(const Simulator& ebosSimulator,
const int perf,
std::vector<Scalar>& mob,
DeferredLogger& deferred_logger) const;
template<class Value>
void getMobility(const Simulator& ebosSimulator,
const int perf,
std::vector<Value>& mob,
DeferredLogger& deferred_logger) const;
void updateWaterMobilityWithPolymer(const Simulator& ebos_simulator,
const int perf,

View File

@@ -618,7 +618,7 @@ namespace Opm
const int cell_idx = this->well_cells_[perf];
const auto& intQuants = ebosSimulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/ 0);
std::vector<EvalWell> mob(this->num_components_, {this->primary_variables_.numWellEq() + Indices::numEq, 0.});
getMobilityEval(ebosSimulator, perf, mob, deferred_logger);
getMobility(ebosSimulator, perf, mob, deferred_logger);
PerforationRates perf_rates;
double trans_mult = ebosSimulator.problem().template rockCompTransMultiplier<double>(intQuants, cell_idx);
@@ -820,15 +820,23 @@ namespace Opm
template<typename TypeTag>
template<class Value>
void
StandardWell<TypeTag>::
getMobilityEval(const Simulator& ebosSimulator,
const int perf,
std::vector<EvalWell>& mob,
DeferredLogger& deferred_logger) const
getMobility(const Simulator& ebosSimulator,
const int perf,
std::vector<Value>& mob,
DeferredLogger& deferred_logger) const
{
auto obtain = [this](const Eval& value)
{
if constexpr (std::is_same_v<Value, Scalar>) {
return getValue(value);
} else {
return this->extendEval(value);
}
};
const int cell_idx = this->well_cells_[perf];
assert (int(mob.size()) == this->num_components_);
const auto& intQuants = ebosSimulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/0);
@@ -838,21 +846,19 @@ namespace Opm
// based on passing the saturation table index
const int satid = this->saturation_table_number_[perf] - 1;
const int satid_elem = materialLawManager->satnumRegionIdx(cell_idx);
if( satid == satid_elem ) { // the same saturation number is used. i.e. just use the mobilty from the cell
if (satid == satid_elem) { // the same saturation number is used. i.e. just use the mobilty from the cell
for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
if (!FluidSystem::phaseIsActive(phaseIdx)) {
continue;
}
const unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
mob[activeCompIdx] = this->extendEval(intQuants.mobility(phaseIdx));
mob[activeCompIdx] = obtain(intQuants.mobility(phaseIdx));
}
if (has_solvent) {
mob[Indices::contiSolventEqIdx] = this->extendEval(intQuants.solventMobility());
mob[Indices::contiSolventEqIdx] = obtain(intQuants.solventMobility());
}
} else {
const auto& paramsCell = materialLawManager->connectionMaterialLawParams(satid, cell_idx);
std::array<Eval,3> relativePerms = { 0.0, 0.0, 0.0 };
MaterialLaw::relativePermeabilities(relativePerms, paramsCell, intQuants.fluidState());
@@ -867,7 +873,7 @@ namespace Opm
}
const unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
mob[activeCompIdx] = this->extendEval(relativePerms[phaseIdx] / intQuants.fluidState().viscosity(phaseIdx));
mob[activeCompIdx] = obtain(relativePerms[phaseIdx] / intQuants.fluidState().viscosity(phaseIdx));
}
// this may not work if viscosity and relperms has been modified?
@@ -885,81 +891,16 @@ namespace Opm
// for the cases related to polymer molecular weight, we assume fully mixing
// as a result, the polymer and water share the same viscosity
if constexpr (!Base::has_polymermw) {
updateWaterMobilityWithPolymer(ebosSimulator, perf, mob, deferred_logger);
}
}
}
template<typename TypeTag>
void
StandardWell<TypeTag>::
getMobilityScalar(const Simulator& ebosSimulator,
const int perf,
std::vector<Scalar>& mob,
DeferredLogger& deferred_logger) const
{
const int cell_idx = this->well_cells_[perf];
assert (int(mob.size()) == this->num_components_);
const auto& intQuants = ebosSimulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/0);
const auto& materialLawManager = ebosSimulator.problem().materialLawManager();
// either use mobility of the perforation cell or calcualte its own
// based on passing the saturation table index
const int satid = this->saturation_table_number_[perf] - 1;
const int satid_elem = materialLawManager->satnumRegionIdx(cell_idx);
if( satid == satid_elem ) { // the same saturation number is used. i.e. just use the mobilty from the cell
for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
if (!FluidSystem::phaseIsActive(phaseIdx)) {
continue;
}
const unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
mob[activeCompIdx] = getValue(intQuants.mobility(phaseIdx));
}
if (has_solvent) {
mob[Indices::contiSolventEqIdx] = getValue(intQuants.solventMobility());
}
} else {
const auto& paramsCell = materialLawManager->connectionMaterialLawParams(satid, cell_idx);
std::array<Eval,3> relativePerms = { 0.0, 0.0, 0.0 };
MaterialLaw::relativePermeabilities(relativePerms, paramsCell, intQuants.fluidState());
// reset the satnumvalue back to original
materialLawManager->connectionMaterialLawParams(satid_elem, cell_idx);
// compute the mobility
for (unsigned phaseIdx = 0; phaseIdx < FluidSystem::numPhases; ++phaseIdx) {
if (!FluidSystem::phaseIsActive(phaseIdx)) {
continue;
}
const unsigned activeCompIdx = Indices::canonicalToActiveComponentIndex(FluidSystem::solventComponentIndex(phaseIdx));
mob[activeCompIdx] = getValue(relativePerms[phaseIdx]) / getValue(intQuants.fluidState().viscosity(phaseIdx));
}
// this may not work if viscosity and relperms has been modified?
if constexpr (has_solvent) {
OPM_DEFLOG_THROW(std::runtime_error, "individual mobility for wells does not work in combination with solvent", deferred_logger);
}
}
// modify the water mobility if polymer is present
if constexpr (has_polymer) {
if (!FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) {
OPM_DEFLOG_THROW(std::runtime_error, "Water is required when polymer is active", deferred_logger);
}
// for the cases related to polymer molecular weight, we assume fully mixing
// as a result, the polymer and water share the same viscosity
if constexpr (!Base::has_polymermw) {
std::vector<EvalWell> mob_eval(this->num_components_, {this->primary_variables_.numWellEq() + Indices::numEq, 0.});
for (size_t i = 0; i < mob.size(); ++i)
mob_eval[i].setValue(mob[i]);
updateWaterMobilityWithPolymer(ebosSimulator, perf, mob_eval, deferred_logger);
for (size_t i = 0; i < mob.size(); ++i) {
mob[i] = getValue(mob_eval[i]);
if constexpr (std::is_same_v<Value, Scalar>) {
std::vector<EvalWell> mob_eval(this->num_components_, {this->primary_variables_.numWellEq() + Indices::numEq, 0.});
for (size_t i = 0; i < mob.size(); ++i)
mob_eval[i].setValue(mob[i]);
updateWaterMobilityWithPolymer(ebosSimulator, perf, mob_eval, deferred_logger);
for (size_t i = 0; i < mob.size(); ++i) {
mob[i] = getValue(mob_eval[i]);
}
} else {
updateWaterMobilityWithPolymer(ebosSimulator, perf, mob, deferred_logger);
}
}
}
@@ -1043,7 +984,7 @@ namespace Opm
for (int perf = 0; perf < this->number_of_perforations_; ++perf) {
std::vector<Scalar> mob(this->num_components_, 0.0);
getMobilityScalar(ebos_simulator, perf, mob, deferred_logger);
getMobility(ebos_simulator, perf, mob, deferred_logger);
const int cell_idx = this->well_cells_[perf];
const auto& int_quantities = ebos_simulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/ 0);
@@ -1450,7 +1391,7 @@ namespace Opm
};
std::vector<EvalWell> mob(this->num_components_, {this->primary_variables_.numWellEq() + Indices::numEq, 0.0});
getMobilityEval(ebosSimulator, static_cast<int>(subsetPerfID), mob, deferred_logger);
getMobility(ebosSimulator, static_cast<int>(subsetPerfID), mob, deferred_logger);
const auto& fs = fluidState(subsetPerfID);
setToZero(connPI);
@@ -1652,7 +1593,7 @@ namespace Opm
const auto& intQuants = ebosSimulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/ 0);
// flux for each perforation
std::vector<Scalar> mob(this->num_components_, 0.);
getMobilityScalar(ebosSimulator, perf, mob, deferred_logger);
getMobility(ebosSimulator, perf, mob, deferred_logger);
double trans_mult = ebosSimulator.problem().template rockCompTransMultiplier<double>(intQuants, cell_idx);
const double Tw = this->well_index_[perf] * trans_mult;
@@ -2509,7 +2450,7 @@ namespace Opm
const int cell_idx = this->well_cells_[perf];
const auto& intQuants = ebosSimulator.model().intensiveQuantities(cell_idx, /*timeIdx=*/ 0);
std::vector<Scalar> mob(this->num_components_, 0.);
getMobilityScalar(ebosSimulator, perf, mob, deferred_logger);
getMobility(ebosSimulator, perf, mob, deferred_logger);
std::vector<Scalar> cq_s(this->num_components_, 0.);
double trans_mult = ebosSimulator.problem().template rockCompTransMultiplier<double>(intQuants, cell_idx);
const double Tw = this->well_index_[perf] * trans_mult;