Compute the group rates that takes into account limitation on the sublevels

opm/simulators/wells/BlackoilWellModelGeneric.cpp

@@ -1044,6 +1044,7 @@ BlackoilWellModelGeneric::
 gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
                           GLiftProdWells& prod_wells,
                           GLiftOptWells& glift_wells,
+                          GasLiftGroupInfo& group_info,
                           GLiftWellStateMap& glift_well_state_map,
                           const int episodeIndex)
 {
@@ -1056,6 +1057,7 @@ gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
                          this->groupState(),
                          prod_wells,
                          glift_wells,
+                         group_info,
                          glift_well_state_map,
                          this->glift_debug
     };

opm/simulators/wells/BlackoilWellModelGeneric.hpp

@@ -47,6 +47,7 @@ namespace Opm {
     class EclipseState;
     class GasLiftSingleWellGeneric;
     class GasLiftWellState;
+    class GasLiftGroupInfo;
     class Group;
     class GuideRateConfig;
     class ParallelWellInfo;
@@ -305,6 +306,7 @@ protected:
     void gasLiftOptimizationStage2(DeferredLogger& deferred_logger,
                                    GLiftProdWells& prod_wells,
                                    GLiftOptWells& glift_wells,
+                                   GasLiftGroupInfo& group_info,
                                    GLiftWellStateMap& map,
                                    const int episodeIndex);
 

opm/simulators/wells/BlackoilWellModel_impl.hpp

@@ -994,7 +994,7 @@ namespace Opm {
             gasLiftOptimizationStage1(
                 deferred_logger, prod_wells, glift_wells, group_info, state_map);
             gasLiftOptimizationStage2(
-                deferred_logger, prod_wells, glift_wells, state_map,
+                deferred_logger, prod_wells, glift_wells, group_info, state_map,
                 ebosSimulator_.episodeIndex());
             if (this->glift_debug) gliftDebugShowALQ(deferred_logger);
             num_wells_changed = glift_wells.size();

opm/simulators/wells/GasLiftGroupInfo.cpp

@@ -75,6 +75,28 @@ gasRate(const std::string& group_name) const
     return group_rate.gasRate();
 }
 
+double
+GasLiftGroupInfo::
+gasPotential(const std::string& group_name) const
+{
+    auto& group_rate = this->group_rate_map_.at(group_name);
+    return group_rate.gasPotential();
+}
+double
+GasLiftGroupInfo::
+waterPotential(const std::string& group_name) const
+{
+    auto& group_rate = this->group_rate_map_.at(group_name);
+    return group_rate.waterPotential();
+}
+double
+GasLiftGroupInfo::
+oilPotential(const std::string& group_name) const
+{
+    auto& group_rate = this->group_rate_map_.at(group_name);
+    return group_rate.oilPotential();
+}
+
 std::optional<double>
 GasLiftGroupInfo::
 gasTarget(const std::string& group_name) const
@@ -101,6 +123,24 @@ getRate(Rate rate_type, const std::string& group_name) const
         throw std::runtime_error("This should not happen");
     }
 }
+double
+GasLiftGroupInfo::
+getPotential(Rate rate_type, const std::string& group_name) const
+{
+    switch (rate_type) {
+    case Rate::oil:
+        return oilPotential(group_name);
+    case Rate::gas:
+        return gasPotential(group_name);
+    case Rate::water:
+        return waterPotential(group_name);
+    case Rate::liquid:
+        return oilPotential(group_name) + waterPotential(group_name);
+    default:
+        // Need this to avoid compiler warning : control reaches end of non-void function
+        throw std::runtime_error("This should not happen");
+    }
+}
 
 std::tuple<double, double, double, double>
 GasLiftGroupInfo::
@@ -438,36 +478,56 @@ displayDebugMessage_(const std::string &msg, const std::string &well_name)
 }
 
 
-std::tuple<double, double, double>
+std::tuple<double, double, double, double, double, double>
 GasLiftGroupInfo::
-getProducerWellRates_(int well_index)
+getProducerWellRates_(const Well* well, int well_index)
 {
     const auto& pu = this->phase_usage_;
     const auto& ws= this->well_state_.well(well_index);
     const auto& wrate = ws.well_potentials;
 
-    const auto oil_rate = pu.phase_used[Oil]
+    const auto oil_pot = pu.phase_used[Oil]
         ? wrate[pu.phase_pos[Oil]]
         : 0.0;
 
-    const auto gas_rate = pu.phase_used[Gas]
+    const auto gas_pot = pu.phase_used[Gas]
         ? wrate[pu.phase_pos[Gas]]
         : 0.0;
 
-    const auto water_rate = pu.phase_used[Water]
+    const auto water_pot = pu.phase_used[Water]
         ? wrate[pu.phase_pos[Water]]
         : 0.0;
 
-    return {oil_rate, gas_rate, water_rate};
+    const auto controls = well->productionControls(this->summary_state_);
+    double oil_rate = oil_pot;
+    if (controls.hasControl(Well::ProducerCMode::ORAT)) {
+        oil_rate = std::min(controls.oil_rate, oil_rate);
+    }
+    double gas_rate = gas_pot;
+    if (controls.hasControl(Well::ProducerCMode::GRAT)) {
+        gas_rate = std::min(controls.gas_rate, gas_rate);
+    }
+    double water_rate = water_pot;
+    if (controls.hasControl(Well::ProducerCMode::WRAT)) {
+        water_rate = std::min(controls.water_rate, water_rate);
+    }
+    if (controls.hasControl(Well::ProducerCMode::LRAT)) {
+        double liquid_rate = oil_rate + water_rate;
+        double liquid_rate_lim = std::min(controls.liquid_rate, liquid_rate);
+        water_rate = water_rate / liquid_rate * liquid_rate_lim;
+        oil_rate = oil_rate / liquid_rate * liquid_rate_lim;
     }
 
-std::tuple<double, double, double, double>
+    return {oil_rate, gas_rate, water_rate, oil_pot, gas_pot, water_pot};
+}
+
+std::tuple<double, double, double, double, double, double, double>
 GasLiftGroupInfo::
 initializeGroupRatesRecursive_(const Group &group)
 {
-    std::array<double,4> rates{};
+    std::array<double,7> rates{};
     if (this->debug) debugStartInitializeGroup(group.name());
-    auto& [oil_rate, water_rate, gas_rate, alq] = rates;
+    auto& [oil_rate, water_rate, gas_rate, oil_potential, water_potential, gas_potential, alq] = rates;
     if (group.wellgroup()) {
         for (const std::string& well_name : group.wells()) {
             // NOTE: we cannot simply use:
@@ -482,17 +542,21 @@ initializeGroupRatesRecursive_(const Group &group)
                 assert(well); // Should never be nullptr
                 const int index = (itr->second).second;
                 if (well->isProducer()) {
-                    auto [sw_oil_rate, sw_gas_rate, sw_water_rate] = getProducerWellRates_(index);
+                    auto [sw_oil_rate, sw_gas_rate, sw_water_rate, sw_oil_pot, sw_gas_pot, sw_water_pot] = getProducerWellRates_(well, index);
                     auto sw_alq = this->well_state_.getALQ(well_name);
                     double factor = well->getEfficiencyFactor();
                     oil_rate += (factor * sw_oil_rate);
                     gas_rate += (factor * sw_gas_rate);
                     water_rate += (factor * sw_water_rate);
+                    oil_potential += (factor * sw_oil_pot);
+                    gas_potential += (factor * sw_gas_pot);
+                    water_potential += (factor * sw_water_pot);
+
                     alq += (factor * sw_alq);
                     if (this->debug) {
                         debugDisplayWellContribution_(
                             group.name(), well_name, factor,
-                            sw_oil_rate, sw_gas_rate, sw_water_rate, sw_alq,
+                            sw_oil_pot, sw_gas_pot, sw_water_pot, sw_alq,
                             oil_rate, gas_rate, water_rate, alq
                         );
                     }
@@ -507,12 +571,16 @@ initializeGroupRatesRecursive_(const Group &group)
                 continue;
             const Group& sub_group = this->schedule_.getGroup(
                 group_name, this->report_step_idx_);
-            auto [sg_oil_rate, sg_gas_rate, sg_water_rate, sg_alq]
+            auto [sg_oil_rate, sg_gas_rate, sg_water_rate,
+                  sg_oil_pot, sg_gas_pot, sg_water_pot, sg_alq]
                 = initializeGroupRatesRecursive_(sub_group);
             const auto gefac = sub_group.getGroupEfficiencyFactor();
             oil_rate += (gefac * sg_oil_rate);
             gas_rate += (gefac * sg_gas_rate);
             water_rate += (gefac * sg_water_rate);
+            oil_potential += (gefac * sg_oil_pot);
+            gas_potential += (gefac * sg_gas_pot);
+            water_potential += (gefac * sg_water_pot);
             alq += (gefac * sg_alq);
         }
     }
@@ -538,14 +606,29 @@ initializeGroupRatesRecursive_(const Group &group)
     }
     if (oil_target || liquid_target || water_target || gas_target || max_total_gas || max_alq) {
         updateGroupIdxMap_(group.name());
+        if(oil_target)
+            oil_rate = std::min(oil_rate, *oil_target);
+        if(gas_target)
+            gas_rate = std::min(gas_rate, *gas_target);
+        if(water_target)
+            water_rate = std::min(water_rate, *water_target);
+        if(liquid_target) {
+            double liquid_rate = oil_rate + water_rate;
+            double liquid_rate_limited = std::min(liquid_rate, *liquid_target);
+            oil_rate = oil_rate / liquid_rate * liquid_rate_limited;
+            water_rate = water_rate / liquid_rate * liquid_rate_limited;
+        }
+
         this->group_rate_map_.try_emplace(group.name(),
-            oil_rate, gas_rate, water_rate, alq, oil_target, gas_target, water_target, liquid_target, max_total_gas, max_alq);
+            oil_rate, gas_rate, water_rate, alq,
+            oil_potential, gas_potential, water_potential,
+            oil_target, gas_target, water_target, liquid_target, max_total_gas, max_alq);
         if (this->debug) {
             debugDisplayUpdatedGroupRates(
                 group.name(), oil_rate, gas_rate, water_rate, alq);
         }
     }
-    return std::make_tuple(oil_rate, gas_rate, water_rate, alq);
+    return std::make_tuple(oil_rate, gas_rate, water_rate, oil_potential, gas_potential, water_potential, alq);
 }
 
 void

opm/simulators/wells/GasLiftGroupInfo.hpp

@@ -91,8 +91,12 @@ public:
 
     double alqRate(const std::string& group_name);
     double gasRate(const std::string& group_name) const;
+    double gasPotential(const std::string& group_name) const;
+    double waterPotential(const std::string& group_name) const;
+    double oilPotential(const std::string& group_name) const;
     int getGroupIdx(const std::string& group_name);
     double getRate(Rate rate_type, const std::string& group_name) const;
+    double getPotential(Rate rate_type, const std::string& group_name) const;
     std::tuple<double,double,double,double> getRates(const int group_idx) const;
     std::optional<double> gasTarget(const std::string& group_name) const;
     std::optional<double> getTarget(
@@ -130,8 +134,9 @@ protected:
     void debugStartInitializeGroup(const std::string& name) const;
     void displayDebugMessage_(const std::string& msg) const override;
     void displayDebugMessage_(const std::string& msg, const std::string& well_name);
-    std::tuple<double, double, double> getProducerWellRates_(const int index);
+    std::tuple<double, double, double, double, double, double>
-    std::tuple<double, double, double, double>
+        getProducerWellRates_(const Well* well, const int index);
+    std::tuple<double, double, double, double, double, double, double>
         initializeGroupRatesRecursive_(const Group &group);
     void initializeWell2GroupMapRecursive_(
         const Group& group, std::vector<std::string>& group_names,
@@ -142,6 +147,7 @@ protected:
     class GroupRates {
     public:
         GroupRates( double oil_rate, double gas_rate, double water_rate, double alq,
+            double oil_potential, double gas_potential, double water_potential,
             std::optional<double> oil_target,
             std::optional<double> gas_target,
             std::optional<double> water_target,
@@ -153,6 +159,9 @@ protected:
             gas_rate_{gas_rate},
             water_rate_{water_rate},
             alq_{alq},
+            oil_potential_{oil_potential},
+            gas_potential_{gas_potential},
+            water_potential_{water_potential},
             oil_target_{oil_target},
             gas_target_{gas_target},
             water_target_{water_target},
@@ -177,6 +186,9 @@ protected:
         double oilRate() const { return oil_rate_; }
         std::optional<double> oilTarget() const { return oil_target_; }
         std::optional<double> liquidTarget() const { return liquid_target_; }
+        double oilPotential() const { return oil_potential_; }
+        double gasPotential() const { return gas_potential_; }
+        double waterPotential() const { return water_potential_; }
 
         void update(double delta_oil, double delta_gas, double delta_water, double delta_alq)
         {
@@ -190,6 +202,9 @@ protected:
         double gas_rate_;
         double water_rate_;
         double alq_;
+        double oil_potential_;
+        double gas_potential_;
+        double water_potential_;
         std::optional<double> oil_target_;
         std::optional<double> gas_target_;
         std::optional<double> water_target_;

opm/simulators/wells/GasLiftSingleWellGeneric.cpp

@@ -247,7 +247,7 @@ GasLiftSingleWellGeneric::checkGroupTargetsViolated(const BasicRates& rates, con
             auto target_opt = this->group_info_.getTarget(rate_type, group_name);
             if (target_opt) {
                 auto delta_rate = new_rates[rate_type] - rates[rate_type];
-                auto new_group_rate = this->group_info_.getRate(rate_type, group_name) + efficiency * delta_rate;
+                auto new_group_rate = this->group_info_.getPotential(rate_type, group_name) + efficiency * delta_rate;
                 if (new_group_rate > *target_opt) {
                     if (this->debug) {
                         const std::string msg

opm/simulators/wells/GasLiftStage2.cpp

@@ -28,7 +28,7 @@
 #include <opm/simulators/wells/GasLiftWellState.hpp>
 #include <opm/simulators/wells/WellInterfaceGeneric.hpp>
 #include <opm/simulators/wells/WellState.hpp>
-
+#include <opm/simulators/wells/GasLiftGroupInfo.hpp>
 #include <cmath>
 #include <optional>
 #include <string>
@@ -47,12 +47,14 @@ GasLiftStage2::GasLiftStage2(
     const GroupState &group_state,
     GLiftProdWells &prod_wells,
     GLiftOptWells &glift_wells,
+    GasLiftGroupInfo& group_info,
     GLiftWellStateMap &state_map,
     bool glift_debug
 ) :
     GasLiftCommon(well_state, group_state, deferred_logger, comm, glift_debug)
     , prod_wells_{prod_wells}
     , stage1_wells_{glift_wells}
+    , group_info_{group_info}
     , well_state_map_{state_map}
     , report_step_idx_{report_step_idx}
     , summary_state_{summary_state}
@@ -278,136 +280,10 @@ std::tuple<double, double, double, double>
 GasLiftStage2::
 getCurrentGroupRates_(const Group &group)
 {
-    auto rates = getCurrentGroupRatesRecursive_(group);
+    return {this->group_info_.oilRate(group.name()),
-    this->comm_.sum(rates.data(), rates.size());
+            this->group_info_.gasRate(group.name()),
-    auto [oil_rate, gas_rate, water_rate, alq] = rates;
+            this->group_info_.waterRate(group.name()),
-    if (this->debug) {
+            this->group_info_.alqRate(group.name())};
-        const std::string msg = fmt::format(
-            "Current group rates for {} : oil: {}, gas: {}, water: {}, alq: {}",
-            group.name(), oil_rate, gas_rate, water_rate, alq);
-        displayDebugMessageOnRank0_(msg);
-    }
-
-    return {oil_rate, gas_rate, water_rate, alq};
-}
-
-
-std::array <double, 4>
-GasLiftStage2::
-getCurrentGroupRatesRecursive_(const Group &group)
-{
-    double oil_rate = 0.0;
-    double gas_rate = 0.0;
-    double water_rate = 0.0;
-    double alq = 0.0;
-    // NOTE: A group can either contain wells or groups, but not both
-    if (group.wellgroup()) {
-        for (const std::string& well_name : group.wells()) {
-            auto [sw_oil_rate, sw_gas_rate, sw_water_rate, sw_alq] =
-                getCurrentWellRates_(well_name, group.name());
-            oil_rate += sw_oil_rate;
-            gas_rate += sw_gas_rate;
-            water_rate += sw_water_rate;
-            alq += sw_alq;
-        }
-
-    }
-    else {
-        for (const std::string& group_name : group.groups()) {
-            if(this->schedule_.back().groups.has(group_name)) {
-                const Group& sub_group =
-                    this->schedule_.getGroup(group_name, this->report_step_idx_);
-                // If groups have efficiency factors to model
-                // synchronized downtime of their subordinate wells
-                // (see keyword GEFAC), their lift gas injection rates
-                // are multiplied by their efficiency factors when
-                // they are added to the lift gas supply rate of the
-                // parent group.
-                const auto gefac = sub_group.getGroupEfficiencyFactor();
-                auto rates = getCurrentGroupRatesRecursive_(sub_group);
-                auto [sg_oil_rate, sg_gas_rate, sg_water_rate, sg_alq] = rates;
-                oil_rate += (gefac * sg_oil_rate);
-                gas_rate += (gefac * sg_gas_rate);
-                water_rate += (gefac * sg_water_rate);
-                alq += (gefac * sg_alq);
-            }
-        }
-    }
-    return {oil_rate, gas_rate, water_rate, alq};
-}
-
-std::tuple<double, double, double, double>
-GasLiftStage2::
-getCurrentWellRates_(const std::string &well_name, const std::string &group_name)
-{
-    double oil_rate, gas_rate, water_rate, alq;
-    bool success = false;
-    const WellInterfaceGeneric *well_ptr = nullptr;
-    std::string debug_info;
-    if (this->stage1_wells_.count(well_name) == 1) {
-        GasLiftSingleWell &gs_well = *(this->stage1_wells_.at(well_name).get());
-        const WellInterfaceGeneric &well = gs_well.getWell();
-        well_ptr = &well;
-        GasLiftWellState &state = *(this->well_state_map_.at(well_name).get());
-        std::tie(oil_rate, gas_rate) = state.getRates();
-        water_rate = state.waterRate();
-        success = true;
-        if ( this->debug) debug_info = "(A)";
-    }
-    else if (this->prod_wells_.count(well_name) == 1) {
-        well_ptr = this->prod_wells_.at(well_name);
-        std::tie(oil_rate, gas_rate, water_rate) = getWellRates_(*well_ptr);
-        success = true;
-        if ( this->debug) debug_info = "(B)";
-    }
-
-    if (well_ptr) {
-        // we only want rates from wells owned by the rank
-        if (!well_state_.wellIsOwned(well_ptr->indexOfWell(), well_name)) {
-            success = false;
-        }
-    }
-
-    if (success) {
-        assert(well_ptr);
-        assert(well_ptr->isProducer());
-        alq = this->well_state_.getALQ(well_name);
-        if (this->debug) {
-            const std::string msg = fmt::format(
-                "Rates {} for well {} : oil: {}, gas: {}, water: {}, alq: {}",
-                debug_info, well_name, oil_rate, gas_rate, water_rate, alq);
-            displayDebugMessage_(msg, group_name);
-        }
-        // If wells have efficiency factors to take account of regular
-        // downtime (see keyword WEFAC), their lift gas injection
-        // rates are multiplied by their efficiency factors when they
-        // are added to the group lift gas supply rate. This is
-        // consistent with the summation of flow rates for wells with
-        // downtime, and preserves the ratio of production rate to
-        // lift gas injection rate.
-        const auto &well_ecl = well_ptr->wellEcl();
-        double factor = well_ecl.getEfficiencyFactor();
-        oil_rate *= factor;
-        gas_rate *= factor;
-        water_rate *= factor;
-        alq *= factor;
-        if (this->debug && (factor != 1)) {
-            const std::string msg = fmt::format(
-                "Well {} : efficiency factor {}. New rates : oil: {}, gas: {}, water: {}, alq: {}",
-                well_name, factor, oil_rate, gas_rate, water_rate, alq);
-            displayDebugMessage_(msg, group_name);
-        }
-    }
-    else {
-        // NOTE: This happens for wells that are not producers, or not active.
-        if (this->debug) {
-            const std::string msg = fmt::format("Could not determine current rates for "
-                "well {}: (not active or injector)", well_name);
-            displayDebugMessage_(msg, group_name);
-        }
-        oil_rate = 0.0; gas_rate = 0.0; water_rate = 0.0; alq = 0.0;
-    }
-    return std::make_tuple(oil_rate, gas_rate, water_rate, alq);
 }
 
 std::optional<double>
@@ -430,24 +306,6 @@ GasLiftStage2::getGroupMaxTotalGas_(const Group &group)
     return std::nullopt; // If GLIFTOPT is missing from schedule, assume unlimited alq
 }
 
-std::tuple<double, double, double>
-GasLiftStage2::
-getWellRates_(const WellInterfaceGeneric &well)
-{
-    const int well_index = well.indexOfWell();
-    const auto& ws = this->well_state_.well(well_index);
-    const auto& pu = well.phaseUsage();
-    auto oil_rate = ws.well_potentials[pu.phase_pos[Oil]];
-    auto water_rate = ws.well_potentials[pu.phase_pos[Water]];
-    double gas_rate = 0.0;
-    // See comment for setupPhaseVariables_() in GasLiftSingleWell_impl.hpp
-    //  about the two-phase oil-water case.
-    if (pu.phase_used[BlackoilPhases::Vapour]) {
-        gas_rate = ws.well_potentials[pu.phase_pos[Gas]];
-    }
-    return {oil_rate, gas_rate, water_rate};
-}
-
 // Find all subordinate wells of a given group.
 //
 // NOTE: A group can either contain wells or groups, not both.

opm/simulators/wells/GasLiftStage2.hpp

@@ -70,6 +70,7 @@ public:
         const GroupState& group_state,
         GLiftProdWells& prod_wells,
         GLiftOptWells& glift_wells,
+        GasLiftGroupInfo& group_info,
         GLiftWellStateMap& state_map,
         bool glift_debug
     );
@@ -89,16 +90,12 @@ protected:
     void displayWarning_(const std::string& msg, const std::string& group_name);
     void displayWarning_(const std::string& msg);
     std::tuple<double, double, double, double> getCurrentGroupRates_(const Group& group);
-    std::array<double,4> getCurrentGroupRatesRecursive_(const Group& group);
-    std::tuple<double, double, double, double> getCurrentWellRates_(
-        const std::string& well_name, const std::string& group_name);
     std::optional<double> getGroupMaxALQ_(const Group &group);
     std::optional<double> getGroupMaxTotalGas_(const Group &group);
     std::vector<GasLiftSingleWell *> getGroupGliftWells_(
         const Group& group);
     void getGroupGliftWellsRecursive_(
         const Group& group, std::vector<GasLiftSingleWell *>& wells);
-    std::tuple<double, double, double> getWellRates_(const WellInterfaceGeneric& well);
     void optimizeGroup_(const Group& group);
     void optimizeGroupsRecursive_(const Group& group);
     void recalculateGradientAndUpdateData_(
@@ -126,6 +123,7 @@ protected:
 
     GLiftProdWells& prod_wells_;
     GLiftOptWells& stage1_wells_;
+    GasLiftGroupInfo& group_info_;
     GLiftWellStateMap& well_state_map_;
 
     int report_step_idx_;