From f30ddffdb7826851b16a0dc85d99cc0e9a432f9c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Fri, 27 Mar 2020 13:27:45 +0100 Subject: [PATCH] Split WellGroupHelpers to cpp/hpp. Also changed namespace name. --- CMakeLists_files.cmake | 2 + .../wells/BlackoilWellModel_impl.hpp | 88 +- opm/simulators/wells/TargetCalculator.hpp | 160 +++ opm/simulators/wells/WellGroupHelpers.cpp | 1165 ++++++++++++++++ opm/simulators/wells/WellGroupHelpers.hpp | 1225 +++-------------- opm/simulators/wells/WellInterface_impl.hpp | 33 +- 6 files changed, 1589 insertions(+), 1084 deletions(-) create mode 100644 opm/simulators/wells/TargetCalculator.hpp create mode 100644 opm/simulators/wells/WellGroupHelpers.cpp diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 9a02b13d8..7cf3171ca 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -39,6 +39,7 @@ list (APPEND MAIN_SOURCE_FILES opm/simulators/utils/ParallelRestart.cpp opm/simulators/wells/VFPProdProperties.cpp opm/simulators/wells/VFPInjProperties.cpp + opm/simulators/wells/WellGroupHelpers.cpp ) if(CUDA_FOUND) @@ -188,6 +189,7 @@ list (APPEND PUBLIC_HEADER_FILES opm/simulators/wells/PerforationData.hpp opm/simulators/wells/RateConverter.hpp opm/simulators/wells/SimFIBODetails.hpp + opm/simulators/wells/TargetCalculator.hpp opm/simulators/wells/WellConnectionAuxiliaryModule.hpp opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp opm/simulators/wells/VFPProperties.hpp diff --git a/opm/simulators/wells/BlackoilWellModel_impl.hpp b/opm/simulators/wells/BlackoilWellModel_impl.hpp index 185cd2416..72dc169f7 100644 --- a/opm/simulators/wells/BlackoilWellModel_impl.hpp +++ b/opm/simulators/wells/BlackoilWellModel_impl.hpp @@ -292,7 +292,7 @@ namespace Opm { } } const Group& fieldGroup = schedule().getGroup("FIELD", timeStepIdx); - wellGroupHelpers::setCmodeGroup(fieldGroup, schedule(), summaryState, timeStepIdx, well_state_); + WellGroupHelpers::setCmodeGroup(fieldGroup, schedule(), summaryState, timeStepIdx, well_state_); // Compute reservoir volumes for RESV controls. rateConverter_.reset(new RateConverterType (phase_usage_, @@ -381,7 +381,7 @@ namespace Opm { //compute well guideRates const auto& comm = ebosSimulator_.vanguard().grid().comm(); - wellGroupHelpers::updateGuideRatesForWells(schedule(), phase_usage_, reportStepIdx, simulationTime, well_state_, comm, guideRate_.get()); + WellGroupHelpers::updateGuideRatesForWells(schedule(), phase_usage_, reportStepIdx, simulationTime, well_state_, comm, guideRate_.get()); } @@ -408,7 +408,7 @@ namespace Opm { well->init(&phase_usage_, depth_, gravity_, number_of_cells_); const Well& wellEcl = schedule().getWell(well_name, timeStepIdx); double well_efficiency_factor = wellEcl.getEfficiencyFactor(); - wellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), timeStepIdx), schedule(), timeStepIdx, well_efficiency_factor); + WellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), timeStepIdx), schedule(), timeStepIdx, well_efficiency_factor); well->setWellEfficiencyFactor(well_efficiency_factor); well->setVFPProperties(vfp_properties_.get()); well->setGuideRate(guideRate_.get()); @@ -1231,23 +1231,23 @@ namespace Opm { // 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); - wellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ false, phase_usage_, well_state_nupcol_, well_state_, groupTargetReduction); + WellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ false, phase_usage_, well_state_nupcol_, well_state_, groupTargetReduction); std::vector groupTargetReductionInj(numPhases(), 0.0); - wellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ true, phase_usage_, well_state_nupcol_, well_state_, groupTargetReductionInj); + WellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ true, phase_usage_, well_state_nupcol_, well_state_, groupTargetReductionInj); 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); + 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::updateGuideRateForGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, /*isInjector*/ true, well_state_, comm, guideRate_.get(), potInj); 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_); + WellGroupHelpers::updateREINForGroups(fieldGroup, schedule(), reportStepIdx, phase_usage_, summaryState, well_state_nupcol_, well_state_); + WellGroupHelpers::updateVREPForGroups(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_); - wellGroupHelpers::updateReservoirRatesInjectionGroups(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_); - wellGroupHelpers::updateGroupProductionRates(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_); - wellGroupHelpers::updateWellRates(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_); + WellGroupHelpers::updateReservoirRatesInjectionGroups(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_); + WellGroupHelpers::updateGroupProductionRates(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_); + WellGroupHelpers::updateWellRates(fieldGroup, schedule(), reportStepIdx, well_state_nupcol_, well_state_); well_state_.communicateGroupRates(comm); // compute wsolvent fraction for REIN wells @@ -1396,7 +1396,7 @@ namespace Opm { for (auto& well : well_container_) { const Well& wellEcl = well->wellEcl(); double well_efficiency_factor = wellEcl.getEfficiencyFactor(); - wellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), reportStepIdx), schedule(), reportStepIdx, well_efficiency_factor); + WellGroupHelpers::accumulateGroupEfficiencyFactor(schedule().getGroup(wellEcl.groupName(), reportStepIdx), schedule(), reportStepIdx, well_efficiency_factor); well->setWellEfficiencyFactor(well_efficiency_factor); } } @@ -1846,7 +1846,7 @@ namespace Opm { if (currentControl != Group::ProductionCMode::ORAT) { double current_rate = 0.0; - current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false); + current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false); // sum over all nodes current_rate = comm.sum(current_rate); @@ -1863,7 +1863,7 @@ namespace Opm { { double current_rate = 0.0; - current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false); + current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false); // sum over all nodes current_rate = comm.sum(current_rate); @@ -1878,7 +1878,7 @@ namespace Opm { if (currentControl != Group::ProductionCMode::GRAT) { double current_rate = 0.0; - current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], false); + current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], false); // sum over all nodes current_rate = comm.sum(current_rate); @@ -1892,8 +1892,8 @@ namespace Opm { if (currentControl != Group::ProductionCMode::LRAT) { double current_rate = 0.0; - current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false); - current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false); + current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false); + current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false); // sum over all nodes current_rate = comm.sum(current_rate); @@ -1913,9 +1913,9 @@ namespace Opm { if (currentControl != Group::ProductionCMode::RESV) { double current_rate = 0.0; - current_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], true); - current_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], true); - current_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], true); + current_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], true); + current_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], true); + current_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], true); // sum over all nodes current_rate = comm.sum(current_rate); @@ -1962,7 +1962,7 @@ namespace Opm { if (currentControl != Group::InjectionCMode::RATE) { double current_rate = 0.0; - current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true); + current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true); // sum over all nodes current_rate = comm.sum(current_rate); @@ -1977,7 +1977,7 @@ namespace Opm { if (currentControl != Group::InjectionCMode::RESV) { double current_rate = 0.0; - current_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true); + current_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true); // sum over all nodes current_rate = comm.sum(current_rate); @@ -1992,13 +1992,13 @@ namespace Opm { { double production_Rate = 0.0; const Group& groupRein = schedule().getGroup(controls.reinj_group, reportStepIdx); - production_Rate += wellGroupHelpers::sumWellRates(groupRein, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/false); + production_Rate += WellGroupHelpers::sumWellRates(groupRein, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/false); // sum over all nodes production_Rate = comm.sum(production_Rate); double current_rate = 0.0; - current_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true); + current_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, phasePos, /*isInjector*/true); // sum over all nodes current_rate = comm.sum(current_rate); @@ -2014,17 +2014,17 @@ namespace Opm { { double voidage_rate = 0.0; const Group& groupVoidage = schedule().getGroup(controls.voidage_group, reportStepIdx); - voidage_rate += wellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false); - voidage_rate += wellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false); - voidage_rate += wellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], false); + voidage_rate += WellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], false); + voidage_rate += WellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], false); + voidage_rate += WellGroupHelpers::sumWellResRates(groupVoidage, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], false); // sum over all nodes voidage_rate = comm.sum(voidage_rate); double total_rate = 0.0; - total_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], true); - total_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], true); - total_rate += wellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], true); + total_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Aqua], true); + total_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Liquid], true); + total_rate += WellGroupHelpers::sumWellResRates(group, schedule(), well_state, reportStepIdx, phase_usage_.phase_pos[BlackoilPhases::Vapour], true); // sum over all nodes total_rate = comm.sum(total_rate); @@ -2044,8 +2044,8 @@ namespace Opm { double sales_rate = 0.0; int gasPos = phase_usage_.phase_pos[BlackoilPhases::Vapour]; - sales_rate += wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, gasPos, /*isInjector*/false); - sales_rate -= wellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, gasPos, /*isInjector*/true); + sales_rate += WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, gasPos, /*isInjector*/false); + sales_rate -= WellGroupHelpers::sumWellRates(group, schedule(), well_state, reportStepIdx, gasPos, /*isInjector*/true); // sum over all nodes sales_rate = comm.sum(sales_rate); @@ -2171,13 +2171,17 @@ namespace Opm { BlackoilWellModel:: checkGroupHigherConstraints(const Group& group, Opm::DeferredLogger& deferred_logger, std::set& switched_groups) { + // Set up coefficients for RESV <-> surface rate conversion. // Use the pvtRegionIdx from the top cell of the first well. // TODO fix this! // This is only used for converting RESV rates. // What is the proper approach? + const int fipnum = 0; const int pvtreg = well_perf_data_.empty() ? pvt_region_idx_[0] : pvt_region_idx_[well_perf_data_[0][0].cell_index]; + std::vector resv_coeff(phase_usage_.num_phases, 0.0); + rateConverter_->calcCoeff(fipnum, pvtreg, resv_coeff); const int reportStepIdx = ebosSimulator_.episodeIndex(); const auto& summaryState = ebosSimulator_.vanguard().summaryState(); @@ -2190,7 +2194,7 @@ namespace Opm { if (!skip && group.isInjectionGroup()) { // Obtain rates for group. for (int phasePos = 0; phasePos < phase_usage_.num_phases; ++phasePos) { - const double local_current_rate = wellGroupHelpers::sumWellRates( + const double local_current_rate = WellGroupHelpers::sumWellRates( group, schedule(), well_state_, reportStepIdx, phasePos, /* isInjector */ true); // Sum over all processes rates[phasePos] = comm.sum(local_current_rate); @@ -2201,7 +2205,7 @@ namespace Opm { const Group::InjectionCMode& currentControl = well_state_.currentInjectionGroupControl(phase, group.name()); if (currentControl != Group::InjectionCMode::FLD) { const Group& parentGroup = schedule().getGroup(group.parent(), reportStepIdx); - const std::pair changed = wellGroupHelpers::checkGroupConstraintsInj( + const std::pair changed = WellGroupHelpers::checkGroupConstraintsInj( group.name(), group.parent(), parentGroup, @@ -2214,8 +2218,7 @@ namespace Opm { group.getGroupEfficiencyFactor(), schedule(), summaryState, - *rateConverter_, - pvtreg, + resv_coeff, deferred_logger); if (changed.first) { switched_groups.insert(group.name()); @@ -2228,7 +2231,7 @@ namespace Opm { if (!skip && group.isProductionGroup()) { // Obtain rates for group. for (int phasePos = 0; phasePos < phase_usage_.num_phases; ++phasePos) { - const double local_current_rate = wellGroupHelpers::sumWellRates( + const double local_current_rate = WellGroupHelpers::sumWellRates( group, schedule(), well_state_, reportStepIdx, phasePos, /* isInjector */ false); // Sum over all processes rates[phasePos] = -comm.sum(local_current_rate); @@ -2237,7 +2240,7 @@ namespace Opm { const Group::ProductionCMode& currentControl = well_state_.currentProductionGroupControl(group.name()); if (currentControl != Group::ProductionCMode::FLD) { const Group& parentGroup = schedule().getGroup(group.parent(), reportStepIdx); - const std::pair changed = wellGroupHelpers::checkGroupConstraintsProd( + const std::pair changed = WellGroupHelpers::checkGroupConstraintsProd( group.name(), group.parent(), parentGroup, @@ -2249,8 +2252,7 @@ namespace Opm { group.getGroupEfficiencyFactor(), schedule(), summaryState, - *rateConverter_, - pvtreg, + resv_coeff, deferred_logger); if (changed.first) { switched_groups.insert(group.name()); @@ -2283,8 +2285,8 @@ namespace Opm { const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(Phase::GAS, group.name()); if( currentGroupControl == Group::InjectionCMode::REIN ) { int gasPos = phase_usage_.phase_pos[BlackoilPhases::Vapour]; - double gasProductionRate = wellGroupHelpers::sumWellRates(group, schedule, wellState, reportStepIdx, gasPos, /*isInjector*/false); - double solventProductionRate = wellGroupHelpers::sumSolventRates(group, schedule, wellState, reportStepIdx, /*isInjector*/false); + double gasProductionRate = WellGroupHelpers::sumWellRates(group, schedule, wellState, reportStepIdx, gasPos, /*isInjector*/false); + double solventProductionRate = WellGroupHelpers::sumSolventRates(group, schedule, wellState, reportStepIdx, /*isInjector*/false); const auto& comm = ebosSimulator_.vanguard().grid().comm(); solventProductionRate = comm.sum(solventProductionRate); diff --git a/opm/simulators/wells/TargetCalculator.hpp b/opm/simulators/wells/TargetCalculator.hpp new file mode 100644 index 000000000..dc2cf2f5d --- /dev/null +++ b/opm/simulators/wells/TargetCalculator.hpp @@ -0,0 +1,160 @@ +/* + Copyright 2020 Equinor ASA. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + + +#ifndef OPM_TARGETCALCULATOR_HEADER_INCLUDED +#define OPM_TARGETCALCULATOR_HEADER_INCLUDED + +#include +#include +#include +#include +#include +#include + +#include +#include +#include +#include + +namespace Opm +{ + +namespace WellGroupHelpers +{ + + /// Based on a group control mode, extract or calculate rates, and + /// provide other conveniences. + class TargetCalculator + { + public: + TargetCalculator(const Group::ProductionCMode cmode, + const PhaseUsage& pu, + const std::vector& resv_coeff) + : cmode_(cmode) + , pu_(pu) + , resv_coeff_(resv_coeff) + { + } + + template + auto calcModeRateFromRates(const RateVec& rates) const + { + // ElemType is just the plain element type of the rates container, + // without any reference, const or volatile modifiers. + using ElemType = std::remove_cv_t>; + switch (cmode_) { + case Group::ProductionCMode::ORAT: { + assert(pu_.phase_used[BlackoilPhases::Liquid]); + const int pos = pu_.phase_pos[BlackoilPhases::Liquid]; + return rates[pos]; + } + case Group::ProductionCMode::WRAT: { + assert(pu_.phase_used[BlackoilPhases::Aqua]); + const int pos = pu_.phase_pos[BlackoilPhases::Aqua]; + return rates[pos]; + } + case Group::ProductionCMode::GRAT: { + assert(pu_.phase_used[BlackoilPhases::Vapour]); + const int pos = pu_.phase_pos[BlackoilPhases::Vapour]; + return rates[pos]; + } + case Group::ProductionCMode::LRAT: { + assert(pu_.phase_used[BlackoilPhases::Liquid]); + assert(pu_.phase_used[BlackoilPhases::Aqua]); + const int opos = pu_.phase_pos[BlackoilPhases::Liquid]; + const int wpos = pu_.phase_pos[BlackoilPhases::Aqua]; + return rates[opos] + rates[wpos]; + } + case Group::ProductionCMode::RESV: { + assert(pu_.phase_used[BlackoilPhases::Liquid]); + assert(pu_.phase_used[BlackoilPhases::Aqua]); + assert(pu_.phase_used[BlackoilPhases::Vapour]); + ElemType mode_rate = zero(); + for (int phase = 0; phase < pu_.num_phases; ++phase) { + mode_rate += rates[phase] * resv_coeff_[phase]; + } + return mode_rate; + } + default: + // Should never be here. + assert(false); + return zero(); + } + } + + double groupTarget(const Group::ProductionControls ctrl) const + { + switch (cmode_) { + case Group::ProductionCMode::ORAT: + return ctrl.oil_target; + case Group::ProductionCMode::WRAT: + return ctrl.water_target; + case Group::ProductionCMode::GRAT: + return ctrl.gas_target; + case Group::ProductionCMode::LRAT: + return ctrl.liquid_target; + case Group::ProductionCMode::RESV: + return ctrl.resv_target; + default: + // Should never be here. + assert(false); + return 0.0; + } + } + + GuideRateModel::Target guideTargetMode() const + { + switch (cmode_) { + case Group::ProductionCMode::ORAT: + return GuideRateModel::Target::OIL; + case Group::ProductionCMode::WRAT: + return GuideRateModel::Target::WAT; + case Group::ProductionCMode::GRAT: + return GuideRateModel::Target::GAS; + case Group::ProductionCMode::LRAT: + return GuideRateModel::Target::LIQ; + case Group::ProductionCMode::RESV: + return GuideRateModel::Target::RES; + default: + // Should never be here. + assert(false); + return GuideRateModel::Target::NONE; + } + } + + private: + template + static ElemType zero() + { + // This is for Evaluation types. + ElemType x; + x = 0.0; + return x; + } + Group::ProductionCMode cmode_; + const PhaseUsage& pu_; + const std::vector& resv_coeff_; + }; + +} // namespace WellGroupHelpers + +} // namespace Opm + +#endif diff --git a/opm/simulators/wells/WellGroupHelpers.cpp b/opm/simulators/wells/WellGroupHelpers.cpp new file mode 100644 index 000000000..ad4644f46 --- /dev/null +++ b/opm/simulators/wells/WellGroupHelpers.cpp @@ -0,0 +1,1165 @@ +/* + Copyright 2019 Norce. + Copyright 2020 Equinor ASA. + + This file is part of the Open Porous Media project (OPM). + + OPM is free software: you can redistribute it and/or modify + it under the terms of the GNU General Public License as published by + the Free Software Foundation, either version 3 of the License, or + (at your option) any later version. + + OPM is distributed in the hope that it will be useful, + but WITHOUT ANY WARRANTY; without even the implied warranty of + MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + GNU General Public License for more details. + + You should have received a copy of the GNU General Public License + along with OPM. If not, see . +*/ + +#include +#include + +#include +#include + +namespace Opm +{ + + +namespace WellGroupHelpers +{ + + + void setCmodeGroup(const Group& group, + const Schedule& schedule, + const SummaryState& summaryState, + const int reportStepIdx, + WellStateFullyImplicitBlackoil& wellState) + { + + for (const std::string& groupName : group.groups()) { + setCmodeGroup( + schedule.getGroup(groupName, reportStepIdx), schedule, summaryState, reportStepIdx, wellState); + } + + // use NONE as default control + const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS}; + for (Phase phase : all) { + if (!wellState.hasInjectionGroupControl(phase, group.name())) { + wellState.setCurrentInjectionGroupControl(phase, group.name(), Group::InjectionCMode::NONE); + } + } + if (!wellState.hasProductionGroupControl(group.name())) { + wellState.setCurrentProductionGroupControl(group.name(), Group::ProductionCMode::NONE); + } + + if (group.isInjectionGroup() + && schedule.hasWellGroupEvent(group.name(), ScheduleEvents::GROUP_INJECTION_UPDATE, reportStepIdx)) { + + for (Phase phase : all) { + if (!group.hasInjectionControl(phase)) + continue; + + const auto& controls = group.injectionControls(phase, summaryState); + wellState.setCurrentInjectionGroupControl(phase, group.name(), controls.cmode); + } + } + + if (group.isProductionGroup() + && schedule.hasWellGroupEvent(group.name(), ScheduleEvents::GROUP_PRODUCTION_UPDATE, reportStepIdx)) { + const auto controls = group.productionControls(summaryState); + wellState.setCurrentProductionGroupControl(group.name(), controls.cmode); + } + + if (schedule.gConSale(reportStepIdx).has(group.name())) { + wellState.setCurrentInjectionGroupControl(Phase::GAS, group.name(), Group::InjectionCMode::SALE); + } + } + + void accumulateGroupEfficiencyFactor(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + double& factor) + { + factor *= group.getGroupEfficiencyFactor(); + if (group.parent() != "FIELD") + accumulateGroupEfficiencyFactor( + schedule.getGroup(group.parent(), reportStepIdx), schedule, reportStepIdx, factor); + } + + double sumWellPhaseRates(const std::vector& rates, + const Group& group, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const int phasePos, + const bool injector) + { + + double rate = 0.0; + for (const std::string& groupName : group.groups()) { + const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); + rate += groupTmp.getGroupEfficiencyFactor() + * sumWellPhaseRates(rates, groupTmp, schedule, wellState, reportStepIdx, phasePos, injector); + } + const auto& end = wellState.wellMap().end(); + for (const std::string& wellName : group.wells()) { + const auto& it = wellState.wellMap().find(wellName); + if (it == end) // the well is not found + continue; + + int well_index = it->second[0]; + + const auto& wellEcl = schedule.getWell(wellName, reportStepIdx); + // only count producers or injectors + if ((wellEcl.isProducer() && injector) || (wellEcl.isInjector() && !injector)) + continue; + + if (wellEcl.getStatus() == Well::Status::SHUT) + continue; + + double factor = wellEcl.getEfficiencyFactor(); + const auto wellrate_index = well_index * wellState.numPhases(); + if (injector) + rate += factor * rates[wellrate_index + phasePos]; + else + rate -= factor * rates[wellrate_index + phasePos]; + } + return rate; + } + + double sumWellRates(const Group& group, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const int phasePos, + const bool injector) + { + return sumWellPhaseRates(wellState.wellRates(), group, schedule, wellState, reportStepIdx, phasePos, injector); + } + + double sumWellResRates(const Group& group, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const int phasePos, + const bool injector) + { + return sumWellPhaseRates( + wellState.wellReservoirRates(), group, schedule, wellState, reportStepIdx, phasePos, injector); + } + + double sumSolventRates(const Group& group, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const bool injector) + { + + double rate = 0.0; + for (const std::string& groupName : group.groups()) { + const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); + rate += groupTmp.getGroupEfficiencyFactor() + * sumSolventRates(groupTmp, schedule, wellState, reportStepIdx, injector); + } + const auto& end = wellState.wellMap().end(); + for (const std::string& wellName : group.wells()) { + const auto& it = wellState.wellMap().find(wellName); + if (it == end) // the well is not found + continue; + + int well_index = it->second[0]; + + const auto& wellEcl = schedule.getWell(wellName, reportStepIdx); + // only count producers or injectors + if ((wellEcl.isProducer() && injector) || (wellEcl.isInjector() && !injector)) + continue; + + if (wellEcl.getStatus() == Well::Status::SHUT) + continue; + + double factor = wellEcl.getEfficiencyFactor(); + if (injector) + rate += factor * wellState.solventWellRate(well_index); + else + rate -= factor * wellState.solventWellRate(well_index); + } + return rate; + } + + void updateGroupTargetReduction(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const bool isInjector, + const PhaseUsage& pu, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState, + std::vector& groupTargetReduction) + { + const int np = wellState.numPhases(); + for (const std::string& subGroupName : group.groups()) { + std::vector subGroupTargetReduction(np, 0.0); + const Group& subGroup = schedule.getGroup(subGroupName, reportStepIdx); + updateGroupTargetReduction( + subGroup, schedule, reportStepIdx, isInjector, pu, wellStateNupcol, wellState, subGroupTargetReduction); + + // accumulate group contribution from sub group + if (isInjector) { + const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS}; + for (Phase phase : all) { + const Group::InjectionCMode& currentGroupControl + = wellState.currentInjectionGroupControl(phase, subGroupName); + 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; + + if (currentGroupControl != Group::InjectionCMode::FLD + && currentGroupControl != Group::InjectionCMode::NONE) { + // Subgroup is under individual control. + groupTargetReduction[phasePos] + += sumWellRates(subGroup, schedule, wellStateNupcol, reportStepIdx, phasePos, isInjector); + } else { + groupTargetReduction[phasePos] += subGroupTargetReduction[phasePos]; + } + } + } else { + const Group::ProductionCMode& currentGroupControl + = wellState.currentProductionGroupControl(subGroupName); + const bool individual_control = (currentGroupControl != Group::ProductionCMode::FLD + && currentGroupControl != Group::ProductionCMode::NONE); + const int num_group_controlled_wells + = groupControlledWells(schedule, wellStateNupcol, reportStepIdx, subGroupName, ""); + if (individual_control || num_group_controlled_wells == 0) { + for (int phase = 0; phase < np; phase++) { + groupTargetReduction[phase] + += sumWellRates(subGroup, schedule, wellStateNupcol, reportStepIdx, phase, isInjector); + } + } else { + // or accumulate directly from the wells if controled from its parents + for (int phase = 0; phase < np; phase++) { + // groupTargetReduction[phase] += subGroupTargetReduction[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]; + const auto wellrate_index = well_index * wellState.numPhases(); + const double efficiency = wellTmp.getEfficiencyFactor(); + // add contributino from wells not under group control + if (isInjector) { + if (wellState.currentInjectionControls()[well_index] != Well::InjectorCMode::GRUP) + for (int phase = 0; phase < np; phase++) { + groupTargetReduction[phase] += wellStateNupcol.wellRates()[wellrate_index + phase] * efficiency; + } + } else { + if (wellState.currentProductionControls()[well_index] != Well::ProducerCMode::GRUP) + for (int phase = 0; phase < np; phase++) { + groupTargetReduction[phase] -= wellStateNupcol.wellRates()[wellrate_index + phase] * efficiency; + } + } + } + const double groupEfficiency = group.getGroupEfficiencyFactor(); + for (double& elem : groupTargetReduction) { + elem *= groupEfficiency; + } + if (isInjector) + wellState.setCurrentInjectionGroupReductionRates(group.name(), groupTargetReduction); + else + wellState.setCurrentProductionGroupReductionRates(group.name(), 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) + { + 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]; + 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) { // the well is found + + int well_index = it->second[0]; + + 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, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState) + { + for (const std::string& groupName : group.groups()) { + const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); + updateVREPForGroups(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState); + } + const int np = wellState.numPhases(); + double resv = 0.0; + for (int phase = 0; phase < np; ++phase) { + resv += sumWellPhaseRates(wellStateNupcol.wellReservoirRates(), + group, + schedule, + wellState, + reportStepIdx, + phase, + /*isInjector*/ false); + } + wellState.setCurrentInjectionVREPRates(group.name(), resv); + } + + void updateReservoirRatesInjectionGroups(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState) + { + for (const std::string& groupName : group.groups()) { + const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); + updateReservoirRatesInjectionGroups(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState); + } + const int np = wellState.numPhases(); + std::vector resv(np, 0.0); + for (int phase = 0; phase < np; ++phase) { + resv[phase] = sumWellPhaseRates(wellStateNupcol.wellReservoirRates(), + group, + schedule, + wellState, + reportStepIdx, + phase, + /*isInjector*/ true); + } + wellState.setCurrentInjectionGroupReservoirRates(group.name(), resv); + } + + void updateWellRates(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState) + { + for (const std::string& groupName : group.groups()) { + const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); + updateWellRates(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState); + } + const int np = wellState.numPhases(); + const auto& end = wellState.wellMap().end(); + for (const std::string& wellName : group.wells()) { + std::vector rates(np, 0.0); + const auto& it = wellState.wellMap().find(wellName); + if (it != end) { // the well is found on this node + int well_index = it->second[0]; + const auto& wellTmp = schedule.getWell(wellName, reportStepIdx); + int sign = 1; + // production wellRates are negative. The users of currentWellRates uses the convention in + // opm-common that production and injection rates are positive. + if (!wellTmp.isInjector()) + sign = -1; + for (int phase = 0; phase < np; ++phase) { + rates[phase] = sign * wellStateNupcol.wellRates()[well_index * np + phase]; + } + } + wellState.setCurrentWellRates(wellName, rates); + } + } + + void updateGroupProductionRates(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState) + { + for (const std::string& groupName : group.groups()) { + const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); + updateGroupProductionRates(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState); + } + const int np = wellState.numPhases(); + std::vector rates(np, 0.0); + for (int phase = 0; phase < np; ++phase) { + rates[phase] = sumWellPhaseRates( + wellStateNupcol.wellRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false); + } + wellState.setCurrentProductionGroupRates(group.name(), rates); + } + + void updateREINForGroups(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const PhaseUsage& pu, + const SummaryState& st, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState) + { + const int np = wellState.numPhases(); + for (const std::string& groupName : group.groups()) { + const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); + updateREINForGroups(groupTmp, schedule, reportStepIdx, pu, st, wellStateNupcol, wellState); + } + + std::vector rein(np, 0.0); + for (int phase = 0; phase < np; ++phase) { + rein[phase] = sumWellPhaseRates( + wellStateNupcol.wellRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false); + } + + // add import rate and substract consumption rate for group for gas + if (schedule.gConSump(reportStepIdx).has(group.name())) { + const auto& gconsump = schedule.gConSump(reportStepIdx).get(group.name(), st); + if (pu.phase_used[BlackoilPhases::Vapour]) { + rein[pu.phase_pos[BlackoilPhases::Vapour]] += gconsump.import_rate; + rein[pu.phase_pos[BlackoilPhases::Vapour]] -= gconsump.consumption_rate; + } + } + + wellState.setCurrentInjectionREINRates(group.name(), rein); + } + + GuideRate::RateVector + getRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& name) + { + const std::vector& rates = well_state.currentWellRates(name); + double oilRate = 0.0; + if (pu.phase_used[BlackoilPhases::Liquid]) + oilRate = rates[pu.phase_pos[BlackoilPhases::Liquid]]; + + double gasRate = 0.0; + if (pu.phase_used[BlackoilPhases::Vapour]) + gasRate = rates[pu.phase_pos[BlackoilPhases::Vapour]]; + + double waterRate = 0.0; + if (pu.phase_used[BlackoilPhases::Aqua]) + waterRate = rates[pu.phase_pos[BlackoilPhases::Aqua]]; + + return GuideRate::RateVector {oilRate, gasRate, waterRate}; + } + + + double getGuideRate(const std::string& name, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const GuideRate* guideRate, + const GuideRateModel::Target target, + const PhaseUsage& pu) + { + if (schedule.hasWell(name, reportStepIdx) || guideRate->has(name)) { + return guideRate->get(name, target, getRateVector(wellState, pu, name)); + } + + double totalGuideRate = 0.0; + const Group& group = schedule.getGroup(name, reportStepIdx); + + for (const std::string& groupName : group.groups()) { + const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName); + if (currentGroupControl == Group::ProductionCMode::FLD + || currentGroupControl == Group::ProductionCMode::NONE) { + // accumulate from sub wells/groups + totalGuideRate += getGuideRate(groupName, schedule, wellState, reportStepIdx, guideRate, target, pu); + } + } + + for (const std::string& wellName : group.wells()) { + const auto& wellTmp = schedule.getWell(wellName, reportStepIdx); + + if (wellTmp.isInjector()) + continue; + + if (wellTmp.getStatus() == Well::Status::SHUT) + continue; + + // Only count wells under group control or the ru + if (!wellState.isProductionGrup(wellName)) + continue; + + totalGuideRate += guideRate->get(wellName, target, getRateVector(wellState, pu, wellName)); + } + return totalGuideRate; + } + + + double getGuideRateInj(const std::string& name, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const GuideRate* guideRate, + const GuideRateModel::Target target, + const Phase& injectionPhase, + const PhaseUsage& pu) + { + if (schedule.hasWell(name, reportStepIdx)) { + return guideRate->get(name, target, getRateVector(wellState, pu, name)); + } + + double totalGuideRate = 0.0; + const Group& group = schedule.getGroup(name, reportStepIdx); + + for (const std::string& groupName : group.groups()) { + const Group::InjectionCMode& currentGroupControl + = wellState.currentInjectionGroupControl(injectionPhase, groupName); + if (currentGroupControl == Group::InjectionCMode::FLD + || currentGroupControl == Group::InjectionCMode::NONE) { + // accumulate from sub wells/groups + totalGuideRate += getGuideRateInj( + groupName, schedule, wellState, reportStepIdx, guideRate, target, injectionPhase, pu); + } + } + + for (const std::string& wellName : group.wells()) { + const auto& wellTmp = schedule.getWell(wellName, reportStepIdx); + + if (!wellTmp.isInjector()) + continue; + + if (wellTmp.getStatus() == Well::Status::SHUT) + continue; + + // Only count wells under group control or the ru + if (!wellState.isInjectionGrup(wellName)) + continue; + + totalGuideRate += guideRate->get(wellName, target, getRateVector(wellState, pu, wellName)); + } + return totalGuideRate; + } + + + + int groupControlledWells(const Schedule& schedule, + const WellStateFullyImplicitBlackoil& well_state, + const int report_step, + const std::string& group_name, + const std::string& always_included_child) + { + 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); + if (included) { + num_wells + += groupControlledWells(schedule, well_state, report_step, 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); + if (included) { + ++num_wells; + } + } + return num_wells; + } + + + FractionCalculator::FractionCalculator(const Schedule& schedule, + const WellStateFullyImplicitBlackoil& well_state, + const int report_step, + const GuideRate* guide_rate, + const GuideRateModel::Target target, + const PhaseUsage& pu) + : schedule_(schedule) + , well_state_(well_state) + , report_step_(report_step) + , guide_rate_(guide_rate) + , target_(target) + , pu_(pu) + { + } + double FractionCalculator::fraction(const std::string& name, + const std::string& control_group_name, + const bool always_include_this) + { + double fraction = 1.0; + std::string current = name; + while (current != control_group_name) { + fraction *= localFraction(current, always_include_this ? name : ""); + current = parent(current); + } + return fraction; + } + double FractionCalculator::localFraction(const std::string& name, const std::string& always_included_child) + { + const double my_guide_rate = guideRate(name, always_included_child); + const Group& parent_group = schedule_.getGroup(parent(name), report_step_); + const double total_guide_rate = guideRateSum(parent_group, always_included_child); + assert(total_guide_rate >= my_guide_rate); + const double guide_rate_epsilon = 1e-12; + return (total_guide_rate > guide_rate_epsilon) ? my_guide_rate / total_guide_rate : 0.0; + } + std::string FractionCalculator::parent(const std::string& name) + { + if (schedule_.hasWell(name)) { + return schedule_.getWell(name, report_step_).groupName(); + } else { + return schedule_.getGroup(name, report_step_).parent(); + } + } + double FractionCalculator::guideRateSum(const Group& group, const std::string& always_included_child) + { + 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); + 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); + if (included) { + total_guide_rate += guideRate(child_well, always_included_child); + } + } + return total_guide_rate; + } + 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)); + } else { + if (groupControlledWells(name, always_included_child) > 0) { + if (guide_rate_->has(name)) { + return guide_rate_->get(name, target_, getGroupRateVector(name)); + } else { + // We are a group, with default guide rate. + // Compute guide rate by accumulating our children's guide rates. + // (only children not under individual control though). + const Group& group = schedule_.getGroup(name, report_step_); + return guideRateSum(group, always_included_child); + } + } else { + // No group-controlled subordinate wells. + return 0.0; + } + } + } + int FractionCalculator::groupControlledWells(const std::string& group_name, + const std::string& always_included_child) + { + return ::Opm::WellGroupHelpers::groupControlledWells( + schedule_, well_state_, report_step_, group_name, always_included_child); + } + + GuideRate::RateVector FractionCalculator::getGroupRateVector(const std::string& group_name) + { + + std::vector groupRates = well_state_.currentProductionGroupRates(group_name); + double oilRate = 0.0; + if (pu_.phase_used[BlackoilPhases::Liquid]) + oilRate = groupRates[pu_.phase_pos[BlackoilPhases::Liquid]]; + + double gasRate = 0.0; + if (pu_.phase_used[BlackoilPhases::Vapour]) + gasRate = groupRates[pu_.phase_pos[BlackoilPhases::Vapour]]; + + double waterRate = 0.0; + if (pu_.phase_used[BlackoilPhases::Aqua]) + waterRate = groupRates[pu_.phase_pos[BlackoilPhases::Aqua]]; + + return GuideRate::RateVector {oilRate, gasRate, waterRate}; + } + + + 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. + + 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.isAvailableForGroupControl()) { + 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; + 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 * efficiencyFactor)) + / efficiencyFactor; + 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 * efficiencyFactor)) + / efficiencyFactor; + 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 * efficiencyFactor)) + / efficiencyFactor; + 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 * efficiencyFactor)) + / efficiencyFactor; + 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; + double inj_rate = wellState.currentInjectionREINRates(group.name())[phasePos]; + if (schedule.gConSump(reportStepIdx).has(group.name())) { + const auto& gconsump = schedule.gConSump(reportStepIdx).get(group.name(), summaryState); + if (pu.phase_used[BlackoilPhases::Vapour]) { + inj_rate += gconsump.import_rate; + inj_rate -= gconsump.consumption_rate; + } + } + const auto& gconsale = schedule.gConSale(reportStepIdx).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 * efficiencyFactor)) / efficiencyFactor; + 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) + { + // Get initial parent, 'bottom' can be a well or a group. + std::string parent; + if (schedule.hasWell(bottom, report_step)) { + parent = schedule.getWell(bottom, report_step).groupName(); + } else { + parent = schedule.getGroup(bottom, report_step).parent(); + } + + // Build the chain from bottom to top. + std::vector chain; + chain.push_back(bottom); + chain.push_back(parent); + while (parent != top) { + parent = schedule.getGroup(parent, report_step).parent(); + chain.push_back(parent); + } + assert(chain.back() == top); + + // Reverse order and return. + std::reverse(chain.begin(), chain.end()); + return chain; + } + + + + + std::pair checkGroupConstraintsProd(const std::string& name, + const std::string& parent, + const Group& group, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const GuideRate* guideRate, + const double* rates, + 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. + + const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(group.name()); + + if (currentGroupControl == Group::ProductionCMode::FLD || currentGroupControl == Group::ProductionCMode::NONE) { + // Return if we are not available for parent group. + if (!group.isAvailableForGroupControl()) { + return std::make_pair(false, 1); + } + // Otherwise: check production share of parent's control. + const auto& parentGroup = schedule.getGroup(group.parent(), reportStepIdx); + return checkGroupConstraintsProd(name, + parent, + parentGroup, + wellState, + reportStepIdx, + guideRate, + rates, + 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.isProductionGroup()) { + 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. + TargetCalculator tcalc(currentGroupControl, pu, resv_coeff); + FractionCalculator fcalc(schedule, wellState, reportStepIdx, guideRate, tcalc.guideTargetMode(), pu); + + auto localFraction = [&](const std::string& child) { return fcalc.localFraction(child, name); }; + + auto localReduction = [&](const std::string& group_name) { + const std::vector& groupTargetReductions + = wellState.currentProductionGroupReductionRates(group_name); + return tcalc.calcModeRateFromRates(groupTargetReductions); + }; + + const double orig_target = tcalc.groupTarget(group.productionControls(summaryState)); + // 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; + double target = orig_target; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + target -= localReduction(chain[ii]); + if (ii == num_ancestors - 1) { + // Final level. Add my reduction back. + target += current_rate * efficiencyFactor; + } else { + // 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], ""); + if (num_gr_ctrl == 0) { + 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 / efficiencyFactor); + return std::make_pair(current_rate > target_rate, target_rate / current_rate); + } + + +} // namespace WellGroupHelpers + +} // namespace Opm diff --git a/opm/simulators/wells/WellGroupHelpers.hpp b/opm/simulators/wells/WellGroupHelpers.hpp index e68be56cf..7d80b9ad6 100644 --- a/opm/simulators/wells/WellGroupHelpers.hpp +++ b/opm/simulators/wells/WellGroupHelpers.hpp @@ -21,255 +21,104 @@ #ifndef OPM_WELLGROUPHELPERS_HEADER_INCLUDED #define OPM_WELLGROUPHELPERS_HEADER_INCLUDED +#include +#include +#include #include #include -#include +#include #include +#include +#include #include -namespace Opm { +namespace Opm +{ - namespace wellGroupHelpers - { +namespace WellGroupHelpers +{ - int groupControlledWells(const Schedule& schedule, - const WellStateFullyImplicitBlackoil& well_state, - const int report_step, - const std::string& group_name, - const std::string& always_included_child); - inline void setCmodeGroup(const Group& group, const Schedule& schedule, const SummaryState& summaryState, const int reportStepIdx, WellStateFullyImplicitBlackoil& wellState) { + void setCmodeGroup(const Group& group, + const Schedule& schedule, + const SummaryState& summaryState, + const int reportStepIdx, + WellStateFullyImplicitBlackoil& wellState); - for (const std::string& groupName : group.groups()) { - setCmodeGroup( schedule.getGroup(groupName, reportStepIdx), schedule, summaryState, reportStepIdx, wellState); - } + void accumulateGroupEfficiencyFactor(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + double& factor); - // use NONE as default control - const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS}; - for (Phase phase : all) { - if (!wellState.hasInjectionGroupControl(phase, group.name())) { - wellState.setCurrentInjectionGroupControl(phase, group.name(), Group::InjectionCMode::NONE); - } - } - if (!wellState.hasProductionGroupControl(group.name())) { - wellState.setCurrentProductionGroupControl(group.name(), Group::ProductionCMode::NONE); - } + double sumWellPhaseRates(const std::vector& rates, + const Group& group, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const int phasePos, + const bool injector); - if (group.isInjectionGroup() && schedule.hasWellGroupEvent(group.name(), ScheduleEvents::GROUP_INJECTION_UPDATE, reportStepIdx)) { + double sumWellRates(const Group& group, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const int phasePos, + const bool injector); - for (Phase phase : all) { - if (!group.hasInjectionControl(phase)) - continue; + double sumWellResRates(const Group& group, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const int phasePos, + const bool injector); - const auto& controls = group.injectionControls(phase, summaryState); - wellState.setCurrentInjectionGroupControl(phase, group.name(), controls.cmode); - } - } + double sumSolventRates(const Group& group, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const bool injector); - if (group.isProductionGroup() && schedule.hasWellGroupEvent(group.name(), ScheduleEvents::GROUP_PRODUCTION_UPDATE, reportStepIdx)) { - const auto controls = group.productionControls(summaryState); - wellState.setCurrentProductionGroupControl(group.name(), controls.cmode); - } - - if (schedule.gConSale(reportStepIdx).has(group.name())) { - wellState.setCurrentInjectionGroupControl(Phase::GAS, group.name(), Group::InjectionCMode::SALE); - } - } - - inline void accumulateGroupEfficiencyFactor(const Group& group, const Schedule& schedule, const int reportStepIdx, double& factor) { - factor *= group.getGroupEfficiencyFactor(); - if (group.parent() != "FIELD") - accumulateGroupEfficiencyFactor(schedule.getGroup(group.parent(), reportStepIdx), schedule, reportStepIdx, factor); - } - - inline double sumWellPhaseRates(const std::vector& rates, const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const int phasePos, - const bool injector) { - - double rate = 0.0; - for (const std::string& groupName : group.groups()) { - const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); - rate += groupTmp.getGroupEfficiencyFactor()*sumWellPhaseRates(rates, groupTmp, schedule, wellState, reportStepIdx, phasePos, injector); - } - const auto& end = wellState.wellMap().end(); - for (const std::string& wellName : group.wells()) { - const auto& it = wellState.wellMap().find( wellName ); - if (it == end) // the well is not found - continue; - - int well_index = it->second[0]; - - const auto& wellEcl = schedule.getWell(wellName, reportStepIdx); - //only count producers or injectors - if ( (wellEcl.isProducer() && injector) || (wellEcl.isInjector() && !injector)) - continue; - - if (wellEcl.getStatus() == Well::Status::SHUT) - continue; - - double factor = wellEcl.getEfficiencyFactor(); - const auto wellrate_index = well_index * wellState.numPhases(); - if (injector) - rate += factor * rates[ wellrate_index + phasePos]; - else - rate -= factor * rates[ wellrate_index + phasePos]; - } - return rate; - } - - inline double sumWellRates(const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const int phasePos, const bool injector) { - return sumWellPhaseRates(wellState.wellRates(), group, schedule, wellState, reportStepIdx, phasePos, injector); - } - - inline double sumWellResRates(const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const int phasePos, const bool injector) { - return sumWellPhaseRates(wellState.wellReservoirRates(), group, schedule, wellState, reportStepIdx, phasePos, injector); - } - - inline double sumSolventRates(const Group& group, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const bool injector) { - - double rate = 0.0; - for (const std::string& groupName : group.groups()) { - const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); - rate += groupTmp.getGroupEfficiencyFactor()*sumSolventRates(groupTmp, schedule, wellState, reportStepIdx, injector); - } - const auto& end = wellState.wellMap().end(); - for (const std::string& wellName : group.wells()) { - const auto& it = wellState.wellMap().find( wellName ); - if (it == end) // the well is not found - continue; - - int well_index = it->second[0]; - - const auto& wellEcl = schedule.getWell(wellName, reportStepIdx); - //only count producers or injectors - if ( (wellEcl.isProducer() && injector) || (wellEcl.isInjector() && !injector)) - continue; - - if (wellEcl.getStatus() == Well::Status::SHUT) - continue; - - double factor = wellEcl.getEfficiencyFactor(); - if (injector) - rate += factor * wellState.solventWellRate(well_index); - else - rate -= factor * wellState.solventWellRate(well_index); - } - return rate; - } - - inline void updateGroupTargetReduction(const Group& group, const Schedule& schedule, const int reportStepIdx, const bool isInjector, const PhaseUsage& pu, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState, std::vector& groupTargetReduction) - { - const int np = wellState.numPhases(); - for (const std::string& subGroupName : group.groups()) { - std::vector subGroupTargetReduction(np, 0.0); - const Group& subGroup = schedule.getGroup(subGroupName, reportStepIdx); - updateGroupTargetReduction(subGroup, schedule, reportStepIdx, isInjector, pu, wellStateNupcol, wellState, subGroupTargetReduction); - - // accumulate group contribution from sub group - if (isInjector) { - const Phase all[] = {Phase::WATER, Phase::OIL, Phase::GAS}; - for (Phase phase : all) { - const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(phase, subGroupName); - 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; - - if (currentGroupControl != Group::InjectionCMode::FLD && - currentGroupControl != Group::InjectionCMode::NONE) { - // Subgroup is under individual control. - groupTargetReduction[phasePos] += sumWellRates(subGroup, schedule, wellStateNupcol, reportStepIdx, phasePos, isInjector); - } else { - groupTargetReduction[phasePos] += subGroupTargetReduction[phasePos]; - } - } - } else { - const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(subGroupName); - const bool individual_control = (currentGroupControl != Group::ProductionCMode::FLD && - currentGroupControl != Group::ProductionCMode::NONE); - const int num_group_controlled_wells = groupControlledWells(schedule, wellStateNupcol, reportStepIdx, subGroupName, ""); - if (individual_control || num_group_controlled_wells == 0) { - for (int phase = 0; phase < np; phase++) { - groupTargetReduction[phase] += sumWellRates(subGroup, schedule, wellStateNupcol, reportStepIdx, phase, isInjector); - } - } else { - // or accumulate directly from the wells if controled from its parents - for (int phase = 0; phase < np; phase++) { - // groupTargetReduction[phase] += subGroupTargetReduction[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]; - const auto wellrate_index = well_index * wellState.numPhases(); - const double efficiency = wellTmp.getEfficiencyFactor(); - // add contributino from wells not under group control - if (isInjector) { - if (wellState.currentInjectionControls()[well_index] != Well::InjectorCMode::GRUP) - for (int phase = 0; phase < np; phase++) { - groupTargetReduction[phase] += wellStateNupcol.wellRates()[wellrate_index + phase] * efficiency; - } - } else { - if (wellState.currentProductionControls()[well_index] != Well::ProducerCMode::GRUP) - for (int phase = 0; phase < np; phase++) { - groupTargetReduction[phase] -= wellStateNupcol.wellRates()[wellrate_index + phase] * efficiency; - } - } - } - const double groupEfficiency = group.getGroupEfficiencyFactor(); - for (double& elem : groupTargetReduction) { - elem *= groupEfficiency; - } - if (isInjector) - wellState.setCurrentInjectionGroupReductionRates(group.name(), groupTargetReduction); - else - wellState.setCurrentProductionGroupReductionRates(group.name(), groupTargetReduction); - } + void updateGroupTargetReduction(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const bool isInjector, + const PhaseUsage& pu, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState, + std::vector& groupTargetReduction); template - inline 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 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); + 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] ) + 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] ) + else if (phase == Phase::WATER && pu.phase_used[BlackoilPhases::Aqua]) phasePos = pu.phase_pos[BlackoilPhases::Aqua]; else continue; @@ -278,14 +127,14 @@ namespace Opm { } } else { const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName); - if (currentGroupControl != Group::ProductionCMode::FLD && currentGroupControl != Group::ProductionCMode::NONE) { + 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); @@ -299,8 +148,8 @@ namespace Opm { 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 + const auto& it = wellState.wellMap().find(wellName); + if (it == end) // the well is not found continue; int well_index = it->second[0]; @@ -313,15 +162,15 @@ namespace Opm { double oilPot = 0.0; if (pu.phase_used[BlackoilPhases::Liquid]) - oilPot = pot [ pu.phase_pos[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]]; + gasPot = pot[pu.phase_pos[BlackoilPhases::Vapour]]; double waterPot = 0.0; if (pu.phase_used[BlackoilPhases::Aqua]) - waterPot = pot [pu.phase_pos[BlackoilPhases::Aqua]]; + waterPot = pot[pu.phase_pos[BlackoilPhases::Aqua]]; const double gefac = group.getGroupEfficiencyFactor(); @@ -337,7 +186,14 @@ namespace Opm { } template - inline void updateGuideRatesForWells(const Schedule& schedule, const PhaseUsage& pu, const int reportStepIdx, const double& simTime, const WellStateFullyImplicitBlackoil& wellState, const Comm& comm, GuideRate* guideRate) { + 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)) { @@ -345,12 +201,12 @@ namespace Opm { double gaspot = 0.0; double waterpot = 0.0; - const auto& it = wellState.wellMap().find( well.name()); - if (it != end) { // the well is found + const auto& it = wellState.wellMap().find(well.name()); + if (it != end) { // the well is found int well_index = it->second[0]; - const auto wpot = wellState.wellPotentials().data() + well_index*wellState.numPhases(); + const auto wpot = wellState.wellPotentials().data() + well_index * wellState.numPhases(); if (pu.phase_used[BlackoilPhases::Liquid] > 0) oilpot = wpot[pu.phase_pos[BlackoilPhases::Liquid]]; @@ -366,228 +222,69 @@ namespace Opm { waterpot = comm.sum(waterpot) * wefac; guideRate->compute(well.name(), reportStepIdx, simTime, oilpot, gaspot, waterpot); } - } - inline void updateVREPForGroups(const Group& group, const Schedule& schedule, const int reportStepIdx, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState) { - for (const std::string& groupName : group.groups()) { - const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); - updateVREPForGroups(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState); - } - const int np = wellState.numPhases(); - double resv = 0.0; - for (int phase = 0; phase < np; ++phase) { - resv += sumWellPhaseRates(wellStateNupcol.wellReservoirRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false); - } - wellState.setCurrentInjectionVREPRates(group.name(), resv); - } + void updateVREPForGroups(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState); - inline void updateReservoirRatesInjectionGroups(const Group& group, const Schedule& schedule, const int reportStepIdx, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState) { - for (const std::string& groupName : group.groups()) { - const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); - updateReservoirRatesInjectionGroups(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState); - } - const int np = wellState.numPhases(); - std::vector resv(np, 0.0); - for (int phase = 0; phase < np; ++phase) { - resv[phase] = sumWellPhaseRates(wellStateNupcol.wellReservoirRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ true); - } - wellState.setCurrentInjectionGroupReservoirRates(group.name(), resv); - } - - inline void updateWellRates(const Group& group, const Schedule& schedule, const int reportStepIdx, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState) { - for (const std::string& groupName : group.groups()) { - const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); - updateWellRates(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState); - } - const int np = wellState.numPhases(); - const auto& end = wellState.wellMap().end(); - for (const std::string& wellName : group.wells()) { - std::vector rates(np, 0.0); - const auto& it = wellState.wellMap().find(wellName); - if (it != end) { // the well is found on this node - int well_index = it->second[0]; - const auto& wellTmp = schedule.getWell(wellName, reportStepIdx); - int sign = 1; - // production wellRates are negative. The users of currentWellRates uses the convention in - // opm-common that production and injection rates are positive. - if (!wellTmp.isInjector()) - sign = -1; - for (int phase = 0; phase < np; ++phase) { - rates[phase] = sign * wellStateNupcol.wellRates()[well_index*np + phase]; - } - } - wellState.setCurrentWellRates(wellName, rates); - } - } - - inline void updateGroupProductionRates(const Group& group, const Schedule& schedule, const int reportStepIdx, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState) { - for (const std::string& groupName : group.groups()) { - const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); - updateGroupProductionRates(groupTmp, schedule, reportStepIdx, wellStateNupcol, wellState); - } - const int np = wellState.numPhases(); - std::vector rates(np, 0.0); - for (int phase = 0; phase < np; ++phase) { - rates[phase] = sumWellPhaseRates(wellStateNupcol.wellRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false); - } - wellState.setCurrentProductionGroupRates(group.name(), rates); - } + void updateReservoirRatesInjectionGroups(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState); - inline void updateREINForGroups(const Group& group, const Schedule& schedule, const int reportStepIdx, const PhaseUsage& pu, const SummaryState& st, const WellStateFullyImplicitBlackoil& wellStateNupcol, WellStateFullyImplicitBlackoil& wellState) { - const int np = wellState.numPhases(); - for (const std::string& groupName : group.groups()) { - const Group& groupTmp = schedule.getGroup(groupName, reportStepIdx); - updateREINForGroups(groupTmp, schedule, reportStepIdx, pu, st, wellStateNupcol, wellState); - } + void updateWellRates(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState); - std::vector rein(np, 0.0); - for (int phase = 0; phase < np; ++phase) { - rein[phase] = sumWellPhaseRates(wellStateNupcol.wellRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false); - } + void updateGroupProductionRates(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState); - // add import rate and substract consumption rate for group for gas - if (schedule.gConSump(reportStepIdx).has(group.name())) { - const auto& gconsump = schedule.gConSump(reportStepIdx).get(group.name(), st); - if (pu.phase_used[BlackoilPhases::Vapour]) { - rein[pu.phase_pos[BlackoilPhases::Vapour]] += gconsump.import_rate; - rein[pu.phase_pos[BlackoilPhases::Vapour]] -= gconsump.consumption_rate; - } - } + void updateREINForGroups(const Group& group, + const Schedule& schedule, + const int reportStepIdx, + const PhaseUsage& pu, + const SummaryState& st, + const WellStateFullyImplicitBlackoil& wellStateNupcol, + WellStateFullyImplicitBlackoil& wellState); - wellState.setCurrentInjectionREINRates(group.name(), rein); - } - - inline GuideRate::RateVector getRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& name) { - const std::vector& rates = well_state.currentWellRates(name); - double oilRate = 0.0; - if (pu.phase_used[BlackoilPhases::Liquid]) - oilRate = rates[pu.phase_pos[BlackoilPhases::Liquid]]; - - double gasRate = 0.0; - if (pu.phase_used[BlackoilPhases::Vapour]) - gasRate = rates[pu.phase_pos[BlackoilPhases::Vapour]]; - - double waterRate = 0.0; - if (pu.phase_used[BlackoilPhases::Aqua]) - waterRate = rates[pu.phase_pos[BlackoilPhases::Aqua]]; - - return GuideRate::RateVector{oilRate, gasRate, waterRate}; - } + GuideRate::RateVector + getRateVector(const WellStateFullyImplicitBlackoil& well_state, const PhaseUsage& pu, const std::string& name); - inline double getGuideRate(const std::string& name, - const Schedule& schedule, - const WellStateFullyImplicitBlackoil& wellState, - const int reportStepIdx, - const GuideRate* guideRate, - const GuideRateModel::Target target, - const PhaseUsage& pu) - { - if (schedule.hasWell(name, reportStepIdx) || guideRate->has(name)) { - return guideRate->get(name, target, getRateVector(wellState, pu, name)); - } - - double totalGuideRate = 0.0; - const Group& group = schedule.getGroup(name, reportStepIdx); - - for (const std::string& groupName : group.groups()) { - const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName); - if (currentGroupControl == Group::ProductionCMode::FLD || currentGroupControl == Group::ProductionCMode::NONE) { - // accumulate from sub wells/groups - totalGuideRate += getGuideRate(groupName, schedule, wellState, reportStepIdx, guideRate, target, pu); - } - } - - for (const std::string& wellName : group.wells()) { - const auto& wellTmp = schedule.getWell(wellName, reportStepIdx); - - if (wellTmp.isInjector()) - continue; - - if (wellTmp.getStatus() == Well::Status::SHUT) - continue; - - // Only count wells under group control or the ru - if (!wellState.isProductionGrup(wellName)) - continue; - - totalGuideRate += guideRate->get(wellName, target, getRateVector(wellState, pu, wellName)); - } - return totalGuideRate; - } + double getGuideRate(const std::string& name, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const GuideRate* guideRate, + const GuideRateModel::Target target, + const PhaseUsage& pu); - inline double getGuideRateInj(const std::string& name, - const Schedule& schedule, - const WellStateFullyImplicitBlackoil& wellState, - const int reportStepIdx, - const GuideRate* guideRate, - const GuideRateModel::Target target, - const Phase& injectionPhase, - const PhaseUsage& pu) - { - if (schedule.hasWell(name, reportStepIdx)) { - return guideRate->get(name, target, getRateVector(wellState, pu, name)); - } - - double totalGuideRate = 0.0; - const Group& group = schedule.getGroup(name, reportStepIdx); - - for (const std::string& groupName : group.groups()) { - const Group::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(injectionPhase, groupName); - if (currentGroupControl == Group::InjectionCMode::FLD || currentGroupControl == Group::InjectionCMode::NONE) { - // accumulate from sub wells/groups - totalGuideRate += getGuideRateInj(groupName, schedule, wellState, reportStepIdx, guideRate, target, injectionPhase, pu); - } - } - - for (const std::string& wellName : group.wells()) { - const auto& wellTmp = schedule.getWell(wellName, reportStepIdx); - - if (!wellTmp.isInjector()) - continue; - - if (wellTmp.getStatus() == Well::Status::SHUT) - continue; - - // Only count wells under group control or the ru - if (!wellState.isInjectionGrup(wellName)) - continue; - - totalGuideRate += guideRate->get(wellName, target, getRateVector(wellState, pu, wellName)); - } - return totalGuideRate; - } + double getGuideRateInj(const std::string& name, + const Schedule& schedule, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const GuideRate* guideRate, + const GuideRateModel::Target target, + const Phase& injectionPhase, + const PhaseUsage& pu); - - inline int - groupControlledWells(const Schedule& schedule, - const WellStateFullyImplicitBlackoil& well_state, - const int report_step, - const std::string& group_name, - const std::string& always_included_child) - { - 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); - if (included) { - num_wells += groupControlledWells(schedule, well_state, report_step, 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); - if (included) { - ++num_wells; - } - } - return num_wells; - } + int groupControlledWells(const Schedule& schedule, + const WellStateFullyImplicitBlackoil& well_state, + const int report_step, + const std::string& group_name, + const std::string& always_included_child); class FractionCalculator @@ -598,134 +295,16 @@ namespace Opm { const int report_step, const GuideRate* guide_rate, const GuideRateModel::Target target, - const PhaseUsage& pu) - : schedule_(schedule) - , well_state_(well_state) - , report_step_(report_step) - , guide_rate_(guide_rate) - , target_(target) - , pu_(pu) - { - } - double fraction(const std::string& name, - const std::string& control_group_name, - const bool always_include_this) - { - double fraction = 1.0; - std::string current = name; - while (current != control_group_name) { - fraction *= localFraction(current, always_include_this ? name : ""); - current = parent(current); - } - return fraction; - } - double localFraction(const std::string& name, const std::string& always_included_child) - { - const double my_guide_rate = guideRate(name, always_included_child); - const Group& parent_group = schedule_.getGroup(parent(name), report_step_); - const double total_guide_rate = guideRateSum(parent_group, always_included_child); - assert(total_guide_rate >= my_guide_rate); - const double guide_rate_epsilon = 1e-12; - return (total_guide_rate > guide_rate_epsilon) - ? my_guide_rate / total_guide_rate - : 0.0; - } + const PhaseUsage& pu); + 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); + private: - std::string parent(const std::string& name) - { - if (schedule_.hasWell(name)) { - return schedule_.getWell(name, report_step_).groupName(); - } else { - return schedule_.getGroup(name, report_step_).parent(); - } - } - double guideRateSum(const Group& group, const std::string& always_included_child) - { - 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); - 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); - if (included) { - total_guide_rate += guideRate(child_well, always_included_child); - } - } - return total_guide_rate; - } - double 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)); - } else { - if (groupControlledWells(name, always_included_child) > 0) { - if (guide_rate_->has(name)) { - return guide_rate_->get(name, target_, getGroupRateVector(name)); - } else { - // We are a group, with default guide rate. - // Compute guide rate by accumulating our children's guide rates. - // (only children not under individual control though). - const Group& group = schedule_.getGroup(name, report_step_); - return guideRateSum(group, always_included_child); - } - } else { - // No group-controlled subordinate wells. - return 0.0; - } - } - } - int groupControlledWells(const std::string& group_name, const std::string& always_included_child) - { - /* - 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); - if (included) { - num_wells += groupControlledWells(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); - if (included) { - ++num_wells; - } - } - return num_wells; - */ - return ::Opm::wellGroupHelpers::groupControlledWells(schedule_, well_state_, report_step_, group_name, always_included_child); - } - - inline GuideRate::RateVector getGroupRateVector(const std::string& group_name) { - - std::vector groupRates = well_state_.currentProductionGroupRates(group_name); - double oilRate = 0.0; - if (pu_.phase_used[BlackoilPhases::Liquid]) - oilRate = groupRates[pu_.phase_pos[BlackoilPhases::Liquid]]; - - double gasRate = 0.0; - if (pu_.phase_used[BlackoilPhases::Vapour]) - gasRate = groupRates[pu_.phase_pos[BlackoilPhases::Vapour]]; - - double waterRate = 0.0; - if (pu_.phase_used[BlackoilPhases::Aqua]) - waterRate = groupRates[pu_.phase_pos[BlackoilPhases::Aqua]]; - - return GuideRate::RateVector{oilRate, gasRate, waterRate}; - } - - + std::string parent(const std::string& name); + double guideRateSum(const Group& group, const std::string& always_included_child); + double guideRate(const std::string& name, const std::string& always_included_child); + 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 WellStateFullyImplicitBlackoil& well_state_; int report_step_; @@ -734,491 +313,75 @@ namespace Opm { PhaseUsage pu_; }; - - inline 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) - { - FractionCalculator calc(schedule, wellState, reportStepIdx, guideRate, target, pu); - return calc.fraction(name, controlGroupName, alwaysIncludeThis); - } - - inline 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) - { - 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; - } + GuideRate::RateVector getGroupRateVector(const std::string& group_name); - template - inline 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 RateConverterType& rateConverter, - const int pvtRegionIdx, - 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. + 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); - 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.isAvailableForGroupControl()) { - 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, - rateConverter, - pvtRegionIdx, - deferred_logger); - } + 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); - // 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 = wellGroupHelpers::fractionFromInjectionPotentials(name, group.name(), schedule, wellState, reportStepIdx, guideRate, target, pu, injectionPhase, true); - double target_fraction = 1.0; - bool constraint_broken = false; - 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*efficiencyFactor)) / efficiencyFactor; - if (current_rate > target_rate) { - constraint_broken = true; - target_fraction = target_rate / current_rate; - } - break; - } - case Group::InjectionCMode::RESV: - { - std::vector convert_coeff(pu.num_phases, 1.0); - rateConverter.calcCoeff(/*fipreg*/ 0, pvtRegionIdx, convert_coeff); - const double coeff = convert_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*efficiencyFactor)) / efficiencyFactor; - 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*efficiencyFactor)) / efficiencyFactor; - if (current_rate > target_rate) { - constraint_broken = true; - target_fraction = target_rate / current_rate; - } - break; - } - case Group::InjectionCMode::VREP: - { - std::vector convert_coeff(pu.num_phases, 1.0); - rateConverter.calcCoeff(/*fipreg*/ 0, pvtRegionIdx, convert_coeff); - const double coeff = convert_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]]*convert_coeff[pu.phase_pos[BlackoilPhases::Aqua]]; - if (groupcontrols.phase != Phase::OIL) - injReduction += groupInjectionReductions[pu.phase_pos[BlackoilPhases::Liquid]]*convert_coeff[pu.phase_pos[BlackoilPhases::Liquid]]; - if (groupcontrols.phase != Phase::GAS) - injReduction += groupInjectionReductions[pu.phase_pos[BlackoilPhases::Vapour]]*convert_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*efficiencyFactor)) / efficiencyFactor; - 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; - double inj_rate = wellState.currentInjectionREINRates(group.name())[phasePos]; - if (schedule.gConSump(reportStepIdx).has(group.name())) { - const auto& gconsump = schedule.gConSump(reportStepIdx).get(group.name(), summaryState); - if (pu.phase_used[BlackoilPhases::Vapour]) { - inj_rate += gconsump.import_rate; - inj_rate -= gconsump.consumption_rate; - } - } - const auto& gconsale = schedule.gConSale(reportStepIdx).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*efficiencyFactor)) / efficiencyFactor; - 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); - } - - template - class TargetCalculator - { - public: - TargetCalculator(const Group::ProductionCMode cmode, - const PhaseUsage& pu, - const RateConverterType& rate_converter, - const int pvt_region_idx) - : cmode_(cmode) - , pu_(pu) - , rate_converter_(rate_converter) - , pvt_region_idx_(pvt_region_idx) - { - } - template - static ElemType zero() - { - // This is for Evaluation types. - ElemType x; - x = 0.0; - return x; - } - //template <> - //static double zero() - //{ - // return 0.0; - //} - template - auto calcModeRateFromRates(const RateVec& rates) const - { - // ElemType is just the plain element type of the rates container, - // without any reference, const or volatile modifiers. - using ElemType = typename std::remove_cv::type>::type; - switch (cmode_) { - case Group::ProductionCMode::ORAT: { - assert(pu_.phase_used[BlackoilPhases::Liquid]); - const int pos = pu_.phase_pos[BlackoilPhases::Liquid]; - return rates[pos]; - } - case Group::ProductionCMode::WRAT: { - assert(pu_.phase_used[BlackoilPhases::Aqua]); - const int pos = pu_.phase_pos[BlackoilPhases::Aqua]; - return rates[pos]; - } - case Group::ProductionCMode::GRAT: { - assert(pu_.phase_used[BlackoilPhases::Vapour]); - const int pos = pu_.phase_pos[BlackoilPhases::Vapour]; - return rates[pos]; - } - case Group::ProductionCMode::LRAT: { - assert(pu_.phase_used[BlackoilPhases::Liquid]); - assert(pu_.phase_used[BlackoilPhases::Aqua]); - const int opos = pu_.phase_pos[BlackoilPhases::Liquid]; - const int wpos = pu_.phase_pos[BlackoilPhases::Aqua]; - return rates[opos] + rates[wpos]; - } - case Group::ProductionCMode::RESV: { - assert(pu_.phase_used[BlackoilPhases::Liquid]); - assert(pu_.phase_used[BlackoilPhases::Aqua]); - assert(pu_.phase_used[BlackoilPhases::Vapour]); - std::vector convert_coeff(pu_.num_phases, 1.0); - rate_converter_.calcCoeff(/*fipreg*/ 0, pvt_region_idx_, convert_coeff); - ElemType mode_rate = zero(); - for (int phase = 0; phase < pu_.num_phases; ++phase) { - mode_rate += rates[phase] * convert_coeff[phase]; - } - return mode_rate; - } - default: - // Should never be here. - assert(false); - return zero(); - } - } - double groupTarget(const Group::ProductionControls ctrl) const - { - switch (cmode_) { - case Group::ProductionCMode::ORAT: - return ctrl.oil_target; - case Group::ProductionCMode::WRAT: - return ctrl.water_target; - case Group::ProductionCMode::GRAT: - return ctrl.gas_target; - case Group::ProductionCMode::LRAT: - return ctrl.liquid_target; - case Group::ProductionCMode::RESV: - return ctrl.resv_target; - default: - // Should never be here. - assert(false); - return 0.0; - } - } - GuideRateModel::Target guideTargetMode() const - { - switch (cmode_) { - case Group::ProductionCMode::ORAT: - return GuideRateModel::Target::OIL; - case Group::ProductionCMode::WRAT: - return GuideRateModel::Target::WAT; - case Group::ProductionCMode::GRAT: - return GuideRateModel::Target::GAS; - case Group::ProductionCMode::LRAT: - return GuideRateModel::Target::LIQ; - case Group::ProductionCMode::RESV: - return GuideRateModel::Target::RES; - default: - // Should never be here. - assert(false); - return GuideRateModel::Target::NONE; - } - } - private: - Group::ProductionCMode cmode_; - const PhaseUsage& pu_; - const RateConverterType& rate_converter_; - int pvt_region_idx_; - }; + 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); - inline std::vector - groupChainTopBot(const std::string& bottom, const std::string& top, const Schedule& schedule, const int report_step) - { - // Get initial parent, 'bottom' can be a well or a group. - std::string parent; - if (schedule.hasWell(bottom, report_step)) { - parent = schedule.getWell(bottom, report_step).groupName(); - } else { - parent = schedule.getGroup(bottom, report_step).parent(); - } - // Build the chain from bottom to top. - std::vector chain; - chain.push_back(bottom); - chain.push_back(parent); - while (parent != top) { - parent = schedule.getGroup(parent, report_step).parent(); - chain.push_back(parent); - } - assert(chain.back() == top); - // Reverse order and return. - std::reverse(chain.begin(), chain.end()); - return chain; - } + std::vector groupChainTopBot(const std::string& bottom, + const std::string& top, + const Schedule& schedule, + const int report_step); - template - inline std::pair checkGroupConstraintsProd(const std::string& name, - const std::string& parent, - const Group& group, - const WellStateFullyImplicitBlackoil& wellState, - const int reportStepIdx, - const GuideRate* guideRate, - const double* rates, - const PhaseUsage& pu, - const double efficiencyFactor, - const Schedule& schedule, - const SummaryState& summaryState, - const RateConverterType& rateConverter, - const int pvtRegionIdx, - 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. - - const Group::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(group.name()); - - if (currentGroupControl == Group::ProductionCMode::FLD || - currentGroupControl == Group::ProductionCMode::NONE) { - // Return if we are not available for parent group. - if (!group.isAvailableForGroupControl()) { - return std::make_pair(false,1); - } - // Otherwise: check production share of parent's control. - const auto& parentGroup = schedule.getGroup(group.parent(), reportStepIdx); - return checkGroupConstraintsProd(name, - parent, - parentGroup, - wellState, - reportStepIdx, - guideRate, - rates, - pu, - efficiencyFactor * group.getGroupEfficiencyFactor(), - schedule, - summaryState, - rateConverter, - pvtRegionIdx, - deferred_logger); - } - - // This can be false for FLD-controlled groups, we must therefore - // check for FLD first (done above). - if (!group.isProductionGroup()) { - 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. - TargetCalculator tcalc(currentGroupControl, pu, rateConverter, pvtRegionIdx); - FractionCalculator fcalc(schedule, wellState, reportStepIdx, guideRate, tcalc.guideTargetMode(), pu); - - auto localFraction = [&](const std::string& child) { - return fcalc.localFraction(child, name); - }; - - auto localReduction = [&](const std::string& group_name) { - const std::vector& groupTargetReductions = wellState.currentProductionGroupReductionRates(group_name); - return tcalc.calcModeRateFromRates(groupTargetReductions); - }; - - const double orig_target = tcalc.groupTarget(group.productionControls(summaryState)); - // 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; - double target = orig_target; - for (size_t ii = 0; ii < num_ancestors; ++ii) { - target -= localReduction(chain[ii]); - if (ii == num_ancestors - 1) { - // Final level. Add my reduction back. - target += current_rate*efficiencyFactor; - } else { - // 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], ""); - if (num_gr_ctrl == 0) { - 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 / efficiencyFactor); - return std::make_pair(current_rate > target_rate, target_rate / current_rate); - } + std::pair checkGroupConstraintsProd(const std::string& name, + const std::string& parent, + const Group& group, + const WellStateFullyImplicitBlackoil& wellState, + const int reportStepIdx, + const GuideRate* guideRate, + const double* rates, + const PhaseUsage& pu, + const double efficiencyFactor, + const Schedule& schedule, + const SummaryState& summaryState, + const std::vector& resv_coeff, + DeferredLogger& deferred_logger); +} // namespace WellGroupHelpers - } // namespace wellGroupHelpers - -} +} // namespace Opm #endif diff --git a/opm/simulators/wells/WellInterface_impl.hpp b/opm/simulators/wells/WellInterface_impl.hpp index e732a25d9..8dbae9653 100644 --- a/opm/simulators/wells/WellInterface_impl.hpp +++ b/opm/simulators/wells/WellInterface_impl.hpp @@ -21,6 +21,7 @@ #include #include +#include namespace Opm { @@ -1773,8 +1774,13 @@ namespace Opm default: throw("Expected WATER, OIL or GAS as type for injector " + name()); } + + // 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. + // Call check for the well's injection phase. - return wellGroupHelpers::checkGroupConstraintsInj(name(), + return WellGroupHelpers::checkGroupConstraintsInj(name(), well_ecl_.groupName(), group, well_state, @@ -1786,8 +1792,7 @@ namespace Opm efficiencyFactor, schedule, summaryState, - rateConverter_, - pvtRegionIdx_, + resv_coeff, deferred_logger); } @@ -1804,7 +1809,11 @@ namespace Opm const SummaryState& summaryState, DeferredLogger& deferred_logger) const { - return wellGroupHelpers::checkGroupConstraintsProd(name(), + // 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. + + return WellGroupHelpers::checkGroupConstraintsProd(name(), well_ecl_.groupName(), group, well_state, @@ -1815,8 +1824,7 @@ namespace Opm efficiencyFactor, schedule, summaryState, - rateConverter_, - pvtRegionIdx_, + resv_coeff, deferred_logger); } @@ -2096,7 +2104,7 @@ namespace Opm const std::vector& groupInjectionReductions = well_state.currentInjectionGroupReductionRates(group.name()); double groupTargetReduction = groupInjectionReductions[phasePos]; - double fraction = wellGroupHelpers::fractionFromInjectionPotentials(well.name(), + double fraction = WellGroupHelpers::fractionFromInjectionPotentials(well.name(), group.name(), schedule, well_state, @@ -2242,8 +2250,13 @@ namespace Opm // 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. - wellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, rateConverter_, pvtRegionIdx_); - wellGroupHelpers::FractionCalculator fcalc(schedule, well_state, current_step_, guide_rate_, tcalc.guideTargetMode(), pu); + + // 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. + + WellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, resv_coeff); + WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, current_step_, guide_rate_, tcalc.guideTargetMode(), pu); auto localFraction = [&](const std::string& child) { return fcalc.localFraction(child, ""); @@ -2255,7 +2268,7 @@ namespace Opm }; const double orig_target = tcalc.groupTarget(group.productionControls(summaryState)); - const auto chain = wellGroupHelpers::groupChainTopBot(name(), group.name(), schedule, current_step_); + 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;