Merge pull request #5166 from akva2/indices_drop_ebos

Drop ebos in index functions
This commit is contained in:
Kai Bao 2024-03-05 23:31:27 +01:00 committed by GitHub
commit cd9a208757
No known key found for this signature in database
GPG Key ID: B5690EEEBB952194
10 changed files with 38 additions and 38 deletions

View File

@ -523,7 +523,7 @@ volumeFractionScaled(const int seg,
// For reservoir rate control, the distr in well control is used for the
// rate conversion coefficients. For the injection well, only the distr of the injection
// phase is not zero.
const double scale = well_.scalingFactor(well_.ebosCompIdxToFlowCompIdx(comp_idx));
const double scale = well_.scalingFactor(well_.modelCompIdxToFlowCompIdx(comp_idx));
if (scale > 0.) {
return this->volumeFraction(seg, comp_idx) / scale;
}
@ -564,17 +564,17 @@ getSegmentRateUpwinding(const int seg,
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)
&& Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx) == comp_idx
&& phase == InjectorType::WATER)
return evaluation_[seg][WQTotal] / well_.scalingFactor(well_.ebosCompIdxToFlowCompIdx(comp_idx));
return evaluation_[seg][WQTotal] / well_.scalingFactor(well_.modelCompIdxToFlowCompIdx(comp_idx));
if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)
&& Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx) == comp_idx
&& phase == InjectorType::OIL)
return evaluation_[seg][WQTotal] / well_.scalingFactor(well_.ebosCompIdxToFlowCompIdx(comp_idx));
return evaluation_[seg][WQTotal] / well_.scalingFactor(well_.modelCompIdxToFlowCompIdx(comp_idx));
if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)
&& Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx) == comp_idx
&& phase == InjectorType::GAS)
return evaluation_[seg][WQTotal] / well_.scalingFactor(well_.ebosCompIdxToFlowCompIdx(comp_idx));
return evaluation_[seg][WQTotal] / well_.scalingFactor(well_.modelCompIdxToFlowCompIdx(comp_idx));
return 0.0;
}
@ -669,6 +669,7 @@ INSTANCE(BlackOilTwoPhaseIndices<0u,0u,2u,0u,false,false,0u,2u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,false,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<0u,0u,0u,1u,false,true,0u,0u,0u>)
INSTANCE(BlackOilTwoPhaseIndices<1u,0u,0u,0u,false,false,0u,0u,0u>)
// Blackoil
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,false,0u,0u>)
INSTANCE(BlackOilIndices<1u,0u,0u,0u,false,false,0u,0u>)
@ -680,7 +681,6 @@ INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,true,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,1u,false,true,0u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,false,1u,0u>)
INSTANCE(BlackOilIndices<0u,0u,0u,0u,false,true,2u,0u>)
INSTANCE(BlackOilIndices<1u,0u,0u,0u,true,false,0u,0u>)
}

View File

@ -369,7 +369,7 @@ namespace Opm
allow_cf, cq_s, perf_press, perf_rates, deferred_logger);
for(int p = 0; p < np; ++p) {
well_flux[this->ebosCompIdxToFlowCompIdx(p)] += cq_s[p];
well_flux[this->modelCompIdxToFlowCompIdx(p)] += cq_s[p];
}
}
}
@ -442,7 +442,7 @@ namespace Opm
well_flux.resize(np, 0.0);
for (int compIdx = 0; compIdx < this->num_components_; ++compIdx) {
const EvalWell rate = well_copy.primary_variables_.getQs(compIdx);
well_flux[this->ebosCompIdxToFlowCompIdx(compIdx)] = rate.value();
well_flux[this->modelCompIdxToFlowCompIdx(compIdx)] = rate.value();
}
debug_cost_counter_ += well_copy.debug_cost_counter_;
}
@ -562,7 +562,7 @@ namespace Opm
well_potentials.resize(np, 0.0);
for (int compIdx = 0; compIdx < this->num_components_; ++compIdx) {
const EvalWell rate = well_copy.primary_variables_.getQs(compIdx);
well_potentials[this->ebosCompIdxToFlowCompIdx(compIdx)] = rate.value();
well_potentials[this->modelCompIdxToFlowCompIdx(compIdx)] = rate.value();
}
debug_cost_counter_ += well_copy.debug_cost_counter_;
return converged;
@ -1361,7 +1361,7 @@ namespace Opm
constexpr int num_eq = MSWEval::numWellEq;
for (int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx){
const EvalWell comp_rate = this->primary_variables_.getQs(comp_idx);
const int idx = this->ebosCompIdxToFlowCompIdx(comp_idx);
const int idx = this->modelCompIdxToFlowCompIdx(comp_idx);
for (size_t pvIdx = 0; pvIdx < num_eq; ++pvIdx) {
// well primary variable derivatives in EvalWell start at position Indices::numEq
ws.implicit_ipr_b[idx] -= x_well[0][pvIdx]*comp_rate.derivative(pvIdx+Indices::numEq);
@ -1879,7 +1879,7 @@ namespace Opm
// store the perf pressure and rates
for (int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx) {
perf_rates[perf*this->number_of_phases_ + this->ebosCompIdxToFlowCompIdx(comp_idx)] = cq_s[comp_idx].value();
perf_rates[perf*this->number_of_phases_ + this->modelCompIdxToFlowCompIdx(comp_idx)] = cq_s[comp_idx].value();
}
perf_press_state[perf] = perf_press.value();

View File

@ -439,7 +439,7 @@ computeProperties(const WellState& well_state,
for (int perf = 0; perf < nperf; ++perf) {
for (int comp = 0; comp < np; ++comp) {
perfRates[perf * well_.numComponents() + comp] = perf_rates_state[perf * np + well_.ebosCompIdxToFlowCompIdx(comp)];
perfRates[perf * well_.numComponents() + comp] = perf_rates_state[perf * np + well_.modelCompIdxToFlowCompIdx(comp)];
}
}
@ -479,9 +479,9 @@ computeProperties(const WellState& well_state,
double total_mobility = 0.0;
double total_invB = 0.;
for (int p = 0; p < np; ++p) {
int ebosPhaseIdx = well_.flowPhaseToEbosPhaseIdx(p);
total_mobility += invB(cell_idx, ebosPhaseIdx) * mobility(cell_idx, ebosPhaseIdx);
total_invB += invB(cell_idx, ebosPhaseIdx);
int modelPhaseIdx = well_.flowPhaseToModelPhaseIdx(p);
total_mobility += invB(cell_idx, modelPhaseIdx) * mobility(cell_idx, modelPhaseIdx);
total_invB += invB(cell_idx, modelPhaseIdx);
}
if constexpr (Indices::enableSolvent) {
total_mobility += solventInverseFormationVolumeFactor(cell_idx) * solventMobility(cell_idx);
@ -494,11 +494,11 @@ computeProperties(const WellState& well_state,
// ratios for those perforations.
constexpr double small_value = 1.e-10;
for (int p = 0; p < np; ++p) {
const int ebosPhaseIdx = well_.flowPhaseToEbosPhaseIdx(p);
const int modelPhaseIdx = well_.flowPhaseToModelPhaseIdx(p);
const auto mob_ratio = non_zero_total_mobility
? mobility(cell_idx, ebosPhaseIdx) / total_mobility
? mobility(cell_idx, modelPhaseIdx) / total_mobility
: small_value / total_invB;
perfRates[perf * well_.numComponents() + p] = well_tw_fraction * invB(cell_idx, ebosPhaseIdx) * mob_ratio;
perfRates[perf * well_.numComponents() + p] = well_tw_fraction * invB(cell_idx, modelPhaseIdx) * mob_ratio;
}
if constexpr (Indices::enableSolvent) {
const auto mob_ratio = non_zero_total_mobility

View File

@ -485,7 +485,7 @@ typename StandardWellPrimaryVariables<FluidSystem,Indices,Scalar>::EvalWell
StandardWellPrimaryVariables<FluidSystem,Indices,Scalar>::
volumeFractionScaled(const int compIdx) const
{
const int legacyCompIdx = well_.ebosCompIdxToFlowCompIdx(compIdx);
const int legacyCompIdx = well_.modelCompIdxToFlowCompIdx(compIdx);
const double scal = well_.scalingFactor(legacyCompIdx);
if (scal > 0)
return this->volumeFraction(compIdx) / scal;

View File

@ -406,7 +406,7 @@ namespace Opm
auto& perf_rate_solvent = perf_data.solvent_rates;
perf_rate_solvent[perf] = cq_s[componentIdx].value();
} else {
perf_rates[perf*np + this->ebosCompIdxToFlowCompIdx(componentIdx)] = cq_s[componentIdx].value();
perf_rates[perf*np + this->modelCompIdxToFlowCompIdx(componentIdx)] = cq_s[componentIdx].value();
}
}
@ -897,7 +897,7 @@ namespace Opm
for (int comp_idx = 0; comp_idx < this->num_components_; ++comp_idx){
EvalWell comp_rate = this->primary_variables_.getQs(comp_idx);
const int idx = this->ebosCompIdxToFlowCompIdx(comp_idx);
const int idx = this->modelCompIdxToFlowCompIdx(comp_idx);
for (size_t pvIdx = 0; pvIdx < nEq; ++pvIdx) {
// well primary variable derivatives in EvalWell start at position Indices::numEq
ws.implicit_ipr_b[idx] -= x_well[0][pvIdx]*comp_rate.derivative(pvIdx+Indices::numEq);
@ -1459,7 +1459,7 @@ namespace Opm
cq_s, perf_rates, deferred_logger);
for(int p = 0; p < np; ++p) {
well_flux[this->ebosCompIdxToFlowCompIdx(p)] += cq_s[p];
well_flux[this->modelCompIdxToFlowCompIdx(p)] += cq_s[p];
}
// the solvent contribution is added to the gas potentials
@ -1631,7 +1631,7 @@ namespace Opm
well_potentials.resize(np, 0.0);
for (int compIdx = 0; compIdx < this->num_components_; ++compIdx) {
const EvalWell rate = well_copy.primary_variables_.getQs(compIdx);
well_potentials[this->ebosCompIdxToFlowCompIdx(compIdx)] = rate.value();
well_potentials[this->modelCompIdxToFlowCompIdx(compIdx)] = rate.value();
}
return converged;
}

View File

@ -224,7 +224,7 @@ checkConstraints(WellState& well_state,
template<typename FluidSystem>
int
WellInterfaceFluidSystem<FluidSystem>::
flowPhaseToEbosPhaseIdx(const int phaseIdx) const
flowPhaseToModelPhaseIdx(const int phaseIdx) const
{
const auto& pu = this->phaseUsage();
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && pu.phase_pos[Water] == phaseIdx)

View File

@ -52,7 +52,7 @@ protected:
static constexpr int INVALIDCOMPLETION = std::numeric_limits<int>::max();
public:
int flowPhaseToEbosPhaseIdx(const int phaseIdx) const;
int flowPhaseToModelPhaseIdx(const int phaseIdx) const;
static constexpr int Water = BlackoilPhases::Aqua;
static constexpr int Oil = BlackoilPhases::Liquid;

View File

@ -59,7 +59,7 @@ WellInterfaceIndices(const Well& well,
template<class FluidSystem, class Indices, class Scalar>
int
WellInterfaceIndices<FluidSystem,Indices,Scalar>::
flowPhaseToEbosCompIdx(const int phaseIdx) const
flowPhaseToModelCompIdx(const int phaseIdx) const
{
const auto& pu = this->phaseUsage();
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && pu.phase_pos[Water] == phaseIdx)
@ -76,7 +76,7 @@ flowPhaseToEbosCompIdx(const int phaseIdx) const
template<class FluidSystem, class Indices, class Scalar>
int
WellInterfaceIndices<FluidSystem,Indices,Scalar>::
ebosCompIdxToFlowCompIdx(const unsigned compIdx) const
modelCompIdxToFlowCompIdx(const unsigned compIdx) const
{
const auto& pu = this->phaseUsage();
if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx) && Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx) == compIdx)

View File

@ -39,8 +39,8 @@ public:
using WellInterfaceFluidSystem<FluidSystem>::Water;
using Eval = DenseAd::Evaluation<Scalar, /*size=*/Indices::numEq>;
int flowPhaseToEbosCompIdx(const int phaseIdx) const;
int ebosCompIdxToFlowCompIdx(const unsigned compIdx) const;
int flowPhaseToModelCompIdx(const int phaseIdx) const;
int modelCompIdxToFlowCompIdx(const unsigned compIdx) const;
double scalingFactor(const int phaseIdx) const;
template <class EvalWell>

View File

@ -1429,12 +1429,12 @@ namespace Opm
const double well_tw_fraction = this->well_index_[perf] / total_tw;
double total_mobility = 0.0;
for (int p = 0; p < np; ++p) {
int ebosPhaseIdx = this->flowPhaseToEbosPhaseIdx(p);
total_mobility += fs.invB(ebosPhaseIdx).value() * intQuants.mobility(ebosPhaseIdx).value();
int modelPhaseIdx = this->flowPhaseToModelPhaseIdx(p);
total_mobility += fs.invB(modelPhaseIdx).value() * intQuants.mobility(modelPhaseIdx).value();
}
for (int p = 0; p < np; ++p) {
int ebosPhaseIdx = this->flowPhaseToEbosPhaseIdx(p);
scaling_factor[p] += well_tw_fraction * fs.invB(ebosPhaseIdx).value() * intQuants.mobility(ebosPhaseIdx).value() / total_mobility;
int modelPhaseIdx = this->flowPhaseToModelPhaseIdx(p);
scaling_factor[p] += well_tw_fraction * fs.invB(modelPhaseIdx).value() * intQuants.mobility(modelPhaseIdx).value() / total_mobility;
}
}
return scaling_factor;
@ -1472,18 +1472,18 @@ namespace Opm
// No nonzero rates.
// Use the computed rate directly
for (int p = 0; p < this->number_of_phases_; ++p) {
ws.surface_rates[p] = well_q_s[this->flowPhaseToEbosCompIdx(p)];
ws.surface_rates[p] = well_q_s[this->flowPhaseToModelCompIdx(p)];
}
return;
}
// Set the currently-zero phase flows to be nonzero in proportion to well_q_s.
const double initial_nonzero_rate = ws.surface_rates[nonzero_rate_index];
const int comp_idx_nz = this->flowPhaseToEbosCompIdx(nonzero_rate_index);
const int comp_idx_nz = this->flowPhaseToModelCompIdx(nonzero_rate_index);
if (std::abs(well_q_s[comp_idx_nz]) > floating_point_error_epsilon) {
for (int p = 0; p < this->number_of_phases_; ++p) {
if (p != nonzero_rate_index) {
const int comp_idx = this->flowPhaseToEbosCompIdx(p);
const int comp_idx = this->flowPhaseToModelCompIdx(p);
double& rate = ws.surface_rates[p];
rate = (initial_nonzero_rate / well_q_s[comp_idx_nz]) * (well_q_s[comp_idx]);
}
@ -1792,8 +1792,8 @@ namespace Opm
// Note: E100's notion of PI value phase mobility includes
// the reciprocal FVF.
const auto connMob =
mobility[this->flowPhaseToEbosCompIdx(p)]
* fs.invB(this->flowPhaseToEbosPhaseIdx(p)).value();
mobility[this->flowPhaseToModelCompIdx(p)]
* fs.invB(this->flowPhaseToModelPhaseIdx(p)).value();
connPI[p] = connPICalc(connMob);
}
@ -1845,7 +1845,7 @@ namespace Opm
}
const auto mt = std::accumulate(mobility.begin(), mobility.end(), 0.0);
connII[phase_pos] = connIICalc(mt * fs.invB(this->flowPhaseToEbosPhaseIdx(phase_pos)).value());
connII[phase_pos] = connIICalc(mt * fs.invB(this->flowPhaseToModelPhaseIdx(phase_pos)).value());
}
} // namespace Opm