reorganize (Std|MS)Well::computeWellPotentials

put shared parts in WellInterfaceGeneric::computeWellPotentials and
WellInterfaceGeneric::checkNegativeWellPotentials

opm/simulators/wells/MultisegmentWell_impl.hpp

@@ -272,53 +272,12 @@ namespace Opm
                           std::vector<double>& well_potentials,
                           DeferredLogger& deferred_logger)
     {
-        const int np = this->number_of_phases_;
+        const auto [compute_potential, bhp_controlled_well] =
-        well_potentials.resize(np, 0.0);
+            this->WellInterfaceGeneric::computeWellPotentials(well_potentials, well_state);
 
-        // Stopped wells have zero potential.
+        if (!compute_potential) {
-        if (this->wellIsStopped()) {
             return;
         }
-        this->operability_status_.has_negative_potentials = false;
-
-        // If the well is pressure controlled the potential equals the rate.
-        bool thp_controlled_well = false;
-        bool bhp_controlled_well = false;
-        const auto& ws = well_state.well(this->index_of_well_);
-        if (this->isInjector()) {
-            const Well::InjectorCMode& current = ws.injection_cmode;
-            if (current == Well::InjectorCMode::THP) {
-                thp_controlled_well = true;
-            }
-            if (current == Well::InjectorCMode::BHP) {
-                bhp_controlled_well = true;
-            }
-        } else {
-            const Well::ProducerCMode& current = ws.production_cmode;
-            if (current == Well::ProducerCMode::THP) {
-                thp_controlled_well = true;
-            }
-            if (current == Well::ProducerCMode::BHP) {
-                bhp_controlled_well = true;
-            }
-        }
-        if (!this->changed_to_open_this_step_ && (thp_controlled_well || bhp_controlled_well)) {
-
-            double total_rate = 0.0;
-            const double sign = this->isInjector() ? 1.0:-1.0;
-            for (int phase = 0; phase < np; ++phase){
-                total_rate += sign * ws.surface_rates[phase];
-            }
-            // for pressure controlled wells the well rates are the potentials
-            // if the rates are trivial we are most probably looking at the newly
-            // opened well, and we therefore make the effort of computing the potentials anyway.
-            if (total_rate > 0) {
-                for (int phase = 0; phase < np; ++phase){
-                    well_potentials[phase] = sign * ws.surface_rates[phase];
-                }
-                return;
-            }
-        }
 
         debug_cost_counter_ = 0;
         // does the well have a THP related constraint?
@@ -332,18 +291,9 @@ namespace Opm
         deferred_logger.debug("Cost in iterations of finding well potential for well "
                               + this->name() + ": " + std::to_string(debug_cost_counter_));
 
-        const double sign = this->isInjector() ? 1.0:-1.0;
+        this->checkNegativeWellPotentials(well_potentials,
-        double total_potential = 0.0;
+                                          this->param_.check_well_operability_,
-        for (int phase = 0; phase < np; ++phase){
+                                          deferred_logger);
-            well_potentials[phase] *= sign;
-            total_potential += well_potentials[phase];
-        }
-        if (total_potential < 0.0 && this->param_.check_well_operability_) {
-            // wells with negative potentials are not operable
-            this->operability_status_.has_negative_potentials = true;
-            const std::string msg = std::string("well ") + this->name() + std::string(": has non negative potentials is not operable");
-            deferred_logger.warning("NEGATIVE_POTENTIALS_INOPERABLE", msg);
-        }
     }
 
 

opm/simulators/wells/StandardWell_impl.hpp

@@ -1714,53 +1714,13 @@ namespace Opm
                           std::vector<double>& well_potentials,
                           DeferredLogger& deferred_logger) // const
     {
-        const int np = this->number_of_phases_;
+        const auto [compute_potential, bhp_controlled_well] =
-        well_potentials.resize(np, 0.0);
+            this->WellInterfaceGeneric::computeWellPotentials(well_potentials, well_state);
 
-        if (this->wellIsStopped()) {
+        if (!compute_potential) {
             return;
         }
 
-        this->operability_status_.has_negative_potentials = false;
-        // If the well is pressure controlled the potential equals the rate.
-        bool thp_controlled_well = false;
-        bool bhp_controlled_well = false;
-        const auto& ws = well_state.well(this->index_of_well_);
-        if (this->isInjector()) {
-            const Well::InjectorCMode& current = ws.injection_cmode;
-            if (current == Well::InjectorCMode::THP) {
-                thp_controlled_well = true;
-            }
-            if (current == Well::InjectorCMode::BHP) {
-                bhp_controlled_well = true;
-            }
-        } else {
-            const Well::ProducerCMode& current = ws.production_cmode;
-            if (current == Well::ProducerCMode::THP) {
-                thp_controlled_well = true;
-            }
-            if (current == Well::ProducerCMode::BHP) {
-                bhp_controlled_well = true;
-            }
-        }
-        if (!this->changed_to_open_this_step_ && (thp_controlled_well || bhp_controlled_well)) {
-
-            double total_rate = 0.0;
-            const double sign = this->isInjector() ? 1.0:-1.0;
-            for (int phase = 0; phase < np; ++phase){
-                total_rate += sign * ws.surface_rates[phase];
-            }
-            // for pressure controlled wells the well rates are the potentials
-            // if the rates are trivial we are most probably looking at the newly
-            // opened well and we therefore make the affort of computing the potentials anyway.
-            if (total_rate > 0) {
-                for (int phase = 0; phase < np; ++phase){
-                    well_potentials[phase] = sign * ws.surface_rates[phase];
-                }
-                return;
-            }
-        }
-
         // does the well have a THP related constraint?
         const auto& summaryState = ebosSimulator.vanguard().summaryState();
         if (!Base::wellHasTHPConstraints(summaryState) || bhp_controlled_well) {
@@ -1773,6 +1733,7 @@ namespace Opm
             // As a remedy the bhp that gives the largest potential is used.
             // For converged cases, ws.bhp <=bhp for injectors and ws.bhp >= bhp,
             // and the potentials will be computed using the limit as expected.
+            const auto& ws = well_state.well(this->index_of_well_);
             if (this->isInjector())
                 bhp = std::max(ws.bhp, bhp);
             else
@@ -1785,18 +1746,9 @@ namespace Opm
             well_potentials = computeWellPotentialWithTHP(ebosSimulator, deferred_logger, well_state);
         }
 
-        const double sign = this->isInjector() ? 1.0:-1.0;
+        this->checkNegativeWellPotentials(well_potentials,
-        double total_potential = 0.0;
+                                          this->param_.check_well_operability_,
-        for (int phase = 0; phase < np; ++phase){
+                                          deferred_logger);
-            well_potentials[phase] *= sign;
-            total_potential += well_potentials[phase];
-        }
-        if (total_potential < 0.0 && this->param_.check_well_operability_) {
-            // wells with negative potentials are not operable
-            this->operability_status_.has_negative_potentials = true;
-            const std::string msg = std::string("well ") + this->name() + std::string(": has negative potentials and is not operable");
-            deferred_logger.warning("NEGATIVE_POTENTIALS_INOPERABLE", msg);
-        }
     }
 
 

opm/simulators/wells/WellInterfaceGeneric.cpp

@@ -575,4 +575,81 @@ int WellInterfaceGeneric::polymerInjTable_() const
     return this->well_ecl_.getPolymerProperties().m_plymwinjtable;
 }
 
+std::pair<bool,bool> WellInterfaceGeneric::
+computeWellPotentials(std::vector<double>& well_potentials,
+                      const WellState& well_state)
+{
+    const int np = this->number_of_phases_;
+    well_potentials.resize(np, 0.0);
+
+    // Stopped wells have zero potential.
+    if (this->wellIsStopped()) {
+        return {false, false};
+    }
+    this->operability_status_.has_negative_potentials = false;
+
+    // If the well is pressure controlled the potential equals the rate.
+    bool thp_controlled_well = false;
+    bool bhp_controlled_well = false;
+    bool compute_potential = true;
+    const auto& ws = well_state.well(this->index_of_well_);
+    if (this->isInjector()) {
+        const Well::InjectorCMode& current = ws.injection_cmode;
+        if (current == Well::InjectorCMode::THP) {
+            thp_controlled_well = true;
+        }
+        if (current == Well::InjectorCMode::BHP) {
+            bhp_controlled_well = true;
+        }
+    } else {
+        const Well::ProducerCMode& current = ws.production_cmode;
+        if (current == Well::ProducerCMode::THP) {
+            thp_controlled_well = true;
+        }
+        if (current == Well::ProducerCMode::BHP) {
+            bhp_controlled_well = true;
+        }
+    }
+
+    if (!this->changed_to_open_this_step_ &&
+        (thp_controlled_well || bhp_controlled_well)) {
+        double total_rate = 0.0;
+        const double sign = this->isInjector() ? 1.0 : -1.0;
+        for (int phase = 0; phase < np; ++phase){
+            total_rate += sign * ws.surface_rates[phase];
+        }
+        // for pressure controlled wells the well rates are the potentials
+        // if the rates are trivial we are most probably looking at the newly
+        // opened well, and we therefore make the effort of computing the potentials anyway.
+        if (total_rate > 0) {
+            for (int phase = 0; phase < np; ++phase){
+                well_potentials[phase] = sign * ws.surface_rates[phase];
+            }
+            compute_potential = false;
+        }
+    }
+
+    return {compute_potential, bhp_controlled_well};
+}
+
+void WellInterfaceGeneric::
+checkNegativeWellPotentials(std::vector<double>& well_potentials,
+                            const bool checkOperability,
+                            DeferredLogger& deferred_logger)
+{
+    const double sign = this->isInjector() ? 1.0 : -1.0;
+    double total_potential = 0.0;
+    for (int phase = 0; phase < this->number_of_phases_; ++phase) {
+        well_potentials[phase] *= sign;
+        total_potential += well_potentials[phase];
+    }
+    if (total_potential < 0.0 && checkOperability) {
+        // wells with negative potentials are not operable
+        this->operability_status_.has_negative_potentials = true;
+        const std::string msg = std::string("well ") + this->name() +
+                                ": has negative potentials and is not operable";
+        deferred_logger.warning("NEGATIVE_POTENTIALS_INOPERABLE", msg);
+    }
+}
+
 } // namespace Opm

opm/simulators/wells/WellInterfaceGeneric.hpp

@@ -214,6 +214,14 @@ protected:
     bool wellUnderZeroRateTarget(const SummaryState& summary_state,
                                  const WellState& well_state) const;
 
+    std::pair<bool,bool>
+    computeWellPotentials(std::vector<double>& well_potentials,
+                          const WellState& well_state);
+
+    void checkNegativeWellPotentials(std::vector<double>& well_potentials,
+                                     const bool checkOperability,
+                                     DeferredLogger& deferred_logger);
+
     // definition of the struct OperabilityStatus
     struct OperabilityStatus {
         bool isOperableAndSolvable() const {