diff --git a/CMakeLists_files.cmake b/CMakeLists_files.cmake index 0fffd3008..c7f2e1cab 100644 --- a/CMakeLists_files.cmake +++ b/CMakeLists_files.cmake @@ -104,6 +104,7 @@ list (APPEND MAIN_SOURCE_FILES opm/simulators/wells/VFPProdProperties.cpp opm/simulators/wells/VFPInjProperties.cpp opm/simulators/wells/WellConvergence.cpp + opm/simulators/wells/WellGroupControls.cpp opm/simulators/wells/WellGroupHelpers.cpp opm/simulators/wells/WellHelpers.cpp opm/simulators/wells/WellInterfaceEval.cpp @@ -385,6 +386,7 @@ list (APPEND PUBLIC_HEADER_FILES opm/simulators/wells/VFPProperties.hpp opm/simulators/wells/WellConnectionAuxiliaryModule.hpp opm/simulators/wells/WellConvergence.hpp + opm/simulators/wells/WellGroupControls.hpp opm/simulators/wells/WellGroupHelpers.hpp opm/simulators/wells/WellHelpers.hpp opm/simulators/wells/WellInterface.hpp diff --git a/opm/simulators/wells/WellGroupControls.cpp b/opm/simulators/wells/WellGroupControls.cpp new file mode 100644 index 000000000..a151b6907 --- /dev/null +++ b/opm/simulators/wells/WellGroupControls.cpp @@ -0,0 +1,525 @@ +/* + Copyright 2017 SINTEF Digital, Mathematics and Cybernetics. + Copyright 2017 Statoil ASA. + Copyright 2018 IRIS + + 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 +#include +#include + +#include + +#include +#include +#include +#include +#include + +#include +#include + +namespace Opm +{ + +template +void WellGroupControls:: +getGroupInjectionControl(const Group& group, + const WellState& well_state, + const GroupState& group_state, + const Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + const EvalWell& bhp, + const EvalWell& injection_rate, + const RateConvFunc& rateConverter, + double efficiencyFactor, + EvalWell& control_eq, + DeferredLogger& deferred_logger) const +{ + // Setting some defaults to silence warnings below. + // Will be overwritten in the switch statement. + Phase injectionPhase = Phase::WATER; + switch (injectorType) { + case InjectorType::WATER: + { + injectionPhase = Phase::WATER; + break; + } + case InjectorType::OIL: + { + injectionPhase = Phase::OIL; + break; + } + case InjectorType::GAS: + { + injectionPhase = Phase::GAS; + break; + } + default: + // Should not be here. + assert(false); + } + + auto currentGroupControl = group_state.injection_control(group.name(), injectionPhase); + if (currentGroupControl == Group::InjectionCMode::FLD || + currentGroupControl == Group::InjectionCMode::NONE) { + if (!group.injectionGroupControlAvailable(injectionPhase)) { + // We cannot go any further up the hierarchy. This could + // be the FIELD group, or any group for which this has + // been set in GCONINJE or GCONPROD. If we are here + // anyway, it is likely that the deck set inconsistent + // requirements, such as GRUP control mode on a well with + // no appropriate controls defined on any of its + // containing groups. We will therefore use the wells' bhp + // limit equation as a fallback. + const auto& controls = well_.wellEcl().injectionControls(summaryState); + control_eq = bhp - controls.bhp_limit; + return; + } else { + // Inject share of parents control + const auto& parent = schedule.getGroup(group.parent(), well_.currentStep()); + efficiencyFactor *= group.getGroupEfficiencyFactor(); + getGroupInjectionControl(parent, well_state, + group_state, schedule, + summaryState, injectorType, + bhp, injection_rate, + rateConverter, + efficiencyFactor, + control_eq, + deferred_logger); + return; + } + } + + const auto& well = well_.wellEcl(); + const auto pu = well_.phaseUsage(); + + if (!group.isInjectionGroup()) { + // use bhp as control eq and let the updateControl code find a valid control + const auto& controls = well.injectionControls(summaryState); + control_eq = bhp - controls.bhp_limit; + return; + } + + // 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. + + // Make conversion factors for RESV <-> surface rates. + std::vector resv_coeff(pu.num_phases, 1.0); + rateConverter(0, well_.pvtRegionIdx(), resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS. + + double sales_target = 0; + if (schedule[well_.currentStep()].gconsale().has(group.name())) { + const auto& gconsale = schedule[well_.currentStep()].gconsale().get(group.name(), summaryState); + sales_target = gconsale.sales_target; + } + WellGroupHelpers::InjectionTargetCalculator tcalc(currentGroupControl, pu, + resv_coeff, group.name(), + sales_target, group_state, + injectionPhase, + group.has_gpmaint_control(injectionPhase, currentGroupControl), + deferred_logger); + WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, + group_state, well_.currentStep(), + well_.guideRate(), + tcalc.guideTargetMode(), + pu, false, injectionPhase); + + auto localFraction = [&](const std::string& child) { + return fcalc.localFraction(child, child); + }; + + auto localReduction = [&](const std::string& group_name) { + const std::vector& groupTargetReductions = group_state.injection_reduction_rates(group_name); + return tcalc.calcModeRateFromRates(groupTargetReductions); + }; + + const double orig_target = tcalc.groupTarget(group.injectionControls(injectionPhase, + summaryState), + deferred_logger); + const auto chain = WellGroupHelpers::groupChainTopBot(well_.name(), group.name(), + schedule, well_.currentStep()); + // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. + const size_t num_ancestors = chain.size() - 1; + double target = orig_target; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + if ((ii == 0) || well_.guideRate()->has(chain[ii], injectionPhase)) { + // Apply local reductions only at the control level + // (top) and for levels where we have a specified + // group guide rate. + target -= localReduction(chain[ii]); + } + target *= localFraction(chain[ii+1]); + } + // Avoid negative target rates coming from too large local reductions. + const double target_rate = std::max(0.0, target / efficiencyFactor); + const auto current_rate = injection_rate; + + control_eq = current_rate - target_rate; +} + +std::optional +WellGroupControls:: +getGroupInjectionTargetRate(const Group& group, + const WellState& well_state, + const GroupState& group_state, + const Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + const RateConvFunc& rateConverter, + double efficiencyFactor, + DeferredLogger& deferred_logger) const +{ + // Setting some defaults to silence warnings below. + // Will be overwritten in the switch statement. + Phase injectionPhase = Phase::WATER; + switch (injectorType) { + case InjectorType::WATER: + { + injectionPhase = Phase::WATER; + break; + } + case InjectorType::OIL: + { + injectionPhase = Phase::OIL; + break; + } + case InjectorType::GAS: + { + injectionPhase = Phase::GAS; + break; + } + default: + // Should not be here. + assert(false); + } + + auto currentGroupControl = group_state.injection_control(group.name(), injectionPhase); + if (currentGroupControl == Group::InjectionCMode::FLD || + currentGroupControl == Group::InjectionCMode::NONE) { + if (!group.injectionGroupControlAvailable(injectionPhase)) { + // We cannot go any further up the hierarchy. This could + // be the FIELD group, or any group for which this has + // been set in GCONINJE or GCONPROD. If we are here + // anyway, it is likely that the deck set inconsistent + // requirements, such as GRUP control mode on a well with + // no appropriate controls defined on any of its + // containing groups. We will therefore use the wells' bhp + // limit equation as a fallback. + return std::nullopt; + } else { + // Inject share of parents control + const auto& parent = schedule.getGroup( group.parent(), well_.currentStep()); + efficiencyFactor *= group.getGroupEfficiencyFactor(); + return getGroupInjectionTargetRate(parent, well_state, group_state, + schedule, summaryState, injectorType, + rateConverter, efficiencyFactor, deferred_logger); + } + } + + const auto pu = well_.phaseUsage(); + + if (!group.isInjectionGroup()) { + return std::nullopt; + } + + // 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. + + // Make conversion factors for RESV <-> surface rates. + std::vector resv_coeff(pu.num_phases, 1.0); + rateConverter(0, well_.pvtRegionIdx(), resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS. + + double sales_target = 0; + if (schedule[well_.currentStep()].gconsale().has(group.name())) { + const auto& gconsale = schedule[well_.currentStep()].gconsale().get(group.name(), summaryState); + sales_target = gconsale.sales_target; + } + WellGroupHelpers::InjectionTargetCalculator tcalc(currentGroupControl, pu, resv_coeff, + group.name(), sales_target, group_state, + injectionPhase, + group.has_gpmaint_control(injectionPhase, currentGroupControl), + deferred_logger); + WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, group_state, + well_.currentStep(), well_.guideRate(), + tcalc.guideTargetMode(), pu, false, injectionPhase); + + auto localFraction = [&](const std::string& child) { + return fcalc.localFraction(child, child); //Note child needs to be passed to always include since the global isGrup map is not updated yet. + }; + + auto localReduction = [&](const std::string& group_name) { + const std::vector& groupTargetReductions = group_state.injection_reduction_rates(group_name); + return tcalc.calcModeRateFromRates(groupTargetReductions); + }; + + const double orig_target = tcalc.groupTarget(group.injectionControls(injectionPhase, summaryState), deferred_logger); + const auto chain = WellGroupHelpers::groupChainTopBot(well_.name(), group.name(), schedule, well_.currentStep()); + // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. + const size_t num_ancestors = chain.size() - 1; + double target = orig_target; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + if ((ii == 0) || well_.guideRate()->has(chain[ii], injectionPhase)) { + // Apply local reductions only at the control level + // (top) and for levels where we have a specified + // group guide rate. + target -= localReduction(chain[ii]); + } + target *= localFraction(chain[ii+1]); + } + return std::max(0.0, target / efficiencyFactor); +} + +template +void WellGroupControls::getGroupProductionControl(const Group& group, + const WellState& well_state, + const GroupState& group_state, + const Schedule& schedule, + const SummaryState& summaryState, + const EvalWell& bhp, + const std::vector& rates, + const RateConvFunc& rateConverter, + double efficiencyFactor, + EvalWell& control_eq) const +{ + const Group::ProductionCMode& currentGroupControl = group_state.production_control(group.name()); + if (currentGroupControl == Group::ProductionCMode::FLD || + currentGroupControl == Group::ProductionCMode::NONE) { + if (!group.productionGroupControlAvailable()) { + // We cannot go any further up the hierarchy. This could + // be the FIELD group, or any group for which this has + // been set in GCONINJE or GCONPROD. If we are here + // anyway, it is likely that the deck set inconsistent + // requirements, such as GRUP control mode on a well with + // no appropriate controls defined on any of its + // containing groups. We will therefore use the wells' bhp + // limit equation as a fallback. + const auto& controls = well_.wellEcl().productionControls(summaryState); + control_eq = bhp - controls.bhp_limit; + return; + } else { + // Produce share of parents control + const auto& parent = schedule.getGroup(group.parent(), well_.currentStep()); + efficiencyFactor *= group.getGroupEfficiencyFactor(); + getGroupProductionControl(parent, well_state, group_state, + schedule, summaryState, bhp, + rates, rateConverter, + efficiencyFactor, control_eq); + return; + } + } + + const auto& well = well_.wellEcl(); + const auto pu = well_.phaseUsage(); + + if (!group.isProductionGroup()) { + // use bhp as control eq and let the updateControl code find a valid control + const auto& controls = well.productionControls(summaryState); + control_eq = bhp - controls.bhp_limit; + return; + } + + // 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. + + // Make conversion factors for RESV <-> surface rates. + std::vector resv_coeff(well_.phaseUsage().num_phases, 1.0); + rateConverter(0, well_.pvtRegionIdx(), resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS. + + // gconsale may adjust the grat target. + // the adjusted rates is send to the targetCalculator + double gratTargetFromSales = 0.0; + if (group_state.has_grat_sales_target(group.name())) + gratTargetFromSales = group_state.grat_sales_target(group.name()); + + WellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, resv_coeff, + gratTargetFromSales, group.name(), + group_state, + group.has_gpmaint_control(currentGroupControl)); + WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, group_state, + well_.currentStep(), + well_.guideRate(), + tcalc.guideTargetMode(), + pu, true, Phase::OIL); + + auto localFraction = [&](const std::string& child) { + return fcalc.localFraction(child, child); + }; + + auto localReduction = [&](const std::string& group_name) { + const std::vector& groupTargetReductions = group_state.production_reduction_rates(group_name); + return tcalc.calcModeRateFromRates(groupTargetReductions); + }; + + const double orig_target = tcalc.groupTarget(group.productionControls(summaryState)); + const auto chain = WellGroupHelpers::groupChainTopBot(well_.name(), group.name(), + schedule, well_.currentStep()); + // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. + const size_t num_ancestors = chain.size() - 1; + double target = orig_target; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + if ((ii == 0) || well_.guideRate()->has(chain[ii])) { + // Apply local reductions only at the control level + // (top) and for levels where we have a specified + // group guide rate. + target -= localReduction(chain[ii]); + } + target *= localFraction(chain[ii+1]); + } + // Avoid negative target rates coming from too large local reductions. + const double target_rate = std::max(0.0, target / efficiencyFactor); + const auto current_rate = -tcalc.calcModeRateFromRates(rates); // Switch sign since 'rates' are negative for producers. + control_eq = current_rate - target_rate; +} + +double WellGroupControls:: +getGroupProductionTargetRate(const Group& group, + const WellState& well_state, + const GroupState& group_state, + const Schedule& schedule, + const SummaryState& summaryState, + const RateConvFunc& rateConverter, + double efficiencyFactor) const +{ + const Group::ProductionCMode& currentGroupControl = group_state.production_control(group.name()); + if (currentGroupControl == Group::ProductionCMode::FLD || + currentGroupControl == Group::ProductionCMode::NONE) { + if (!group.productionGroupControlAvailable()) { + return 1.0; + } else { + // Produce share of parents control + const auto& parent = schedule.getGroup(group.parent(), well_.currentStep()); + efficiencyFactor *= group.getGroupEfficiencyFactor(); + return getGroupProductionTargetRate(parent, well_state, group_state, + schedule, summaryState, + rateConverter, efficiencyFactor); + } + } + + const auto pu = well_.phaseUsage(); + + if (!group.isProductionGroup()) { + return 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. + + // Make conversion factors for RESV <-> surface rates. + std::vector resv_coeff(well_.phaseUsage().num_phases, 1.0); + rateConverter(0, well_.pvtRegionIdx(), resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS. + + // gconsale may adjust the grat target. + // the adjusted rates is send to the targetCalculator + double gratTargetFromSales = 0.0; + if (group_state.has_grat_sales_target(group.name())) + gratTargetFromSales = group_state.grat_sales_target(group.name()); + + WellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, resv_coeff, gratTargetFromSales, group.name(), group_state, group.has_gpmaint_control(currentGroupControl)); + WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, group_state, + well_.currentStep(), + well_.guideRate(), + tcalc.guideTargetMode(), + pu, true, Phase::OIL); + + auto localFraction = [&](const std::string& child) { + return fcalc.localFraction(child, child); //Note child needs to be passed to always include since the global isGrup map is not updated yet. + }; + + auto localReduction = [&](const std::string& group_name) { + const std::vector& groupTargetReductions = group_state.production_reduction_rates(group_name); + return tcalc.calcModeRateFromRates(groupTargetReductions); + }; + + const double orig_target = tcalc.groupTarget(group.productionControls(summaryState)); + const auto chain = WellGroupHelpers::groupChainTopBot(well_.name(), group.name(), + schedule, well_.currentStep()); + // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. + const size_t num_ancestors = chain.size() - 1; + double target = orig_target; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + if ((ii == 0) || well_.guideRate()->has(chain[ii])) { + // Apply local reductions only at the control level + // (top) and for levels where we have a specified + // group guide rate. + target -= localReduction(chain[ii]); + } + target *= localFraction(chain[ii+1]); + } + // Avoid negative target rates coming from too large local reductions. + const double target_rate = std::max(0.0, target / efficiencyFactor); + const auto& ws = well_state.well(well_.indexOfWell()); + const auto& rates = ws.surface_rates; + const auto current_rate = -tcalc.calcModeRateFromRates(rates); // Switch sign since 'rates' are negative for producers. + double scale = 1.0; + if (target_rate == 0.0) { + return 0.0; + } + + if (current_rate > 1e-14) + scale = target_rate/current_rate; + return scale; +} + +#define INSTANCE(...) \ +template void WellGroupControls:: \ +getGroupInjectionControl<__VA_ARGS__>(const Group&, \ + const WellState&, \ + const GroupState&, \ + const Schedule&, \ + const SummaryState&, \ + const InjectorType&, \ + const __VA_ARGS__& bhp, \ + const __VA_ARGS__& injection_rate, \ + const RateConvFunc& rateConverter, \ + double efficiencyFactor, \ + __VA_ARGS__& control_eq, \ + DeferredLogger& deferred_logger) const; \ +template void WellGroupControls:: \ +getGroupProductionControl<__VA_ARGS__>(const Group&, \ + const WellState&, \ + const GroupState&, \ + const Schedule&, \ + const SummaryState&, \ + const __VA_ARGS__& bhp, \ + const std::vector<__VA_ARGS__>&, \ + const RateConvFunc& rateConverter, \ + double efficiencyFactor, \ + __VA_ARGS__& control_eq) const; + +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) +INSTANCE(DenseAd::Evaluation) + +} // namespace Opm diff --git a/opm/simulators/wells/WellGroupControls.hpp b/opm/simulators/wells/WellGroupControls.hpp new file mode 100644 index 000000000..23a351826 --- /dev/null +++ b/opm/simulators/wells/WellGroupControls.hpp @@ -0,0 +1,103 @@ +/* + Copyright 2017 SINTEF Digital, Mathematics and Cybernetics. + Copyright 2017 Statoil ASA. + Copyright 2017 IRIS + Copyright 2019 Norce + + 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_WELL_GROUP_CONTROLS_HEADER_INCLUDED +#define OPM_WELL_GROUP_CONTROLS_HEADER_INCLUDED + +#include +#include +#include + +namespace Opm +{ + +class DeferredLogger; +class Group; +class GroupState; +enum class InjectorType; +using RegionId = int; +class Schedule; +class SummaryState; +class WellInterfaceGeneric; +class WellState; + +//! \brief Class for computing well group controls. +class WellGroupControls { +public: + //! \brief Constructor sets reference to well. + WellGroupControls(const WellInterfaceGeneric& well) : well_(well) {} + + using RateConvFunc = std::function&)>; + + template + void getGroupInjectionControl(const Group& group, + const WellState& well_state, + const GroupState& group_state, + const Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + const EvalWell& bhp, + const EvalWell& injection_rate, + const RateConvFunc& rateConverter, + double efficiencyFactor, + EvalWell& control_eq, + DeferredLogger& deferred_logger) const; + + std::optional + getGroupInjectionTargetRate(const Group& group, + const WellState& well_state, + const GroupState& group_state, + const Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + const RateConvFunc& rateConverter, + double efficiencyFactor, + DeferredLogger& deferred_logger) const; + + template + void getGroupProductionControl(const Group& group, + const WellState& well_state, + const GroupState& group_state, + const Schedule& schedule, + const SummaryState& summaryState, + const EvalWell& bhp, + const std::vector& rates, + const RateConvFunc& rateConverter, + double efficiencyFactor, + EvalWell& control_eq) const; + + double getGroupProductionTargetRate(const Group& group, + const WellState& well_state, + const GroupState& group_state, + const Schedule& schedule, + const SummaryState& summaryState, + const RateConvFunc& rateConverter, + double efficiencyFactor) const; + +private: + const WellInterfaceGeneric& well_; //!< Reference to well interface +}; + +} + +#endif // OPM_WELL_GROUP_CONTROLS_HEADER_INCLUDED diff --git a/opm/simulators/wells/WellInterfaceEval.cpp b/opm/simulators/wells/WellInterfaceEval.cpp index 0117fc8a8..7be74d96a 100644 --- a/opm/simulators/wells/WellInterfaceEval.cpp +++ b/opm/simulators/wells/WellInterfaceEval.cpp @@ -30,9 +30,8 @@ #include #include #include -#include #include -#include +#include #include #include #include @@ -50,218 +49,6 @@ WellInterfaceEval(const WellInterfaceFluidSystem& baseif) : baseif_(baseif) {} -template -template -void -WellInterfaceEval:: -getGroupInjectionControl(const Group& group, - const WellState& well_state, - const GroupState& group_state, - const Schedule& schedule, - const SummaryState& summaryState, - const InjectorType& injectorType, - const EvalWell& bhp, - const EvalWell& injection_rate, - EvalWell& control_eq, - double efficiencyFactor, - DeferredLogger& deferred_logger) const -{ - // Setting some defaults to silence warnings below. - // Will be overwritten in the switch statement. - Phase injectionPhase = Phase::WATER; - switch (injectorType) { - case InjectorType::WATER: - { - injectionPhase = Phase::WATER; - break; - } - case InjectorType::OIL: - { - injectionPhase = Phase::OIL; - break; - } - case InjectorType::GAS: - { - injectionPhase = Phase::GAS; - break; - } - default: - // Should not be here. - assert(false); - } - - auto currentGroupControl = group_state.injection_control(group.name(), injectionPhase); - if (currentGroupControl == Group::InjectionCMode::FLD || - currentGroupControl == Group::InjectionCMode::NONE) { - if (!group.injectionGroupControlAvailable(injectionPhase)) { - // We cannot go any further up the hierarchy. This could - // be the FIELD group, or any group for which this has - // been set in GCONINJE or GCONPROD. If we are here - // anyway, it is likely that the deck set inconsistent - // requirements, such as GRUP control mode on a well with - // no appropriate controls defined on any of its - // containing groups. We will therefore use the wells' bhp - // limit equation as a fallback. - const auto& controls = baseif_.wellEcl().injectionControls(summaryState); - control_eq = bhp - controls.bhp_limit; - return; - } else { - // Inject share of parents control - const auto& parent = schedule.getGroup( group.parent(), baseif_.currentStep()); - efficiencyFactor *= group.getGroupEfficiencyFactor(); - getGroupInjectionControl(parent, well_state, group_state, schedule, summaryState, injectorType, bhp, injection_rate, control_eq, efficiencyFactor, deferred_logger); - return; - } - } - - const auto& well = baseif_.wellEcl(); - const auto pu = baseif_.phaseUsage(); - - if (!group.isInjectionGroup()) { - // use bhp as control eq and let the updateControl code find a valid control - const auto& controls = well.injectionControls(summaryState); - control_eq = bhp - controls.bhp_limit; - return; - } - - // 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. - - // Make conversion factors for RESV <-> surface rates. - std::vector resv_coeff(pu.num_phases, 1.0); - baseif_.rateConverter().calcInjCoeff(0, baseif_.pvtRegionIdx(), resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS. - - double sales_target = 0; - if (schedule[baseif_.currentStep()].gconsale().has(group.name())) { - const auto& gconsale = schedule[baseif_.currentStep()].gconsale().get(group.name(), summaryState); - sales_target = gconsale.sales_target; - } - WellGroupHelpers::InjectionTargetCalculator tcalc(currentGroupControl, pu, resv_coeff, group.name(), sales_target, group_state, injectionPhase, group.has_gpmaint_control(injectionPhase, currentGroupControl), deferred_logger); - WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, group_state, baseif_.currentStep(), baseif_.guideRate(), tcalc.guideTargetMode(), pu, false, injectionPhase); - - auto localFraction = [&](const std::string& child) { - return fcalc.localFraction(child, child); - }; - - auto localReduction = [&](const std::string& group_name) { - const std::vector& groupTargetReductions = group_state.injection_reduction_rates(group_name); - return tcalc.calcModeRateFromRates(groupTargetReductions); - }; - - const double orig_target = tcalc.groupTarget(group.injectionControls(injectionPhase, summaryState), deferred_logger); - const auto chain = WellGroupHelpers::groupChainTopBot(baseif_.name(), group.name(), schedule, baseif_.currentStep()); - // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. - const size_t num_ancestors = chain.size() - 1; - double target = orig_target; - for (size_t ii = 0; ii < num_ancestors; ++ii) { - if ((ii == 0) || baseif_.guideRate()->has(chain[ii], injectionPhase)) { - // Apply local reductions only at the control level - // (top) and for levels where we have a specified - // group guide rate. - target -= localReduction(chain[ii]); - } - target *= localFraction(chain[ii+1]); - } - // Avoid negative target rates coming from too large local reductions. - const double target_rate = std::max(0.0, target / efficiencyFactor); - const auto current_rate = injection_rate; - - control_eq = current_rate - target_rate; -} - -template -template -void -WellInterfaceEval:: -getGroupProductionControl(const Group& group, - const WellState& well_state, - const GroupState& group_state, - const Schedule& schedule, - const SummaryState& summaryState, - const EvalWell& bhp, - const std::vector& rates, - EvalWell& control_eq, - double efficiencyFactor) const -{ - const Group::ProductionCMode& currentGroupControl = group_state.production_control(group.name()); - if (currentGroupControl == Group::ProductionCMode::FLD || - currentGroupControl == Group::ProductionCMode::NONE) { - if (!group.productionGroupControlAvailable()) { - // We cannot go any further up the hierarchy. This could - // be the FIELD group, or any group for which this has - // been set in GCONINJE or GCONPROD. If we are here - // anyway, it is likely that the deck set inconsistent - // requirements, such as GRUP control mode on a well with - // no appropriate controls defined on any of its - // containing groups. We will therefore use the wells' bhp - // limit equation as a fallback. - const auto& controls = baseif_.wellEcl().productionControls(summaryState); - control_eq = bhp - controls.bhp_limit; - return; - } else { - // Produce share of parents control - const auto& parent = schedule.getGroup(group.parent(), baseif_.currentStep()); - efficiencyFactor *= group.getGroupEfficiencyFactor(); - getGroupProductionControl(parent, well_state, group_state, schedule, summaryState, bhp, rates, control_eq, efficiencyFactor); - return; - } - } - - const auto& well = baseif_.wellEcl(); - const auto pu = baseif_.phaseUsage(); - - if (!group.isProductionGroup()) { - // use bhp as control eq and let the updateControl code find a valid control - const auto& controls = well.productionControls(summaryState); - control_eq = bhp - controls.bhp_limit; - return; - } - - // 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. - - // Make conversion factors for RESV <-> surface rates. - std::vector resv_coeff(baseif_.phaseUsage().num_phases, 1.0); - baseif_.rateConverter().calcCoeff(0, baseif_.pvtRegionIdx(), resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS. - - // gconsale may adjust the grat target. - // the adjusted rates is send to the targetCalculator - double gratTargetFromSales = 0.0; - if (group_state.has_grat_sales_target(group.name())) - gratTargetFromSales = group_state.grat_sales_target(group.name()); - - WellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, resv_coeff, gratTargetFromSales, group.name(), group_state, group.has_gpmaint_control(currentGroupControl)); - WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, group_state, baseif_.currentStep(), baseif_.guideRate(), tcalc.guideTargetMode(), pu, true, Phase::OIL); - - auto localFraction = [&](const std::string& child) { - return fcalc.localFraction(child, child); - }; - - auto localReduction = [&](const std::string& group_name) { - const std::vector& groupTargetReductions = group_state.production_reduction_rates(group_name); - return tcalc.calcModeRateFromRates(groupTargetReductions); - }; - - const double orig_target = tcalc.groupTarget(group.productionControls(summaryState)); - const auto chain = WellGroupHelpers::groupChainTopBot(baseif_.name(), group.name(), schedule, baseif_.currentStep()); - // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. - const size_t num_ancestors = chain.size() - 1; - double target = orig_target; - for (size_t ii = 0; ii < num_ancestors; ++ii) { - if ((ii == 0) || baseif_.guideRate()->has(chain[ii])) { - // Apply local reductions only at the control level - // (top) and for levels where we have a specified - // group guide rate. - target -= localReduction(chain[ii]); - } - target *= localFraction(chain[ii+1]); - } - // Avoid negative target rates coming from too large local reductions. - const double target_rate = std::max(0.0, target / efficiencyFactor); - const auto current_rate = -tcalc.calcModeRateFromRates(rates); // Switch sign since 'rates' are negative for producers. - control_eq = current_rate - target_rate; -} - template template void @@ -364,7 +151,18 @@ assembleControlEqProd_(const WellState& well_state, active_rates[pu.phase_pos[canonical_phase]] = rates[canonical_phase]; } } - getGroupProductionControl(group, well_state, group_state, schedule, summaryState, bhp, active_rates, control_eq, efficiencyFactor); + auto rCoeff = [this](const RegionId id, const int region, std::vector& coeff) + { + baseif_.rateConverter().calcCoeff(id, region, coeff); + }; + WellGroupControls(baseif_).getGroupProductionControl(group, well_state, + group_state, + schedule, + summaryState, + bhp, active_rates, + rCoeff, + efficiencyFactor, + control_eq); break; } case Well::ProducerCMode::CMODE_UNDEFINED: { @@ -438,17 +236,22 @@ assembleControlEqInj_(const WellState& well_state, case Well::InjectorCMode::GRUP: { assert(baseif_.wellEcl().isAvailableForGroupControl()); const auto& group = schedule.getGroup(baseif_.wellEcl().groupName(), baseif_.currentStep()); - this->getGroupInjectionControl(group, - well_state, - group_state, - schedule, - summaryState, - injectorType, - bhp, - injection_rate, - control_eq, - efficiencyFactor, - deferred_logger); + auto rCoeff = [this](const RegionId id, const int region, std::vector& coeff) + { + baseif_.rateConverter().calcInjCoeff(id, region, coeff); + }; + WellGroupControls(baseif_).getGroupInjectionControl(group, + well_state, + group_state, + schedule, + summaryState, + injectorType, + bhp, + injection_rate, + rCoeff, + efficiencyFactor, + control_eq, + deferred_logger); break; } case Well::InjectorCMode::CMODE_UNDEFINED: { diff --git a/opm/simulators/wells/WellInterfaceEval.hpp b/opm/simulators/wells/WellInterfaceEval.hpp index 754e6567b..33159024f 100644 --- a/opm/simulators/wells/WellInterfaceEval.hpp +++ b/opm/simulators/wells/WellInterfaceEval.hpp @@ -56,30 +56,6 @@ public: const SummaryState& summaryState, const double rho, DeferredLogger& deferred_logger) const; - template - void getGroupInjectionControl(const Group& group, - const WellState& well_state, - const GroupState& group_state, - const Schedule& schedule, - const SummaryState& summaryState, - const InjectorType& injectorType, - const EvalWell& bhp, - const EvalWell& injection_rate, - EvalWell& control_eq, - double efficiencyFactor, - DeferredLogger& deferred_logger) const; - - - template - void getGroupProductionControl(const Group& group, - const WellState& well_state, - const GroupState& group_state, - const Schedule& schedule, - const SummaryState& summaryState, - const EvalWell& bhp, - const std::vector& rates, - EvalWell& control_eq, - double efficiencyFactor) const; template void assembleControlEqProd(const WellState& well_state, diff --git a/opm/simulators/wells/WellInterfaceFluidSystem.cpp b/opm/simulators/wells/WellInterfaceFluidSystem.cpp index 36f468d2f..0ba58228f 100644 --- a/opm/simulators/wells/WellInterfaceFluidSystem.cpp +++ b/opm/simulators/wells/WellInterfaceFluidSystem.cpp @@ -32,6 +32,7 @@ #include #include #include +#include #include #include @@ -531,97 +532,18 @@ getGroupInjectionTargetRate(const Group& group, double efficiencyFactor, DeferredLogger& deferred_logger) const { - // Setting some defaults to silence warnings below. - // Will be overwritten in the switch statement. - Phase injectionPhase = Phase::WATER; - switch (injectorType) { - case InjectorType::WATER: + auto rCoeff = [this](const int id, const int region, std::vector& coeff) { - injectionPhase = Phase::WATER; - break; - } - case InjectorType::OIL: - { - injectionPhase = Phase::OIL; - break; - } - case InjectorType::GAS: - { - injectionPhase = Phase::GAS; - break; - } - default: - // Should not be here. - assert(false); - } - - auto currentGroupControl = group_state.injection_control(group.name(), injectionPhase); - if (currentGroupControl == Group::InjectionCMode::FLD || - currentGroupControl == Group::InjectionCMode::NONE) { - if (!group.injectionGroupControlAvailable(injectionPhase)) { - // We cannot go any further up the hierarchy. This could - // be the FIELD group, or any group for which this has - // been set in GCONINJE or GCONPROD. If we are here - // anyway, it is likely that the deck set inconsistent - // requirements, such as GRUP control mode on a well with - // no appropriate controls defined on any of its - // containing groups. We will therefore use the wells' bhp - // limit equation as a fallback. - return std::nullopt; - } else { - // Inject share of parents control - const auto& parent = schedule.getGroup( group.parent(), currentStep()); - efficiencyFactor *= group.getGroupEfficiencyFactor(); - return getGroupInjectionTargetRate(parent, well_state, group_state, schedule, summaryState, injectorType, efficiencyFactor, deferred_logger); - } - } - - const auto pu = phaseUsage(); - - if (!group.isInjectionGroup()) { - return std::nullopt; - } - - // 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. - - // Make conversion factors for RESV <-> surface rates. - std::vector resv_coeff(pu.num_phases, 1.0); - rateConverter_.calcCoeff(0, pvtRegionIdx(), resv_coeff); // FIPNUM region 0 here, should use FIPNUM from WELSPECS. - - double sales_target = 0; - if (schedule[currentStep()].gconsale().has(group.name())) { - const auto& gconsale = schedule[currentStep()].gconsale().get(group.name(), summaryState); - sales_target = gconsale.sales_target; - } - WellGroupHelpers::InjectionTargetCalculator tcalc(currentGroupControl, pu, resv_coeff, group.name(), sales_target, group_state, injectionPhase, group.has_gpmaint_control(injectionPhase, currentGroupControl), deferred_logger); - WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, group_state, currentStep(), guideRate(), tcalc.guideTargetMode(), pu, false, injectionPhase); - - auto localFraction = [&](const std::string& child) { - return fcalc.localFraction(child, child); //Note child needs to be passed to always include since the global isGrup map is not updated yet. + this->rateConverter().calcCoeff(id, region, coeff); }; - auto localReduction = [&](const std::string& group_name) { - const std::vector& groupTargetReductions = group_state.injection_reduction_rates(group_name); - return tcalc.calcModeRateFromRates(groupTargetReductions); - }; - - const double orig_target = tcalc.groupTarget(group.injectionControls(injectionPhase, summaryState), deferred_logger); - const auto chain = WellGroupHelpers::groupChainTopBot(name(), group.name(), schedule, currentStep()); - // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. - const size_t num_ancestors = chain.size() - 1; - double target = orig_target; - for (size_t ii = 0; ii < num_ancestors; ++ii) { - if ((ii == 0) || guideRate()->has(chain[ii], injectionPhase)) { - // Apply local reductions only at the control level - // (top) and for levels where we have a specified - // group guide rate. - target -= localReduction(chain[ii]); - } - target *= localFraction(chain[ii+1]); - } - return std::max(0.0, target / efficiencyFactor); + return WellGroupControls(*this).getGroupInjectionTargetRate(group, well_state, + group_state, schedule, + summaryState, injectorType, + rCoeff, efficiencyFactor, + deferred_logger); } + template double WellInterfaceFluidSystem:: @@ -632,78 +554,17 @@ getGroupProductionTargetRate(const Group& group, const SummaryState& summaryState, double efficiencyFactor) const { - const Group::ProductionCMode& currentGroupControl = group_state.production_control(group.name()); - if (currentGroupControl == Group::ProductionCMode::FLD || - currentGroupControl == Group::ProductionCMode::NONE) { - if (!group.productionGroupControlAvailable()) { - return 1.0; - } else { - // Produce share of parents control - const auto& parent = schedule.getGroup(group.parent(), currentStep()); - efficiencyFactor *= group.getGroupEfficiencyFactor(); - return getGroupProductionTargetRate(parent, well_state, group_state, schedule, summaryState, efficiencyFactor); - } - } - - const auto pu = phaseUsage(); - - if (!group.isProductionGroup()) { - return 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. - - // 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. - - // gconsale may adjust the grat target. - // the adjusted rates is send to the targetCalculator - double gratTargetFromSales = 0.0; - if (group_state.has_grat_sales_target(group.name())) - gratTargetFromSales = group_state.grat_sales_target(group.name()); - - WellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, resv_coeff, gratTargetFromSales, group.name(), group_state, group.has_gpmaint_control(currentGroupControl)); - WellGroupHelpers::FractionCalculator fcalc(schedule, well_state, group_state, currentStep(), guideRate(), tcalc.guideTargetMode(), pu, true, Phase::OIL); - - auto localFraction = [&](const std::string& child) { - return fcalc.localFraction(child, child); //Note child needs to be passed to always include since the global isGrup map is not updated yet. + auto rCoeff = [this](const int id, const int region, std::vector& coeff) + { + this->rateConverter().calcCoeff(id, region, coeff); }; - auto localReduction = [&](const std::string& group_name) { - const std::vector& groupTargetReductions = group_state.production_reduction_rates(group_name); - return tcalc.calcModeRateFromRates(groupTargetReductions); - }; - - const double orig_target = tcalc.groupTarget(group.productionControls(summaryState)); - const auto chain = WellGroupHelpers::groupChainTopBot(name(), group.name(), schedule, currentStep()); - // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. - const size_t num_ancestors = chain.size() - 1; - double target = orig_target; - for (size_t ii = 0; ii < num_ancestors; ++ii) { - if ((ii == 0) || guideRate()->has(chain[ii])) { - // Apply local reductions only at the control level - // (top) and for levels where we have a specified - // group guide rate. - target -= localReduction(chain[ii]); - } - target *= localFraction(chain[ii+1]); - } - // Avoid negative target rates coming from too large local reductions. - const double target_rate = std::max(0.0, target / efficiencyFactor); - const auto& ws = well_state.well(this->index_of_well_); - const auto& rates = ws.surface_rates; - const auto current_rate = -tcalc.calcModeRateFromRates(rates); // Switch sign since 'rates' are negative for producers. - double scale = 1.0; - if (target_rate == 0.0) { - return 0.0; - } - - if (current_rate > 1e-14) - scale = target_rate/current_rate; - return scale; + return WellGroupControls(*this).getGroupProductionTargetRate(group, well_state, + group_state, schedule, + summaryState, + rCoeff, efficiencyFactor); } + template class WellInterfaceFluidSystem>; } // namespace Opm