Merge pull request #3748 from hakonhagland/glift_swg_fx2

Refactor some methods in `GasLiftSingleWellGeneric.cpp`

opm/simulators/wells/GasLiftGroupInfo.cpp

@@ -84,15 +84,54 @@ gasTarget(const std::string& group_name)
     return group_rate.gasTarget();
 }
 
+double
+GasLiftGroupInfo::
+getRate(Rate rate_type, const std::string& group_name)
+{
+    switch (rate_type) {
+    case Rate::oil:
+        return oilRate(group_name);
+    case Rate::gas:
+        return gasRate(group_name);
+    case Rate::water:
+        return waterRate(group_name);
+    case Rate::liquid:
+        return oilRate(group_name) + waterRate(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::
-getRates(int group_idx)
+getRates(const int group_idx)
 {
     const auto& group_name = groupIdxToName(group_idx);
     auto& rates = this->group_rate_map_.at(group_name);
     return std::make_tuple(rates.oilRate(), rates.gasRate(), rates.waterRate(), rates.alq());
 }
 
+std::optional<double>
+GasLiftGroupInfo::
+getTarget(Rate rate_type, const std::string& group_name)
+{
+    switch (rate_type) {
+    case Rate::oil:
+        return oilTarget(group_name);
+    case Rate::gas:
+        return gasTarget(group_name);
+    case Rate::water:
+        return waterTarget(group_name);
+    case Rate::liquid:
+        return liquidTarget(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::vector<std::pair<std::string,double>>&
 GasLiftGroupInfo::
 getWellGroups(const std::string& well_name)
@@ -189,6 +228,23 @@ liquidTarget(const std::string &group_name)
     return group_rate.liquidTarget();
 }
 
+const std::string
+GasLiftGroupInfo::
+rateToString(Rate rate) {
+    switch (rate) {
+    case Rate::oil:
+        return "oil";
+    case Rate::gas:
+        return "gas";
+    case Rate::water:
+        return "water";
+    case Rate::liquid:
+        return "liquid";
+    default:
+        throw std::runtime_error("This should not happen");
+    }
+}
+
 void
 GasLiftGroupInfo::
 update(

opm/simulators/wells/GasLiftGroupInfo.hpp

@@ -67,6 +67,8 @@ class GasLiftGroupInfo
     static const int Oil = BlackoilPhases::Liquid;
     static const int Gas = BlackoilPhases::Vapour;
 public:
+    enum class Rate {oil, gas, water, liquid};
+
     using GLiftEclWells = std::map<std::string,std::pair<const Well *,int>>;
     GasLiftGroupInfo(
         GLiftEclWells& ecl_wells,
@@ -84,13 +86,16 @@ public:
     double alqRate(const std::string& group_name);
     double gasRate(const std::string& group_name);
     int getGroupIdx(const std::string& group_name);
-    std::tuple<double,double,double,double> getRates(int group_idx);
+    double getRate(Rate rate_type, const std::string& group_name);
+    std::tuple<double,double,double,double> getRates(const int group_idx);
     std::optional<double> gasTarget(const std::string& group_name);
+    std::optional<double> getTarget(Rate rate_type, const std::string& group_name);
     const std::string& groupIdxToName(int group_idx);
     bool hasWell(const std::string& well_name);
     void initialize();
     std::optional<double> maxAlq(const std::string& group_name);
     double oilRate(const std::string& group_name);
+    static const std::string rateToString(Rate rate);
     double waterRate(const std::string& group_name);
     std::optional<double> oilTarget(const std::string& group_name);
     std::optional<double> waterTarget(const std::string& group_name);

opm/simulators/wells/GasLiftSingleWellGeneric.cpp

@@ -601,105 +601,34 @@ getOilRateWithLimit_(const std::vector<double>& potentials) const
     return { new_rate, limited};
 }
 
-
 std::pair<double, bool>
 GasLiftSingleWellGeneric::
-getOilRateWithGroupLimit_(const double new_oil_rate, const double oil_rate) const
+getOilRateWithGroupLimit_(double new_oil_rate, double oil_rate) const
 {
-    const double delta_oil = new_oil_rate - oil_rate;
+    [[maybe_unused]] auto [rate, gr_name, efficiency]
-    if (delta_oil > 0) {
+        = getRateWithGroupLimit_(Rate::oil, new_oil_rate, oil_rate);
-      // It is required that the production rate for a given group is
+    bool limited = gr_name != nullptr;
-      // is less than or equal to its target rate, see assert() below.
+    return {rate, limited};
-      // Then it only makes sense to check if the group target is exceeded
-      //  if delta_oil > 0
-      const auto &pairs =
-          this->group_info_.getWellGroups(this->well_name_);
-      bool limit = false;
-      double limited_oil_rate = new_oil_rate;
-      double gr_oil_target, new_gr_oil_rate;
-      const std::string *group_name = nullptr;
-      for (const auto& [group_name_temp, efficiency] : pairs) {
-          auto gr_oil_target_opt = this->group_info_.oilTarget(group_name_temp);
-          if (gr_oil_target_opt) {
-            double gr_oil_target_temp = *gr_oil_target_opt;
-            double gr_oil_rate_temp =
-                this->group_info_.oilRate(group_name_temp);
-            assert(gr_oil_rate_temp <= gr_oil_target_temp);
-            double new_gr_oil_rate_temp = gr_oil_rate_temp + efficiency * delta_oil;
-            if (new_gr_oil_rate_temp > gr_oil_target_temp) {
-                double limited_oil_rate_temp =
-                    oil_rate + (gr_oil_target_temp - gr_oil_rate_temp) / efficiency;
-                if (limited_oil_rate_temp < limited_oil_rate) {
-                    limit = true;
-                    group_name = &group_name_temp;
-                    limited_oil_rate = limited_oil_rate_temp;
-                    gr_oil_target = gr_oil_target_temp;
-                    new_gr_oil_rate = new_gr_oil_rate_temp;
-                }
-            }
-         }
-      }
-      if (this->debug_ && limit) {
-          const std::string msg = fmt::format(
-              "limiting oil rate from {} to {} to meet group target {} "
-              "for group {}. Computed group rate was: {}",
-              new_oil_rate, limited_oil_rate, gr_oil_target,
-              *group_name, new_gr_oil_rate);
-          displayDebugMessage_(msg);
-          return { limited_oil_rate, /*limit=*/true};
-      }
-    }
-    return { new_oil_rate, /*limit=*/false};
 }
 
 std::pair<double, bool>
 GasLiftSingleWellGeneric::
-getGasRateWithGroupLimit_(const double new_gas_rate, const double gas_rate) const
+getGasRateWithGroupLimit_(double new_gas_rate, double gas_rate) const
 {
-    const double delta_gas = new_gas_rate - gas_rate;
+    [[maybe_unused]] auto [rate, gr_name, efficiency]
-    const auto &pairs =
+        = getRateWithGroupLimit_(Rate::gas, new_gas_rate, gas_rate);
-        this->group_info_.getWellGroups(this->well_name_);
+    bool limited = gr_name != nullptr;
-    for (const auto &[group_name, efficiency] : pairs) {
+    return {rate, limited};
-        auto gr_gas_target_opt = this->group_info_.gasTarget(group_name);
-        if (gr_gas_target_opt) {
-            double gr_gas_rate =
-                this->group_info_.gasRate(group_name);
-            double new_gr_gas_rate = gr_gas_rate + efficiency * delta_gas;
-            if (new_gr_gas_rate > *gr_gas_target_opt) {
-                const std::string msg = fmt::format("limiting gas rate to group target: "
-                    "computed group rate: {}, target: {}", new_gr_gas_rate, *gr_gas_target_opt);
-                displayDebugMessage_(msg);
-                double new_rate = gas_rate + (*gr_gas_target_opt - gr_gas_rate) / efficiency;
-                return { std::min(new_rate, new_gas_rate), /*limit=*/true};
-            }
-        }
-    }
-    return { new_gas_rate, /*limit=*/false};
 }
 
 std::pair<double, bool>
 GasLiftSingleWellGeneric::
-getWaterRateWithGroupLimit_(const double new_water_rate, const double water_rate) const
+getWaterRateWithGroupLimit_(double new_water_rate, double water_rate) const
 {
-    const double delta_water = new_water_rate - water_rate;
+    [[maybe_unused]] auto [rate, gr_name, efficiency] = getRateWithGroupLimit_(
-    const auto &pairs =
+                                            Rate::water, new_water_rate, water_rate);
-        this->group_info_.getWellGroups(this->well_name_);
+    bool limited = gr_name != nullptr;
-    for (const auto &[group_name, efficiency] : pairs) {
+    return {rate, limited};
-        auto gr_water_target_opt = this->group_info_.waterTarget(group_name);
-        if (gr_water_target_opt) {
-            double gr_water_rate =
-                this->group_info_.waterRate(group_name);
-            double new_gr_water_rate = gr_water_rate + efficiency * delta_water;
-            if (new_gr_water_rate > *gr_water_target_opt) {
-                const std::string msg = fmt::format("limiting water rate to group target: "
-                    "computed group rate: {}, target: {}", new_gr_water_rate, *gr_water_target_opt);
-                displayDebugMessage_(msg);
-                double new_rate = water_rate + (*gr_water_target_opt - gr_water_rate) / efficiency;
-                return { std::min(new_rate, new_water_rate), /*limit=*/true};
-            }
-        }
-    }
-    return { new_water_rate, /*limit=*/false};
 }
 
 std::tuple<double, double, bool, bool>
@@ -707,33 +636,90 @@ GasLiftSingleWellGeneric::
 getLiquidRateWithGroupLimit_(const double new_oil_rate, const double oil_rate,
                              const double new_water_rate, const double water_rate) const
 {
-    const double delta_water = new_water_rate - water_rate;
+    auto liquid_rate = oil_rate + water_rate;
-    const double delta_oil = new_oil_rate - oil_rate;
+    auto new_liquid_rate = new_oil_rate + new_water_rate;
-    const auto &pairs =
+    auto [liquid_rate_limited, group_name, efficiency]
-        this->group_info_.getWellGroups(this->well_name_);
+        = getRateWithGroupLimit_(Rate::liquid, new_liquid_rate, liquid_rate);
-    for (const auto &[group_name, efficiency] : pairs) {
+    bool limited = group_name != nullptr;
-        auto gr_liquid_target_opt = this->group_info_.liquidTarget(group_name);
+    if (limited) {
-        if (gr_liquid_target_opt) {
+        // the oil, gas, and water cases can be handled directly by
-            double gr_water_rate =
+        //  getRateWithGroupLimit_() above. However, for the liquid case
-                this->group_info_.waterRate(group_name);
+        //  we must do some postprocessing. I chose to include it here
-            double gr_oil_rate =
+        //  instead of cluttering up getRateWithGroupLimit_() with this
-                this->group_info_.oilRate(group_name);
+        //  special case.
+        double delta_water = new_water_rate - water_rate;
+        double delta_oil = new_oil_rate - oil_rate;
+
+        double gr_water_rate = this->group_info_.waterRate(*group_name);
+        double gr_oil_rate = this->group_info_.oilRate(*group_name);
+
+        // NOTE: these rates are too large according to the limited liquid rate
+        //  but it does not matter since we are only using them to calculate
+        //  the fraction of the liquid corresponding to the oil phase
         double new_gr_water_rate = gr_water_rate + efficiency * delta_water;
         double new_gr_oil_rate = gr_oil_rate + efficiency * delta_oil;
         double new_gr_liquid_rate = new_gr_water_rate + new_gr_oil_rate;
-            if (new_gr_liquid_rate > *gr_liquid_target_opt) {
+
-                const std::string msg = fmt::format("limiting liquid rate to group target: "
-                    "computed group rate: {}, target: {}", new_gr_liquid_rate, *gr_liquid_target_opt);
-                displayDebugMessage_(msg);
         double oil_fraction = new_gr_oil_rate / new_gr_liquid_rate;
-                double water_rate_limited = water_rate + (1.0 - oil_fraction) * (new_gr_liquid_rate - *gr_liquid_target_opt) / efficiency;
+        double delta_liquid = liquid_rate_limited - liquid_rate;
-                double oil_rate_limited = oil_rate + oil_fraction * (new_gr_liquid_rate - *gr_liquid_target_opt) / efficiency;
+        auto limited_oil_rate = oil_rate + oil_fraction * delta_liquid;
-                return { std::min(oil_rate_limited, new_oil_rate), std::min(water_rate_limited, new_water_rate), /*limit=*/true, /*limit=*/true};
+        auto limited_water_rate = water_rate + (1.0 - oil_fraction) * delta_liquid;
+        return {limited_oil_rate, limited_water_rate, limited, limited};
+    }
+    return {new_oil_rate, new_water_rate, limited, limited};
+}
+
+std::tuple<double, const std::string*, double>
+GasLiftSingleWellGeneric::
+getRateWithGroupLimit_(
+    Rate rate_type, const double new_rate, const double old_rate) const
+{
+    const double delta_rate = new_rate - old_rate;
+    if (delta_rate > 0) {
+      // It is required that the production rate for a given group is
+      // is less than or equal to its target rate, see assert() below.
+      // Then it only makes sense to check if the group target is exceeded
+      //  if delta_rate > 0
+      const auto &pairs =
+          this->group_info_.getWellGroups(this->well_name_);
+      double limited_rate = new_rate;
+      double gr_target, new_gr_rate, efficiency;
+      const std::string *group_name = nullptr;
+      for (const auto& [group_name_temp, efficiency_temp] : pairs) {
+          auto gr_target_opt = this->group_info_.getTarget(rate_type, group_name_temp);
+          if (gr_target_opt) {
+            double gr_target_temp = *gr_target_opt;
+            double gr_rate_temp =
+                this->group_info_.getRate(rate_type, group_name_temp);
+            assert(gr_rate_temp <= gr_target_temp);
+            double new_gr_rate_temp = gr_rate_temp + efficiency_temp * delta_rate;
+            if (new_gr_rate_temp > gr_target_temp) {
+                double limited_rate_temp =
+                    old_rate + (gr_target_temp - gr_rate_temp) / efficiency_temp;
+                if (limited_rate_temp < limited_rate) {
+                    group_name = &group_name_temp;
+                    efficiency = efficiency_temp;
+                    limited_rate = limited_rate_temp;
+                    gr_target = gr_target_temp;
+                    new_gr_rate = new_gr_rate_temp;
                 }
             }
          }
-    return { new_oil_rate, new_water_rate, /*limit=*/false, /*limit=*/false};
       }
+      if (this->debug_ && group_name) {
+          const std::string msg = fmt::format(
+              "limiting {} rate from {} to {} to meet group target {} "
+              "for group {}. Computed group rate was: {}",
+              GasLiftGroupInfo::rateToString(rate_type),
+              new_rate, limited_rate, gr_target,
+              *group_name, new_gr_rate);
+          displayDebugMessage_(msg);
+          return { limited_rate, group_name, efficiency };
+      }
+    }
+    return { new_rate, /*group_name =*/nullptr, /*efficiency dummy value*/0.0 };
+}
+
 
 std::tuple<double,double,double, bool, bool,bool>
 GasLiftSingleWellGeneric::

opm/simulators/wells/GasLiftSingleWellGeneric.hpp

@@ -57,6 +57,7 @@ protected:
 
 public:
     using GLiftSyncGroups = std::set<int>;
+    using Rate = GasLiftGroupInfo::Rate;
     struct GradInfo
     {
         GradInfo() { }
@@ -180,14 +181,22 @@ protected:
     std::pair<double, bool> getGasRateWithLimit_(const std::vector<double>& potentials) const;
     std::tuple<double,double,double,bool,bool,bool>
     getInitialRatesWithLimit_(const std::vector<double>& potentials);
-    std::pair<double, bool> getOilRateWithLimit_(const std::vector<double>& potentials) const;
+    std::tuple<double, const std::string*, double> getRateWithGroupLimit_(
-    std::pair<double, bool> getWaterRateWithLimit_(const std::vector<double>& potentials) const;
+                  Rate rate_type, const double new_rate, const double old_rate) const;
+    std::pair<double, bool> getOilRateWithLimit_(
+                  const std::vector<double>& potentials) const;
+    std::pair<double, bool> getWaterRateWithLimit_(
+                  const std::vector<double>& potentials) const;
 
-    std::pair<double, bool> getOilRateWithGroupLimit_(const double new_oil_rate, const double oil_rate) const;
+    std::pair<double, bool> getOilRateWithGroupLimit_(
-    std::pair<double, bool> getGasRateWithGroupLimit_(const double new_gas_rate, const double gas_rate) const;
+                           double new_oil_rate, double oil_rate) const;
-    std::pair<double, bool> getWaterRateWithGroupLimit_(const double new_water_rate, const double water_rate) const;
+    std::pair<double, bool> getGasRateWithGroupLimit_(
-    std::tuple<double,double,bool,bool> getLiquidRateWithGroupLimit_(double new_oil_rate, const double oil_rate,
+                           double new_gas_rate, double gas_rate) const;
-                                                                     double new_water_rate, const double water_rate) const;
+    std::pair<double, bool> getWaterRateWithGroupLimit_(
+                           double new_water_rate, double water_rate) const;
+    std::tuple<double,double,bool,bool> getLiquidRateWithGroupLimit_(
+                           const double new_oil_rate, const double oil_rate,
+                           const double new_water_rate, const double water_rate) const;
 
     std::tuple<double,double,bool,bool,double>
     increaseALQtoPositiveOilRate_(double alq,