From 766d02cacce5ed87feb7fba59a458ff69d08493f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Atgeirr=20Fl=C3=B8=20Rasmussen?= Date: Tue, 24 Mar 2020 11:07:14 +0100 Subject: [PATCH] Unify group control equation codes. --- opm/simulators/wells/MultisegmentWell.hpp | 25 +- .../wells/MultisegmentWell_impl.hpp | 373 +++++------------- opm/simulators/wells/StandardWell.hpp | 21 +- opm/simulators/wells/StandardWell_impl.hpp | 355 ++--------------- opm/simulators/wells/WellInterface.hpp | 22 ++ opm/simulators/wells/WellInterface_impl.hpp | 249 ++++++++++++ 6 files changed, 435 insertions(+), 610 deletions(-) diff --git a/opm/simulators/wells/MultisegmentWell.hpp b/opm/simulators/wells/MultisegmentWell.hpp index 6a1ad001b..a253c9147 100644 --- a/opm/simulators/wells/MultisegmentWell.hpp +++ b/opm/simulators/wells/MultisegmentWell.hpp @@ -283,11 +283,10 @@ namespace Opm void initMatrixAndVectors(const int num_cells) const; - // protected functions - // EvalWell getBhp(); this one should be something similar to getSegmentPressure(); - // EvalWell getQs(); this one should be something similar to getSegmentRates() - // EValWell wellVolumeFractionScaled, wellVolumeFraction, wellSurfaceVolumeFraction ... these should have different names, and probably will be needed. - // bool crossFlowAllowed(const Simulator& ebosSimulator) const; probably will be needed + EvalWell getBhp() const; + EvalWell getQs(const int comp_idx) const; + EvalWell getWQTotal() const; + // xw = inv(D)*(rw - C*x) void recoverSolutionWell(const BVector& x, BVectorWell& xw) const; @@ -381,8 +380,20 @@ namespace Opm const Well::ProductionControls& prod_controls, Opm::DeferredLogger& deferred_logger); - void assembleGroupProductionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, EvalWell& control_eq, double efficincyFactor, Opm::DeferredLogger& deferred_logger); - void assembleGroupInjectionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, const InjectorType& injectorType, EvalWell& control_eq, double efficincyFactor, Opm::DeferredLogger& deferred_logger); + void assembleGroupProductionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + EvalWell& control_eq, + double efficiencyFactor); + void assembleGroupInjectionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + EvalWell& control_eq, + double efficiencyFactor, + Opm::DeferredLogger& deferred_logger); void assemblePressureEq(const int seg) const; diff --git a/opm/simulators/wells/MultisegmentWell_impl.hpp b/opm/simulators/wells/MultisegmentWell_impl.hpp index 15bba590f..34c7e30d1 100644 --- a/opm/simulators/wells/MultisegmentWell_impl.hpp +++ b/opm/simulators/wells/MultisegmentWell_impl.hpp @@ -690,7 +690,7 @@ namespace Opm } if (pressure_controlled_well) { for (int compIdx = 0; compIdx < num_components_; ++compIdx) { - const EvalWell rate = this->getSegmentRate(0, compIdx); + const EvalWell rate = this->getQs(compIdx); well_potentials[ebosCompIdxToFlowCompIdx(compIdx)] = rate.value(); } return; @@ -790,7 +790,7 @@ namespace Opm const int np = number_of_phases_; well_flux.resize(np, 0.0); for (int compIdx = 0; compIdx < num_components_; ++compIdx) { - const EvalWell rate = well_copy.getSegmentRate(0, compIdx); + const EvalWell rate = well_copy.getQs(compIdx); well_flux[ebosCompIdxToFlowCompIdx(compIdx)] = rate.value(); } debug_cost_counter_ += well_copy.debug_cost_counter_; @@ -1621,6 +1621,18 @@ namespace Opm + template + typename MultisegmentWell::EvalWell + MultisegmentWell:: + getBhp() const + { + return getSegmentPressure(0); + } + + + + + template typename MultisegmentWell::EvalWell MultisegmentWell:: @@ -1634,6 +1646,18 @@ namespace Opm + template + typename MultisegmentWell::EvalWell + MultisegmentWell:: + getQs(const int comp_idx) const + { + return getSegmentRate(0, comp_idx); + } + + + + + template typename MultisegmentWell::EvalWell MultisegmentWell:: @@ -1690,6 +1714,18 @@ namespace Opm + template + typename MultisegmentWell::EvalWell + MultisegmentWell:: + getWQTotal() const + { + return getSegmentGTotal(0); + } + + + + + template void MultisegmentWell:: @@ -1763,7 +1799,7 @@ namespace Opm double efficiencyFactor = well.getEfficiencyFactor(); if (wellIsStopped_) { - control_eq = getSegmentGTotal(0); + control_eq = getWQTotal(); } else if (this->isInjector() ) { const Opm::Well::InjectorCMode& current = well_state.currentInjectionControls()[well_index]; const auto& controls = inj_controls; @@ -1796,7 +1832,7 @@ namespace Opm switch(current) { case Well::InjectorCMode::RATE: { - control_eq = getSegmentGTotal(0) / scaling - controls.surface_rate; + control_eq = getWQTotal() / scaling - controls.surface_rate; break; } @@ -1828,7 +1864,7 @@ namespace Opm } - control_eq = coeff*getSegmentGTotal(0) / scaling - controls.reservoir_rate; + control_eq = coeff*getWQTotal() / scaling - controls.reservoir_rate; break; } @@ -1836,22 +1872,22 @@ namespace Opm { std::vector rates(3, 0.); if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { - rates[ Water ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)); + rates[ Water ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)); } if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) { - rates[ Oil ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); + rates[ Oil ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); } if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { - rates[ Gas ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)); + rates[ Gas ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)); } - control_eq = getSegmentPressure(0) - calculateBhpFromThp(rates, well, summaryState, deferred_logger); + control_eq = getBhp() - calculateBhpFromThp(rates, well, summaryState, deferred_logger); break; } case Well::InjectorCMode::BHP: { const auto& bhp = controls.bhp_limit; - control_eq = getSegmentPressure(0) - bhp; + control_eq = getBhp() - bhp; break; } @@ -1881,21 +1917,21 @@ namespace Opm case Well::ProducerCMode::ORAT: { assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)); - const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); + const EvalWell& rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); control_eq = rate - controls.oil_rate; break; } case Well::ProducerCMode::WRAT: { assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)); - const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)); + const EvalWell& rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)); control_eq = rate - controls.water_rate; break; } case Well::ProducerCMode::GRAT: { assert(FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)); - const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)); + const EvalWell& rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)); control_eq = rate - controls.gas_rate; break; @@ -1904,8 +1940,8 @@ namespace Opm { assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)); assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)); - EvalWell rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)) - -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); + EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)) + -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); control_eq = rate - controls.liquid_rate; break; } @@ -1919,7 +1955,7 @@ namespace Opm std::vector convert_coeff(number_of_phases_, 1.0); Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff); for (int phase = 0; phase < number_of_phases_; ++phase) { - total_rate += getSegmentRate(0, flowPhaseToEbosCompIdx(phase) ) * convert_coeff[phase]; + total_rate += getQs(flowPhaseToEbosCompIdx(phase) ) * convert_coeff[phase]; } if (controls.prediction_mode) { @@ -1945,22 +1981,22 @@ namespace Opm } case Well::ProducerCMode::BHP: { - control_eq = getSegmentPressure(0) - controls.bhp_limit; + control_eq = getBhp() - controls.bhp_limit; break; } case Well::ProducerCMode::THP: { std::vector rates(3, 0.); if (FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)) { - rates[ Water ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)); + rates[ Water ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)); } if (FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)) { - rates[ Oil ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); + rates[ Oil ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); } if (FluidSystem::phaseIsActive(FluidSystem::gasPhaseIdx)) { - rates[ Gas ] = getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)); + rates[ Gas ] = getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)); } - control_eq = getSegmentPressure(0) - calculateBhpFromThp(rates, well, summaryState, deferred_logger); + control_eq = getBhp() - calculateBhpFromThp(rates, well, summaryState, deferred_logger); break; } case Well::ProducerCMode::GRUP: @@ -1968,7 +2004,7 @@ namespace Opm assert(well.isAvailableForGroupControl()); const auto& group = schedule.getGroup( well.groupName(), current_step_ ); - assembleGroupProductionControl(group, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger); + assembleGroupProductionControl(group, well_state, schedule, summaryState, control_eq, efficiencyFactor); break; } case Well::ProducerCMode::CMODE_UNDEFINED: @@ -2112,168 +2148,54 @@ namespace Opm template void - MultisegmentWell:: - assembleGroupInjectionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, const InjectorType& injectorType, EvalWell& control_eq, double efficiencyFactor, Opm::DeferredLogger& deferred_logger) + MultisegmentWell::assembleGroupInjectionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + EvalWell& control_eq, + double efficiencyFactor, + Opm::DeferredLogger& deferred_logger) { - - if (!group.isAvailableForGroupControl()) { - // 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_ecl_.injectionControls(summaryState); - control_eq = getSegmentPressure(0) - controls.bhp_limit; - return; - } - - const auto& well = well_ecl_; + // The total rate primary variable is the scaled surface rate + // sum, also for injectors, so we must scale it to get a + // proper surface injection rate. This is different from + // standard wells, where the primary variable is unscaled if + // the well is an injector. const auto& pu = phaseUsage(); - - int phasePos; - Well::GuideRateTarget wellTarget; - Phase injectionPhase; double scaling = 1.0; - switch (injectorType) { case InjectorType::WATER: { - phasePos = pu.phase_pos[BlackoilPhases::Aqua]; - wellTarget = Well::GuideRateTarget::WAT; - injectionPhase = Phase::WATER; scaling = scalingFactor(pu.phase_pos[BlackoilPhases::Aqua]); break; } case InjectorType::OIL: { - phasePos = pu.phase_pos[BlackoilPhases::Liquid]; - wellTarget = Well::GuideRateTarget::OIL; - injectionPhase = Phase::OIL; scaling = scalingFactor(pu.phase_pos[BlackoilPhases::Liquid]); break; } case InjectorType::GAS: { - phasePos = pu.phase_pos[BlackoilPhases::Vapour]; - wellTarget = Well::GuideRateTarget::GAS; - injectionPhase = Phase::GAS; scaling = scalingFactor(pu.phase_pos[BlackoilPhases::Vapour]); break; } default: - throw("Expected WATER, OIL or GAS as type for injectors " + well.name()); - } - - const Group::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(injectionPhase, group.name()); - if (currentGroupControl == Group::InjectionCMode::FLD || - currentGroupControl == Group::InjectionCMode::NONE) { - // Inject share of parents control - const auto& parent = schedule.getGroup( group.parent(), current_step_ ); - - efficiencyFactor *= group.getGroupEfficiencyFactor(); - - assembleGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, control_eq, efficiencyFactor, deferred_logger); - return; - } - - assert(group.hasInjectionControl(injectionPhase)); - const auto& groupcontrols = group.injectionControls(injectionPhase, summaryState); - - const std::vector& groupInjectionReductions = well_state.currentInjectionGroupReductionRates(group.name()); - double groupTargetReduction = groupInjectionReductions[phasePos]; - double fraction = wellGroupHelpers::fractionFromInjectionPotentials(well.name(), group.name(), schedule,well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(wellTarget), pu, injectionPhase,false); - - switch(currentGroupControl) { - case Group::InjectionCMode::NONE: - { - // The NONE case is handled earlier + // Should not be here. assert(false); - break; } - case Group::InjectionCMode::RATE: - { - double target = std::max(0.0, (groupcontrols.surface_max_rate - groupTargetReduction)) / efficiencyFactor; - control_eq = getSegmentGTotal(0) / scaling - fraction * target; - break; - } - case Group::InjectionCMode::RESV: - { - std::vector convert_coeff(number_of_phases_, 1.0); - Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff); - double coeff = convert_coeff[phasePos]; - double target = std::max(0.0, (groupcontrols.resv_max_rate/coeff - groupTargetReduction)) / efficiencyFactor; - control_eq = getSegmentGTotal(0) / scaling - fraction * target; - break; - } - case Group::InjectionCMode::REIN: - { - double productionRate = well_state.currentInjectionREINRates(groupcontrols.reinj_group)[phasePos]; - productionRate /= efficiencyFactor; - double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction)) / efficiencyFactor; - control_eq = getSegmentGTotal(0) / scaling - fraction * target; - break; - } - case Group::InjectionCMode::VREP: - { - std::vector convert_coeff(number_of_phases_, 1.0); - Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff); - double coeff = convert_coeff[phasePos]; - double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group)*groupcontrols.target_void_fraction; - - double injReduction = 0.0; - std::vector groupInjectionReservoirRates = well_state.currentInjectionGroupReservoirRates(group.name()); - if (groupcontrols.phase != Phase::WATER) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Aqua]]; - - if (groupcontrols.phase != Phase::OIL) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Liquid]]; - - if (groupcontrols.phase != Phase::GAS) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Vapour]]; - - voidageRate -= injReduction; - - double target = std::max(0.0, ( voidageRate/coeff - groupTargetReduction)) / efficiencyFactor; - control_eq = getSegmentGTotal(0) / scaling - fraction * target; - break; - } - case Group::InjectionCMode::FLD: - { - // The FLD case is handled earlier - assert(false); - break; - } - case Group::InjectionCMode::SALE: - { - // only for gas injectors - assert (phasePos == pu.phase_pos[BlackoilPhases::Vapour]); - - // Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate; - double inj_rate = well_state.currentInjectionREINRates(group.name())[phasePos]; - if (schedule.gConSump(current_step_).has(group.name())) { - const auto& gconsump = schedule.gConSump(current_step_).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(current_step_).get(group.name(), summaryState); - inj_rate -= gconsale.sales_target; - - inj_rate /= efficiencyFactor; - double target = std::max(0.0, (inj_rate - groupTargetReduction)) / efficiencyFactor; - control_eq = getSegmentGTotal(0) /scaling - fraction * target; - break; - } - - default: - OPM_DEFLOG_THROW(std::runtime_error, "Unvalid group control specified for group " + well.groupName(), deferred_logger ); - } - + const EvalWell injection_rate = getWQTotal() / scaling; + // Call the generic implementation. + Base::getGroupInjectionControl(group, + well_state, + schedule, + summaryState, + injectorType, + getBhp(), + injection_rate, + control_eq, + efficiencyFactor); } @@ -2281,122 +2203,26 @@ namespace Opm template void - MultisegmentWell:: - assembleGroupProductionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, EvalWell& control_eq, double efficiencyFactor, Opm::DeferredLogger& deferred_logger) + MultisegmentWell::assembleGroupProductionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + EvalWell& control_eq, + double efficiencyFactor) { - - const auto& well = well_ecl_; const auto pu = phaseUsage(); - - const Group::ProductionCMode& currentGroupControl = well_state.currentProductionGroupControl(group.name()); - if (currentGroupControl == Group::ProductionCMode::FLD || - currentGroupControl == Group::ProductionCMode::NONE) { - if (!group.isAvailableForGroupControl()) { - // 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_ecl_.productionControls(summaryState); - control_eq = getSegmentPressure(0) - controls.bhp_limit; - return; - } else { - // Produce share of parents control - const auto& parent = schedule.getGroup( group.parent(), current_step_ ); - efficiencyFactor *= group.getGroupEfficiencyFactor(); - assembleGroupProductionControl(parent, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger); - return; + std::vector rates(pu.num_phases); + const int compIndices[3] = { FluidSystem::waterCompIdx, FluidSystem::oilCompIdx, FluidSystem::gasCompIdx }; + const BlackoilPhases::PhaseIndex phases[3] = { BlackoilPhases::Aqua, BlackoilPhases::Liquid, BlackoilPhases::Vapour }; + for (int canonical_phase = 0; canonical_phase < 3; ++canonical_phase) { + const auto phase = phases[canonical_phase]; + if (pu.phase_used[phase]) { + const auto compIdx = compIndices[canonical_phase]; + rates[pu.phase_pos[phase]] = getQs(Indices::canonicalToActiveComponentIndex(compIdx)); } } - if (!group.isProductionGroup()) { - // use bhp as control eq and let the updateControl code find a vallied control - const auto& controls = well.productionControls(summaryState); - control_eq = getSegmentPressure(0) - controls.bhp_limit; - return; - } - - const auto& groupcontrols = group.productionControls(summaryState); - const std::vector& groupTargetReductions = well_state.currentProductionGroupReductionRates(group.name()); - - switch(currentGroupControl) { - case Group::ProductionCMode::NONE: - { - // The NONE case is handled earlier - assert(false); - break; - } - case Group::ProductionCMode::ORAT: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]]; - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::OIL), pu); - - const double rate_target = std::max(0.0, groupcontrols.oil_target - groupTargetReduction) / efficiencyFactor; - assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)); - const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); - control_eq = rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::WRAT: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Water]]; - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::WAT), pu); - - const double rate_target = std::max(0.0, groupcontrols.gas_target - groupTargetReduction) / efficiencyFactor; - assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)); - const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)); - control_eq = rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::GRAT: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Gas]]; - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::GAS), pu); - - const double rate_target = std::max(0.0, groupcontrols.gas_target - groupTargetReduction) / efficiencyFactor; - assert(FluidSystem::phaseIsActive(FluidSystem::gasCompIdx)); - const EvalWell& rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)); - control_eq = rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::LRAT: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]] + groupTargetReductions[pu.phase_pos[Water]]; - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::LIQ), pu); - - const double rate_target = std::max(0.0, groupcontrols.liquid_target - groupTargetReduction) / efficiencyFactor; assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)); - EvalWell rate = -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)) - -getSegmentRate(0, Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); - control_eq = rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::CRAT: - { - OPM_DEFLOG_THROW(std::runtime_error, "CRAT group control not implemented for producers", deferred_logger ); - break; - } - case Group::ProductionCMode::RESV: - { - OPM_DEFLOG_THROW(std::runtime_error, "RESV group control not implemented for producers", deferred_logger ); - break; - } - case Group::ProductionCMode::PRBL: - { - OPM_DEFLOG_THROW(std::runtime_error, "PRBL group control not implemented for producers", deferred_logger ); - break; - } - case Group::ProductionCMode::FLD: - { - // The FLD case is handled earlier - assert(false); - break; - } - default: - OPM_DEFLOG_THROW(std::runtime_error, "Unvallied group control specified for group " + well.groupName(), deferred_logger ); - } + Base::getGroupProductionControl(group, well_state, schedule, summaryState, getBhp(), rates, control_eq, efficiencyFactor); } @@ -4110,5 +3936,4 @@ namespace Opm const double sign = mass_rate <= 0. ? 1.0 : -1.0; return sign * (friction_pressure_loss + constriction_pressure_loss); } - -} + } diff --git a/opm/simulators/wells/StandardWell.hpp b/opm/simulators/wells/StandardWell.hpp index 96042e2a3..0de8e9419 100644 --- a/opm/simulators/wells/StandardWell.hpp +++ b/opm/simulators/wells/StandardWell.hpp @@ -410,10 +410,25 @@ namespace Opm void updateThp(WellState& well_state, Opm::DeferredLogger& deferred_logger) const; - void assembleControlEq(const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, Opm::DeferredLogger& deferred_logger); + void assembleControlEq(const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + Opm::DeferredLogger& deferred_logger); - void assembleGroupProductionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, EvalWell& control_eq, double efficincyFactor, Opm::DeferredLogger& deferred_logger); - void assembleGroupInjectionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, const InjectorType& injectorType, EvalWell& control_eq, double efficincyFactor, Opm::DeferredLogger& deferred_logger); + void assembleGroupProductionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + EvalWell& control_eq, + double efficiencyFactor); + void assembleGroupInjectionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + EvalWell& control_eq, + double efficiencyFactor, + Opm::DeferredLogger& deferred_logger); // handle the non reasonable fractions due to numerical overshoot void processFractions() const; diff --git a/opm/simulators/wells/StandardWell_impl.hpp b/opm/simulators/wells/StandardWell_impl.hpp index 3b3172268..7087a419c 100644 --- a/opm/simulators/wells/StandardWell_impl.hpp +++ b/opm/simulators/wells/StandardWell_impl.hpp @@ -764,8 +764,10 @@ namespace Opm template void - StandardWell:: - assembleControlEq(const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, Opm::DeferredLogger& deferred_logger) + StandardWell::assembleControlEq(const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + Opm::DeferredLogger& deferred_logger) { EvalWell control_eq(numWellEq_ + numEq, 0.); @@ -952,7 +954,7 @@ namespace Opm assert(well.isAvailableForGroupControl()); const auto& group = schedule.getGroup( well.groupName(), current_step_ ); - assembleGroupProductionControl(group, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger); + assembleGroupProductionControl(group, well_state, schedule, summaryState, control_eq, efficiencyFactor); break; } case Well::ProducerCMode::CMODE_UNDEFINED: @@ -978,157 +980,24 @@ namespace Opm template void - StandardWell:: - assembleGroupInjectionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, const InjectorType& injectorType, EvalWell& control_eq, double efficiencyFactor, Opm::DeferredLogger& deferred_logger) + StandardWell::assembleGroupInjectionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + EvalWell& control_eq, + double efficiencyFactor, + Opm::DeferredLogger& deferred_logger) { - if (!group.isAvailableForGroupControl()) { - // 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_ecl_.injectionControls(summaryState); - control_eq = getBhp() - controls.bhp_limit; - return; - } - - const auto& well = well_ecl_; - const auto pu = phaseUsage(); - - int phasePos; - Well::GuideRateTarget wellTarget; - Phase injectionPhase; - - switch (injectorType) { - case InjectorType::WATER: - { - phasePos = pu.phase_pos[BlackoilPhases::Aqua]; - wellTarget = Well::GuideRateTarget::WAT; - injectionPhase = Phase::WATER; - break; - } - case InjectorType::OIL: - { - phasePos = pu.phase_pos[BlackoilPhases::Liquid]; - wellTarget = Well::GuideRateTarget::OIL; - injectionPhase = Phase::OIL; - break; - } - case InjectorType::GAS: - { - phasePos = pu.phase_pos[BlackoilPhases::Vapour]; - wellTarget = Well::GuideRateTarget::GAS; - injectionPhase = Phase::GAS; - break; - } - default: - throw("Expected WATER, OIL or GAS as type for injectors " + well.name()); - } - const Group::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(injectionPhase, group.name()); - - if (currentGroupControl == Group::InjectionCMode::FLD || - currentGroupControl == Group::InjectionCMode::NONE) { - // Inject share of parents control - const auto& parent = schedule.getGroup( group.parent(), current_step_ ); - efficiencyFactor *= group.getGroupEfficiencyFactor(); - assembleGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, control_eq, efficiencyFactor, deferred_logger); - return; - } - - assert(group.hasInjectionControl(injectionPhase)); - const auto& groupcontrols = group.injectionControls(injectionPhase, summaryState); - - const std::vector& groupInjectionReductions = well_state.currentInjectionGroupReductionRates(group.name()); - double groupTargetReduction = groupInjectionReductions[phasePos]; - double fraction = wellGroupHelpers::fractionFromInjectionPotentials(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(wellTarget), pu, injectionPhase,false); - switch(currentGroupControl) { - case Group::InjectionCMode::NONE: - { - // The NONE case is handled earlier - assert(false); - break; - } - case Group::InjectionCMode::RATE: - { - double target = std::max(0.0, (groupcontrols.surface_max_rate - groupTargetReduction)) / efficiencyFactor; - control_eq = getWQTotal() - fraction * target; - break; - } - case Group::InjectionCMode::RESV: - { - std::vector convert_coeff(number_of_phases_, 1.0); - Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff); - double coeff = convert_coeff[phasePos]; - double target = std::max(0.0, (groupcontrols.resv_max_rate/coeff - groupTargetReduction)) / efficiencyFactor; - control_eq = getWQTotal() - fraction * target; - break; - } - case Group::InjectionCMode::REIN: - { - double productionRate = well_state.currentInjectionREINRates(groupcontrols.reinj_group)[phasePos]; - double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction)) / efficiencyFactor; - control_eq = getWQTotal() - fraction * target; - break; - } - case Group::InjectionCMode::VREP: - { - std::vector convert_coeff(number_of_phases_, 1.0); - Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff); - double coeff = convert_coeff[phasePos]; - double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group)*groupcontrols.target_void_fraction; - - double injReduction = 0.0; - std::vector groupInjectionReservoirRates = well_state.currentInjectionGroupReservoirRates(group.name()); - if (groupcontrols.phase != Phase::WATER) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Aqua]]; - - if (groupcontrols.phase != Phase::OIL) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Liquid]]; - - if (groupcontrols.phase != Phase::GAS) - injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Vapour]]; - - voidageRate -= injReduction; - - double target = std::max(0.0, ( voidageRate/coeff - groupTargetReduction)) / efficiencyFactor; - control_eq = getWQTotal() - fraction * target; - break; - } - case Group::InjectionCMode::FLD: - { - // The FLD case is handled earlier - assert(false); - break; - } - case Group::InjectionCMode::SALE: - { - // only for gas injectors - assert (phasePos == pu.phase_pos[BlackoilPhases::Vapour]); - - // Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate; - double inj_rate = well_state.currentInjectionREINRates(group.name())[phasePos]; - if (schedule.gConSump(current_step_).has(group.name())) { - const auto& gconsump = schedule.gConSump(current_step_).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(current_step_).get(group.name(), summaryState); - inj_rate -= gconsale.sales_target; - - double target = std::max(0.0, (inj_rate - groupTargetReduction)) / efficiencyFactor; - control_eq = getWQTotal() - fraction * target; - break; - } - default: - OPM_DEFLOG_THROW(std::runtime_error, "Unvalid group control specified for group " + well.groupName(), deferred_logger ); - } - - + Base::getGroupInjectionControl(group, + well_state, + schedule, + summaryState, + injectorType, + getBhp(), + getWQTotal(), + control_eq, + efficiencyFactor); } @@ -1136,78 +1005,14 @@ namespace Opm template void - StandardWell:: - assembleGroupProductionControl(const Group& group, const WellState& well_state, const Opm::Schedule& schedule, const SummaryState& summaryState, EvalWell& control_eq, double efficiencyFactor, Opm::DeferredLogger& deferred_logger) + StandardWell::assembleGroupProductionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + EvalWell& control_eq, + double efficiencyFactor) { - const auto& well = well_ecl_; const auto pu = phaseUsage(); - - const Group::ProductionCMode& currentGroupControl = well_state.currentProductionGroupControl(group.name()); - if (currentGroupControl == Group::ProductionCMode::FLD || - currentGroupControl == Group::ProductionCMode::NONE) { - if (!group.isAvailableForGroupControl()) { - // 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_ecl_.productionControls(summaryState); - control_eq = getBhp() - controls.bhp_limit; - return; - } else { - // Produce share of parents control - const auto& parent = schedule.getGroup( group.parent(), current_step_ ); - efficiencyFactor *= group.getGroupEfficiencyFactor(); - assembleGroupProductionControl(parent, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger); - return; - } - } - - if (!group.isProductionGroup()) { - // use bhp as control eq and let the updateControl code find a vallied control - const auto& controls = well.productionControls(summaryState); - control_eq = getBhp() - controls.bhp_limit; - return; - } - - // ------------------------------- New code start -------------------------------- - - // 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, Base::rateConverter_, Base::pvtRegionIdx_); - wellGroupHelpers::FractionCalculator fcalc(schedule, well_state, current_step_, Base::guide_rate_, tcalc.guideTargetMode(), pu); - - auto localFraction = [&](const std::string& child) { - return fcalc.localFraction(child, ""); - }; - - auto localReduction = [&](const std::string& group_name) { - const std::vector& groupTargetReductions = well_state.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 auto chain = wellGroupHelpers::groupChainTopBot(name(), group.name(), schedule, current_step_); - // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. - const size_t num_ancestors = chain.size() - 1; - double target = orig_target; - for (size_t ii = 0; ii < num_ancestors; ++ii) { - target -= localReduction(chain[ii]); - // Next lines: different from in WellGroupHelpers.hpp - // if (ii == num_ancestors - 1) { - // // Final level. Add my reduction back. - // target += current_rate*efficiencyFactor; - // } - target *= localFraction(chain[ii+1]); - } - const double target_rate = target / efficiencyFactor; - std::vector rates(pu.num_phases); const int compIndices[3] = { FluidSystem::waterCompIdx, FluidSystem::oilCompIdx, FluidSystem::gasCompIdx }; const BlackoilPhases::PhaseIndex phases[3] = { BlackoilPhases::Aqua, BlackoilPhases::Liquid, BlackoilPhases::Vapour }; @@ -1218,110 +1023,8 @@ namespace Opm rates[pu.phase_pos[phase]] = getQs(Indices::canonicalToActiveComponentIndex(compIdx)); } } - const auto current_rate = -tcalc.calcModeRateFromRates(rates); // Switch sign since 'rates' are negative for producers. - control_eq = current_rate - target_rate; - - // ------------------------------- New code end -------------------------------- - - -#if 0 - const auto& groupcontrols = group.productionControls(summaryState); - const std::vector& groupTargetReductions = well_state.currentProductionGroupReductionRates(group.name()); - - switch(currentGroupControl) { - case Group::ProductionCMode::NONE: - { - // The NONE case is handled earlier - assert(false); - } - case Group::ProductionCMode::ORAT: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]]; - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::OIL), pu); - - const double rate_target = std::max(0.0, groupcontrols.oil_target - groupTargetReduction) / efficiencyFactor; - assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)); - EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); - control_eq = rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::WRAT: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Water]]; - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::WAT), pu); - - const double rate_target = std::max(0.0, groupcontrols.water_target - groupTargetReduction) / efficiencyFactor; - assert(FluidSystem::phaseIsActive(FluidSystem::waterPhaseIdx)); - EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)); - control_eq = rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::GRAT: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Gas]]; - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::GAS), pu); - - const double rate_target = std::max(0.0, groupcontrols.gas_target - groupTargetReduction) / efficiencyFactor; - assert(FluidSystem::phaseIsActive(FluidSystem::gasCompIdx)); - EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::gasCompIdx)); - control_eq = rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::LRAT: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]] + groupTargetReductions[pu.phase_pos[Water]]; - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::LIQ), pu); - - const double rate_target = std::max(0.0, groupcontrols.liquid_target - groupTargetReduction) / efficiencyFactor; - assert(FluidSystem::phaseIsActive(FluidSystem::oilPhaseIdx)); - - EvalWell rate = -getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::waterCompIdx)) - - getQs(Indices::canonicalToActiveComponentIndex(FluidSystem::oilCompIdx)); - control_eq = rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::CRAT: - { - OPM_DEFLOG_THROW(std::runtime_error, "CRAT group control not implemented for producers", deferred_logger ); - break; - } - case Group::ProductionCMode::RESV: - { - double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]] - + groupTargetReductions[pu.phase_pos[Gas]] - + groupTargetReductions[pu.phase_pos[Water]]; - - double fraction = wellGroupHelpers::fractionFromGuideRates(well.name(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, GuideRateModel::convert_target(Well::GuideRateTarget::RES), pu); - - EvalWell total_rate(numWellEq_ + numEq, 0.); // reservoir rate - std::vector convert_coeff(number_of_phases_, 1.0); - Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff); - for (int phase = 0; phase < number_of_phases_; ++phase) { - total_rate -= getQs( flowPhaseToEbosCompIdx(phase) ) * convert_coeff[phase]; - } - - const double rate_target = std::max(0.0, groupcontrols.resv_target - groupTargetReduction) / efficiencyFactor; - assert(FluidSystem::phaseIsActive(FluidSystem::gasCompIdx)); - control_eq = total_rate - fraction * rate_target; - break; - } - case Group::ProductionCMode::PRBL: - { - OPM_DEFLOG_THROW(std::runtime_error, "PRBL group control not implemented for producers", deferred_logger ); - break; - } - case Group::ProductionCMode::FLD: - { - // The FLD case is handled earlier - assert(false); - break; - } - - default: - OPM_DEFLOG_THROW(std::runtime_error, "Unvallied group control specified for group " + well.groupName(), deferred_logger ); - } -#endif + Base::getGroupProductionControl(group, well_state, schedule, summaryState, getBhp(), rates, control_eq, efficiencyFactor); } diff --git a/opm/simulators/wells/WellInterface.hpp b/opm/simulators/wells/WellInterface.hpp index eabe290b9..2076df47a 100644 --- a/opm/simulators/wells/WellInterface.hpp +++ b/opm/simulators/wells/WellInterface.hpp @@ -508,6 +508,28 @@ namespace Opm const Schedule& schedule, const SummaryState& summaryState, DeferredLogger& deferred_logger) const; + + template + void getGroupInjectionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + const EvalWell& bhp, + const EvalWell& injection_rate, + EvalWell& control_eq, + double efficiencyFactor); + + template + void getGroupProductionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const EvalWell& bhp, + const std::vector& rates, + EvalWell& control_eq, + double efficiencyFactor); + }; diff --git a/opm/simulators/wells/WellInterface_impl.hpp b/opm/simulators/wells/WellInterface_impl.hpp index 935b67bfc..e4a670b24 100644 --- a/opm/simulators/wells/WellInterface_impl.hpp +++ b/opm/simulators/wells/WellInterface_impl.hpp @@ -1814,6 +1814,255 @@ namespace Opm + template + template + void + WellInterface::getGroupInjectionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const InjectorType& injectorType, + const EvalWell& bhp, + const EvalWell& injection_rate, + EvalWell& control_eq, + double efficiencyFactor) + { + if (!group.isAvailableForGroupControl()) { + // 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_ecl_.injectionControls(summaryState); + control_eq = bhp - controls.bhp_limit; + return; + } + + const auto& well = well_ecl_; + const auto pu = phaseUsage(); + + int phasePos = -1; + Well::GuideRateTarget wellTarget; + Phase injectionPhase; + + switch (injectorType) { + case InjectorType::WATER: + { + phasePos = pu.phase_pos[BlackoilPhases::Aqua]; + wellTarget = Well::GuideRateTarget::WAT; + injectionPhase = Phase::WATER; + break; + } + case InjectorType::OIL: + { + phasePos = pu.phase_pos[BlackoilPhases::Liquid]; + wellTarget = Well::GuideRateTarget::OIL; + injectionPhase = Phase::OIL; + break; + } + case InjectorType::GAS: + { + phasePos = pu.phase_pos[BlackoilPhases::Vapour]; + wellTarget = Well::GuideRateTarget::GAS; + injectionPhase = Phase::GAS; + break; + } + default: + // Should not be here. + assert(false); + } + const Group::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(injectionPhase, group.name()); + + if (currentGroupControl == Group::InjectionCMode::FLD || + currentGroupControl == Group::InjectionCMode::NONE) { + // Inject share of parents control + const auto& parent = schedule.getGroup( group.parent(), current_step_ ); + efficiencyFactor *= group.getGroupEfficiencyFactor(); + getGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, bhp, injection_rate, control_eq, efficiencyFactor); + return; + } + + assert(group.hasInjectionControl(injectionPhase)); + const auto& groupcontrols = group.injectionControls(injectionPhase, summaryState); + + const std::vector& groupInjectionReductions = well_state.currentInjectionGroupReductionRates(group.name()); + double groupTargetReduction = groupInjectionReductions[phasePos]; + double fraction = wellGroupHelpers::fractionFromInjectionPotentials(well.name(), + group.name(), + schedule, + well_state, + current_step_, + guide_rate_, + GuideRateModel::convert_target(wellTarget), + pu, + injectionPhase, + false); + switch (currentGroupControl) { + case Group::InjectionCMode::NONE: + { + // The NONE case is handled earlier + assert(false); + break; + } + case Group::InjectionCMode::RATE: + { + double target = std::max(0.0, (groupcontrols.surface_max_rate - groupTargetReduction)) / efficiencyFactor; + control_eq = injection_rate - fraction * target; + break; + } + case Group::InjectionCMode::RESV: + { + std::vector convert_coeff(number_of_phases_, 1.0); + rateConverter_.calcCoeff(/*fipreg*/ 0, pvtRegionIdx_, convert_coeff); + double coeff = convert_coeff[phasePos]; + double target = std::max(0.0, (groupcontrols.resv_max_rate/coeff - groupTargetReduction)) / efficiencyFactor; + control_eq = injection_rate - fraction * target; + break; + } + case Group::InjectionCMode::REIN: + { + double productionRate = well_state.currentInjectionREINRates(groupcontrols.reinj_group)[phasePos]; + double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction)) / efficiencyFactor; + control_eq = injection_rate - fraction * target; + break; + } + case Group::InjectionCMode::VREP: + { + std::vector convert_coeff(number_of_phases_, 1.0); + rateConverter_.calcCoeff(/*fipreg*/ 0, pvtRegionIdx_, convert_coeff); + double coeff = convert_coeff[phasePos]; + double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group)*groupcontrols.target_void_fraction; + + double injReduction = 0.0; + std::vector groupInjectionReservoirRates = well_state.currentInjectionGroupReservoirRates(group.name()); + if (groupcontrols.phase != Phase::WATER) + injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Aqua]]; + + if (groupcontrols.phase != Phase::OIL) + injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Liquid]]; + + if (groupcontrols.phase != Phase::GAS) + injReduction += groupInjectionReservoirRates[pu.phase_pos[BlackoilPhases::Vapour]]; + + voidageRate -= injReduction; + + double target = std::max(0.0, ( voidageRate/coeff - groupTargetReduction)) / efficiencyFactor; + control_eq = injection_rate - fraction * target; + break; + } + case Group::InjectionCMode::FLD: + { + // The FLD case is handled earlier + assert(false); + break; + } + case Group::InjectionCMode::SALE: + { + // only for gas injectors + assert (phasePos == pu.phase_pos[BlackoilPhases::Vapour]); + + // Gas injection rate = Total gas production rate + gas import rate - gas consumption rate - sales rate; + double inj_rate = well_state.currentInjectionREINRates(group.name())[phasePos]; + if (schedule.gConSump(current_step_).has(group.name())) { + const auto& gconsump = schedule.gConSump(current_step_).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(current_step_).get(group.name(), summaryState); + inj_rate -= gconsale.sales_target; + + double target = std::max(0.0, (inj_rate - groupTargetReduction)) / efficiencyFactor; + control_eq = injection_rate - fraction * target; + break; + } + // default: + // OPM_DEFLOG_THROW(std::runtime_error, "Unvalid group control specified for group " + well.groupName(), deferred_logger ); + } + } + + + + + template + template + void + WellInterface::getGroupProductionControl(const Group& group, + const WellState& well_state, + const Opm::Schedule& schedule, + const SummaryState& summaryState, + const EvalWell& bhp, + const std::vector& rates, + EvalWell& control_eq, + double efficiencyFactor) + { + const auto& well = well_ecl_; + const auto pu = phaseUsage(); + + const Group::ProductionCMode& currentGroupControl = well_state.currentProductionGroupControl(group.name()); + if (currentGroupControl == Group::ProductionCMode::FLD || + currentGroupControl == Group::ProductionCMode::NONE) { + if (!group.isAvailableForGroupControl()) { + // 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_ecl_.productionControls(summaryState); + control_eq = bhp - controls.bhp_limit; + return; + } else { + // Produce share of parents control + const auto& parent = schedule.getGroup( group.parent(), current_step_ ); + efficiencyFactor *= group.getGroupEfficiencyFactor(); + getGroupProductionControl(parent, well_state, schedule, summaryState, bhp, rates, control_eq, efficiencyFactor); + return; + } + } + + if (!group.isProductionGroup()) { + // use bhp as control eq and let the updateControl code find a vallied 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. + wellGroupHelpers::TargetCalculator tcalc(currentGroupControl, pu, rateConverter_, pvtRegionIdx_); + wellGroupHelpers::FractionCalculator fcalc(schedule, well_state, current_step_, guide_rate_, tcalc.guideTargetMode(), pu); + + auto localFraction = [&](const std::string& child) { + return fcalc.localFraction(child, ""); + }; + + auto localReduction = [&](const std::string& group_name) { + const std::vector& groupTargetReductions = well_state.currentProductionGroupReductionRates(group_name); + return tcalc.calcModeRateFromRates(groupTargetReductions); + }; + + const double orig_target = tcalc.groupTarget(group.productionControls(summaryState)); + const auto chain = wellGroupHelpers::groupChainTopBot(name(), group.name(), schedule, current_step_); + // Because 'name' is the last of the elements, and not an ancestor, we subtract one below. + const size_t num_ancestors = chain.size() - 1; + double target = orig_target; + for (size_t ii = 0; ii < num_ancestors; ++ii) { + target -= localReduction(chain[ii]); + target *= localFraction(chain[ii+1]); + } + const double target_rate = target / efficiencyFactor; + + const auto current_rate = -tcalc.calcModeRateFromRates(rates); // Switch sign since 'rates' are negative for producers. + + control_eq = current_rate - target_rate; + } } // namespace Opm