diff --git a/ebos/eclmpiserializer.hh b/ebos/eclmpiserializer.hh index 4c137d521..c3d08cb0b 100644 --- a/ebos/eclmpiserializer.hh +++ b/ebos/eclmpiserializer.hh @@ -260,16 +260,24 @@ public: (*this)(d); }; + auto keyHandle = [&](auto& d) + { + if constexpr (is_pair::value) + pair(d); + else + (*this)(d); + }; + if (m_op == Operation::PACKSIZE) { m_packSize += Mpi::packSize(data.size(), m_comm); for (auto& it : data) { - m_packSize += Mpi::packSize(it.first, m_comm); + keyHandle(it.first); handle(it.second); } } else if (m_op == Operation::PACK) { Mpi::pack(data.size(), m_buffer, m_position, m_comm); for (auto& it : data) { - Mpi::pack(it.first, m_buffer, m_position, m_comm); + keyHandle(it.first); handle(it.second); } } else if (m_op == Operation::UNPACK) { @@ -277,7 +285,7 @@ public: Mpi::unpack(size, m_buffer, m_position, m_comm); for (size_t i = 0; i < size; ++i) { Key key; - Mpi::unpack(key, m_buffer, m_position, m_comm); + keyHandle(key); Data entry; handle(entry); data.insert(std::make_pair(key, entry)); diff --git a/opm/simulators/wells/BlackoilWellModel.hpp b/opm/simulators/wells/BlackoilWellModel.hpp index 865a36641..f168ab290 100644 --- a/opm/simulators/wells/BlackoilWellModel.hpp +++ b/opm/simulators/wells/BlackoilWellModel.hpp @@ -391,7 +391,7 @@ namespace Opm { void updateWellControls(Opm::DeferredLogger& deferred_logger, const bool checkGroupControls); - void updateAndCommunicateGroupData(); + void updateAndCommunicateGroupData(Opm::DeferredLogger& deferred_logger); void updateNetworkPressures(); // setting the well_solutions_ based on well_state. @@ -502,6 +502,7 @@ namespace Opm { void assignGroupControl(const Group& group, data::GroupData& gdata) const; data::GuideRateValue getGuideRateValues(const Well& well) const; data::GuideRateValue getGuideRateValues(const Group& group) const; + data::GuideRateValue getGuideRateInjectionGroupValues(const Group& group) const; void getGuideRateValues(const GuideRate::RateVector& qs, const bool is_inj, const std::string& wgname, diff --git a/opm/simulators/wells/BlackoilWellModel_impl.hpp b/opm/simulators/wells/BlackoilWellModel_impl.hpp index ebd3b22d2..0e2124d07 100644 --- a/opm/simulators/wells/BlackoilWellModel_impl.hpp +++ b/opm/simulators/wells/BlackoilWellModel_impl.hpp @@ -306,7 +306,7 @@ namespace Opm { well_state_ = previous_well_state_; well_state_.disableGliftOptimization(); - updateAndCommunicateGroupData(); + updateAndCommunicateGroupData(local_deferredLogger); const int reportStepIdx = ebosSimulator_.episodeIndex(); const double simulationTime = ebosSimulator_.time(); int exception_thrown = 0; @@ -1227,7 +1227,7 @@ namespace Opm { // For no well active globally we simply return. if( !wellsActive() ) return ; - updateAndCommunicateGroupData(); + updateAndCommunicateGroupData(deferred_logger); updateNetworkPressures(); @@ -1241,7 +1241,7 @@ namespace Opm { // Check group's constraints from higher levels. updateGroupHigherControls(deferred_logger, switched_groups); - updateAndCommunicateGroupData(); + updateAndCommunicateGroupData(deferred_logger); // Check wells' group constraints and communicate. for (const auto& well : well_container_) { @@ -1251,7 +1251,7 @@ namespace Opm { switched_wells.insert(well->name()); } } - updateAndCommunicateGroupData(); + updateAndCommunicateGroupData(deferred_logger); } // Check individual well constraints and communicate. @@ -1262,7 +1262,7 @@ namespace Opm { const auto mode = WellInterface::IndividualOrGroup::Individual; well->updateWellControl(ebosSimulator_, mode, well_state_, deferred_logger); } - updateAndCommunicateGroupData(); + updateAndCommunicateGroupData(deferred_logger); } @@ -1305,7 +1305,7 @@ namespace Opm { template void BlackoilWellModel:: - updateAndCommunicateGroupData() + updateAndCommunicateGroupData(Opm::DeferredLogger& deferred_logger) { const int reportStepIdx = ebosSimulator_.episodeIndex(); const Group& fieldGroup = schedule().getGroup("FIELD", reportStepIdx); @@ -1321,6 +1321,7 @@ namespace Opm { well_state_nupcol_ = well_state_; } + const auto& summaryState = ebosSimulator_.vanguard().summaryState(); // the group target reduction rates needs to be update since wells may have swicthed to/from GRUP control // Currently the group target reduction does not honor NUPCOL. TODO: is that true? std::vector groupTargetReduction(numPhases(), 0.0); @@ -1328,13 +1329,13 @@ namespace Opm { std::vector groupTargetReductionInj(numPhases(), 0.0); WellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ true, phase_usage_, *guideRate_, well_state_nupcol_, well_state_, groupTargetReductionInj); + // the guiderate update should not be part of this function + // remove in seperate PR since it affects existing functionality const double simulationTime = ebosSimulator_.time(); std::vector pot(numPhases(), 0.0); - WellGroupHelpers::updateGuideRateForGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, /*isInjector*/ false, well_state_, comm, guideRate_.get(), pot); - std::vector potInj(numPhases(), 0.0); - WellGroupHelpers::updateGuideRateForGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, /*isInjector*/ true, well_state_, comm, guideRate_.get(), potInj); + WellGroupHelpers::updateGuideRateForProductionGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, well_state_, comm, guideRate_.get(), pot); + WellGroupHelpers::updateGuideRatesForInjectionGroups(fieldGroup, schedule(), summaryState, phase_usage_, reportStepIdx, well_state_, guideRate_.get(), deferred_logger); - const auto& summaryState = ebosSimulator_.vanguard().summaryState(); WellGroupHelpers::updateREINForGroups(fieldGroup, schedule(), reportStepIdx, phase_usage_, summaryState, well_state_nupcol_, well_state_); WellGroupHelpers::updateVREPForGroups(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_); @@ -2873,6 +2874,11 @@ namespace Opm { gr[gname].production = this->getGuideRateValues(group); } + if (this->guideRate_->has(gname, Opm::Phase::WATER) + || this->guideRate_->has(gname, Opm::Phase::GAS)) { + gr[gname].injection = this->getGuideRateInjectionGroupValues(group); + } + const auto parent = group.parent(); if (parent == "FIELD") { continue; } @@ -2951,12 +2957,32 @@ namespace Opm { } const auto qs = WellGroupHelpers:: - getRateVector(this->well_state_, this->phase_usage_, wname); + getWellRateVector(this->well_state_, this->phase_usage_, wname); this->getGuideRateValues(qs, well.isInjector(), wname, grval); return grval; } + template + data::GuideRateValue + BlackoilWellModel:: + getGuideRateInjectionGroupValues(const Group& group) const + { + auto grval = data::GuideRateValue{}; + + assert (this->guideRate_ != nullptr); + + const auto& gname = group.name(); + if (this->guideRate_->has(gname, Opm::Phase::GAS)) { + grval.set(data::GuideRateValue::Item::Gas, + this->guideRate_->get(gname, Opm::Phase::GAS)); + } + if (this->guideRate_->has(gname, Opm::Phase::WATER)) { + grval.set(data::GuideRateValue::Item::Water, + this->guideRate_->get(gname, Opm::Phase::WATER)); + } + return grval; + } template data::GuideRateValue @@ -2968,8 +2994,9 @@ namespace Opm { assert (this->guideRate_ != nullptr); const auto& gname = group.name(); - if (! this->well_state_.hasProductionGroupRates(gname)) { - // No flow rates for 'gname' -- might be before group comes + + if ( ! this->well_state_.hasProductionGroupRates(gname)) { + // No flow rates for production group 'gname' -- might be before group comes // online (e.g., for the initial condition before simulation // starts). return grval; @@ -2981,7 +3008,7 @@ namespace Opm { } const auto qs = WellGroupHelpers:: - getProductionGroupRateVector(this->well_state_, this->phase_usage_, gname); + getProductionGroupRateVector(this->well_state_, this->phase_usage_, gname); const auto is_inj = false; // This procedure only applies to G*PGR. this->getGuideRateValues(qs, is_inj, gname, grval); diff --git a/opm/simulators/wells/TargetCalculator.hpp b/opm/simulators/wells/TargetCalculator.hpp index 903e2aba5..2d13f17ba 100644 --- a/opm/simulators/wells/TargetCalculator.hpp +++ b/opm/simulators/wells/TargetCalculator.hpp @@ -163,6 +163,130 @@ namespace WellGroupHelpers const double group_grat_target_from_sales_; }; + + /// Based on a group control mode, extract or calculate rates, and + /// provide other conveniences. + class InjectionTargetCalculator + { + public: + InjectionTargetCalculator(const Group::InjectionCMode& cmode, + const PhaseUsage& pu, + const std::vector& resv_coeff, + const std::string& group_name, + const double sales_target, + const WellStateFullyImplicitBlackoil& well_state, + const Phase& injection_phase, + DeferredLogger& deferred_logger) + : cmode_(cmode) + , pu_(pu) + , resv_coeff_(resv_coeff) + , group_name_(group_name) + , sales_target_(sales_target) + , well_state_(well_state) + { + // initialize to avoid warning + pos_ = pu.phase_pos[BlackoilPhases::Aqua]; + target_ = GuideRateModel::Target::WAT; + + switch (injection_phase) { + case Phase::WATER: { + pos_ = pu.phase_pos[BlackoilPhases::Aqua]; + target_ = GuideRateModel::Target::WAT; + break; + } + case Phase::OIL: { + pos_ = pu.phase_pos[BlackoilPhases::Liquid]; + target_ = GuideRateModel::Target::OIL; + break; + } + case Phase::GAS: { + pos_ = pu.phase_pos[BlackoilPhases::Vapour]; + target_ = GuideRateModel::Target::GAS; + break; + } + default: + OPM_DEFLOG_THROW(std::logic_error, + "Invalid injection phase in InjectionTargetCalculator", + deferred_logger); + } + } + + + template + auto calcModeRateFromRates(const RateVec& rates) const + { + return rates[pos_]; + } + + double groupTarget(const Group::InjectionControls& ctrl, Opm::DeferredLogger& deferred_logger) const + { + switch (cmode_) { + case Group::InjectionCMode::RATE: + return ctrl.surface_max_rate; + case Group::InjectionCMode::RESV: + return ctrl.resv_max_rate; + case Group::InjectionCMode::REIN: { + double production_rate = well_state_.currentInjectionREINRates(ctrl.reinj_group)[pos_]; + return ctrl.target_reinj_fraction * production_rate; + } + case Group::InjectionCMode::VREP: { + const std::vector& group_injection_reductions + = well_state_.currentInjectionGroupReductionRates(group_name_); + double voidage_rate + = well_state_.currentInjectionVREPRates(ctrl.voidage_group) * ctrl.target_void_fraction; + double inj_reduction = 0.0; + if (ctrl.phase != Phase::WATER) + inj_reduction += group_injection_reductions[pu_.phase_pos[BlackoilPhases::Aqua]] + * resv_coeff_[pu_.phase_pos[BlackoilPhases::Aqua]]; + if (ctrl.phase != Phase::OIL) + inj_reduction += group_injection_reductions[pu_.phase_pos[BlackoilPhases::Liquid]] + * resv_coeff_[pu_.phase_pos[BlackoilPhases::Liquid]]; + if (ctrl.phase != Phase::GAS) + inj_reduction += group_injection_reductions[pu_.phase_pos[BlackoilPhases::Vapour]] + * resv_coeff_[pu_.phase_pos[BlackoilPhases::Vapour]]; + voidage_rate -= inj_reduction; + return voidage_rate / resv_coeff_[pos_]; + } + case Group::InjectionCMode::SALE: { + assert(pos_ == pu_.phase_pos[BlackoilPhases::Vapour]); + // Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate; + // Gas import and consumption is already included in the REIN rates + double inj_rate = well_state_.currentInjectionREINRates(group_name_)[pos_]; + inj_rate -= sales_target_; + return inj_rate; + } + default: + OPM_DEFLOG_THROW(std::logic_error, + "Invalid Group::InjectionCMode in InjectionTargetCalculator", + deferred_logger); + return 0.0; + } + } + + GuideRateModel::Target guideTargetMode() const + { + return target_; + } + + private: + template + static ElemType zero() + { + // This is for Evaluation types. + ElemType x; + x = 0.0; + return x; + } + Group::InjectionCMode cmode_; + const PhaseUsage& pu_; + const std::vector& resv_coeff_; + const std::string& group_name_; + double sales_target_; + const WellStateFullyImplicitBlackoil& well_state_; + int pos_; + GuideRateModel::Target target_; + + }; } // namespace WellGroupHelpers } // namespace Opm diff --git a/opm/simulators/wells/WellGroupHelpers.cpp b/opm/simulators/wells/WellGroupHelpers.cpp index 24ead3e59..5120a89c7 100644 --- a/opm/simulators/wells/WellGroupHelpers.cpp +++ b/opm/simulators/wells/WellGroupHelpers.cpp @@ -23,7 +23,6 @@ #include #include #include - #include #include #include @@ -227,6 +226,63 @@ namespace WellGroupHelpers return gefac * rate; } + void updateGuideRatesForInjectionGroups(const Group& group, + const Schedule& schedule, + const SummaryState& summaryState, + const Opm::PhaseUsage& pu, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellState, + GuideRate* guideRate, + Opm::DeferredLogger& deferred_logger) + { + for (const std::string& groupName : group.groups()) { + const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); + updateGuideRatesForInjectionGroups( + groupTmp, schedule, summaryState, pu, reportStepIdx, wellState, guideRate, deferred_logger); + } + const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS}; + for (Phase phase : all) { + + if(!group.hasInjectionControl(phase)) + continue; + + double guideRateValue = 0.0; + const auto& controls = group.injectionControls(phase, summaryState); + switch (controls.guide_rate_def){ + case Group::GuideRateInjTarget::RATE: + break; + case Group::GuideRateInjTarget::VOID: + { + guideRateValue = wellState.currentInjectionVREPRates(group.name()); + break; + } + case Group::GuideRateInjTarget::NETV: + { + guideRateValue = wellState.currentInjectionVREPRates(group.name()); + const std::vector& injRES = wellState.currentInjectionGroupReservoirRates(group.name()); + if (phase != Phase::OIL && pu.phase_used[BlackoilPhases::Liquid]) + guideRateValue -= injRES[pu.phase_pos[BlackoilPhases::Liquid]]; + if (phase != Phase::GAS && pu.phase_used[BlackoilPhases::Vapour]) + guideRateValue -= injRES[pu.phase_pos[BlackoilPhases::Vapour]]; + if (phase != Phase::WATER && pu.phase_used[BlackoilPhases::Aqua]) + guideRateValue -= injRES[pu.phase_pos[BlackoilPhases::Aqua]]; + break; + } + case Group::GuideRateInjTarget::RESV: + OPM_DEFLOG_THROW(std::runtime_error, "GUIDE PHASE RESV not implemented. Group " + group.name(), deferred_logger); + case Group::GuideRateInjTarget::POTN: + break; + case Group::GuideRateInjTarget::NO_GUIDE_RATE: + break; + default: + OPM_DEFLOG_THROW(std::logic_error, + "Invalid GuideRateInjTarget in updateGuideRatesForInjectionGroups", + deferred_logger); + } + guideRate->compute(group.name(), phase, reportStepIdx, guideRateValue); + } + } + void updateGroupTargetReduction(const Group& group, const Schedule& schedule, const int reportStepIdx, @@ -282,7 +338,7 @@ namespace WellGroupHelpers const bool individual_control = (currentGroupControl != Group::ProductionCMode::FLD && currentGroupControl != Group::ProductionCMode::NONE); const int num_group_controlled_wells - = groupControlledWells(schedule, wellStateNupcol, reportStepIdx, subGroupName, ""); + = groupControlledWells(schedule, wellStateNupcol, reportStepIdx, subGroupName, "", !isInjector, /*injectionPhaseNotUsed*/Phase::OIL); if (individual_control || num_group_controlled_wells == 0) { for (int phase = 0; phase < np; phase++) { groupTargetReduction[phase] @@ -350,140 +406,6 @@ namespace WellGroupHelpers } - /* - template - void updateGuideRateForGroups(const Group& group, const Schedule& schedule, const PhaseUsage& pu, const int - reportStepIdx, const double& simTime, const bool isInjector, WellStateFullyImplicitBlackoil& wellState, const - Comm& comm, GuideRate* guideRate, std::vector& pot) - { - const int np = pu.num_phases; - for (const std::string& groupName : group.groups()) { - std::vector thisPot(np, 0.0); - const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); - updateGuideRateForGroups(groupTmp, schedule, pu, reportStepIdx, simTime, isInjector, wellState, comm, - guideRate, thisPot); - - // accumulate group contribution from sub group unconditionally - if (isInjector) { - const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS}; - for (Phase phase : all) { - const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(phase, - groupName); int phasePos; if (phase == Phase::GAS && pu.phase_used[BlackoilPhases::Vapour] ) phasePos = - pu.phase_pos[BlackoilPhases::Vapour]; else if (phase == Phase::OIL && pu.phase_used[BlackoilPhases::Liquid]) - phasePos = pu.phase_pos[BlackoilPhases::Liquid]; - else if (phase == Phase::WATER && pu.phase_used[BlackoilPhases::Aqua] ) - phasePos = pu.phase_pos[BlackoilPhases::Aqua]; - else - continue; - - pot[phasePos] += thisPot[phasePos]; - } - } else { - const Group::ProductionCMode& currentGroupControl = - wellState.currentProductionGroupControl(groupName); if (currentGroupControl != Group::ProductionCMode::FLD && - currentGroupControl != Group::ProductionCMode::NONE) { continue; - } - for (int phase = 0; phase < np; phase++) { - pot[phase] += thisPot[phase]; - } - } - - } - for (const std::string& wellName : group.wells()) { - const auto& wellTmp = schedule.getWell(wellName, reportStepIdx); - - if (wellTmp.isProducer() && isInjector) - continue; - - if (wellTmp.isInjector() && !isInjector) - continue; - - if (wellTmp.getStatus() == Well::Status::SHUT) - continue; - const auto& end = wellState.wellMap().end(); - const auto& it = wellState.wellMap().find( wellName ); - if (it == end) // the well is not found - continue; - - int well_index = it->second[0]; - - if (! wellIsOwned(well_index, wellName, wellState) ) // Only sum once - { - continue; - } - - const auto wellrate_index = well_index * wellState.numPhases(); - // add contribution from wells unconditionally - for (int phase = 0; phase < np; phase++) { - pot[phase] += wellState.wellPotentials()[wellrate_index + phase]; - } - } - - double oilPot = 0.0; - if (pu.phase_used[BlackoilPhases::Liquid]) - oilPot = pot [ pu.phase_pos[BlackoilPhases::Liquid]]; - - double gasPot = 0.0; - if (pu.phase_used[BlackoilPhases::Vapour]) - gasPot = pot [ pu.phase_pos[BlackoilPhases::Vapour]]; - - double waterPot = 0.0; - if (pu.phase_used[BlackoilPhases::Aqua]) - waterPot = pot [pu.phase_pos[BlackoilPhases::Aqua]]; - - const double gefac = group.getGroupEfficiencyFactor(); - - oilPot = comm.sum(oilPot) * gefac; - gasPot = comm.sum(gasPot) * gefac; - waterPot = comm.sum(waterPot) * gefac; - - if (isInjector) { - wellState.setCurrentGroupInjectionPotentials(group.name(), pot); - } else { - guideRate->compute(group.name(), reportStepIdx, simTime, oilPot, gasPot, waterPot); - } - } - */ - - - /* - template - void updateGuideRatesForWells(const Schedule& schedule, const PhaseUsage& pu, const int reportStepIdx, const - double& simTime, const WellStateFullyImplicitBlackoil& wellState, const Comm& comm, GuideRate* guideRate) { - - const auto& end = wellState.wellMap().end(); - for (const auto& well : schedule.getWells(reportStepIdx)) { - double oilpot = 0.0; - double gaspot = 0.0; - double waterpot = 0.0; - - const auto& it = wellState.wellMap().find( well.name()); - if (it != end && wellState.wellIsOwned(it->second[0], well.name())) - { - // the well is found and owned - - const auto wpot = wellState.wellPotentials().data() + well_index*wellState.numPhases(); - if (pu.phase_used[BlackoilPhases::Liquid] > 0) - oilpot = wpot[pu.phase_pos[BlackoilPhases::Liquid]]; - - if (pu.phase_used[BlackoilPhases::Vapour] > 0) - gaspot = wpot[pu.phase_pos[BlackoilPhases::Vapour]]; - - if (pu.phase_used[BlackoilPhases::Aqua] > 0) - waterpot = wpot[pu.phase_pos[BlackoilPhases::Aqua]]; - - const double wefac = well.getEfficiencyFactor(); - } - oilpot = comm.sum(oilpot) * wefac; - gaspot = comm.sum(gaspot) * wefac; - waterpot = comm.sum(waterpot) * wefac; - guideRate->compute(well.name(), reportStepIdx, simTime, oilpot, gaspot, waterpot); - } - - } - */ - - void updateVREPForGroups(const Group& group, const Schedule& schedule, const int reportStepIdx, @@ -582,6 +504,7 @@ namespace WellGroupHelpers wellState.setCurrentProductionGroupRates(group.name(), rates); } + void updateREINForGroups(const Group& group, const Schedule& schedule, const int reportStepIdx, @@ -738,7 +661,7 @@ namespace WellGroupHelpers GuideRate::RateVector - getRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& name) + getWellRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& name) { return getGuideRateVector(well_state.currentWellRates(name), pu); } @@ -757,8 +680,12 @@ namespace WellGroupHelpers const GuideRateModel::Target target, const PhaseUsage& pu) { - if (schedule.hasWell(name, reportStepIdx) || guideRate->has(name)) { - return guideRate->get(name, target, getRateVector(wellState, pu, name)); + if (schedule.hasWell(name, reportStepIdx)) { + return guideRate->get(name, target, getWellRateVector(wellState, pu, name)); + } + + if (guideRate->has(name)) { + return guideRate->get(name, target, getProductionGroupRateVector(wellState, pu, name)); } double totalGuideRate = 0.0; @@ -786,7 +713,7 @@ namespace WellGroupHelpers if (!wellState.isProductionGrup(wellName)) continue; - totalGuideRate += guideRate->get(wellName, target, getRateVector(wellState, pu, wellName)); + totalGuideRate += guideRate->get(wellName, target, getWellRateVector(wellState, pu, wellName)); } return totalGuideRate; } @@ -802,7 +729,11 @@ namespace WellGroupHelpers const PhaseUsage& pu) { if (schedule.hasWell(name, reportStepIdx)) { - return guideRate->get(name, target, getRateVector(wellState, pu, name)); + return guideRate->get(name, target, getWellRateVector(wellState, pu, name)); + } + + if (guideRate->has(name, injectionPhase)) { + return guideRate->get(name, injectionPhase); } double totalGuideRate = 0.0; @@ -832,7 +763,7 @@ namespace WellGroupHelpers if (!wellState.isInjectionGrup(wellName)) continue; - totalGuideRate += guideRate->get(wellName, target, getRateVector(wellState, pu, wellName)); + totalGuideRate += guideRate->get(wellName, target, getWellRateVector(wellState, pu, wellName)); } return totalGuideRate; } @@ -843,21 +774,35 @@ namespace WellGroupHelpers const WellStateFullyImplicitBlackoil& well_state, const int report_step, const std::string& group_name, - const std::string& always_included_child) + const std::string& always_included_child, + const bool is_production_group, + const Phase injection_phase) { const Group& group = schedule.getGroup(group_name, report_step); int num_wells = 0; for (const std::string& child_group : group.groups()) { - const auto ctrl = well_state.currentProductionGroupControl(child_group); - const bool included = (ctrl == Group::ProductionCMode::FLD) || (ctrl == Group::ProductionCMode::NONE) - || (child_group == always_included_child); + + bool included = (child_group == always_included_child); + if (is_production_group) { + const auto ctrl = well_state.currentProductionGroupControl(child_group); + included = included || (ctrl == Group::ProductionCMode::FLD) || (ctrl == Group::ProductionCMode::NONE); + } else { + const auto ctrl = well_state.currentInjectionGroupControl(injection_phase, child_group); + included = included || (ctrl == Group::InjectionCMode::FLD) || (ctrl == Group::InjectionCMode::NONE); + } + if (included) { num_wells - += groupControlledWells(schedule, well_state, report_step, child_group, always_included_child); + += groupControlledWells(schedule, well_state, report_step, child_group, always_included_child, is_production_group, injection_phase); } } for (const std::string& child_well : group.wells()) { - const bool included = (well_state.isProductionGrup(child_well)) || (child_well == always_included_child); + bool included = (child_well == always_included_child); + if (is_production_group) { + included = included || well_state.isProductionGrup(child_well); + } else { + included = included || well_state.isInjectionGrup(child_well); + } if (included) { ++num_wells; } @@ -866,17 +811,23 @@ namespace WellGroupHelpers } FractionCalculator::FractionCalculator(const Schedule& schedule, + const SummaryState& summary_state, const WellStateFullyImplicitBlackoil& well_state, const int report_step, const GuideRate* guide_rate, const GuideRateModel::Target target, - const PhaseUsage& pu) + const PhaseUsage& pu, + const bool is_producer, + const Phase injection_phase) : schedule_(schedule) + , summary_state_(summary_state) , well_state_(well_state) , report_step_(report_step) , guide_rate_(guide_rate) , target_(target) , pu_(pu) + , is_producer_(is_producer) + , injection_phase_(injection_phase) { } double FractionCalculator::fraction(const std::string& name, @@ -912,15 +863,26 @@ namespace WellGroupHelpers { double total_guide_rate = 0.0; for (const std::string& child_group : group.groups()) { - const auto ctrl = well_state_.currentProductionGroupControl(child_group); - const bool included = (ctrl == Group::ProductionCMode::FLD) || (ctrl == Group::ProductionCMode::NONE) - || (child_group == always_included_child); + bool included = (child_group == always_included_child); + if (is_producer_) { + const auto ctrl = well_state_.currentProductionGroupControl(child_group); + included = included || (ctrl == Group::ProductionCMode::FLD) || (ctrl == Group::ProductionCMode::NONE); + } else { + const auto ctrl = well_state_.currentInjectionGroupControl(injection_phase_, child_group); + included = included || (ctrl == Group::InjectionCMode::FLD) || (ctrl == Group::InjectionCMode::NONE); + } if (included) { total_guide_rate += guideRate(child_group, always_included_child); } } for (const std::string& child_well : group.wells()) { - const bool included = (well_state_.isProductionGrup(child_well)) || (child_well == always_included_child); + bool included = (child_well == always_included_child); + if (is_producer_) { + included = included || well_state_.isProductionGrup(child_well); + } else { + included = included || well_state_.isInjectionGrup(child_well); + } + if (included) { total_guide_rate += guideRate(child_well, always_included_child); } @@ -930,11 +892,13 @@ namespace WellGroupHelpers double FractionCalculator::guideRate(const std::string& name, const std::string& always_included_child) { if (schedule_.hasWell(name, report_step_)) { - return guide_rate_->get(name, target_, getRateVector(well_state_, pu_, name)); + return guide_rate_->get(name, target_, getWellRateVector(well_state_, pu_, name)); } else { if (groupControlledWells(name, always_included_child) > 0) { - if (guide_rate_->has(name)) { + if (is_producer_ && guide_rate_->has(name)) { return guide_rate_->get(name, target_, getGroupRateVector(name)); + } else if (!is_producer_ && guide_rate_->has(name, injection_phase_)) { + return guide_rate_->get(name, injection_phase_); } else { // We are a group, with default guide rate. // Compute guide rate by accumulating our children's guide rates. @@ -951,250 +915,16 @@ namespace WellGroupHelpers const std::string& always_included_child) { return ::Opm::WellGroupHelpers::groupControlledWells( - schedule_, well_state_, report_step_, group_name, always_included_child); + schedule_, well_state_, report_step_, group_name, always_included_child, is_producer_, injection_phase_); } GuideRate::RateVector FractionCalculator::getGroupRateVector(const std::string& group_name) { + assert(is_producer_); return getProductionGroupRateVector(this->well_state_, this->pu_, group_name); } - double fractionFromGuideRates(const std::string& name, - const std::string& controlGroupName, - const Schedule& schedule, - const WellStateFullyImplicitBlackoil& wellState, - const int reportStepIdx, - const GuideRate* guideRate, - const GuideRateModel::Target target, - const PhaseUsage& pu, - const bool alwaysIncludeThis) - { - FractionCalculator calc(schedule, wellState, reportStepIdx, guideRate, target, pu); - return calc.fraction(name, controlGroupName, alwaysIncludeThis); - } - - double fractionFromInjectionPotentials(const std::string& name, - const std::string& controlGroupName, - const Schedule& schedule, - const WellStateFullyImplicitBlackoil& wellState, - const int reportStepIdx, - const GuideRate* guideRate, - const GuideRateModel::Target target, - const PhaseUsage& pu, - const Phase& injectionPhase, - const bool alwaysIncludeThis) - { - double thisGuideRate - = getGuideRateInj(name, schedule, wellState, reportStepIdx, guideRate, target, injectionPhase, pu); - double controlGroupGuideRate = getGuideRateInj( - controlGroupName, schedule, wellState, reportStepIdx, guideRate, target, injectionPhase, pu); - if (alwaysIncludeThis) - controlGroupGuideRate += thisGuideRate; - - assert(controlGroupGuideRate >= thisGuideRate); - const double guideRateEpsilon = 1e-12; - return (controlGroupGuideRate > guideRateEpsilon) ? thisGuideRate / controlGroupGuideRate : 0.0; - } - - - std::pair checkGroupConstraintsInj(const std::string& name, - const std::string& parent, - const Group& group, - const WellStateFullyImplicitBlackoil& wellState, - const int reportStepIdx, - const GuideRate* guideRate, - const double* rates, - Phase injectionPhase, - const PhaseUsage& pu, - const double efficiencyFactor, - const Schedule& schedule, - const SummaryState& summaryState, - const std::vector& resv_coeff, - DeferredLogger& deferred_logger) - { - // When called for a well ('name' is a well name), 'parent' - // will be the name of 'group'. But if we recurse, 'name' and - // 'parent' will stay fixed while 'group' will be higher up - // in the group tree. - // efficiency factor is the well efficiency factor for the first group the well is - // part of. Later it is the accumulated factor including the group efficiency factor - // of the child of group. - - const Group::InjectionCMode& currentGroupControl - = wellState.currentInjectionGroupControl(injectionPhase, group.name()); - if (currentGroupControl == Group::InjectionCMode::FLD || currentGroupControl == Group::InjectionCMode::NONE) { - // Return if we are not available for parent group. - if (!group.injectionGroupControlAvailable(injectionPhase)) { - return std::make_pair(false, 1.0); - } - // Otherwise: check injection share of parent's control. - const auto& parentGroup = schedule.getGroup(group.parent(), reportStepIdx); - return checkGroupConstraintsInj(name, - parent, - parentGroup, - wellState, - reportStepIdx, - guideRate, - rates, - injectionPhase, - pu, - efficiencyFactor * group.getGroupEfficiencyFactor(), - schedule, - summaryState, - resv_coeff, - deferred_logger); - } - - // If we are here, we are at the topmost group to be visited in the recursion. - // This is the group containing the control we will check against. - - // This can be false for FLD-controlled groups, we must therefore - // check for FLD first (done above). - if (!group.isInjectionGroup()) { - return std::make_pair(false, 1.0); - } - - int phasePos; - GuideRateModel::Target target; - - switch (injectionPhase) { - case Phase::WATER: { - phasePos = pu.phase_pos[BlackoilPhases::Aqua]; - target = GuideRateModel::Target::WAT; - break; - } - case Phase::OIL: { - phasePos = pu.phase_pos[BlackoilPhases::Liquid]; - target = GuideRateModel::Target::OIL; - break; - } - case Phase::GAS: { - phasePos = pu.phase_pos[BlackoilPhases::Vapour]; - target = GuideRateModel::Target::GAS; - break; - } - default: - OPM_DEFLOG_THROW( - std::logic_error, "Expected WATER, OIL or GAS as injecting type for " + name, deferred_logger); - } - - assert(group.hasInjectionControl(injectionPhase)); - const auto& groupcontrols = group.injectionControls(injectionPhase, summaryState); - - const std::vector& groupInjectionReductions - = wellState.currentInjectionGroupReductionRates(group.name()); - const double groupTargetReduction = groupInjectionReductions[phasePos]; - double fraction = fractionFromInjectionPotentials( - name, group.name(), schedule, wellState, reportStepIdx, guideRate, target, pu, injectionPhase, true); - double target_fraction = 1.0; - bool constraint_broken = false; - double efficiencyFactorInclGroup = efficiencyFactor * group.getGroupEfficiencyFactor(); - - switch (currentGroupControl) { - case Group::InjectionCMode::RATE: { - const double current_rate = rates[phasePos]; - const double target_rate = fraction - * std::max(0.0, - (groupcontrols.surface_max_rate - groupTargetReduction + current_rate * efficiencyFactorInclGroup)) - / efficiencyFactorInclGroup; - if (current_rate > target_rate) { - constraint_broken = true; - target_fraction = target_rate / current_rate; - } - break; - } - case Group::InjectionCMode::RESV: { - const double coeff = resv_coeff[phasePos]; - const double current_rate = rates[phasePos]; - const double target_rate = fraction - * std::max(0.0, - (groupcontrols.resv_max_rate / coeff - groupTargetReduction - + current_rate * efficiencyFactorInclGroup)) - / efficiencyFactorInclGroup; - if (current_rate > target_rate) { - constraint_broken = true; - target_fraction = target_rate / current_rate; - } - break; - } - case Group::InjectionCMode::REIN: { - double productionRate = wellState.currentInjectionREINRates(groupcontrols.reinj_group)[phasePos]; - const double current_rate = rates[phasePos]; - const double target_rate = fraction - * std::max(0.0, - (groupcontrols.target_reinj_fraction * productionRate - groupTargetReduction - + current_rate * efficiencyFactorInclGroup)) - / efficiencyFactorInclGroup; - if (current_rate > target_rate) { - constraint_broken = true; - target_fraction = target_rate / current_rate; - } - break; - } - case Group::InjectionCMode::VREP: { - const double coeff = resv_coeff[phasePos]; - double voidageRate - = wellState.currentInjectionVREPRates(groupcontrols.voidage_group) * groupcontrols.target_void_fraction; - - double injReduction = 0.0; - if (groupcontrols.phase != Phase::WATER) - injReduction += groupInjectionReductions[pu.phase_pos[BlackoilPhases::Aqua]] - * resv_coeff[pu.phase_pos[BlackoilPhases::Aqua]]; - if (groupcontrols.phase != Phase::OIL) - injReduction += groupInjectionReductions[pu.phase_pos[BlackoilPhases::Liquid]] - * resv_coeff[pu.phase_pos[BlackoilPhases::Liquid]]; - if (groupcontrols.phase != Phase::GAS) - injReduction += groupInjectionReductions[pu.phase_pos[BlackoilPhases::Vapour]] - * resv_coeff[pu.phase_pos[BlackoilPhases::Vapour]]; - voidageRate -= injReduction; - - const double current_rate = rates[phasePos]; - const double target_rate = fraction - * std::max(0.0, (voidageRate / coeff - groupTargetReduction + current_rate * efficiencyFactorInclGroup)) - / efficiencyFactorInclGroup; - if (current_rate > target_rate) { - constraint_broken = true; - target_fraction = target_rate / current_rate; - } - break; - } - case Group::InjectionCMode::SALE: { - // only for gas injectors - assert(phasePos == pu.phase_pos[BlackoilPhases::Vapour]); - - // Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate; - // Gas import and consumption is already included in the REIN rates - double inj_rate = wellState.currentInjectionREINRates(group.name())[phasePos]; - const auto& gconsale = schedule[reportStepIdx].gconsale().get(group.name(), summaryState); - inj_rate -= gconsale.sales_target; - - const double current_rate = rates[phasePos]; - const double target_rate = fraction - * std::max(0.0, (inj_rate - groupTargetReduction + current_rate * efficiencyFactorInclGroup)) / efficiencyFactorInclGroup; - if (current_rate > target_rate) { - constraint_broken = true; - target_fraction = target_rate / current_rate; - } - break; - } - case Group::InjectionCMode::NONE: { - assert(false); // Should already be handled at the top. - } - case Group::InjectionCMode::FLD: { - assert(false); // Should already be handled at the top. - } - default: - OPM_DEFLOG_THROW( - std::runtime_error, "Invalid group control specified for group " + group.name(), deferred_logger); - } - - return std::make_pair(constraint_broken, target_fraction); - } - - - - std::vector groupChainTopBot(const std::string& bottom, const std::string& top, const Schedule& schedule, const int report_step) { @@ -1286,7 +1016,7 @@ namespace WellGroupHelpers gratTargetFromSales = wellState.currentGroupGratTargetFromSales(group.name()); TargetCalculator tcalc(currentGroupControl, pu, resv_coeff, gratTargetFromSales); - FractionCalculator fcalc(schedule, wellState, reportStepIdx, guideRate, tcalc.guideTargetMode(), pu); + FractionCalculator fcalc(schedule, summaryState, wellState, reportStepIdx, guideRate, tcalc.guideTargetMode(), pu, true, Phase::OIL); auto localFraction = [&](const std::string& child) { return fcalc.localFraction(child, name); }; @@ -1333,7 +1063,7 @@ namespace WellGroupHelpers // the current well to be always included, because we // want to know the situation that applied to the // calculation of reductions. - const int num_gr_ctrl = groupControlledWells(schedule, wellState, reportStepIdx, chain[ii + 1], ""); + const int num_gr_ctrl = groupControlledWells(schedule, wellState, reportStepIdx, chain[ii + 1], "", /*is_producer*/true, /*injectionPhaseNotUsed*/Phase::OIL); if (num_gr_ctrl == 0) { if (guideRate->has(chain[ii + 1])) { target += localReduction(chain[ii + 1]); @@ -1347,6 +1077,130 @@ namespace WellGroupHelpers return std::make_pair(current_rate > target_rate, target_rate / current_rate); } + std::pair checkGroupConstraintsInj(const std::string& name, + const std::string& parent, + const Group& group, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const GuideRate* guideRate, + const double* rates, + Phase injectionPhase, + const PhaseUsage& pu, + const double efficiencyFactor, + const Schedule& schedule, + const SummaryState& summaryState, + const std::vector& resv_coeff, + DeferredLogger& deferred_logger) + { + // When called for a well ('name' is a well name), 'parent' + // will be the name of 'group'. But if we recurse, 'name' and + // 'parent' will stay fixed while 'group' will be higher up + // in the group tree. + // efficiencyfactor is the well efficiency factor for the first group the well is + // part of. Later it is the accumulated factor including the group efficiency factor + // of the child of group. + + const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(injectionPhase, group.name()); + + if (currentGroupControl == Group::InjectionCMode::FLD || currentGroupControl == Group::InjectionCMode::NONE) { + // Return if we are not available for parent group. + if (!group.injectionGroupControlAvailable(injectionPhase)) { + return std::make_pair(false, 1); + } + // Otherwise: check production share of parent's control. + const auto& parentGroup = schedule.getGroup(group.parent(), reportStepIdx); + return checkGroupConstraintsInj(name, + parent, + parentGroup, + wellState, + reportStepIdx, + guideRate, + rates, + injectionPhase, + pu, + efficiencyFactor * group.getGroupEfficiencyFactor(), + schedule, + summaryState, + resv_coeff, + deferred_logger); + } + + // This can be false for FLD-controlled groups, we must therefore + // check for FLD first (done above). + if (!group.isInjectionGroup()) { + return std::make_pair(false, 1.0); + } + + // If we are here, we are at the topmost group to be visited in the recursion. + // This is the group containing the control we will check against. + double sales_target = 0; + if (schedule[reportStepIdx].gconsale().has(group.name())) { + const auto& gconsale = schedule[reportStepIdx].gconsale().get(group.name(), summaryState); + sales_target = gconsale.sales_target; + } + InjectionTargetCalculator tcalc(currentGroupControl, pu, resv_coeff, group.name(), sales_target, wellState, injectionPhase, deferred_logger); + FractionCalculator fcalc(schedule, summaryState, wellState, reportStepIdx, guideRate, tcalc.guideTargetMode(), pu, false, injectionPhase); + + auto localFraction = [&](const std::string& child) { return fcalc.localFraction(child, name); }; + + auto localReduction = [&](const std::string& group_name) { + const std::vector& groupTargetReductions + = wellState.currentInjectionGroupReductionRates(group_name); + return tcalc.calcModeRateFromRates(groupTargetReductions); + }; + + const double orig_target = tcalc.groupTarget(group.injectionControls(injectionPhase, summaryState), deferred_logger); + // Assume we have a chain of groups as follows: BOTTOM -> MIDDLE -> TOP. + // Then ... + // TODO finish explanation. + const double current_rate + = tcalc.calcModeRateFromRates(rates); // Switch sign since 'rates' are negative for producers. + const auto chain = groupChainTopBot(name, group.name(), schedule, reportStepIdx); + // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. + const size_t num_ancestors = chain.size() - 1; + // we need to find out the level where the current well is applied to the local reduction + size_t local_reduction_level = 0; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + if ((ii == 0) || guideRate->has(chain[ii], injectionPhase)) { + local_reduction_level = ii; + } + } + + double efficiencyFactorInclGroup = efficiencyFactor * group.getGroupEfficiencyFactor(); + double target = orig_target; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + if ((ii == 0) || guideRate->has(chain[ii], injectionPhase)) { + // Apply local reductions only at the control level + // (top) and for levels where we have a specified + // group guide rate. + target -= localReduction(chain[ii]); + + // Add my reduction back at the level where it is included in the local reduction + if (local_reduction_level == ii ) + target += current_rate * efficiencyFactorInclGroup; + } + if (ii < num_ancestors - 1) { + // Not final level. Add sub-level reduction back, if + // it was nonzero due to having no group-controlled + // wells. Note that we make this call without setting + // the current well to be always included, because we + // want to know the situation that applied to the + // calculation of reductions. + const int num_gr_ctrl = groupControlledWells(schedule, wellState, reportStepIdx, chain[ii + 1], "", /*is_producer*/false, injectionPhase); + if (num_gr_ctrl == 0) { + if (guideRate->has(chain[ii + 1], injectionPhase)) { + target += localReduction(chain[ii + 1]); + } + } + } + target *= localFraction(chain[ii + 1]); + } + // Avoid negative target rates comming from too large local reductions. + const double target_rate = std::max(1e-12, target / efficiencyFactorInclGroup); + return std::make_pair(current_rate > target_rate, target_rate / current_rate); + } + + } // namespace WellGroupHelpers diff --git a/opm/simulators/wells/WellGroupHelpers.hpp b/opm/simulators/wells/WellGroupHelpers.hpp index 658119f1e..fdaa68a29 100644 --- a/opm/simulators/wells/WellGroupHelpers.hpp +++ b/opm/simulators/wells/WellGroupHelpers.hpp @@ -95,16 +95,15 @@ namespace WellGroupHelpers std::vector& groupTargetReduction); template - void updateGuideRateForGroups(const Group& group, - const Schedule& schedule, - const PhaseUsage& pu, - const int reportStepIdx, - const double& simTime, - const bool isInjector, - WellStateFullyImplicitBlackoil& wellState, - const Comm& comm, - GuideRate* guideRate, - std::vector& pot) + void updateGuideRateForProductionGroups(const Group& group, + const Schedule& schedule, + const PhaseUsage& pu, + const int reportStepIdx, + const double& simTime, + WellStateFullyImplicitBlackoil& wellState, + const Comm& comm, + GuideRate* guideRate, + std::vector& pot) { const int np = pu.num_phases; for (const std::string& groupName : group.groups()) { @@ -112,44 +111,25 @@ namespace WellGroupHelpers const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); // Note that group effiency factors for groupTmp are applied in updateGuideRateForGroups - updateGuideRateForGroups( - groupTmp, schedule, pu, reportStepIdx, simTime, isInjector, wellState, comm, guideRate, thisPot); + updateGuideRateForProductionGroups( + groupTmp, schedule, pu, reportStepIdx, simTime, wellState, comm, guideRate, thisPot); // accumulate group contribution from sub group unconditionally - if (isInjector) { - const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS}; - for (Phase phase : all) { - int phasePos; - if (phase == Phase::GAS && pu.phase_used[BlackoilPhases::Vapour]) - phasePos = pu.phase_pos[BlackoilPhases::Vapour]; - else if (phase == Phase::OIL && pu.phase_used[BlackoilPhases::Liquid]) - phasePos = pu.phase_pos[BlackoilPhases::Liquid]; - else if (phase == Phase::WATER && pu.phase_used[BlackoilPhases::Aqua]) - phasePos = pu.phase_pos[BlackoilPhases::Aqua]; - else - continue; - - pot[phasePos] += thisPot[phasePos]; - } - } else { - const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName); - if (currentGroupControl != Group::ProductionCMode::FLD + const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName); + if (currentGroupControl != Group::ProductionCMode::FLD && currentGroupControl != Group::ProductionCMode::NONE) { - continue; - } - for (int phase = 0; phase < np; phase++) { - pot[phase] += thisPot[phase]; - } + continue; } + for (int phase = 0; phase < np; phase++) { + pot[phase] += thisPot[phase]; + } + } for (const std::string& wellName : group.wells()) { const auto& wellTmp = schedule.getWell(wellName, reportStepIdx); const auto wefac = wellTmp.getEfficiencyFactor(); - if (wellTmp.isProducer() && isInjector) - continue; - - if (wellTmp.isInjector() && !isInjector) + if (wellTmp.isInjector()) continue; if (wellTmp.getStatus() == Well::Status::SHUT) @@ -195,12 +175,7 @@ namespace WellGroupHelpers oilPot = comm.sum(oilPot); gasPot = comm.sum(gasPot); waterPot = comm.sum(waterPot); - - if (isInjector) { - wellState.setCurrentGroupInjectionPotentials(group.name(), pot); - } else { - guideRate->compute(group.name(), reportStepIdx, simTime, oilPot, gasPot, waterPot); - } + guideRate->compute(group.name(), reportStepIdx, simTime, oilPot, gasPot, waterPot); } template @@ -242,6 +217,14 @@ namespace WellGroupHelpers } } + void updateGuideRatesForInjectionGroups(const Group& group, + const Schedule& schedule, + const SummaryState& summaryState, + const Opm::PhaseUsage& pu, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellState, + GuideRate* guideRate, + Opm::DeferredLogger& deferred_logger); void updateVREPForGroups(const Group& group, const Schedule& schedule, @@ -283,7 +266,7 @@ namespace WellGroupHelpers const int report_time_step); GuideRate::RateVector - getRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& name); + getWellRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& name); GuideRate::RateVector getProductionGroupRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& group_name); @@ -310,18 +293,23 @@ namespace WellGroupHelpers const WellStateFullyImplicitBlackoil& well_state, const int report_step, const std::string& group_name, - const std::string& always_included_child); + const std::string& always_included_child, + const bool is_production_group, + const Phase injection_phase); class FractionCalculator { public: FractionCalculator(const Schedule& schedule, + const SummaryState& summary_state, const WellStateFullyImplicitBlackoil& well_state, const int report_step, const GuideRate* guide_rate, const GuideRateModel::Target target, - const PhaseUsage& pu); + const PhaseUsage& pu, + const bool is_producer, + const Phase injection_phase); double fraction(const std::string& name, const std::string& control_group_name, const bool always_include_this); double localFraction(const std::string& name, const std::string& always_included_child); @@ -332,36 +320,17 @@ namespace WellGroupHelpers int groupControlledWells(const std::string& group_name, const std::string& always_included_child); GuideRate::RateVector getGroupRateVector(const std::string& group_name); const Schedule& schedule_; + const SummaryState& summary_state_; const WellStateFullyImplicitBlackoil& well_state_; int report_step_; const GuideRate* guide_rate_; GuideRateModel::Target target_; PhaseUsage pu_; + bool is_producer_; + Phase injection_phase_; }; - - double fractionFromGuideRates(const std::string& name, - const std::string& controlGroupName, - const Schedule& schedule, - const WellStateFullyImplicitBlackoil& wellState, - const int reportStepIdx, - const GuideRate* guideRate, - const GuideRateModel::Target target, - const PhaseUsage& pu, - const bool alwaysIncludeThis = false); - - double fractionFromInjectionPotentials(const std::string& name, - const std::string& controlGroupName, - const Schedule& schedule, - const WellStateFullyImplicitBlackoil& wellState, - const int reportStepIdx, - const GuideRate* guideRate, - const GuideRateModel::Target target, - const PhaseUsage& pu, - const Phase& injectionPhase, - const bool alwaysIncludeThis = false); - std::pair checkGroupConstraintsInj(const std::string& name, const std::string& parent, const Group& group, diff --git a/opm/simulators/wells/WellInterface.hpp b/opm/simulators/wells/WellInterface.hpp index 1a47aeb4d..d020cc99e 100644 --- a/opm/simulators/wells/WellInterface.hpp +++ b/opm/simulators/wells/WellInterface.hpp @@ -597,7 +597,8 @@ namespace Opm const EvalWell& bhp, const EvalWell& injection_rate, EvalWell& control_eq, - double efficiencyFactor); + double efficiencyFactor, + Opm::DeferredLogger& deferred_logger); template void getGroupProductionControl(const Group& group, diff --git a/opm/simulators/wells/WellInterface_impl.hpp b/opm/simulators/wells/WellInterface_impl.hpp index 372a9dc1d..8b7e200b6 100644 --- a/opm/simulators/wells/WellInterface_impl.hpp +++ b/opm/simulators/wells/WellInterface_impl.hpp @@ -2016,7 +2016,8 @@ namespace Opm bhp, injection_rate, control_eq, - efficiencyFactor); + efficiencyFactor, + deferred_logger); break; } case Well::InjectorCMode::CMODE_UNDEFINED: { @@ -2142,47 +2143,36 @@ namespace Opm - template template void WellInterface::getGroupInjectionControl(const Group& group, - const WellState& well_state, - const Opm::Schedule& schedule, - const SummaryState& summaryState, - const InjectorType& injectorType, - const EvalWell& bhp, - const EvalWell& injection_rate, - EvalWell& control_eq, - double efficiencyFactor) + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + const EvalWell& bhp, + const EvalWell& injection_rate, + EvalWell& control_eq, + double efficiencyFactor, + Opm::DeferredLogger& deferred_logger) { - const auto& well = well_ecl_; - const auto pu = phaseUsage(); - // Setting some defaults to silence warnings below. // Will be overwritten in the switch statement. - int phasePos = -1; - Well::GuideRateTarget wellTarget = Well::GuideRateTarget::UNDEFINED; Phase injectionPhase = Phase::WATER; switch (injectorType) { case InjectorType::WATER: { - phasePos = pu.phase_pos[BlackoilPhases::Aqua]; - wellTarget = Well::GuideRateTarget::WAT; injectionPhase = Phase::WATER; break; } case InjectorType::OIL: { - phasePos = pu.phase_pos[BlackoilPhases::Liquid]; - wellTarget = Well::GuideRateTarget::OIL; injectionPhase = Phase::OIL; break; } case InjectorType::GAS: { - phasePos = pu.phase_pos[BlackoilPhases::Vapour]; - wellTarget = Well::GuideRateTarget::GAS; injectionPhase = Phase::GAS; break; } @@ -2192,7 +2182,6 @@ namespace Opm } const Group::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(injectionPhase, group.name()); - if (currentGroupControl == Group::InjectionCMode::FLD || currentGroupControl == Group::InjectionCMode::NONE) { if (!group.injectionGroupControlAvailable(injectionPhase)) { @@ -2211,109 +2200,68 @@ namespace Opm // Inject share of parents control const auto& parent = schedule.getGroup( group.parent(), current_step_ ); efficiencyFactor *= group.getGroupEfficiencyFactor(); - getGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, bhp, injection_rate, control_eq, efficiencyFactor); + getGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, bhp, injection_rate, control_eq, efficiencyFactor, deferred_logger); return; } } efficiencyFactor *= group.getGroupEfficiencyFactor(); - assert(group.hasInjectionControl(injectionPhase)); - const auto& groupcontrols = group.injectionControls(injectionPhase, summaryState); + const auto& well = well_ecl_; + const auto pu = phaseUsage(); - const std::vector& groupInjectionReductions = well_state.currentInjectionGroupReductionRates(group.name()); - double groupTargetReduction = groupInjectionReductions[phasePos]; - double fraction = WellGroupHelpers::fractionFromInjectionPotentials(well.name(), - group.name(), - schedule, - well_state, - current_step_, - guide_rate_, - GuideRateModel::convert_target(wellTarget), - pu, - injectionPhase, - false); - switch (currentGroupControl) { - case Group::InjectionCMode::NONE: - { - // The NONE case is handled earlier - assert(false); - break; + if (!group.isInjectionGroup()) { + // use bhp as control eq and let the updateControl code find a valid control + const auto& controls = well.injectionControls(summaryState); + control_eq = bhp - controls.bhp_limit; + return; } - case Group::InjectionCMode::RATE: - { - double target = std::max(0.0, (groupcontrols.surface_max_rate - groupTargetReduction)) / efficiencyFactor; - control_eq = injection_rate - fraction * target; - break; - } - case Group::InjectionCMode::RESV: - { - std::vector convert_coeff(number_of_phases_, 1.0); - rateConverter_.calcCoeff(/*fipreg*/ 0, pvtRegionIdx_, convert_coeff); - double coeff = convert_coeff[phasePos]; - double target = std::max(0.0, (groupcontrols.resv_max_rate/coeff - groupTargetReduction)) / efficiencyFactor; - control_eq = injection_rate - fraction * target; - break; - } - case Group::InjectionCMode::REIN: - { - double productionRate = well_state.currentInjectionREINRates(groupcontrols.reinj_group)[phasePos]; - double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction)) / efficiencyFactor; - control_eq = injection_rate - fraction * target; - break; - } - case Group::InjectionCMode::VREP: - { - std::vector convert_coeff(number_of_phases_, 1.0); - rateConverter_.calcCoeff(/*fipreg*/ 0, pvtRegionIdx_, convert_coeff); - double coeff = convert_coeff[phasePos]; - double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group)*groupcontrols.target_void_fraction; - double injReduction = 0.0; - std::vector groupInjectionReservoirRates = well_state.currentInjectionGroupReservoirRates(group.name()); - if (groupcontrols.phase != Phase::WATER) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Aqua]]; + // If we are here, we are at the topmost group to be visited in the recursion. + // This is the group containing the control we will check against. - if (groupcontrols.phase != Phase::OIL) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Liquid]]; + // Make conversion factors for RESV <-> surface rates. + std::vector resv_coeff(phaseUsage().num_phases, 1.0); + rateConverter_.calcCoeff(0, pvtRegionIdx_, resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS. - if (groupcontrols.phase != Phase::GAS) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Vapour]]; - - voidageRate -= injReduction; - - double target = std::max(0.0, ( voidageRate/coeff - groupTargetReduction)) / efficiencyFactor; - control_eq = injection_rate - fraction * target; - break; - } - case Group::InjectionCMode::FLD: - { - // The FLD case is handled earlier - assert(false); - break; - } - case Group::InjectionCMode::SALE: - { - // only for gas injectors - assert (phasePos == pu.phase_pos[BlackoilPhases::Vapour]); - - // Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate; - // The import and consumption is already included in the REIN rates. - double inj_rate = well_state.currentInjectionREINRates(group.name())[phasePos]; + double sales_target = 0; + if (schedule[current_step_].gconsale().has(group.name())) { const auto& gconsale = schedule[current_step_].gconsale().get(group.name(), summaryState); - inj_rate -= gconsale.sales_target; + sales_target = gconsale.sales_target; + } + WellGroupHelpers::InjectionTargetCalculator tcalc(currentGroupControl, pu, resv_coeff, group.name(), sales_target, well_state, injectionPhase, deferred_logger); + WellGroupHelpers::FractionCalculator fcalc(schedule, summaryState, well_state, current_step_, guide_rate_, tcalc.guideTargetMode(), pu, false, injectionPhase); - double target = std::max(0.0, (inj_rate - groupTargetReduction)) / efficiencyFactor; - control_eq = injection_rate - fraction * target; - break; - } - // default: - // OPM_DEFLOG_THROW(std::runtime_error, "Unvalid group control specified for group " + well.groupName(), deferred_logger ); + auto localFraction = [&](const std::string& child) { + return fcalc.localFraction(child, ""); + }; + + auto localReduction = [&](const std::string& group_name) { + const std::vector& groupTargetReductions = well_state.currentInjectionGroupReductionRates(group_name); + return tcalc.calcModeRateFromRates(groupTargetReductions); + }; + + const double orig_target = tcalc.groupTarget(group.injectionControls(injectionPhase, summaryState), deferred_logger); + const auto chain = WellGroupHelpers::groupChainTopBot(name(), group.name(), schedule, current_step_); + // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. + const size_t num_ancestors = chain.size() - 1; + double target = orig_target; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + if ((ii == 0) || guide_rate_->has(chain[ii], injectionPhase)) { + // Apply local reductions only at the control level + // (top) and for levels where we have a specified + // group guide rate. + target -= localReduction(chain[ii]); + } + target *= localFraction(chain[ii+1]); } + // Avoid negative target rates coming from too large local reductions. + const double target_rate = std::max(0.0, target / efficiencyFactor); + const auto current_rate = injection_rate; // Switch sign since 'rates' are negative for producers. + control_eq = current_rate - target_rate; } - template template void @@ -2375,7 +2323,7 @@ namespace Opm gratTargetFromSales = well_state.currentGroupGratTargetFromSales(group.name()); WellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, resv_coeff, gratTargetFromSales); - WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, current_step_, guide_rate_, tcalc.guideTargetMode(), pu); + WellGroupHelpers::FractionCalculator fcalc(schedule, summaryState, well_state, current_step_, guide_rate_, tcalc.guideTargetMode(), pu, true, Phase::OIL); auto localFraction = [&](const std::string& child) { return fcalc.localFraction(child, ""); diff --git a/opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp b/opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp index e9ede1547..8a19d89ce 100644 --- a/opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp +++ b/opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp @@ -449,6 +449,7 @@ namespace Opm return this->production_group_rates.find(groupName) != this->production_group_rates.end(); } + void setCurrentProductionGroupReductionRates(const std::string& groupName, const std::vector& target ) { production_group_reduction_rates[groupName] = target; }