[test] different approach

opm/simulators/wells/BlackoilWellModelGeneric.cpp

@@ -964,16 +964,35 @@ updateEclWellsConstraints(const int              timeStepIdx,
     this->wellUpdateLoop(sim_update.affected_wells.begin(),
                          sim_update.affected_wells.end(),
                          timeStepIdx,
-                         [this, &st]
-                         (const auto wellIdx, const auto& well)
-    {
-        auto& ws = this->wellState().well(wellIdx);
-        // whether the well was SHUT before applying the action
-        ws.was_shut_before_action_applied = (ws.status == WellStatus::SHUT);
-
-        ws.updateStatus(well.getStatus());
-        ws.update_type_and_targets(well, st);
-    });
+                         [this, &st, timeStepIdx](const auto wellIdx, const auto& well) {
+                             auto& ws = this->wellState().well(wellIdx);
+                             // whether the well was SHUT before applying the action
+                             ws.was_shut_before_action_applied = (ws.status == WellStatus::SHUT);
+
+                             ws.updateStatus(well.getStatus());
+                             ws.update_type_and_targets(well, st);
+
+                             // New production controls also invalidate any THP limit the
+                             // network has imposed on the well, matching the treatment of
+                             // PRODUCTION_UPDATE events at report step initialization. Updates
+                             // that leave the production controls untouched (e.g. a pure
+                             // WELOPEN) keep the retained limit. A well that is still attached
+                             // to a network node receives a fresh limit at the next network
+                             // balance.
+                             if (well.isProducer()
+                                 && this->schedule_[timeStepIdx].wellgroup_events().hasEvent(
+                                     well.name(), ScheduleEvents::PRODUCTION_UPDATE)) {
+                                 this->genNetwork_.eraseImposedThpLimit(well.name());
+                                 const auto genWell
+                                     = std::find_if(this->well_container_generic_.begin(),
+                                                    this->well_container_generic_.end(),
+                                                    [&wname = well.name()](const auto* w) { return w->name() == wname;
+                                                    });
+                                 if (genWell != this->well_container_generic_.end()) {
+                                     (*genWell)->setDynamicThpLimit(std::optional<Scalar> {});
+                                 }
+                             }
+                         });
 }
 
 template<typename Scalar, typename IndexTraits>

opm/simulators/wells/BlackoilWellModelNetworkGeneric.cpp

@@ -27,8 +27,9 @@
 
 #include <opm/material/fluidsystems/BlackOilDefaultFluidSystemIndices.hpp>
 
-#include <opm/input/eclipse/Schedule/Schedule.hpp>
+#include <opm/input/eclipse/Schedule/Events.hpp>
 #include <opm/input/eclipse/Schedule/Network/Balance.hpp>
+#include <opm/input/eclipse/Schedule/Schedule.hpp>
 
 #include <opm/simulators/wells/BlackoilWellModelGeneric.hpp>
 #include <opm/simulators/wells/GroupStateHelper.hpp>
@@ -238,13 +239,18 @@ updatePressures(const int reportStepIdx,
         // Producers only, since we so far only support the
         // "extended" network model (properties defined by
         // BRANPROP and NODEPROP) which only applies to producers.
-        if (well->isProducer() && well->wellEcl().predictionMode()) {
+        // A well without a VFP table cannot translate the network nodal
+        // pressure into a BHP limit, so it cannot be put under THP control,
+        // while its rates still contribute to the network flows.
+        if (well->isProducer() && well->wellEcl().predictionMode()
+            && well->wellEcl().vfp_table_number() > 0) {
             const auto it = node_pressures_.find(well->wellEcl().groupName());
             if (it != node_pressures_.end()) {
                 // The well belongs to a group with has a network pressure constraint,
                 // set the dynamic THP constraint of the well accordingly.
                 const Scalar new_limit = it->second;
                 well->setDynamicThpLimit(new_limit);
+                network_imposed_thp_limits_[well->name()] = new_limit;
                 SingleWellState<Scalar, IndexTraits>& ws = well_model_.wellState()[well->indexOfWell()];
                 const bool thp_is_limit = ws.production_cmode == Well::ProducerCMode::THP;
                 // TODO: not sure why the thp is NOT updated properly elsewhere
@@ -318,6 +324,14 @@ template<typename Scalar, typename IndexTraits>
 void BlackoilWellModelNetworkGeneric<Scalar, IndexTraits>::
 initialize(const int report_step)
 {
+    // Wells with new production controls specified at this report step do not
+    // keep any THP limit imposed earlier by the network.
+    if (!network_imposed_thp_limits_.empty()) {
+        const auto& events = well_model_.schedule()[report_step].wellgroup_events();
+        std::erase_if(network_imposed_thp_limits_, [&events](const auto& item) {
+            return events.hasEvent(item.first, ScheduleEvents::PRODUCTION_UPDATE);
+        });
+    }
     const auto& network = well_model_.schedule()[report_step].network();
     if (network.active() && !node_pressures_.empty()) {
         for (auto& well : well_model_.genericWells()) {
@@ -333,12 +347,24 @@ initializeWell(WellInterfaceGeneric<Scalar,IndexTraits>& well)
     // Producers only, since we so far only support the
     // "extended" network model (properties defined by
     // BRANPROP and NODEPROP) which only applies to producers.
-    if (well.isProducer() && !node_pressures_.empty()) {
+    // A well without a VFP table cannot translate the network nodal
+    // pressure into a BHP limit, so it cannot be put under THP control,
+    // while its rates still contribute to the network flows.
+    if (well.isProducer() && well.wellEcl().vfp_table_number() > 0) {
         const auto it = this->node_pressures_.find(well.wellEcl().groupName());
         if (it != this->node_pressures_.end()) {
             // The well belongs to a group which has a network nodal pressure,
             // set the dynamic THP constraint based on the network nodal pressure
             well.setDynamicThpLimit(it->second);
+            network_imposed_thp_limits_[well.name()] = it->second;
+        } else {
+            // The well is not (or no longer) attached to a network node. A THP
+            // limit imposed earlier by the network remains in force until new
+            // production controls are specified for the well (Eclipse behavior).
+            const auto imposed = network_imposed_thp_limits_.find(well.name());
+            if (imposed != network_imposed_thp_limits_.end()) {
+                well.setDynamicThpLimit(imposed->second);
+            }
         }
     }
 }
@@ -374,7 +400,9 @@ operator==(const BlackoilWellModelNetworkGeneric<Scalar,IndexTraits>& rhs) const
         && this->node_pressures_ == rhs.node_pressures_
         && this->last_valid_node_pressures_ == rhs.last_valid_node_pressures_
         && this->branch_data_ == rhs.branch_data_
-        && this->last_valid_branch_data_ == rhs.last_valid_branch_data_;
+        && this->last_valid_branch_data_ == rhs.last_valid_branch_data_
+        && this->network_imposed_thp_limits_ == rhs.network_imposed_thp_limits_
+        && this->last_valid_network_imposed_thp_limits_ == rhs.last_valid_network_imposed_thp_limits_;
 }
 
 template class BlackoilWellModelNetworkGeneric<double, BlackOilDefaultFluidSystemIndices>;

opm/simulators/wells/BlackoilWellModelNetworkGeneric.hpp

@@ -64,6 +64,14 @@ class BlackoilWellModelNetworkGeneric
     void setNodePressures(const std::map<std::string, Scalar>& values)
     { node_pressures_ = values; }
 
+    // do not use, only needed for serialization testing
+    void setNetworkImposedThpLimits(const std::map<std::string, Scalar>& current,
+                                    const std::map<std::string, Scalar>& last_valid)
+    {
+        network_imposed_thp_limits_ = current;
+        last_valid_network_imposed_thp_limits_ = last_valid;
+    }
+
     void setFromRestart(const std::optional<std::map<std::string, double>>& restart_pressures);
 
     //! \brief Initialize wells according to network configuration.
@@ -72,6 +80,19 @@ class BlackoilWellModelNetworkGeneric
     //! \brief Initialize a single well according to network configuration.
     void initializeWell(WellInterfaceGeneric<Scalar,IndexTraits>& well);
 
+    //! \brief Forget any THP limit the network has imposed on well \p wname.
+    //!
+    //! Needed when new production controls are specified for the well
+    //! outside the regular report-step processing (e.g. through ACTIONX),
+    //! so that a well detached from the network does not keep a stale
+    //! network-imposed THP limit. The entry is erased from the last-valid
+    //! state as well, since schedule changes survive retried time steps.
+    void eraseImposedThpLimit(const std::string& wname)
+    {
+        network_imposed_thp_limits_.erase(wname);
+        last_valid_network_imposed_thp_limits_.erase(wname);
+    }
+
     /// Checks if network is active (at least one network well on prediction).
     void updateActiveState(const int report_step);
 
@@ -99,12 +120,14 @@ class BlackoilWellModelNetworkGeneric
     {
         this->last_valid_node_pressures_ = this->node_pressures_;
         this->last_valid_branch_data_ = this->branch_data_;
+        this->last_valid_network_imposed_thp_limits_ = this->network_imposed_thp_limits_;
     }
 
     void resetState()
     {
         this->node_pressures_ = this->last_valid_node_pressures_;
         this->branch_data_ = this->last_valid_branch_data_;
+        this->network_imposed_thp_limits_ = this->last_valid_network_imposed_thp_limits_;
     }
 
     template<class Serializer>
@@ -114,6 +137,8 @@ class BlackoilWellModelNetworkGeneric
         serializer(last_valid_node_pressures_);
         serializer(branch_data_);
         serializer(last_valid_branch_data_);
+        serializer(network_imposed_thp_limits_);
+        serializer(last_valid_network_imposed_thp_limits_);
     }
 
     bool operator==(const BlackoilWellModelNetworkGeneric<Scalar,IndexTraits>& rhs) const;
@@ -138,6 +163,12 @@ class BlackoilWellModelNetworkGeneric
     std::map<std::string, Scalar> last_valid_node_pressures_;
     // Valid network branch pressure drops and flow rates for output (outlet branch for production network, inlet branch for injection network) for safe restart after failed iterations
     std::map<std::string, data::BranchData> last_valid_branch_data_;
+    // The latest THP limit the network imposed on each well. A well that is
+    // detached from the network keeps its network-imposed THP limit until new
+    // production controls are specified for it (matching Eclipse behavior).
+    std::map<std::string, Scalar> network_imposed_thp_limits_;
+    // Valid network-imposed THP limits for safe restart after failed iterations
+    std::map<std::string, Scalar> last_valid_network_imposed_thp_limits_;
 };
 
 } // namespace Opm

opm/simulators/wells/SingleWellState.cpp

@@ -53,7 +53,6 @@ SingleWellState(const std::string& name_,
     , reservoir_rates(pu.numActivePhases())
     , prev_surface_rates(pu.numActivePhases())
     , perf_data(perf_input.size(), !is_producer, pu.numActivePhases())
-    , trivial_group_target(false)
     , use_group_target_fallback(false)
 {
     for (std::size_t perf = 0; perf < perf_input.size(); perf++) {
@@ -407,7 +406,6 @@ bool SingleWellState<Scalar, IndexTraits>::operator==(const SingleWellState& rhs
            this->frac_rate == rhs.frac_rate &&
            this->perf_data == rhs.perf_data &&
            this->filtrate_conc == rhs.filtrate_conc &&
-           this->trivial_group_target == rhs.trivial_group_target &&
            this->segments == rhs.segments &&
            this->events == rhs.events &&
            this->injection_cmode == rhs.injection_cmode &&

opm/simulators/wells/SingleWellState.hpp

@@ -78,7 +78,6 @@ class SingleWellState {
         serializer(reservoir_rates);
         serializer(prev_surface_rates);
         serializer(frac_rate);
-        serializer(trivial_group_target);
         serializer(segments);
         serializer(events);
         serializer(injection_cmode);
@@ -156,7 +155,6 @@ class SingleWellState {
     std::vector<Scalar> prev_surface_rates;
     Scalar frac_rate{0.0};
     PerfData<Scalar> perf_data;
-    bool trivial_group_target;
     std::optional<GroupTarget> group_target;
     std::optional<GroupTarget> group_target_fallback;
     bool use_group_target_fallback;

opm/simulators/wells/WellConstraints.cpp

@@ -292,9 +292,12 @@ activeProductionConstraint(const SingleWellState<Scalar, IndexTraits>& ws,
     if (well_.wellHasTHPConstraints(summaryState) && currentControl != Well::ProducerCMode::THP) {
         const auto& thp = well_.getTHPConstraint(summaryState);
         Scalar current_thp = ws.thp;
-        // For trivial group targets (for instance caused by NETV) we dont want to flip to THP control.
-        const bool dont_check = (currentControl == Well::ProducerCMode::GRUP && ws.trivial_group_target);
-        if (thp > current_thp && !dont_check) {
+        // Wells under a zero group rate target are handled before this check
+        // is reached: stoppedOrZeroRateTarget() in updateWellControl() and
+        // wellUnderZeroRateTarget() in the local-iteration switching both
+        // return early for them. The THP limit is checked unconditionally
+        // here.
+        if (thp > current_thp) {
             // If WVFPEXP item 4 is set to YES1 or YES2
             // switching to THP is prevented if the well will
             // produce at a higher rate with THP control

opm/simulators/wells/WellInterface_impl.hpp

@@ -1733,11 +1733,6 @@ namespace Opm
                     for (int p = 0; p<np; ++p) {
                         ws.surface_rates[p] *= scale;
                     }
-                    ws.trivial_group_target = false;
-                } else {
-                    // If group target is trivial we dont want to flip to other controls. To avoid oscillation we store
-                    // this information in the well state and explicitly check for this condition when evaluating well controls.
-                    ws.trivial_group_target = true;
                 }
                 break;
             }

tests/test_RestartSerialization.cpp

@@ -330,6 +330,7 @@ class BlackoilWellModelGenericTest : public BlackoilWellModelGeneric<double, Ind
         closed_this_step_ = {"test1", "test2"};
         guideRate_.setSerializationTestData();
         genNetwork_.setNodePressures({{"test3", 4.0}});
+        genNetwork_.setNetworkImposedThpLimits({{"test4", 5.0}}, {{"test5", 6.0}});
         active_wgstate_ = WGState<double, IndexTraits>::serializationTestObject(dummy);
         last_valid_wgstate_ = WGState<double, IndexTraits>::serializationTestObject(dummy);
         nupcol_wgstate_ = WGState<double, IndexTraits>::serializationTestObject(dummy);