From 12fa7a4ac8a8607e80ef01a7d28a05e9c7bb489c Mon Sep 17 00:00:00 2001 From: Kai Bao Date: Mon, 28 Jun 2021 23:12:18 +0200 Subject: [PATCH] putting all the guide rate update function in one single function to make sure we only do once checking of the expiration of the guide rates. --- .../wells/BlackoilWellModel_impl.hpp | 6 ++-- opm/simulators/wells/WellGroupHelpers.cpp | 36 ++++++++++++++++++- opm/simulators/wells/WellGroupHelpers.hpp | 14 ++++++++ 3 files changed, 51 insertions(+), 5 deletions(-) diff --git a/opm/simulators/wells/BlackoilWellModel_impl.hpp b/opm/simulators/wells/BlackoilWellModel_impl.hpp index 79f97cc52..f317627f1 100644 --- a/opm/simulators/wells/BlackoilWellModel_impl.hpp +++ b/opm/simulators/wells/BlackoilWellModel_impl.hpp @@ -322,10 +322,8 @@ namespace Opm { const auto& summaryState = ebosSimulator_.vanguard().summaryState(); std::vector pot(numPhases(), 0.0); const Group& fieldGroup = schedule().getGroup("FIELD", reportStepIdx); - WellGroupHelpers::updateGuideRateForProductionGroups(fieldGroup, schedule(), phase_usage_, reportStepIdx, simulationTime, this->wellState(), this->groupState(), comm, &guideRate_, pot); - WellGroupHelpers::updateGuideRatesForInjectionGroups(fieldGroup, schedule(), summaryState, phase_usage_, reportStepIdx, this->wellState(), this->groupState(), &guideRate_, local_deferredLogger); - WellGroupHelpers::updateGuideRatesForWells(schedule(), phase_usage_, reportStepIdx, simulationTime, this->wellState(), comm, &guideRate_); - + WellGroupHelpers::updateGuideRates(fieldGroup, schedule(), summaryState, this->phase_usage_, reportStepIdx, simulationTime, + this->wellState(), this->groupState(), comm, &this->guideRate_, pot, local_deferredLogger); try { // Compute initial well solution for new wells and injectors that change injection type i.e. WAG. for (auto& well : well_container_) { diff --git a/opm/simulators/wells/WellGroupHelpers.cpp b/opm/simulators/wells/WellGroupHelpers.cpp index 481f5251f..fbf7c06f5 100644 --- a/opm/simulators/wells/WellGroupHelpers.cpp +++ b/opm/simulators/wells/WellGroupHelpers.cpp @@ -1306,6 +1306,27 @@ namespace WellGroupHelpers return std::make_pair(current_rate > target_rate, scale); } + + template + void updateGuideRates(const Group& group, + const Schedule& schedule, + const SummaryState& summary_state, + const PhaseUsage& pu, + const int report_step, + const double sim_time, + WellState& well_state, + const GroupState& group_state, + const Comm& comm, + GuideRate* guide_rate, + std::vector& pot, + Opm::DeferredLogger& deferred_logger) + { + guide_rate->updateGuideRateExpiration(sim_time, report_step); + updateGuideRateForProductionGroups(group, schedule, pu, report_step, sim_time, well_state, group_state, comm, guide_rate, pot); + updateGuideRatesForInjectionGroups(group, schedule, summary_state, pu, report_step, well_state, group_state, guide_rate, deferred_logger); + updateGuideRatesForWells(schedule, pu, report_step, sim_time, well_state, comm, guide_rate); + } + template void updateGuideRateForProductionGroups(const Group& group, const Schedule& schedule, @@ -1454,7 +1475,20 @@ namespace WellGroupHelpers const double& simTime, \ const WellState& wellState, \ const Dune::CollectiveCommunication<__VA_ARGS__>& comm, \ - GuideRate* guideRate); + GuideRate* guideRate); \ + template \ + void updateGuideRates>(const Group& group, \ + const Schedule& schedule, \ + const SummaryState& summary_state, \ + const PhaseUsage& pu, \ + const int report_step, \ + const double sim_time, \ + WellState& well_state, \ + const GroupState& group_state, \ + const Dune::CollectiveCommunication<__VA_ARGS__>& comm,\ + GuideRate* guide_rate, \ + std::vector& pot,\ + Opm::DeferredLogger& deferred_logger); #if HAVE_MPI INSTANCE_WELLGROUP_HELPERS(MPI_Comm) diff --git a/opm/simulators/wells/WellGroupHelpers.hpp b/opm/simulators/wells/WellGroupHelpers.hpp index a8f958857..bbe9a3f3b 100644 --- a/opm/simulators/wells/WellGroupHelpers.hpp +++ b/opm/simulators/wells/WellGroupHelpers.hpp @@ -100,6 +100,20 @@ namespace WellGroupHelpers GroupState& group_state, std::vector& groupTargetReduction); + template + void updateGuideRates(const Group& group, + const Schedule& schedule, + const SummaryState& summary_state, + const PhaseUsage& pu, + int report_step, + double sim_time, + WellState& well_state, + const GroupState& group_state, + const Comm& comm, + GuideRate* guide_rate, + std::vector& pot, + Opm::DeferredLogger& deferred_logge); + template void updateGuideRateForProductionGroups(const Group& group, const Schedule& schedule,