Support nupcol

Only update groupTargetReduction, group REIN rates and group VREP rates when iteration < nupcol

opm/simulators/wells/BlackoilWellModel_impl.hpp

@@ -788,7 +788,17 @@ namespace Opm {
             // only check group controls for iterationIdx smaller then nupcol
             const int reportStepIdx = ebosSimulator_.episodeIndex();
             const int nupcol = schedule().getNupcol(reportStepIdx);
-            updateWellControls(local_deferredLogger, iterationIdx < nupcol);
+            if (iterationIdx < nupcol) {
+                const Group2& fieldGroup = schedule().getGroup2("FIELD", reportStepIdx);
+                std::vector<double> groupTargetReduction(numPhases(), 0.0);
+                wellGroupHelpers::updateGroupTargetReduction(fieldGroup, schedule(), reportStepIdx, /*isInjector*/ false, well_state_, groupTargetReduction);
+                std::vector<double> rein(numPhases(), 0.0);
+                wellGroupHelpers::updateREINForGroups(fieldGroup, schedule(), reportStepIdx, well_state_, rein);
+                double resv = 0.0;
+                wellGroupHelpers::updateVREPForGroups(fieldGroup, schedule(), reportStepIdx, well_state_, resv);
+            }
+
+            updateWellControls(local_deferredLogger, true);
             // Set the well primary variables based on the value of well solutions
             initPrimaryVariablesEvaluation();
 

opm/simulators/wells/MultisegmentWell_impl.hpp

@@ -1980,7 +1980,17 @@ namespace Opm
         const auto& well = well_ecl_;
         const auto pu = phaseUsage();
         const Group2::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(group.name());
-        if (!group.isInjectionGroup() && currentGroupControl!=Group2::InjectionCMode::FLD )
+        if (currentGroupControl == Group2::InjectionCMode::FLD) {
+            // Inject share of parents control
+            const auto& parent = schedule.getGroup2( group.parent(), current_step_ );
+            if (group.getTransferGroupEfficiencyFactor())
+                efficiencyFactor *= group.getGroupEfficiencyFactor();
+
+            assembleGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, control_eq, efficiencyFactor, deferred_logger);
+            return;
+        }
+
+        if (!group.isInjectionGroup())
             return;
 
         const auto& groupcontrols = group.injectionControls(summaryState);
@@ -2019,8 +2029,8 @@ namespace Opm
             throw("Expected WATER, OIL or GAS as type for injectors " + well.name());
         }
 
-        double groupTargetReduction = 0.0;
-        wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/true, groupTargetReduction);
+        const std::vector<double>& groupTargetReductions = well_state.currentProductionGroupReductionRates(group.name());
+        double groupTargetReduction = groupTargetReductions[phasePos];
         double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, wellTarget, /*isInjector*/true);
         wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, groupTarget, /*isInjector*/true, fraction);
 
@@ -2047,9 +2057,7 @@ namespace Opm
         }
         case Group2::InjectionCMode::REIN:
         {
-            double productionRate = 0.0;
-            const Group2& reinGroup = schedule.getGroup2( groupcontrols.reinj_group, current_step_ );
-            productionRate += wellGroupHelpers::sumWellRates(reinGroup, schedule, well_state, current_step_, phasePos, /*isInjector*/false);
+            double productionRate = well_state.currentInjectionVREPRates(groupcontrols.reinj_group);
             productionRate /= efficiencyFactor;
             double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction));
             control_eq = getSegmentGTotal(0) / scaling - fraction * target;
@@ -2060,11 +2068,7 @@ namespace Opm
             std::vector<double> convert_coeff(number_of_phases_, 1.0);
             Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
             double coeff = convert_coeff[phasePos];
-            double voidageRate = 0.0;
-            const Group2& voidageGroup = schedule.getGroup2( groupcontrols.voidage_group, current_step_ );
-            voidageRate += wellGroupHelpers::sumWellResRates(voidageGroup, schedule, well_state, current_step_, pu.phase_pos[BlackoilPhases::Aqua], /*injector*/false);
-            voidageRate += wellGroupHelpers::sumWellResRates(voidageGroup, schedule, well_state, current_step_, pu.phase_pos[BlackoilPhases::Liquid], /*injector*/false);
-            voidageRate += wellGroupHelpers::sumWellResRates(voidageGroup, schedule, well_state, current_step_, pu.phase_pos[BlackoilPhases::Vapour], /*injector*/false);
+            double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group);
             voidageRate /= efficiencyFactor;
             double target = std::max(0.0, ( groupcontrols.target_void_fraction*voidageRate/coeff - groupTargetReduction));
             control_eq = getSegmentGTotal(0) / scaling - fraction * target ;
@@ -2072,11 +2076,8 @@ namespace Opm
         }
         case Group2::InjectionCMode::FLD:
         {
-            // Inject share of parents control
-            const auto& parent = schedule.getGroup2( group.parent(), current_step_ );
-            if (group.getTransferGroupEfficiencyFactor())
-                efficiencyFactor *= group.getGroupEfficiencyFactor();
-            assembleGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, control_eq, efficiencyFactor, deferred_logger);
+            // The FLD case is handled earlier
+            assert(false);
             break;
         }
 
@@ -2100,10 +2101,21 @@ namespace Opm
         const auto pu = phaseUsage();
 
         const Group2::ProductionCMode& currentGroupControl = well_state.currentProductionGroupControl(group.name());
-        if (!group.isProductionGroup() && currentGroupControl != Group2::ProductionCMode::FLD)
+
+        if (currentGroupControl == Group2::ProductionCMode::FLD ) {
+            // Produce share of parents control
+            const auto& parent = schedule.getGroup2( group.parent(), current_step_ );
+            if (group.getTransferGroupEfficiencyFactor())
+                efficiencyFactor *= group.getGroupEfficiencyFactor();
+
+            assembleGroupProductionControl(parent, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger);
+            return;
+        }
+        if (!group.isProductionGroup())
             return;
 
         const auto& groupcontrols = group.productionControls(summaryState);
+        const std::vector<double>& groupTargetReductions = well_state.currentProductionGroupReductionRates(group.name());
 
         switch(currentGroupControl) {
         case Group2::ProductionCMode::NONE:
@@ -2113,9 +2125,7 @@ namespace Opm
         }
         case Group2::ProductionCMode::ORAT:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Oil];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+            double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]];
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::OIL, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::OIL, /*isInjector*/false, fraction);
 
@@ -2127,9 +2137,7 @@ namespace Opm
         }
         case Group2::ProductionCMode::WRAT:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Water];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+            double groupTargetReduction = groupTargetReductions[pu.phase_pos[Water]];
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::WAT, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::WAT, /*isInjector*/false, fraction);
 
@@ -2141,9 +2149,7 @@ namespace Opm
         }
         case Group2::ProductionCMode::GRAT:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Gas];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+           double groupTargetReduction = groupTargetReductions[pu.phase_pos[Gas]];
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::GAS, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::GAS, /*isInjector*/false, fraction);
             const double rate_target = std::max(0.0, groupcontrols.gas_target / efficiencyFactor - groupTargetReduction);
@@ -2154,11 +2160,7 @@ namespace Opm
         }
         case Group2::ProductionCMode::LRAT:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Water];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
-            phasePos = pu.phase_pos[Oil];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+            double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]] + groupTargetReductions[pu.phase_pos[Water]];
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::LIQ, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::LIQ, /*isInjector*/false, fraction);
 
@@ -2187,15 +2189,10 @@ namespace Opm
         }
         case Group2::ProductionCMode::FLD:
         {
-            // Produce share of parents control
-            const auto& parent = schedule.getGroup2( group.parent(), current_step_ );
-            if (group.getTransferGroupEfficiencyFactor())
-                efficiencyFactor *= group.getGroupEfficiencyFactor();
-
-            assembleGroupProductionControl(parent, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger);
+            // 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 );
         }

opm/simulators/wells/StandardWell_impl.hpp

@@ -961,15 +961,24 @@ namespace Opm
         const auto pu = phaseUsage();
         const auto& groupcontrols = group.injectionControls(summaryState);
         const Group2::InjectionCMode& currentGroupControl = well_state.currentInjectionGroupControl(group.name());
-        // it may have a injection group higher up in tree
-        if (!group.isInjectionGroup() && currentGroupControl != Group2::InjectionCMode::FLD)
+
+        if (currentGroupControl == Group2::InjectionCMode::FLD) {
+            // Inject share of parents control
+            const auto& parent = schedule.getGroup2( group.parent(), current_step_ );
+            if (group.getTransferGroupEfficiencyFactor())
+                efficiencyFactor *= group.getGroupEfficiencyFactor();
+
+            assembleGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, control_eq, efficiencyFactor, deferred_logger);
+            return;
+        }
+
+        if (!group.isInjectionGroup())
             return;
 
         int phasePos;
         Well2::GuideRateTarget wellTarget;
         Group2::GuideRateTarget groupTarget;
 
-
         switch (injectorType) {
         case Well2::InjectorType::WATER:
         {
@@ -996,8 +1005,8 @@ namespace Opm
             throw("Expected WATER, OIL or GAS as type for injectors " + well.name());
         }
 
-        double groupTargetReduction = 0.0;
-        wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/true, groupTargetReduction);
+        const std::vector<double>& groupTargetReductions = well_state.currentProductionGroupReductionRates(group.name());
+        double groupTargetReduction = groupTargetReductions[phasePos];
         double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, wellTarget, /*isInjector*/true);
         wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, groupTarget, /*isInjector*/true, fraction);
 
@@ -1024,9 +1033,7 @@ namespace Opm
         }
         case Group2::InjectionCMode::REIN:
         {
-            double productionRate = 0.0;
-            const Group2& reinGroup = schedule.getGroup2( groupcontrols.reinj_group, current_step_ );
-            productionRate += wellGroupHelpers::sumWellRates(reinGroup, schedule, well_state, current_step_, phasePos, /*isInjector*/false);
+            double productionRate = well_state.currentInjectionVREPRates(groupcontrols.reinj_group);
             productionRate /= efficiencyFactor;
             double target = std::max(0.0, (groupcontrols.target_reinj_fraction*productionRate - groupTargetReduction));
             control_eq = getWQTotal() - fraction * target;            
@@ -1037,11 +1044,7 @@ namespace Opm
             std::vector<double> convert_coeff(number_of_phases_, 1.0);
             Base::rateConverter_.calcCoeff(/*fipreg*/ 0, Base::pvtRegionIdx_, convert_coeff);
             double coeff = convert_coeff[phasePos];
-            double voidageRate = 0.0;
-            const Group2& voidageGroup = schedule.getGroup2( groupcontrols.voidage_group, current_step_ );
-            voidageRate += wellGroupHelpers::sumWellResRates(voidageGroup, schedule, well_state, current_step_, pu.phase_pos[BlackoilPhases::Aqua], /*injector*/false);
-            voidageRate += wellGroupHelpers::sumWellResRates(voidageGroup, schedule, well_state, current_step_, pu.phase_pos[BlackoilPhases::Liquid], /*injector*/false);
-            voidageRate += wellGroupHelpers::sumWellResRates(voidageGroup, schedule, well_state, current_step_, pu.phase_pos[BlackoilPhases::Vapour], /*injector*/false);
+            double voidageRate = well_state.currentInjectionVREPRates(groupcontrols.voidage_group);
             voidageRate /= efficiencyFactor;
             double target = std::max(0.0, ( groupcontrols.target_void_fraction*voidageRate/coeff - groupTargetReduction));
             control_eq = getWQTotal() - fraction * target;
@@ -1049,15 +1052,10 @@ namespace Opm
         }
         case Group2::InjectionCMode::FLD:
         {
-            // Inject share of parents control
-            const auto& parent = schedule.getGroup2( group.parent(), current_step_ );
-            if (group.getTransferGroupEfficiencyFactor())
-                efficiencyFactor *= group.getGroupEfficiencyFactor();
-
-            assembleGroupInjectionControl(parent, well_state, schedule, summaryState, injectorType, control_eq, efficiencyFactor, deferred_logger);
+            // The FLD case is handled earlier
+            assert(false);
             break;
         }
-
         default:
             OPM_DEFLOG_THROW(std::runtime_error, "Unvalid group control specified for group "  + well.groupName(), deferred_logger );
         }
@@ -1078,11 +1076,23 @@ namespace Opm
         const auto pu = phaseUsage();
 
         const Group2::ProductionCMode& currentGroupControl = well_state.currentProductionGroupControl(group.name());
-        if (!group.isProductionGroup() && currentGroupControl != Group2::ProductionCMode::FLD)
+        if (currentGroupControl == Group2::ProductionCMode::FLD) {
+            // Produce share of parents control
+            const auto& parent = schedule.getGroup2( group.parent(), current_step_ );
+            if (group.getTransferGroupEfficiencyFactor())
+                efficiencyFactor *= group.getGroupEfficiencyFactor();
+
+            assembleGroupProductionControl(parent, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger);
+            return;
+        }
+
+        if (!group.isProductionGroup())
             return;
 
         const auto& groupcontrols = group.productionControls(summaryState);
 
+        const std::vector<double>& groupTargetReductions = well_state.currentProductionGroupReductionRates(group.name());
+
         switch(currentGroupControl) {
         case Group2::ProductionCMode::NONE:
         {
@@ -1091,9 +1101,7 @@ namespace Opm
         }
         case Group2::ProductionCMode::ORAT:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Oil];                    
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+            double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]];
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::OIL, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::OIL, /*isInjector*/false, fraction);
 
@@ -1105,9 +1113,7 @@ namespace Opm
         }
         case Group2::ProductionCMode::WRAT:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Water];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+            double groupTargetReduction = groupTargetReductions[pu.phase_pos[Water]];
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::WAT, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::WAT, /*isInjector*/false, fraction);
 
@@ -1119,9 +1125,7 @@ namespace Opm
         }
         case Group2::ProductionCMode::GRAT:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Gas];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+            double groupTargetReduction = groupTargetReductions[pu.phase_pos[Gas]];
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::GAS, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::GAS, /*isInjector*/false, fraction);
 
@@ -1133,11 +1137,7 @@ namespace Opm
         }
         case Group2::ProductionCMode::LRAT:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Water];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
-            phasePos = pu.phase_pos[Oil];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+            double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]] + groupTargetReductions[pu.phase_pos[Water]];
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::LIQ, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::LIQ, /*isInjector*/false, fraction);
 
@@ -1156,13 +1156,9 @@ namespace Opm
         }
         case Group2::ProductionCMode::RESV:
         {
-            double groupTargetReduction = 0.0;
-            int phasePos = pu.phase_pos[Water];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
-            phasePos = pu.phase_pos[Oil];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
-            phasePos = pu.phase_pos[Gas];
-            wellGroupHelpers::computeGroupTargetReduction(group, well_state, schedule, current_step_, phasePos, /*isInjector*/false, groupTargetReduction);
+            double groupTargetReduction = groupTargetReductions[pu.phase_pos[Oil]]
+                    + groupTargetReductions[pu.phase_pos[Gas]]
+                    + groupTargetReductions[pu.phase_pos[Water]];
 
             double fraction = wellGroupHelpers::wellFractionFromGuideRates(well, schedule, well_state, current_step_, Base::guide_rate_, Well2::GuideRateTarget::RES, /*isInjector*/false);
             wellGroupHelpers::accumulateGroupFractions(well.groupName(), group.name(), schedule, well_state, current_step_, Base::guide_rate_, Group2::GuideRateTarget::RES, /*isInjector*/false, fraction);
@@ -1186,12 +1182,8 @@ namespace Opm
         }
         case Group2::ProductionCMode::FLD:
         {
-            // Produce share of parents control
-            const auto& parent = schedule.getGroup2( group.parent(), current_step_ );
-            if (group.getTransferGroupEfficiencyFactor())
-                efficiencyFactor *= group.getGroupEfficiencyFactor();
-
-            assembleGroupProductionControl(parent, well_state, schedule, summaryState, control_eq, efficiencyFactor, deferred_logger);
+            // The FLD case is handled earlier
+            assert(false);
             break;
         }
 

opm/simulators/wells/WellGroupHelpers.hpp

@@ -211,6 +211,77 @@ namespace Opm {
         return rate;
     }
 
+    inline void updateGroupTargetReduction(const Group2& group, const Schedule& schedule, const int reportStepIdx, const bool isInjector, WellStateFullyImplicitBlackoil& wellState, std::vector<double>& groupTargetReduction)
+    {
+        const int np = wellState.numPhases();
+        for (const std::string& groupName : group.groups()) {
+            std::vector<double> thisGroupTargetReduction(np, 0.0);
+            const Group2& groupTmp = schedule.getGroup2(groupName, reportStepIdx);
+            updateGroupTargetReduction(groupTmp, schedule, reportStepIdx, isInjector, wellState, thisGroupTargetReduction);
+
+            // accumulate group contribution from sub group
+            if (isInjector) {
+                const Group2::InjectionCMode& currentGroupControl = wellState.currentInjectionGroupControl(groupName);
+                if (currentGroupControl != Group2::InjectionCMode::FLD) {
+                    for (int phase = 0; phase < np; phase++) {
+                        groupTargetReduction[phase] += sumWellRates(group, schedule, wellState, reportStepIdx, phase, isInjector);
+                    }
+                    continue;
+                }
+            } else {
+                const Group2::ProductionCMode& currentGroupControl = wellState.currentProductionGroupControl(groupName);
+                if (currentGroupControl != Group2::ProductionCMode::FLD) {
+                    for (int phase = 0; phase < np; phase++) {
+                        groupTargetReduction[phase] += sumWellRates(group, schedule, wellState, reportStepIdx, phase, isInjector);
+                    }
+                    continue;
+                }
+            }
+            // or accumulate directly from the wells if controled from its parents
+            for (int phase = 0; phase < np; phase++) {
+                groupTargetReduction[phase] += thisGroupTargetReduction[phase];
+            }
+        }
+        //std::vector<double> groupTargetReduction(np, 0.0);
+        for (const std::string& wellName : group.wells()) {
+            const auto& wellTmp = schedule.getWell2(wellName, reportStepIdx);
+
+            if (wellTmp.isProducer() && isInjector)
+                continue;
+
+            if (wellTmp.isInjector() && !isInjector)
+                continue;
+
+            if (wellTmp.getStatus() == Well2::Status::SHUT)
+                continue;
+
+            const auto& end = wellState.wellMap().end();
+            const auto& it = wellState.wellMap().find( wellName );
+            if (it == end)  // the well is not found
+                continue;
+
+            int well_index = it->second[0];
+            const auto wellrate_index = well_index * wellState.numPhases();
+            // add contributino from wells not under group control
+            if (isInjector) {
+                if (wellState.currentInjectionControls()[well_index] != Well2::InjectorCMode::GRUP)
+                    for (int phase = 0; phase < np; phase++) {
+                        groupTargetReduction[phase] += wellState.wellRates()[wellrate_index + phase];
+                    }
+            } else {
+                if (wellState.currentProductionControls()[well_index] !=  Well2::ProducerCMode::GRUP)
+                    for (int phase = 0; phase < np; phase++) {
+                        groupTargetReduction[phase] -= wellState.wellRates()[wellrate_index + phase];
+                    }
+            }
+        }
+        if (isInjector)
+            wellState.setCurrentInjectionGroupReductionRates(group.name(), groupTargetReduction);
+        else
+            wellState.setCurrentProductionGroupReductionRates(group.name(), groupTargetReduction);
+    }
+
+
 
     inline void updateGuideRateForGroups(const Group2& group, const Schedule& schedule, const PhaseUsage& pu, const int reportStepIdx, const double& simTime, GuideRate* guideRate, WellStateFullyImplicitBlackoil& wellState) {
         for (const std::string& groupName : group.groups()) {
@@ -238,6 +309,37 @@ namespace Opm {
         guideRate->compute(group.name(), reportStepIdx, simTime, oilPot, gasPot, waterPot);
     }
 
+    inline void updateVREPForGroups(const Group2& group, const Schedule& schedule, const int reportStepIdx, WellStateFullyImplicitBlackoil& wellState, double& resv) {
+        for (const std::string& groupName : group.groups()) {
+            const Group2& groupTmp = schedule.getGroup2(groupName, reportStepIdx);
+            double thisResv = 0.0;
+            updateVREPForGroups(groupTmp, schedule, reportStepIdx, wellState, thisResv);
+            resv += thisResv;
+        }
+        const int np = wellState.numPhases();
+        for (int phase = 0; phase < np; ++phase) {
+            resv += sumWellPhaseRates(wellState.wellReservoirRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false);
+        }
+
+        wellState.setCurrentInjectionVREPRates(group.name(), resv);
+    }
+
+    inline void updateREINForGroups(const Group2& group, const Schedule& schedule, const int reportStepIdx, WellStateFullyImplicitBlackoil& wellState, std::vector<double>& rein) {
+        const int np = wellState.numPhases();
+        for (const std::string& groupName : group.groups()) {
+            const Group2& groupTmp = schedule.getGroup2(groupName, reportStepIdx);
+            std::vector<double> thisRein(np, 0.0);
+            updateREINForGroups(groupTmp, schedule, reportStepIdx, wellState, thisRein);
+            for (int phase = 0; phase < np; ++phase) {
+                rein[phase] = thisRein[phase];
+            }
+        }
+        for (int phase = 0; phase < np; ++phase) {
+            rein[phase] = sumWellPhaseRates(wellState.wellRates(), group, schedule, wellState, reportStepIdx, phase, /*isInjector*/ false);
+        }
+        wellState.setCurrentInjectionREINRates(group.name(), rein);
+    }
+
     inline double wellFractionFromGuideRates(const Well2& well, const Schedule& schedule, const WellStateFullyImplicitBlackoil& wellState, const int reportStepIdx, const GuideRate* guideRate, const Well2::GuideRateTarget& wellTarget, const bool isInjector) {
         double groupTotalGuideRate = 0.0;
         const Group2& groupTmp = schedule.getGroup2(well.groupName(), reportStepIdx);

opm/simulators/wells/WellStateFullyImplicitBlackoil.hpp

@@ -334,7 +334,7 @@ namespace Opm
             auto it = current_production_group_controls_.find(groupName);
 
             if (it == current_production_group_controls_.end())
-                OPM_THROW(std::logic_error, "Could not find any well control for " << groupName);
+                OPM_THROW(std::logic_error, "Could not find any control for production group " << groupName);
 
             return it->second;
         }
@@ -347,11 +347,61 @@ namespace Opm
             auto it = current_injection_group_controls_.find(groupName);
 
             if (it == current_injection_group_controls_.end())
-                OPM_THROW(std::logic_error, "Could not find any well control for " << groupName);
+                OPM_THROW(std::logic_error, "Could not find any control for injection group " << groupName);
 
             return it->second;
         }
 
+
+        void setCurrentProductionGroupReductionRates(const std::string& groupName, const std::vector<double>& target ) {
+            production_group_reduction_rates[groupName] = target;
+        }
+        const std::vector<double>& currentProductionGroupReductionRates(const std::string& groupName) const {
+            auto it = production_group_reduction_rates.find(groupName);
+
+            if (it == production_group_reduction_rates.end())
+                OPM_THROW(std::logic_error, "Could not find any reduction rates for production group  " << groupName);
+
+            return it->second;
+        }
+
+        void setCurrentInjectionGroupReductionRates(const std::string& groupName, const std::vector<double>& target ) {
+            injection_group_reduction_rates[groupName] = target;
+        }
+        const std::vector<double>& currentInjectionGroupReductionRates(const std::string& groupName) const {
+            auto it = injection_group_reduction_rates.find(groupName);
+
+            if (it == injection_group_reduction_rates.end())
+                OPM_THROW(std::logic_error, "Could not find any reduction rates for injection group " << groupName);
+
+            return it->second;
+        }
+
+        void setCurrentInjectionVREPRates(const std::string& groupName, const double& target ) {
+            injection_group_vrep_rates[groupName] = target;
+        }
+        const double& currentInjectionVREPRates(const std::string& groupName) const {
+            auto it = injection_group_vrep_rates.find(groupName);
+
+            if (it == injection_group_vrep_rates.end())
+                OPM_THROW(std::logic_error, "Could not find any VREP rates for group " << groupName);
+
+            return it->second;
+        }
+
+        void setCurrentInjectionREINRates(const std::string& groupName, const std::vector<double>& target ) {
+            injection_group_rein_rates[groupName] = target;
+        }
+        const std::vector<double>& currentInjectionREINRates(const std::string& groupName) const {
+            auto it = injection_group_rein_rates.find(groupName);
+
+            if (it == injection_group_rein_rates.end())
+                OPM_THROW(std::logic_error, "Could not find any REIN rates for group " << groupName);
+
+            return it->second;
+        }
+
+
         data::Wells report(const PhaseUsage &pu, const int* globalCellIdxMap) const override
         {
             data::Wells res = WellState::report(pu, globalCellIdxMap);
@@ -783,6 +833,10 @@ namespace Opm
         std::map<std::string, Group2::ProductionCMode> current_production_group_controls_;
         std::map<std::string, Group2::InjectionCMode> current_injection_group_controls_;
 
+        std::map<std::string, std::vector<double>> production_group_reduction_rates;
+        std::map<std::string, std::vector<double>> injection_group_reduction_rates;
+        std::map<std::string, double> injection_group_vrep_rates;
+        std::map<std::string, std::vector<double>> injection_group_rein_rates;
 
         std::vector<double> perfRateSolvent_;