Merge pull request #5693 from bska/private-fetch-well-pi

Make WELPI Infrastructure for ACTIONX Fully Private

opm/simulators/flow/ActionHandler.cpp

@@ -22,6 +22,7 @@
 */
 
 #include <config.h>
+
 #include <opm/simulators/flow/ActionHandler.hpp>
 
 #include <opm/common/OpmLog/OpmLog.hpp>
@@ -29,6 +30,7 @@
 #include <opm/common/TimingMacros.hpp>
 
 #include <opm/input/eclipse/EclipseState/EclipseState.hpp>
+
 #include <opm/input/eclipse/Schedule/Action/ActionContext.hpp>
 #include <opm/input/eclipse/Schedule/Action/Actions.hpp>
 #include <opm/input/eclipse/Schedule/Action/ActionX.hpp>
@@ -40,6 +42,8 @@
 #include <opm/input/eclipse/Schedule/Well/WellMatcher.hpp>
 
 #include <opm/simulators/wells/BlackoilWellModelGeneric.hpp>
+
+#include <opm/simulators/utils/ParallelCommunication.hpp>
 #include <opm/simulators/utils/ParallelSerialization.hpp>
 
 #include <chrono>
@@ -128,6 +132,61 @@ namespace {
 
         Opm::OpmLog::debug("ACTION_NOT_TRIGGERED", message);
     }
+
+    template <typename Scalar, class WellModel>
+    std::unordered_map<std::string, Scalar>
+    fetchWellPI(const int                          reportStep,
+                const Opm::Schedule&               schedule,
+                const WellModel&                   wellModel,
+                const Opm::Action::ActionX&        action,
+                const std::vector<std::string>&    matching_wells,
+                const Opm::Parallel::Communication comm)
+    {
+        auto wellpi = std::unordered_map<std::string, Scalar> {};
+
+        const auto wellpi_wells = action.wellpi_wells
+            (schedule.wellMatcher(reportStep), matching_wells);
+
+        if (wellpi_wells.empty()) {
+            return wellpi;
+        }
+
+        const auto num_wells = schedule[reportStep].well_order().size();
+
+        std::vector<Scalar> wellpi_vector(num_wells);
+        for (const auto& wname : wellpi_wells) {
+            if (wellModel.hasWell(wname)) {
+                const auto& well = schedule.getWell(wname, reportStep);
+                wellpi_vector[well.seqIndex()] = wellModel.wellPI(wname);
+            }
+        }
+
+        if (comm.size() > 1) {
+            std::vector<Scalar> wellpi_buffer(num_wells * comm.size());
+            comm.gather(wellpi_vector.data(), wellpi_buffer.data(), num_wells, 0);
+
+            if (comm.rank() == 0) {
+                for (int rank = 1; rank < comm.size(); ++rank) {
+                    for (std::size_t well_index = 0; well_index < num_wells; ++well_index) {
+                        const auto global_index = rank*num_wells + well_index;
+                        const auto value = wellpi_buffer[global_index];
+                        if (value != Scalar{0}) {
+                            wellpi_vector[well_index] = value;
+                        }
+                    }
+                }
+            }
+
+            comm.broadcast(wellpi_vector.data(), wellpi_vector.size(), 0);
+        }
+
+        for (const auto& wname : wellpi_wells) {
+            const auto& well = schedule.getWell(wname, reportStep);
+            wellpi.insert_or_assign(wname, wellpi_vector[well.seqIndex()]);
+        }
+
+        return wellpi;
+    }
 } // Anonymous namespace
 
 namespace Opm {
@@ -180,13 +239,13 @@ applyActions(const int reportStep,
             logActivePyAction(pyaction->name(), ts);
         }
 
-        this->applySimulatorUpdate(reportStep, sim_update, commit_wellstate, transUp);
+        this->applySimulatorUpdate(reportStep, sim_update, transUp, commit_wellstate);
     }
 
     auto non_triggered = 0;
     const auto simTime = asTimeT(now);
-    for (const auto& action : actions.pending(actionState_, simTime)) {
+    for (const auto& action : actions.pending(this->actionState_, simTime)) {
-        const auto actionResult = action->eval(context);
+        auto actionResult = action->eval(context);
         if (! actionResult) {
             ++non_triggered;
             logInactiveAction(action->name(), ts);
@@ -197,13 +256,15 @@ applyActions(const int reportStep,
 
         logActiveAction(action->name(), matching_wells, ts);
 
-        const auto wellpi = this->fetchWellPI(reportStep, *action, matching_wells);
+        const auto wellpi = fetchWellPI<Scalar>
+            (reportStep, this->schedule_, this->wellModel_,
+             *action, matching_wells, this->comm_);
 
         const auto sim_update = this->schedule_
             .applyAction(reportStep, *action, matching_wells, wellpi);
 
-        this->applySimulatorUpdate(reportStep, sim_update, commit_wellstate, transUp);
+        this->applySimulatorUpdate(reportStep, sim_update, transUp, commit_wellstate);
-        actionState_.add_run(*action, simTime, std::move(actionResult));
+        this->actionState_.add_run(*action, simTime, std::move(actionResult));
     }
 
     if (non_triggered > 0) {
@@ -223,8 +284,8 @@ template<class Scalar>
 void ActionHandler<Scalar>::
 applySimulatorUpdate(const int report_step,
                      const SimulatorUpdate& sim_update,
-                     bool& commit_wellstate,
+                     const TransFunc& updateTrans,
-                     const TransFunc& updateTrans)
+                     bool& commit_wellstate)
 {
     OPM_TIMEBLOCK(applySimulatorUpdate);
 
@@ -244,53 +305,6 @@ applySimulatorUpdate(const int report_step,
     }
 }
 
-template<class Scalar>
-std::unordered_map<std::string, Scalar>
-ActionHandler<Scalar>::
-fetchWellPI(const int reportStep,
-            const Action::ActionX& action,
-            const std::vector<std::string>& matching_wells) const
-{
-  auto wellpi_wells = action.wellpi_wells
-      (this->schedule_.wellMatcher(reportStep), matching_wells);
-
-  if (wellpi_wells.empty()) {
-      return {};
-  }
-
-  const auto num_wells = schedule_[reportStep].well_order().size();
-  std::vector<Scalar> wellpi_vector(num_wells);
-  for (const auto& wname : wellpi_wells) {
-      if (this->wellModel_.hasWell(wname)) {
-          const auto& well = schedule_.getWell( wname, reportStep );
-          wellpi_vector[well.seqIndex()] = this->wellModel_.wellPI(wname);
-      }
-  }
-
-  if (comm_.size() > 1) {
-      std::vector<Scalar> wellpi_buffer(num_wells * comm_.size());
-      comm_.gather( wellpi_vector.data(), wellpi_buffer.data(), num_wells, 0 );
-      if (comm_.rank() == 0) {
-          for (int rank = 1; rank < comm_.size(); rank++) {
-              for (std::size_t well_index=0; well_index < num_wells; well_index++) {
-                  const auto global_index = rank*num_wells + well_index;
-                  const auto value = wellpi_buffer[global_index];
-                  if (value != 0)
-                      wellpi_vector[well_index] = value;
-              }
-          }
-      }
-      comm_.broadcast(wellpi_vector.data(), wellpi_vector.size(), 0);
-  }
-
-  std::unordered_map<std::string, Scalar> wellpi;
-  for (const auto& wname : wellpi_wells) {
-      const auto& well = schedule_.getWell( wname, reportStep );
-      wellpi[wname] = wellpi_vector[ well.seqIndex() ];
-  }
-  return wellpi;
-}
-
 template<class Scalar>
 void ActionHandler<Scalar>::
 evalUDQAssignments(const unsigned episodeIdx,

opm/simulators/flow/ActionHandler.hpp

@@ -31,19 +31,20 @@
 #include <unordered_map>
 #include <vector>
 
+namespace Opm::Action {
+    class State;
+} // namespace Opm::Action
+
 namespace Opm {
-
-namespace Action {
-class ActionX;
-class State;
-}
-
 template<class Scalar> class BlackoilWellModelGeneric;
 class EclipseState;
 class Schedule;
 struct SimulatorUpdate;
 class SummaryState;
 class UDQState;
+} // namespace Opm
+
+namespace Opm {
 
 //! \brief Class handling Action support in simulator
 template<class Scalar>
@@ -53,6 +54,22 @@ public:
     //! \brief Function handle to update transmissiblities.
     using TransFunc = std::function<void(bool)>;
 
+    /// Constructor.
+    ///
+    /// \param[in,out] ecl_state Container of static properties such as
+    /// permeability and transmissibility.
+    ///
+    /// \param[in,out] schedule Container of dynamic objects, such as wells.
+    ///
+    /// \param[in,out] actionState Dynamic state object for all actions.
+    ///
+    /// \param[in,out] summaryState Dynamic state object for all summary
+    /// vectors.
+    ///
+    /// \param[in,out] wellModel Simulation wells on this rank.
+    ///
+    /// \param[in] comm MPI communicator object linking all simulation
+    /// ranks.
     ActionHandler(EclipseState& ecl_state,
                   Schedule& schedule,
                   Action::State& actionState,
@@ -60,36 +77,70 @@ public:
                   BlackoilWellModelGeneric<Scalar>& wellModel,
                   Parallel::Communication comm);
 
+    /// Run all pending actions.
+    ///
+    /// \param[in] reportStep Zero-based report step index.
+    ///
+    /// \param[in] sim_time Elapsed time since simulation start.
+    ///
+    /// \param[in] updateTrans Call-back for affecting transmissibility
+    /// updates.  Typically invoked if the action triggers a keyword like
+    /// MULTZ.
     void applyActions(int reportStep,
                       double sim_time,
                       const TransFunc& updateTrans);
 
-    //! \brief Evaluates UDQ assign statements.
+    /// \brief Evaluates UDQ assign statements.
+    ///
+    /// \param[in] episodeIdx Zero-based report step index.
+    ///
+    /// \param[in,out] udq_state Dynamic state of all user-defined
+    /// quantities.
     void evalUDQAssignments(const unsigned episodeIdx,
                             UDQState& udq_state);
 
-  private:
+private:
-    /*
+    /// Convey dynamic updates triggered by an action block back to the
-       This function is run after applyAction has been completed in the Schedule
+    /// running simulator.
-       implementation. The sim_update argument should have members & flags for
+    ///
-       the simulator properties which need to be updated. This functionality is
+    /// This function is run after applyAction has been completed in the
-       probably not complete.
+    /// Schedule implementation.  The sim_update argument should have
-    */
+    /// members & flags for the simulator properties which need to be
+    /// updated.  This functionality is probably not complete.
+    ///
+    /// \param[in] report_step Zero-based report step index.
+    ///
+    /// \param[in] sim_update Action's resulting simulator update.
+    ///
+    /// \param[in] updateTrans Call-back for affecting transmissibility
+    /// updates.  Typically invoked if the action triggers a keyword like
+    /// MULTZ.
+    ///
+    /// \param[out] commit_wellstate Whether or not the action affected any
+    /// simulation wells which, in turn, may require rebuilding internal
+    /// data structures in the simulator and therefore would require
+    /// preserving the dynamic well and group states prior to doing so.
     void applySimulatorUpdate(int report_step,
                               const SimulatorUpdate& sim_update,
-                              bool& commit_wellstate,
+                              const TransFunc& updateTrans,
-                              const TransFunc& updateTrans);
+                              bool& commit_wellstate);
-
-    std::unordered_map<std::string, Scalar>
-    fetchWellPI(int reportStep,
-                const Action::ActionX& action,
-                const std::vector<std::string>& matching_wells) const;
 
+    /// Static properties such as permeability and transmissibility.
     EclipseState& ecl_state_;
+
+    /// Dynamic objects such as wells.
     Schedule& schedule_;
+
+    /// Dynamic state for all actions--e.g., their run count.
     Action::State& actionState_;
+
+    /// Dynamic state for all user-defined quantities.
     SummaryState& summaryState_;
+
+    /// Simulation wells on this rank.
     BlackoilWellModelGeneric<Scalar>& wellModel_;
+
+    /// MPI communicator object linking all simulation ranks.
     Parallel::Communication comm_;
 };