Better separation of concerns in the core (#3971)

Fixes #3907

Move callbacks from `communication.cpp` to the relevant source files. This hides implementation details of many submodules, which no longer need to expose local functions in the public interface with doxygen blocks stating for example "Used by \ref mpi_place_particle, should not be used elsewhere". This also reduces the visibility of `MpiCallbacks.hpp`, `communication.hpp`, `Particle.hpp` and `cuda_init.hpp`. When the head function and worker function have the same code, remove code duplication using `mpi_call_all()`. Remove the `CALLBACK_LIST` macro and split `mpi_gather_stats()` into independent callbacks.

src/core/CMakeLists.txt

@@ -17,6 +17,7 @@ set(EspressoCore_SRC
     global.cpp
     grid.cpp
     immersed_boundaries.cpp
+    interactions.cpp
     event.cpp
     integrate.cpp
     npt.cpp

src/core/RuntimeErrorCollector.cpp

@@ -24,9 +24,9 @@
 #include <boost/mpi/collectives.hpp>
 
 #include <iostream>
+#include <sstream>
 #include <utility>
 
-using namespace std;
 using boost::mpi::all_reduce;
 using boost::mpi::communicator;
 
@@ -56,39 +56,40 @@ void RuntimeErrorCollector::message(RuntimeError::ErrorLevel level,
                                     const std::string &msg,
                                     const char *function, const char *file,
                                     const int line) {
-  m_errors.emplace_back(level, m_comm.rank(), msg, string(function),
-                        string(file), line);
+  m_errors.emplace_back(level, m_comm.rank(), msg, std::string(function),
+                        std::string(file), line);
 }
 
-void RuntimeErrorCollector::warning(const string &msg, const char *function,
-                                    const char *file, const int line) {
+void RuntimeErrorCollector::warning(const std::string &msg,
+                                    const char *function, const char *file,
+                                    const int line) {
   m_errors.emplace_back(RuntimeError::ErrorLevel::WARNING, m_comm.rank(), msg,
-                        string(function), string(file), line);
+                        std::string(function), std::string(file), line);
 }
 
 void RuntimeErrorCollector::warning(const char *msg, const char *function,
                                     const char *file, const int line) {
-  warning(string(msg), function, file, line);
+  warning(std::string(msg), function, file, line);
 }
 
-void RuntimeErrorCollector::warning(const ostringstream &mstr,
+void RuntimeErrorCollector::warning(const std::ostringstream &mstr,
                                     const char *function, const char *file,
                                     const int line) {
   warning(mstr.str(), function, file, line);
 }
 
-void RuntimeErrorCollector::error(const string &msg, const char *function,
+void RuntimeErrorCollector::error(const std::string &msg, const char *function,
                                   const char *file, const int line) {
   m_errors.emplace_back(RuntimeError::ErrorLevel::ERROR, m_comm.rank(), msg,
-                        string(function), string(file), line);
+                        std::string(function), std::string(file), line);
 }
 
 void RuntimeErrorCollector::error(const char *msg, const char *function,
                                   const char *file, const int line) {
-  error(string(msg), function, file, line);
+  error(std::string(msg), function, file, line);
 }
 
-void RuntimeErrorCollector::error(const ostringstream &mstr,
+void RuntimeErrorCollector::error(const std::ostringstream &mstr,
                                   const char *function, const char *file,
                                   const int line) {
   error(mstr.str(), function, file, line);
@@ -111,8 +112,8 @@ int RuntimeErrorCollector::count(RuntimeError::ErrorLevel level) {
 
 void RuntimeErrorCollector::clear() { m_errors.clear(); }
 
-vector<RuntimeError> RuntimeErrorCollector::gather() {
-  vector<RuntimeError> all_errors{};
+std::vector<RuntimeError> RuntimeErrorCollector::gather() {
+  std::vector<RuntimeError> all_errors{};
   std::swap(all_errors, m_errors);
 
   Utils::Mpi::gather_buffer(all_errors, m_comm);

src/core/RuntimeErrorCollector.hpp

@@ -20,6 +20,7 @@
 #ifndef ERROR_HANDLING_RUNTIMEERRORCOLLECTOR_HPP
 #define ERROR_HANDLING_RUNTIMEERRORCOLLECTOR_HPP
 
+#include <sstream>
 #include <string>
 #include <vector>
 

src/core/actor/DipolarBarnesHut.cpp

@@ -19,8 +19,8 @@
 
 #include "DipolarBarnesHut.hpp"
 #include "EspressoSystemInterface.hpp"
-#include "communication.hpp"
 #include "config.hpp"
+#include "electrostatics_magnetostatics/common.hpp"
 #include "energy.hpp"
 #include "forces.hpp"
 #include "grid.hpp"

src/core/actor/DipolarDirectSum.cpp

@@ -21,7 +21,7 @@
 
 #include "DipolarDirectSum.hpp"
 #include "EspressoSystemInterface.hpp"
-#include "communication.hpp"
+#include "electrostatics_magnetostatics/common.hpp"
 #include "energy.hpp"
 #include "forces.hpp"
 

src/core/bonded_interactions/angle_cosine.cpp

@@ -24,7 +24,7 @@
  */
 #include "angle_cosine.hpp"
 
-#include "communication.hpp"
+#include "interactions.hpp"
 
 #include <utils/constants.hpp>
 

src/core/bonded_interactions/angle_cossquare.cpp

@@ -24,7 +24,7 @@
  */
 #include "angle_cossquare.hpp"
 
-#include "communication.hpp"
+#include "interactions.hpp"
 
 int angle_cossquare_set_params(int bond_type, double bend, double phi0) {
   if (bond_type < 0)

src/core/bonded_interactions/angle_harmonic.cpp

@@ -24,7 +24,7 @@
  */
 #include "angle_harmonic.hpp"
 
-#include "communication.hpp"
+#include "interactions.hpp"
 
 #include <utils/constants.hpp>
 

src/core/bonded_interactions/bonded_coulomb.cpp

@@ -23,7 +23,7 @@
  *  Implementation of \ref bonded_coulomb.hpp
  */
 #include "bonded_coulomb.hpp"
-#include "communication.hpp"
+#include "interactions.hpp"
 
 #include <utils/constants.hpp>
 

src/core/bonded_interactions/bonded_coulomb_sr.cpp

@@ -25,7 +25,7 @@
 #include "bonded_coulomb_sr.hpp"
 
 #ifdef ELECTROSTATICS
-#include "communication.hpp"
+#include "interactions.hpp"
 
 #include <utils/constants.hpp>
 

src/core/bonded_interactions/bonded_interaction_data.cpp

@@ -17,7 +17,7 @@
  * along with this program.  If not, see <http://www.gnu.org/licenses/>.
  */
 #include "bonded_interaction_data.hpp"
-#include "communication.hpp"
+#include "interactions.hpp"
 
 #include <boost/algorithm/cxx11/any_of.hpp>
 #include <boost/range/algorithm/copy.hpp>

src/core/bonded_interactions/bonded_tab.cpp

@@ -20,8 +20,8 @@
  */
 #include "bonded_interactions/bonded_tab.hpp"
 
-#include "communication.hpp"
 #include "errorhandling.hpp"
+#include "interactions.hpp"
 
 int tabulated_bonded_set_params(int bond_type,
                                 TabulatedBondedInteraction tab_type, double min,

src/core/bonded_interactions/dihedral.cpp

@@ -23,7 +23,7 @@
  *  Implementation of \ref dihedral.hpp
  */
 #include "dihedral.hpp"
-#include "communication.hpp"
+#include "interactions.hpp"
 
 int dihedral_set_params(int bond_type, int mult, double bend, double phase) {
   if (bond_type < 0)

src/core/bonded_interactions/fene.cpp

@@ -24,7 +24,7 @@
  */
 
 #include "fene.hpp"
-#include "communication.hpp"
+#include "interactions.hpp"
 
 #include <utils/constants.hpp>
 #include <utils/math/sqr.hpp>

src/core/bonded_interactions/harmonic.cpp

@@ -23,7 +23,7 @@
  *  Implementation of \ref harmonic.hpp
  */
 #include "harmonic.hpp"
-#include "communication.hpp"
+#include "interactions.hpp"
 
 #include <utils/constants.hpp>
 

src/core/bonded_interactions/quartic.cpp

@@ -23,7 +23,7 @@
  *  Implementation of \ref quartic.hpp
  */
 #include "quartic.hpp"
-#include "communication.hpp"
+#include "interactions.hpp"
 
 int quartic_set_params(int bond_type, double k0, double k1, double r,
                        double r_cut) {

src/core/bonded_interactions/thermalized_bond.cpp

@@ -25,9 +25,9 @@
 
 #include "thermalized_bond.hpp"
 #include "bonded_interaction_data.hpp"
-#include "communication.hpp"
 #include "global.hpp"
 #include "integrate.hpp"
+#include "interactions.hpp"
 
 int n_thermalized_bonds = 0;
 

src/core/bonded_interactions/umbrella.cpp

@@ -23,7 +23,7 @@
  *  Implementation of \ref umbrella.hpp
  */
 #include "umbrella.hpp"
-#include "communication.hpp"
+#include "interactions.hpp"
 
 #ifdef UMBRELLA
 

src/core/cells.cpp

@@ -30,6 +30,7 @@
 #include "event.hpp"
 #include "grid.hpp"
 #include "integrate.hpp"
+#include "particle_data.hpp"
 
 #include "DomainDecomposition.hpp"
 #include "ParticleDecomposition.hpp"
@@ -75,7 +76,7 @@ std::vector<std::pair<int, int>> get_pairs(double distance) {
   return ret;
 }
 
-void mpi_get_pairs_slave(int, int) {
+void mpi_get_pairs_local(int, int) {
   on_observable_calc();
 
   double distance;
@@ -86,11 +87,10 @@ void mpi_get_pairs_slave(int, int) {
   Utils::Mpi::gather_buffer(local_pairs, comm_cart);
 }
 
-/**
- * @brief Collect pairs from all nodes.
- */
+REGISTER_CALLBACK(mpi_get_pairs_local)
+
 std::vector<std::pair<int, int>> mpi_get_pairs(double distance) {
-  mpi_call(mpi_get_pairs_slave, 0, 0);
+  mpi_call(mpi_get_pairs_local, 0, 0);
   on_observable_calc();
 
   boost::mpi::broadcast(comm_cart, distance, 0);
@@ -102,11 +102,28 @@ std::vector<std::pair<int, int>> mpi_get_pairs(double distance) {
   return pairs;
 }
 
-/************************************************************
- *            Exported Functions                            *
- ************************************************************/
+void mpi_resort_particles_local(int global_flag, int) {
+  cell_structure.resort_particles(global_flag);
+
+  boost::mpi::gather(
+      comm_cart, static_cast<int>(cell_structure.local_particles().size()), 0);
+}
+
+REGISTER_CALLBACK(mpi_resort_particles_local)
+
+std::vector<int> mpi_resort_particles(int global_flag) {
+  mpi_call(mpi_resort_particles_local, global_flag, 0);
+  cell_structure.resort_particles(global_flag);
 
-/************************************************************/
+  clear_particle_node();
+
+  std::vector<int> n_parts;
+  boost::mpi::gather(comm_cart,
+                     static_cast<int>(cell_structure.local_particles().size()),
+                     n_parts, 0);
+
+  return n_parts;
+}
 
 void cells_re_init(int new_cs) {
   switch (new_cs) {
@@ -124,7 +141,9 @@ void cells_re_init(int new_cs) {
   on_cell_structure_change();
 }
 
-/*************************************************/
+REGISTER_CALLBACK(cells_re_init)
+
+void mpi_bcast_cell_structure(int cs) { mpi_call_all(cells_re_init, cs); }
 
 void check_resort_particles() {
   const double skin2 = Utils::sqr(skin / 2.0);
@@ -140,7 +159,6 @@ void check_resort_particles() {
   cell_structure.set_resort_particles(level);
 }
 
-/*************************************************/
 void cells_update_ghosts(unsigned data_parts) {
   /* data parts that are only updated on resort */
   auto constexpr resort_only_parts =
@@ -184,6 +202,12 @@ const DomainDecomposition *get_domain_decomposition() {
       Utils::as_const(cell_structure).decomposition());
 }
 
-void cells_set_use_verlet_lists(bool use_verlet_lists) {
+void mpi_set_use_verlet_lists_local(bool use_verlet_lists) {
   cell_structure.use_verlet_list = use_verlet_lists;
 }
+
+REGISTER_CALLBACK(mpi_set_use_verlet_lists_local)
+
+void mpi_set_use_verlet_lists(bool use_verlet_lists) {
+  mpi_call_all(mpi_set_use_verlet_lists_local, use_verlet_lists);
+}

src/core/cells.hpp

@@ -46,21 +46,17 @@
 #include <utility>
 #include <vector>
 
-/** \name Flags for exchange_and_sort_particles: whether to do a global
+/** Flags for particle exchange and resorting: whether to do a global
  *  exchange or assume that particles did not move much (faster, used
  *  during integration, where moving far is a catastrophe anyways).
  */
-/*@{*/
-
 enum {
-  /** Flag for exchange_and_sort_particles : Do neighbor exchange. */
+  /** Do neighbor exchange. */
   CELL_NEIGHBOR_EXCHANGE = 0,
-  /** Flag for exchange_and_sort_particles : Do global exchange. */
+  /** Do global exchange. */
   CELL_GLOBAL_EXCHANGE = 1
 };
 
-/*@}*/
-
 /** Type of cell structure in use. */
 extern CellStructure cell_structure;
 
@@ -69,12 +65,16 @@ extern CellStructure cell_structure;
  */
 void cells_re_init(int new_cs);
 
+/** Change the cell structure on all nodes. */
+void mpi_bcast_cell_structure(int cs);
+
 /**
- * @brief Set use_verlet_lists
+ * @brief Set @ref CellStructure::use_verlet_list
+ * "cell_structure::use_verlet_list"
  *
- * @param use_verlet_lists Should verlet lists be used?
+ * @param use_verlet_lists Should Verlet lists be used?
  */
-void cells_set_use_verlet_lists(bool use_verlet_lists);
+void mpi_set_use_verlet_lists(bool use_verlet_lists);
 
 /** Update ghost information. If needed,
  *  the particles are also resorted.
@@ -91,6 +91,17 @@ std::vector<std::pair<int, int>> mpi_get_pairs(double distance);
 /** Check if a particle resorting is required. */
 void check_resort_particles();
 
+/**
+ * @brief Resort the particles.
+ *
+ * This function resorts the particles on the nodes.
+ *
+ * @param global_flag If true a global resort is done,
+ *        if false particles are only exchanges between neighbors.
+ * @return The number of particles on the nodes after the resort.
+ */
+std::vector<int> mpi_resort_particles(int global_flag);
+
 /**
  * @brief Find the cell in which a particle is stored.
  *
@@ -104,8 +115,7 @@ Cell *find_current_cell(const Particle &p);
 /**
  * @brief Return a pointer to the global DomainDecomposition.
  *
- * @return Pointer to the decomposition if it is set and is
- * DomainDecomposition, nullptr otherwise.
+ * @return Pointer to the decomposition if it is set, nullptr otherwise.
  */
 const DomainDecomposition *get_domain_decomposition();