diff --git a/src/core/CMakeLists.txt b/src/core/CMakeLists.txt index a2a807c89e9..6e04e7906bc 100644 --- a/src/core/CMakeLists.txt +++ b/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 diff --git a/src/core/RuntimeErrorCollector.cpp b/src/core/RuntimeErrorCollector.cpp index 45e8758d902..fc9cd790fc6 100644 --- a/src/core/RuntimeErrorCollector.cpp +++ b/src/core/RuntimeErrorCollector.cpp @@ -24,9 +24,9 @@ #include #include +#include #include -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 RuntimeErrorCollector::gather() { - vector all_errors{}; +std::vector RuntimeErrorCollector::gather() { + std::vector all_errors{}; std::swap(all_errors, m_errors); Utils::Mpi::gather_buffer(all_errors, m_comm); diff --git a/src/core/RuntimeErrorCollector.hpp b/src/core/RuntimeErrorCollector.hpp index 2079ac06659..8ea7b3f047f 100644 --- a/src/core/RuntimeErrorCollector.hpp +++ b/src/core/RuntimeErrorCollector.hpp @@ -20,6 +20,7 @@ #ifndef ERROR_HANDLING_RUNTIMEERRORCOLLECTOR_HPP #define ERROR_HANDLING_RUNTIMEERRORCOLLECTOR_HPP +#include #include #include diff --git a/src/core/actor/DipolarBarnesHut.cpp b/src/core/actor/DipolarBarnesHut.cpp index 2dbb406deae..4872470cb49 100644 --- a/src/core/actor/DipolarBarnesHut.cpp +++ b/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" diff --git a/src/core/actor/DipolarDirectSum.cpp b/src/core/actor/DipolarDirectSum.cpp index fb2334a9d28..bfc47ee4451 100644 --- a/src/core/actor/DipolarDirectSum.cpp +++ b/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" diff --git a/src/core/bonded_interactions/angle_cosine.cpp b/src/core/bonded_interactions/angle_cosine.cpp index cda46f777e8..3cefd9a0e11 100644 --- a/src/core/bonded_interactions/angle_cosine.cpp +++ b/src/core/bonded_interactions/angle_cosine.cpp @@ -24,7 +24,7 @@ */ #include "angle_cosine.hpp" -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/bonded_interactions/angle_cossquare.cpp b/src/core/bonded_interactions/angle_cossquare.cpp index 0362b807830..8f3619d8ac4 100644 --- a/src/core/bonded_interactions/angle_cossquare.cpp +++ b/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) diff --git a/src/core/bonded_interactions/angle_harmonic.cpp b/src/core/bonded_interactions/angle_harmonic.cpp index 8e4557de0c7..076af06d7fe 100644 --- a/src/core/bonded_interactions/angle_harmonic.cpp +++ b/src/core/bonded_interactions/angle_harmonic.cpp @@ -24,7 +24,7 @@ */ #include "angle_harmonic.hpp" -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/bonded_interactions/bonded_coulomb.cpp b/src/core/bonded_interactions/bonded_coulomb.cpp index dd50c510391..81cd94fe295 100644 --- a/src/core/bonded_interactions/bonded_coulomb.cpp +++ b/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 diff --git a/src/core/bonded_interactions/bonded_coulomb_sr.cpp b/src/core/bonded_interactions/bonded_coulomb_sr.cpp index d338e42a333..edfbe3ccc36 100644 --- a/src/core/bonded_interactions/bonded_coulomb_sr.cpp +++ b/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 diff --git a/src/core/bonded_interactions/bonded_interaction_data.cpp b/src/core/bonded_interactions/bonded_interaction_data.cpp index 56f47a662f4..1b8fb8483c9 100644 --- a/src/core/bonded_interactions/bonded_interaction_data.cpp +++ b/src/core/bonded_interactions/bonded_interaction_data.cpp @@ -17,7 +17,7 @@ * along with this program. If not, see . */ #include "bonded_interaction_data.hpp" -#include "communication.hpp" +#include "interactions.hpp" #include #include diff --git a/src/core/bonded_interactions/bonded_tab.cpp b/src/core/bonded_interactions/bonded_tab.cpp index 1a2e056ab8c..a5824eba939 100644 --- a/src/core/bonded_interactions/bonded_tab.cpp +++ b/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, diff --git a/src/core/bonded_interactions/dihedral.cpp b/src/core/bonded_interactions/dihedral.cpp index 93c9a9333e8..a38bb976ebe 100644 --- a/src/core/bonded_interactions/dihedral.cpp +++ b/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) diff --git a/src/core/bonded_interactions/fene.cpp b/src/core/bonded_interactions/fene.cpp index 4cdec115b06..fb96121a44b 100644 --- a/src/core/bonded_interactions/fene.cpp +++ b/src/core/bonded_interactions/fene.cpp @@ -24,7 +24,7 @@ */ #include "fene.hpp" -#include "communication.hpp" +#include "interactions.hpp" #include #include diff --git a/src/core/bonded_interactions/harmonic.cpp b/src/core/bonded_interactions/harmonic.cpp index 1d23b6f47d9..eaf47b5c46b 100644 --- a/src/core/bonded_interactions/harmonic.cpp +++ b/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 diff --git a/src/core/bonded_interactions/quartic.cpp b/src/core/bonded_interactions/quartic.cpp index 31284ca40b6..2106179edab 100644 --- a/src/core/bonded_interactions/quartic.cpp +++ b/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) { diff --git a/src/core/bonded_interactions/thermalized_bond.cpp b/src/core/bonded_interactions/thermalized_bond.cpp index 3776c34362d..5de659968d7 100644 --- a/src/core/bonded_interactions/thermalized_bond.cpp +++ b/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; diff --git a/src/core/bonded_interactions/umbrella.cpp b/src/core/bonded_interactions/umbrella.cpp index 8baf49107a2..fea29021b13 100644 --- a/src/core/bonded_interactions/umbrella.cpp +++ b/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 diff --git a/src/core/cells.cpp b/src/core/cells.cpp index ceaa1dee2f2..558a2f45193 100644 --- a/src/core/cells.cpp +++ b/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> 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> 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> 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(cell_structure.local_particles().size()), 0); +} + +REGISTER_CALLBACK(mpi_resort_particles_local) + +std::vector 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 n_parts; + boost::mpi::gather(comm_cart, + static_cast(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); +} diff --git a/src/core/cells.hpp b/src/core/cells.hpp index eb2a6dca973..e0c2ad5c52c 100644 --- a/src/core/cells.hpp +++ b/src/core/cells.hpp @@ -46,21 +46,17 @@ #include #include -/** \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> 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 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(); diff --git a/src/core/communication.cpp b/src/core/communication.cpp index f8800e24d4a..4296399c643 100644 --- a/src/core/communication.cpp +++ b/src/core/communication.cpp @@ -31,40 +31,12 @@ #include "errorhandling.hpp" -#include "CellStructure.hpp" -#include "EspressoSystemInterface.hpp" -#include "bonded_interactions/bonded_interaction_data.hpp" -#include "cells.hpp" -#include "cuda_interface.hpp" -#include "energy.hpp" #include "event.hpp" -#include "galilei.hpp" -#include "global.hpp" #include "grid.hpp" -#include "grid_based_algorithms/lb.hpp" -#include "grid_based_algorithms/lb_interface.hpp" -#include "integrate.hpp" -#include "integrators/steepest_descent.hpp" -#include "npt.hpp" -#include "particle_data.hpp" -#include "pressure.hpp" - -#include "electrostatics_magnetostatics/coulomb.hpp" -#include "electrostatics_magnetostatics/dipole.hpp" -#include "electrostatics_magnetostatics/icc.hpp" -#include "electrostatics_magnetostatics/mdlc_correction.hpp" - -#include "serialization/IA_parameters.hpp" #include -#include -#include -#include -#include #include -using namespace std; - namespace Communication { auto const &mpi_datatype_cache = boost::mpi::detail::mpi_datatype_cache(); std::shared_ptr mpi_env; @@ -88,38 +60,6 @@ using Communication::mpiCallbacks; int this_node = -1; int n_nodes = -1; -// if you want to add a callback, add it here, and here only -#define CALLBACK_LIST \ - CB(mpi_who_has_slave) \ - CB(mpi_place_particle_slave) \ - CB(mpi_bcast_ia_params_slave) \ - CB(mpi_gather_stats_slave) \ - CB(mpi_bcast_coulomb_params_slave) \ - CB(mpi_remove_particle_slave) \ - CB(mpi_rescale_particles_slave) \ - CB(mpi_bcast_steepest_descent_worker) \ - CB(mpi_bcast_cuda_global_part_vars_slave) \ - CB(mpi_resort_particles_slave) \ - CB(mpi_get_pairs_slave) \ - CB(mpi_get_particles_slave) \ - CB(mpi_rotate_system_slave) \ - CB(mpi_update_particle_slave) - -// create the forward declarations -#define CB(name) void name(int node, int param); -#ifndef DOXYGEN -/* this conditional on DOXYGEN prevents an interaction in Doxygen between - * CALLBACK_LIST and whatever follows next, e.g. a function "int foo();" - * would otherwise become "CALLBACK_LIST int foo();" */ -CALLBACK_LIST -#endif - -#undef CB - -/********************************************** - * procedures - **********************************************/ - #if defined(OPEN_MPI) namespace { /** Workaround for "Read -1, expected XXXXXXX, errno = 14" that sometimes @@ -191,10 +131,6 @@ void init(std::shared_ptr mpi_env) { Communication::m_callbacks = std::make_unique(comm_cart); -#define CB(name) Communication::m_callbacks->add(&(name)); - CALLBACK_LIST -#undef CB - ErrorHandling::init_error_handling(mpiCallbacks()); on_program_start(); @@ -210,441 +146,7 @@ std::shared_ptr mpi_init() { return std::make_shared(); } -/****************** PLACE/PLACE NEW PARTICLE ************/ - -void mpi_place_particle(int node, int id, const Utils::Vector3d &pos) { - mpi_call(mpi_place_particle_slave, node, id); - - if (node == this_node) - local_place_particle(id, pos, 0); - else { - comm_cart.send(node, SOME_TAG, pos); - } - - cell_structure.set_resort_particles(Cells::RESORT_GLOBAL); - on_particle_change(); -} - -void mpi_place_particle_slave(int pnode, int part) { - if (pnode == this_node) { - Utils::Vector3d pos; - comm_cart.recv(0, SOME_TAG, pos); - local_place_particle(part, pos, 0); - } - - cell_structure.set_resort_particles(Cells::RESORT_GLOBAL); - on_particle_change(); -} - -boost::optional mpi_place_new_particle_slave(int part, - Utils::Vector3d const &pos) { - auto p = local_place_particle(part, pos, 1); - - on_particle_change(); - - if (p) { - return comm_cart.rank(); - } - - return {}; -} - -REGISTER_CALLBACK_ONE_RANK(mpi_place_new_particle_slave) - -int mpi_place_new_particle(int id, const Utils::Vector3d &pos) { - return mpi_call(Communication::Result::one_rank, mpi_place_new_particle_slave, - id, pos); -} - -/****************** REMOVE PARTICLE ************/ -void mpi_remove_particle(int pnode, int part) { - mpi_call_all(mpi_remove_particle_slave, pnode, part); -} - -void mpi_remove_particle_slave(int pnode, int part) { - if (part != -1) { - cell_structure.remove_particle(part); - } else { - cell_structure.remove_all_particles(); - } - - on_particle_change(); -} - -/********************* STEEPEST DESCENT ********/ -static int mpi_steepest_descent_slave(int steps, int) { - return integrate(steps, -1); -} -REGISTER_CALLBACK_MASTER_RANK(mpi_steepest_descent_slave) - -int mpi_steepest_descent(int steps) { - return mpi_call(Communication::Result::master_rank, - mpi_steepest_descent_slave, steps, 0); -} - -/********************* INTEGRATE ********/ -static int mpi_integrate_slave(int n_steps, int reuse_forces) { - integrate(n_steps, reuse_forces); - - return check_runtime_errors_local(); -} -REGISTER_CALLBACK_REDUCTION(mpi_integrate_slave, std::plus()) - -int mpi_integrate(int n_steps, int reuse_forces) { - return mpi_call(Communication::Result::reduction, std::plus(), - mpi_integrate_slave, n_steps, reuse_forces); -} - -/*************** BCAST IA ************/ -static void mpi_bcast_all_ia_params_slave() { - boost::mpi::broadcast(comm_cart, ia_params, 0); -} - -REGISTER_CALLBACK(mpi_bcast_all_ia_params_slave) - -void mpi_bcast_all_ia_params() { mpi_call_all(mpi_bcast_all_ia_params_slave); } - -inline bool is_tabulated_bond(BondedInteraction const type) { - return (type == BONDED_IA_TABULATED_DISTANCE or - type == BONDED_IA_TABULATED_ANGLE or - type == BONDED_IA_TABULATED_DIHEDRAL); -} - -REGISTER_CALLBACK(mpi_bcast_ia_params_slave) - -void mpi_bcast_ia_params(int i, int j) { - mpi_call_all(mpi_bcast_ia_params_slave, i, j); -} - -void mpi_bcast_ia_params_slave(int i, int j) { - if (j >= 0) { - // non-bonded interaction parameters - boost::mpi::broadcast(comm_cart, *get_ia_param(i, j), 0); - } else { - // bonded interaction parameters - if (this_node) { - make_bond_type_exist(i); // realloc bonded_ia_params on slave nodes! - if (is_tabulated_bond(bonded_ia_params[i].type)) { - delete bonded_ia_params[i].p.tab.pot; - } - } - MPI_Bcast(&(bonded_ia_params[i]), sizeof(Bonded_ia_parameters), MPI_BYTE, 0, - comm_cart); - // for tabulated potentials we have to send the tables extra - if (is_tabulated_bond(bonded_ia_params[i].type)) { - if (this_node) { - bonded_ia_params[i].p.tab.pot = new TabulatedPotential(); - } - boost::mpi::broadcast(comm_cart, *bonded_ia_params[i].p.tab.pot, 0); - } - } - - on_short_range_ia_change(); -} - -/*************** BCAST IA SIZE ************/ - -REGISTER_CALLBACK(realloc_ia_params) -void mpi_realloc_ia_params(int ns) { mpi_call_all(realloc_ia_params, ns); } - -/*************** GATHER ************/ -void mpi_gather_stats(GatherStats job, double *result) { - auto job_slave = static_cast(job); - switch (job) { - case GatherStats::energy: - mpi_call(mpi_gather_stats_slave, -1, job_slave); - energy_calc(sim_time); - break; - case GatherStats::pressure: - mpi_call(mpi_gather_stats_slave, -1, job_slave); - pressure_calc(); - break; - case GatherStats::lb_fluid_momentum: - mpi_call(mpi_gather_stats_slave, -1, job_slave); - lb_calc_fluid_momentum(result, lbpar, lbfields, lblattice); - break; -#ifdef LB_BOUNDARIES - case GatherStats::lb_boundary_forces: - mpi_call(mpi_gather_stats_slave, -1, job_slave); - lb_collect_boundary_forces(result); - break; -#endif - default: - fprintf( - stderr, - "%d: INTERNAL ERROR: illegal request %d for mpi_gather_stats_slave\n", - this_node, job_slave); - errexit(); - } -} - -void mpi_gather_stats_slave(int, int job_slave) { - auto job = static_cast(job_slave); - switch (job) { - case GatherStats::energy: - energy_calc(sim_time); - break; - case GatherStats::pressure: - pressure_calc(); - break; - case GatherStats::lb_fluid_momentum: - lb_calc_fluid_momentum(nullptr, lbpar, lbfields, lblattice); - break; -#ifdef LB_BOUNDARIES - case GatherStats::lb_boundary_forces: - lb_collect_boundary_forces(nullptr); - break; -#endif - default: - fprintf( - stderr, - "%d: INTERNAL ERROR: illegal request %d for mpi_gather_stats_slave\n", - this_node, job_slave); - errexit(); - } -} - -/*************** TIME STEP ************/ -void mpi_set_time_step_slave(double dt) { - time_step = dt; - - on_parameter_change(FIELD_TIMESTEP); -} -REGISTER_CALLBACK(mpi_set_time_step_slave) - -void mpi_set_time_step(double time_s) { - if (time_s <= 0.) - throw std::invalid_argument("time_step must be > 0."); - if (lb_lbfluid_get_lattice_switch() != ActiveLB::NONE) - check_tau_time_step_consistency(lb_lbfluid_get_tau(), time_s); - mpi_call_all(mpi_set_time_step_slave, time_s); -} - -/*************** BCAST COULOMB ************/ -void mpi_bcast_coulomb_params() { -#if defined(ELECTROSTATICS) || defined(DIPOLES) - mpi_call(mpi_bcast_coulomb_params_slave, 1, 0); - mpi_bcast_coulomb_params_slave(-1, 0); -#endif -} - -void mpi_bcast_coulomb_params_slave(int, int) { - -#ifdef ELECTROSTATICS - MPI_Bcast(&coulomb, sizeof(Coulomb_parameters), MPI_BYTE, 0, comm_cart); - - Coulomb::bcast_coulomb_params(); -#endif - -#ifdef DIPOLES - MPI_Bcast(&dipole, sizeof(Dipole_parameters), MPI_BYTE, 0, comm_cart); - - Dipole::set_method_local(dipole.method); - - Dipole::bcast_params(comm_cart); -#endif - -#if defined(ELECTROSTATICS) || defined(DIPOLES) - on_coulomb_change(); -#endif -} - -/****************** RESCALE PARTICLES ************/ - -void mpi_rescale_particles(int dir, double scale) { - int pnode; - - mpi_call(mpi_rescale_particles_slave, -1, dir); - for (pnode = 0; pnode < n_nodes; pnode++) { - if (pnode == this_node) { - local_rescale_particles(dir, scale); - } else { - MPI_Send(&scale, 1, MPI_DOUBLE, pnode, SOME_TAG, comm_cart); - } - } - on_particle_change(); -} - -void mpi_rescale_particles_slave(int, int dir) { - double scale = 0.0; - MPI_Recv(&scale, 1, MPI_DOUBLE, 0, SOME_TAG, comm_cart, MPI_STATUS_IGNORE); - local_rescale_particles(dir, scale); - on_particle_change(); -} - -/*************** BCAST CELL STRUCTURE *****************/ - -REGISTER_CALLBACK(cells_re_init) - -void mpi_bcast_cell_structure(int cs) { mpi_call_all(cells_re_init, cs); } - -REGISTER_CALLBACK(cells_set_use_verlet_lists) - -void mpi_set_use_verlet_lists(bool use_verlet_lists) { - mpi_call_all(cells_set_use_verlet_lists, use_verlet_lists); -} - -/*************** BCAST NPTISO GEOM *****************/ - -#ifdef NPT -REGISTER_CALLBACK(mpi_bcast_nptiso_geom_worker) - -void mpi_bcast_nptiso_geom() { - mpi_call(mpi_bcast_nptiso_geom_worker, -1, 0); - mpi_bcast_nptiso_geom_worker(-1, 0); -} -#endif - -/*************** BCAST STEEPEST DESCENT *****************/ - -void mpi_bcast_steepest_descent() { - mpi_call(mpi_bcast_steepest_descent_worker, -1, 0); - mpi_bcast_steepest_descent_worker(-1, 0); -} - -/******************* BCAST CUDA GLOBAL PART VARS ********************/ - -void mpi_bcast_cuda_global_part_vars() { -#ifdef CUDA - mpi_call(mpi_bcast_cuda_global_part_vars_slave, 1, - 0); // third parameter is meaningless - mpi_bcast_cuda_global_part_vars_slave(-1, 0); -#endif -} - -void mpi_bcast_cuda_global_part_vars_slave(int, int) { -#ifdef CUDA - MPI_Bcast(gpu_get_global_particle_vars_pointer_host(), - sizeof(CUDA_global_part_vars), MPI_BYTE, 0, comm_cart); - espressoSystemInterface.requestParticleStructGpu(); -#endif -} - -/********************* SET EXCLUSION ********/ -#ifdef EXCLUSIONS -void mpi_send_exclusion_slave(int part1, int part2, int _delete) { - local_change_exclusion(part1, part2, _delete); - on_particle_change(); -} - -REGISTER_CALLBACK(mpi_send_exclusion_slave) - -void mpi_send_exclusion(int part1, int part2, int _delete) { - mpi_call(mpi_send_exclusion_slave, part1, part2, _delete); - mpi_send_exclusion_slave(part1, part2, _delete); -} -#endif - -/********************* ICCP3M INIT ********/ -#ifdef ELECTROSTATICS -void mpi_iccp3m_init_slave(const iccp3m_struct &iccp3m_cfg_) { - iccp3m_cfg = iccp3m_cfg_; - - on_particle_charge_change(); - check_runtime_errors(comm_cart); -} - -REGISTER_CALLBACK(mpi_iccp3m_init_slave) - -int mpi_iccp3m_init() { - mpi_call(mpi_iccp3m_init_slave, iccp3m_cfg); - - on_particle_charge_change(); - return check_runtime_errors(comm_cart); -} -#endif - -/********************* CALC MU MAX ********/ - -#ifdef DP3M -REGISTER_CALLBACK(calc_mu_max) - -void mpi_bcast_max_mu() { mpi_call_all(calc_mu_max); } -#endif - -/***** GALILEI TRANSFORM AND ASSOCIATED FUNCTIONS ****/ -void mpi_kill_particle_motion_slave(int rotation) { - local_kill_particle_motion(rotation, cell_structure.local_particles()); - on_particle_change(); -} - -REGISTER_CALLBACK(mpi_kill_particle_motion_slave) - -void mpi_kill_particle_motion(int rotation) { - mpi_call_all(mpi_kill_particle_motion_slave, rotation); -} - -void mpi_kill_particle_forces_slave(int torque) { - local_kill_particle_forces(torque, cell_structure.local_particles()); - on_particle_change(); -} - -REGISTER_CALLBACK(mpi_kill_particle_forces_slave) - -void mpi_kill_particle_forces(int torque) { - mpi_call_all(mpi_kill_particle_forces_slave, torque); -} - -struct pair_sum { - template - auto operator()(std::pair l, std::pair r) const { - return std::pair{l.first + r.first, l.second + r.second}; - } -}; - -Utils::Vector3d mpi_system_CMS() { - auto const data = - mpi_call(Communication::Result::reduction, pair_sum{}, local_system_CMS); - return data.first / data.second; -} - -REGISTER_CALLBACK_REDUCTION(local_system_CMS_velocity, pair_sum{}) - -Utils::Vector3d mpi_system_CMS_velocity() { - auto const data = mpi_call(Communication::Result::reduction, pair_sum{}, - local_system_CMS_velocity); - return data.first / data.second; -} - -REGISTER_CALLBACK_REDUCTION(local_system_CMS, pair_sum{}) - -void mpi_galilei_transform_slave(Utils::Vector3d const &cmsvel) { - local_galilei_transform(cmsvel); - on_particle_change(); -} - -REGISTER_CALLBACK(mpi_galilei_transform_slave) - -void mpi_galilei_transform() { - auto const cmsvel = mpi_system_CMS_velocity(); - - mpi_call_all(mpi_galilei_transform_slave, cmsvel); -} - -/*********************** MAIN LOOP for slaves ****************/ - void mpi_loop() { if (this_node != 0) mpiCallbacks().loop(); } - -std::vector mpi_resort_particles(int global_flag) { - mpi_call(mpi_resort_particles_slave, global_flag, 0); - cell_structure.resort_particles(global_flag); - - clear_particle_node(); - - std::vector n_parts; - boost::mpi::gather(comm_cart, - static_cast(cell_structure.local_particles().size()), - n_parts, 0); - - return n_parts; -} - -void mpi_resort_particles_slave(int global_flag, int) { - cell_structure.resort_particles(global_flag); - - boost::mpi::gather( - comm_cart, static_cast(cell_structure.local_particles().size()), 0); -} diff --git a/src/core/communication.hpp b/src/core/communication.hpp index 32b5b008074..4288c70d01a 100644 --- a/src/core/communication.hpp +++ b/src/core/communication.hpp @@ -39,23 +39,16 @@ * following: * - write the @c mpi_* function that is executed on the master * - write the @c mpi_*_slave function - * - in communication.cpp add your slave function to \ref CALLBACK_LIST or - * register it with one of the @c REGISTER_CALLBACK macros + * - register your slave function with one of the @c REGISTER_CALLBACK macros * * After this, your procedure is free to do anything. However, it has * to be in (MPI) sync with what your new @c mpi_*_slave does. This * procedure is called immediately after the broadcast with the - * arbitrary integer as parameter. To this aim it has also to be added - * to \ref CALLBACK_LIST. + * arbitrary integer as parameter. */ #include "MpiCallbacks.hpp" -/* Includes needed by callbacks. */ -#include "Particle.hpp" -#include "cuda_init.hpp" -#include "grid_based_algorithms/lb_constants.hpp" - #include #include @@ -65,13 +58,6 @@ extern int this_node; extern int n_nodes; /** The communicator */ extern boost::mpi::communicator comm_cart; -/** Statistics to calculate */ -enum class GatherStats : int { - energy, - pressure, - lb_fluid_momentum, - lb_boundary_forces -}; /** * Default MPI tag used by callbacks. @@ -153,147 +139,6 @@ auto mpi_call(Tag tag, TagArg &&tag_arg, R (*fp)(Args...), ArgRef &&... args) { /** Process requests from master node. Slave nodes main loop. */ void mpi_loop(); -/** Move particle to a position on a node. - * Also calls \ref on_particle_change. - * \param id the particle to move. - * \param node the node to attach it to. - * \param pos the particles position. - */ -void mpi_place_particle(int node, int id, const Utils::Vector3d &pos); - -/** Create particle at a position on a node. - * Also calls \ref on_particle_change. - * \param id the particle to create. - * \param pos the particles position. - */ -int mpi_place_new_particle(int id, const Utils::Vector3d &pos); - -/** Send exclusions. - * Also calls \ref on_particle_change. - * \param part identity of first particle of the exclusion. - * \param part2 identity of second particle of the exclusion. - * \param _delete if true, do not add the exclusion, rather delete it if found - */ -void mpi_send_exclusion(int part, int part2, int _delete); - -/** Remove a particle. - * Also calls \ref on_particle_change. - * \param id the particle to remove. - * \param node the node it is attached to. - */ -void mpi_remove_particle(int node, int id); - -/** Start integrator. - * @param n_steps how many steps to do. - * @param reuse_forces whether to trust the old forces for the first half step - * @return nonzero on error - */ -int mpi_integrate(int n_steps, int reuse_forces); - -/** Steepest descent main integration loop - * - * Integration stops when the maximal force is lower than the user limit - * @ref SteepestDescentParameters::f_max "f_max" or when the maximal number - * of steps @p steps is reached. - * - * @param steps Maximal number of integration steps - * @return number of integrated steps - */ -int mpi_steepest_descent(int steps); - -/** Broadcast steepest descent parameters */ -void mpi_bcast_steepest_descent(); - -void mpi_bcast_all_ia_params(); - -/** Send new IA params. - * Also calls \ref on_short_range_ia_change. - * - * Used for both bonded and non-bonded interaction parameters. Therefore - * @p i and @p j are used depending on their value: - * - * \param i particle type for non-bonded interaction parameters / - * bonded interaction type number. - * \param j if not negative: particle type for non-bonded interaction - * parameters / if negative: flag for bonded interaction - */ -void mpi_bcast_ia_params(int i, int j); - -/** Resize \ref ia_params. - * \param s the new size for \ref ia_params. - */ -void mpi_realloc_ia_params(int s); - -/** Gather data for analysis. - * \param[in] job what to do: - * \arg for \ref GatherStats::energy, calculate and reduce (sum up) - * energies, using \ref energy_calc. - * \arg for \ref GatherStats::pressure, calculate and reduce (sum up) - * pressure, using \ref pressure_calc. - * \arg for \ref GatherStats::lb_fluid_momentum, use - * \ref lb_calc_fluid_momentum. - * \arg for \ref GatherStats::lb_boundary_forces, use - * \ref lb_collect_boundary_forces. - * \param[out] result where to store values gathered by - * \ref GatherStats::lb_fluid_momentum, - * \ref GatherStats::lb_boundary_forces - */ -void mpi_gather_stats(GatherStats job, double *result = nullptr); - -/** Send new \ref time_step and rescale the velocities accordingly. */ -void mpi_set_time_step(double time_step); - -/** Send new Coulomb parameters. */ -void mpi_bcast_coulomb_params(); - -/** Rescale all particle positions in direction @p dir by a factor @p scale. */ -void mpi_rescale_particles(int dir, double scale); - -/** Change the cell structure on all nodes. */ -void mpi_bcast_cell_structure(int cs); - -void mpi_set_use_verlet_lists(bool use_verlet_lists); - -/** Broadcast nptiso geometry parameters to all nodes. */ -void mpi_bcast_nptiso_geom(); - -/** Broadcast @ref CUDA_global_part_vars structure */ -void mpi_bcast_cuda_global_part_vars(); - -/** Perform iccp3m initialization. - * @return non-zero value on error - */ -int mpi_iccp3m_init(); - -/** Calculate the maximal dipole moment in the system (part of MDLC) */ -void mpi_bcast_max_mu(); - -/** @name Galilei and other - * - set all particle velocities and rotational inertias to zero - * - set all forces and torques on the particles to zero - * - calculate the centre of mass (CMS) - * - calculate the velocity of the CMS - * - remove the CMS velocity from the system - */ -/*@{*/ -void mpi_kill_particle_motion(int rotation); -void mpi_kill_particle_forces(int torque); -Utils::Vector3d mpi_system_CMS(); -Utils::Vector3d mpi_system_CMS_velocity(); -void mpi_galilei_transform(); -/*@}*/ - -/** - * @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 mpi_resort_particles(int global_flag); - namespace Communication { /** * @brief Init globals for communication. diff --git a/src/core/cuda_interface.cpp b/src/core/cuda_interface.cpp index 909a3da4c0a..3b6c07a8a27 100644 --- a/src/core/cuda_interface.cpp +++ b/src/core/cuda_interface.cpp @@ -21,6 +21,7 @@ #ifdef CUDA +#include "EspressoSystemInterface.hpp" #include "communication.hpp" #include "config.hpp" #include "grid.hpp" @@ -160,6 +161,16 @@ void cuda_mpi_send_forces(const ParticleRange &particles, } } -void cuda_bcast_global_part_params() { mpi_bcast_cuda_global_part_vars(); } +void cuda_bcast_global_part_params_local(int, int) { + MPI_Bcast(gpu_get_global_particle_vars_pointer_host(), + sizeof(CUDA_global_part_vars), MPI_BYTE, 0, comm_cart); + espressoSystemInterface.requestParticleStructGpu(); +} + +REGISTER_CALLBACK(cuda_bcast_global_part_params_local) + +void cuda_bcast_global_part_params() { + mpi_call_all(cuda_bcast_global_part_params_local, -1, -1); +} #endif /* ifdef CUDA */ diff --git a/src/core/dpd.cpp b/src/core/dpd.cpp index cee30551c44..78ca46960c6 100644 --- a/src/core/dpd.cpp +++ b/src/core/dpd.cpp @@ -24,16 +24,19 @@ #include "dpd.hpp" #ifdef DPD +#include "MpiCallbacks.hpp" #include "cells.hpp" #include "communication.hpp" #include "event.hpp" #include "grid.hpp" #include "integrate.hpp" +#include "interactions.hpp" #include "nonbonded_interactions/nonbonded_interaction_data.hpp" #include "random.hpp" #include "thermostat.hpp" #include +#include #include #include diff --git a/src/core/electrostatics_magnetostatics/CMakeLists.txt b/src/core/electrostatics_magnetostatics/CMakeLists.txt index 1a473eabbf0..498dcf79e30 100644 --- a/src/core/electrostatics_magnetostatics/CMakeLists.txt +++ b/src/core/electrostatics_magnetostatics/CMakeLists.txt @@ -1,6 +1,7 @@ target_sources( EspressoCore - PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/debye_hueckel.cpp + PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/common.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/debye_hueckel.cpp ${CMAKE_CURRENT_SOURCE_DIR}/elc.cpp ${CMAKE_CURRENT_SOURCE_DIR}/icc.cpp ${CMAKE_CURRENT_SOURCE_DIR}/magnetic_non_p3m_methods.cpp diff --git a/src/core/electrostatics_magnetostatics/common.cpp b/src/core/electrostatics_magnetostatics/common.cpp new file mode 100644 index 00000000000..ac9808a88d4 --- /dev/null +++ b/src/core/electrostatics_magnetostatics/common.cpp @@ -0,0 +1,55 @@ +/* + * Copyright (C) 2010-2020 The ESPResSo project + * Copyright (C) 2002,2003,2004,2005,2006,2007,2008,2009,2010 + * Max-Planck-Institute for Polymer Research, Theory Group + * + * This file is part of ESPResSo. + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ + +#include "communication.hpp" +#include "event.hpp" + +#include "electrostatics_magnetostatics/common.hpp" + +#include "electrostatics_magnetostatics/coulomb.hpp" +#include "electrostatics_magnetostatics/dipole.hpp" + +#include + +void mpi_bcast_coulomb_params_local(int, int) { +#ifdef ELECTROSTATICS + MPI_Bcast(&coulomb, sizeof(Coulomb_parameters), MPI_BYTE, 0, comm_cart); + Coulomb::bcast_coulomb_params(); +#endif + +#ifdef DIPOLES + MPI_Bcast(&dipole, sizeof(Dipole_parameters), MPI_BYTE, 0, comm_cart); + Dipole::set_method_local(dipole.method); + Dipole::bcast_params(comm_cart); +#endif + +#if defined(ELECTROSTATICS) || defined(DIPOLES) + on_coulomb_change(); +#endif +} + +REGISTER_CALLBACK(mpi_bcast_coulomb_params_local) + +void mpi_bcast_coulomb_params() { +#if defined(ELECTROSTATICS) || defined(DIPOLES) + mpi_call_all(mpi_bcast_coulomb_params_local, -1, -1); +#endif +} diff --git a/src/core/electrostatics_magnetostatics/common.hpp b/src/core/electrostatics_magnetostatics/common.hpp new file mode 100644 index 00000000000..58609d58755 --- /dev/null +++ b/src/core/electrostatics_magnetostatics/common.hpp @@ -0,0 +1,25 @@ +/* + * Copyright (C) 2020 The ESPResSo project + * + * This file is part of ESPResSo. + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ +#ifndef ELECTROSTATICS_MAGNETOSTATICS_COMMON_HPP +#define ELECTROSTATICS_MAGNETOSTATICS_COMMON_HPP + +/** Send new Coulomb/Dipole parameters. */ +void mpi_bcast_coulomb_params(); + +#endif diff --git a/src/core/electrostatics_magnetostatics/coulomb.cpp b/src/core/electrostatics_magnetostatics/coulomb.cpp index 41f5f3122dd..524472e2553 100644 --- a/src/core/electrostatics_magnetostatics/coulomb.cpp +++ b/src/core/electrostatics_magnetostatics/coulomb.cpp @@ -22,6 +22,7 @@ #ifdef ELECTROSTATICS #include "cells.hpp" #include "communication.hpp" +#include "electrostatics_magnetostatics/common.hpp" #include "electrostatics_magnetostatics/debye_hueckel.hpp" #include "electrostatics_magnetostatics/elc.hpp" #include "electrostatics_magnetostatics/icc.hpp" diff --git a/src/core/electrostatics_magnetostatics/debye_hueckel.cpp b/src/core/electrostatics_magnetostatics/debye_hueckel.cpp index 0770dc833ae..8be9afee689 100644 --- a/src/core/electrostatics_magnetostatics/debye_hueckel.cpp +++ b/src/core/electrostatics_magnetostatics/debye_hueckel.cpp @@ -26,7 +26,7 @@ #include "debye_hueckel.hpp" #ifdef ELECTROSTATICS -#include "communication.hpp" +#include "electrostatics_magnetostatics/common.hpp" Debye_hueckel_params dh_params{}; diff --git a/src/core/electrostatics_magnetostatics/dipole.cpp b/src/core/electrostatics_magnetostatics/dipole.cpp index d52d70251b4..9c28dec04bd 100644 --- a/src/core/electrostatics_magnetostatics/dipole.cpp +++ b/src/core/electrostatics_magnetostatics/dipole.cpp @@ -23,6 +23,7 @@ #include "actor/DipolarBarnesHut.hpp" #include "actor/DipolarDirectSum.hpp" +#include "electrostatics_magnetostatics/common.hpp" #include "electrostatics_magnetostatics/magnetic_non_p3m_methods.hpp" #include "electrostatics_magnetostatics/mdlc_correction.hpp" #include "electrostatics_magnetostatics/p3m-dipolar.hpp" diff --git a/src/core/electrostatics_magnetostatics/elc.cpp b/src/core/electrostatics_magnetostatics/elc.cpp index a938f692349..cbfd73d283c 100644 --- a/src/core/electrostatics_magnetostatics/elc.cpp +++ b/src/core/electrostatics_magnetostatics/elc.cpp @@ -30,6 +30,7 @@ #include #include +#include "electrostatics_magnetostatics/common.hpp" #include "electrostatics_magnetostatics/coulomb.hpp" #include "electrostatics_magnetostatics/elc.hpp" #include "electrostatics_magnetostatics/p3m.hpp" diff --git a/src/core/electrostatics_magnetostatics/icc.cpp b/src/core/electrostatics_magnetostatics/icc.cpp index 7400223d9f3..e8a4411ffda 100644 --- a/src/core/electrostatics_magnetostatics/icc.cpp +++ b/src/core/electrostatics_magnetostatics/icc.cpp @@ -210,4 +210,20 @@ void init_forces_iccp3m(const ParticleRange &particles, } } +void mpi_iccp3m_init_local(const iccp3m_struct &iccp3m_cfg_) { + iccp3m_cfg = iccp3m_cfg_; + + on_particle_charge_change(); + check_runtime_errors(comm_cart); +} + +REGISTER_CALLBACK(mpi_iccp3m_init_local) + +int mpi_iccp3m_init() { + mpi_call(mpi_iccp3m_init_local, iccp3m_cfg); + + on_particle_charge_change(); + return check_runtime_errors(comm_cart); +} + #endif diff --git a/src/core/electrostatics_magnetostatics/icc.hpp b/src/core/electrostatics_magnetostatics/icc.hpp index befec352658..8f0debf4ab0 100644 --- a/src/core/electrostatics_magnetostatics/icc.hpp +++ b/src/core/electrostatics_magnetostatics/icc.hpp @@ -103,6 +103,11 @@ void iccp3m_alloc_lists(); */ int iccp3m_sanity_check(); +/** Perform iccp3m initialization. + * @return non-zero value on error + */ +int mpi_iccp3m_init(); + #endif /* ELECTROSTATICS */ #endif /* ICCP3M_H */ diff --git a/src/core/electrostatics_magnetostatics/magnetic_non_p3m_methods.cpp b/src/core/electrostatics_magnetostatics/magnetic_non_p3m_methods.cpp index 321196110a5..c7d4ab8dbf1 100644 --- a/src/core/electrostatics_magnetostatics/magnetic_non_p3m_methods.cpp +++ b/src/core/electrostatics_magnetostatics/magnetic_non_p3m_methods.cpp @@ -24,7 +24,7 @@ #ifdef DIPOLES #include "cells.hpp" #include "communication.hpp" -#include "dipole.hpp" +#include "electrostatics_magnetostatics/common.hpp" #include "electrostatics_magnetostatics/dipole.hpp" #include "errorhandling.hpp" #include "grid.hpp" diff --git a/src/core/electrostatics_magnetostatics/mdlc_correction.cpp b/src/core/electrostatics_magnetostatics/mdlc_correction.cpp index 62b6b3bf69f..df8e3063112 100644 --- a/src/core/electrostatics_magnetostatics/mdlc_correction.cpp +++ b/src/core/electrostatics_magnetostatics/mdlc_correction.cpp @@ -27,6 +27,7 @@ #include "Particle.hpp" #include "cells.hpp" #include "communication.hpp" +#include "electrostatics_magnetostatics/common.hpp" #include "electrostatics_magnetostatics/dipole.hpp" #include "electrostatics_magnetostatics/p3m-dipolar.hpp" #include "errorhandling.hpp" @@ -36,19 +37,26 @@ #include -DLC_struct dlc_params = {1e100, 0, 0, 0, 0}; +#include +#include -static double mu_max; +DLC_struct dlc_params = {1e100, 0, 0, 0, 0}; -void calc_mu_max() { - auto local_particles = cell_structure.local_particles(); - mu_max = std::accumulate( +/** Calculate the maximal dipole moment in the system */ +double calc_mu_max() { + auto const local_particles = cell_structure.local_particles(); + auto const mu_max_local = std::accumulate( local_particles.begin(), local_particles.end(), 0.0, [](double mu, Particle const &p) { return std::max(mu, p.p.dipm); }); - MPI_Allreduce(MPI_IN_PLACE, &mu_max, 1, MPI_DOUBLE, MPI_MAX, comm_cart); + double mu_max; + boost::mpi::reduce(comm_cart, mu_max_local, mu_max, + boost::mpi::maximum(), 0); + return mu_max; } +REGISTER_CALLBACK_MASTER_RANK(calc_mu_max) + inline double g1_DLC_dip(double g, double x) { auto const c = g / x; auto const cc2 = c * c; @@ -409,8 +417,8 @@ double add_mdlc_energy_corrections(const ParticleRange &particles) { /** Compute the cut-off in the DLC dipolar part to get a certain accuracy. * We assume particles to have all the same value of the dipolar momentum - * modulus (@ref mu_max). @ref mu_max is taken as the largest value of mu - * inside the system. If we assume the gap has a width @c gap_size (within + * modulus, which is taken as the largest value of mu inside the system. + * If we assume the gap has a width @c gap_size (within * which there is no particles): Lz = h + gap_size * * BE CAREFUL: @@ -419,9 +427,10 @@ double add_mdlc_energy_corrections(const ParticleRange &particles) { * it makes no sense to have an accurate result for DLC-dipolar. */ int mdlc_tune(double error) { - mpi_bcast_max_mu(); /* we take the maximum dipole in the system, to be sure - that the errors in the other case - will be equal or less than for this one */ + /* we take the maximum dipole in the system, to be sure that the errors + * in the other case will be equal or less than for this one */ + auto const mu_max = mpi_call(Communication::Result::master_rank, calc_mu_max); + auto const mu_max_sq = mu_max * mu_max; const double n = get_n_part(); auto const lx = box_geo.length()[0]; @@ -455,8 +464,7 @@ int mdlc_tune(double error) { 9.0 * exp(-2.0 * gc * h) * g1_DLC_dip(gc, lz + h)); auto const fa1 = 0.5 * sqrt(Utils::pi() / (2.0 * a)) * fa0; auto const fa2 = g2_DLC_dip(gc, lz); - auto const de = - n * (mu_max * mu_max) / (4.0 * (exp(gc * lz) - 1.0)) * (fa1 + fa2); + auto const de = n * mu_max_sq / (4.0 * (exp(gc * lz) - 1.0)) * (fa1 + fa2); if (de < error) { flag = true; break; diff --git a/src/core/electrostatics_magnetostatics/mdlc_correction.hpp b/src/core/electrostatics_magnetostatics/mdlc_correction.hpp index 0221dde72a5..9ee47f92e0d 100644 --- a/src/core/electrostatics_magnetostatics/mdlc_correction.hpp +++ b/src/core/electrostatics_magnetostatics/mdlc_correction.hpp @@ -81,8 +81,6 @@ int mdlc_set_params(double maxPWerror, double gap_size, double far_cut); int mdlc_sanity_checks(); void add_mdlc_force_corrections(const ParticleRange &particles); double add_mdlc_energy_corrections(const ParticleRange &particles); -/** Calculate the maximal dipole moment in the system */ -void calc_mu_max(); #endif // DP3M #endif diff --git a/src/core/electrostatics_magnetostatics/mmm1d.cpp b/src/core/electrostatics_magnetostatics/mmm1d.cpp index 070103da13c..dafd3df6021 100644 --- a/src/core/electrostatics_magnetostatics/mmm1d.cpp +++ b/src/core/electrostatics_magnetostatics/mmm1d.cpp @@ -28,13 +28,13 @@ #ifdef ELECTROSTATICS #include "cells.hpp" -#include "communication.hpp" #include "errorhandling.hpp" #include "grid.hpp" #include "mmm-common.hpp" #include "specfunc.hpp" #include "tuning.hpp" +#include "electrostatics_magnetostatics/common.hpp" #include "electrostatics_magnetostatics/coulomb.hpp" #include diff --git a/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp b/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp index 757d66c7adc..a2f660aaa96 100644 --- a/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp +++ b/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp @@ -33,6 +33,7 @@ */ #include "electrostatics_magnetostatics/p3m-dipolar.hpp" +#include "electrostatics_magnetostatics/common.hpp" #include "electrostatics_magnetostatics/dp3m_influence_function.hpp" #ifdef DP3M diff --git a/src/core/electrostatics_magnetostatics/p3m.cpp b/src/core/electrostatics_magnetostatics/p3m.cpp index c90d98226dc..c32a6f35db2 100644 --- a/src/core/electrostatics_magnetostatics/p3m.cpp +++ b/src/core/electrostatics_magnetostatics/p3m.cpp @@ -29,6 +29,7 @@ #include "Particle.hpp" #include "cells.hpp" #include "communication.hpp" +#include "electrostatics_magnetostatics/common.hpp" #include "electrostatics_magnetostatics/coulomb.hpp" #include "electrostatics_magnetostatics/elc.hpp" #include "electrostatics_magnetostatics/p3m_influence_function.hpp" diff --git a/src/core/electrostatics_magnetostatics/reaction_field.cpp b/src/core/electrostatics_magnetostatics/reaction_field.cpp index 3b504f176f9..30271a11f6c 100644 --- a/src/core/electrostatics_magnetostatics/reaction_field.cpp +++ b/src/core/electrostatics_magnetostatics/reaction_field.cpp @@ -25,7 +25,7 @@ #include "reaction_field.hpp" #ifdef ELECTROSTATICS -#include "communication.hpp" +#include "common.hpp" Reaction_field_params rf_params{}; diff --git a/src/core/energy.cpp b/src/core/energy.cpp index aa487a7c9ad..22d54a9e95f 100644 --- a/src/core/energy.cpp +++ b/src/core/energy.cpp @@ -31,6 +31,7 @@ #include "energy_inline.hpp" #include "event.hpp" #include "forces.hpp" +#include "integrate.hpp" #include "nonbonded_interactions/nonbonded_interaction_data.hpp" #include "reduce_observable_stat.hpp" @@ -46,9 +47,6 @@ Observable_stat obs_energy{1}; Observable_stat const &get_obs_energy() { return obs_energy; } -/** Reduce the system energy from all MPI ranks. */ -void master_energy_calc() { mpi_gather_stats(GatherStats::energy); } - void energy_calc(const double time) { if (!interactions_sanity_checks()) return; @@ -107,7 +105,11 @@ void energy_calc(const double time) { } } -void update_energy() { master_energy_calc(); } +void update_energy_local(int, int) { energy_calc(sim_time); } + +REGISTER_CALLBACK(update_energy_local) + +void update_energy() { mpi_call_all(update_energy_local, -1, -1); } void calc_long_range_energies(const ParticleRange &particles) { #ifdef ELECTROSTATICS diff --git a/src/core/galilei.cpp b/src/core/galilei.cpp index b261a864612..58828b785b4 100644 --- a/src/core/galilei.cpp +++ b/src/core/galilei.cpp @@ -20,11 +20,18 @@ */ #include "galilei.hpp" +#include "ParticleRange.hpp" #include "cells.hpp" +#include "communication.hpp" #include "event.hpp" #include "grid.hpp" #include +#include + +#include + +#include /** Stop particle motion by setting the velocity of each particle to zero. */ void local_kill_particle_motion(int omega, const ParticleRange &particles) { @@ -82,3 +89,60 @@ void local_galilei_transform(const Utils::Vector3d &cmsvel) { p.m.v -= cmsvel; } } + +void mpi_kill_particle_motion_local(int rotation) { + local_kill_particle_motion(rotation, cell_structure.local_particles()); + on_particle_change(); +} + +REGISTER_CALLBACK(mpi_kill_particle_motion_local) + +void mpi_kill_particle_motion(int rotation) { + mpi_call_all(mpi_kill_particle_motion_local, rotation); +} + +void mpi_kill_particle_forces_local(int torque) { + local_kill_particle_forces(torque, cell_structure.local_particles()); + on_particle_change(); +} + +REGISTER_CALLBACK(mpi_kill_particle_forces_local) + +void mpi_kill_particle_forces(int torque) { + mpi_call_all(mpi_kill_particle_forces_local, torque); +} + +struct pair_sum { + template + auto operator()(std::pair l, std::pair r) const { + return std::pair{l.first + r.first, l.second + r.second}; + } +}; + +Utils::Vector3d mpi_system_CMS() { + auto const data = + mpi_call(Communication::Result::reduction, pair_sum{}, local_system_CMS); + return data.first / data.second; +} + +REGISTER_CALLBACK_REDUCTION(local_system_CMS_velocity, pair_sum{}) + +Utils::Vector3d mpi_system_CMS_velocity() { + auto const data = mpi_call(Communication::Result::reduction, pair_sum{}, + local_system_CMS_velocity); + return data.first / data.second; +} + +REGISTER_CALLBACK_REDUCTION(local_system_CMS, pair_sum{}) + +void mpi_galilei_transform_local(Utils::Vector3d const &cmsvel) { + local_galilei_transform(cmsvel); + on_particle_change(); +} + +REGISTER_CALLBACK(mpi_galilei_transform_local) + +void mpi_galilei_transform() { + auto const cmsvel = mpi_system_CMS_velocity(); + mpi_call_all(mpi_galilei_transform_local, cmsvel); +} diff --git a/src/core/galilei.hpp b/src/core/galilei.hpp index 18b3d1ec735..a0b30a1d3b5 100644 --- a/src/core/galilei.hpp +++ b/src/core/galilei.hpp @@ -23,13 +23,23 @@ #include -#include "ParticleRange.hpp" -#include +/** Stop particle motion by setting the velocity of each particle to zero. + * @param rotation if true, also set particle angular velocities to zero + */ +void mpi_kill_particle_motion(int rotation); + +/** Set all the forces acting on the particles to zero. + * @param torque if true, also set particle torques to zero + */ +void mpi_kill_particle_forces(int torque); + +/** Calculate the CMS of the system */ +Utils::Vector3d mpi_system_CMS(); + +/** Calculate the CMS velocity of the system */ +Utils::Vector3d mpi_system_CMS_velocity(); -void local_kill_particle_motion(int, const ParticleRange &particles); -void local_kill_particle_forces(int, const ParticleRange &particles); -std::pair local_system_CMS(); -std::pair local_system_CMS_velocity(); -void local_galilei_transform(const Utils::Vector3d &cmsvel); +/** Remove the CMS velocity */ +void mpi_galilei_transform(); #endif diff --git a/src/core/grid.cpp b/src/core/grid.cpp index d71a53ea235..4cf5841120b 100644 --- a/src/core/grid.cpp +++ b/src/core/grid.cpp @@ -28,6 +28,7 @@ #include "cells.hpp" #include "communication.hpp" #include "global.hpp" +#include "particle_data.hpp" #include #include diff --git a/src/core/grid_based_algorithms/lb_boundaries.cpp b/src/core/grid_based_algorithms/lb_boundaries.cpp index 90b15909bc8..74639c0b8d3 100644 --- a/src/core/grid_based_algorithms/lb_boundaries.cpp +++ b/src/core/grid_based_algorithms/lb_boundaries.cpp @@ -253,6 +253,14 @@ void lb_init_boundaries() { } } +#if defined(LB_BOUNDARIES) +void lb_collect_boundary_forces_local(int, int) { + lb_collect_boundary_forces(nullptr); +} + +REGISTER_CALLBACK(lb_collect_boundary_forces_local) +#endif + Utils::Vector3d lbboundary_get_force(LBBoundary const *lbb) { Utils::Vector3d force{}; #if defined(LB_BOUNDARIES) || defined(LB_BOUNDARIES_GPU) @@ -271,7 +279,8 @@ Utils::Vector3d lbboundary_get_force(LBBoundary const *lbb) { #endif } else if (lattice_switch == ActiveLB::CPU) { #if defined(LB_BOUNDARIES) - mpi_gather_stats(GatherStats::lb_boundary_forces, forces.data()); + mpi_call(lb_collect_boundary_forces_local, -1, -1); + lb_collect_boundary_forces(forces.data()); #endif } auto const container_index = std::distance(lbboundaries.begin(), it); diff --git a/src/core/grid_based_algorithms/lb_interface.cpp b/src/core/grid_based_algorithms/lb_interface.cpp index e67b3361c70..55824f9534c 100644 --- a/src/core/grid_based_algorithms/lb_interface.cpp +++ b/src/core/grid_based_algorithms/lb_interface.cpp @@ -1221,6 +1221,12 @@ const Lattice &lb_lbfluid_get_lattice() { return lblattice; } ActiveLB lb_lbfluid_get_lattice_switch() { return lattice_switch; } +void mpi_lb_lbfluid_calc_fluid_momentum_local(int, int) { + lb_calc_fluid_momentum(nullptr, lbpar, lbfields, lblattice); +} + +REGISTER_CALLBACK(mpi_lb_lbfluid_calc_fluid_momentum_local) + Utils::Vector3d lb_lbfluid_calc_fluid_momentum() { Utils::Vector3d fluid_momentum{}; if (lattice_switch == ActiveLB::GPU) { @@ -1228,7 +1234,8 @@ Utils::Vector3d lb_lbfluid_calc_fluid_momentum() { lb_calc_fluid_momentum_GPU(fluid_momentum.data()); #endif } else if (lattice_switch == ActiveLB::CPU) { - mpi_gather_stats(GatherStats::lb_fluid_momentum, fluid_momentum.data()); + mpi_call(mpi_lb_lbfluid_calc_fluid_momentum_local, -1, -1); + lb_calc_fluid_momentum(fluid_momentum.data(), lbpar, lbfields, lblattice); } return fluid_momentum; } diff --git a/src/core/immersed_boundary/ImmersedBoundaries.cpp b/src/core/immersed_boundary/ImmersedBoundaries.cpp index fab9f2063e2..094df214855 100644 --- a/src/core/immersed_boundary/ImmersedBoundaries.cpp +++ b/src/core/immersed_boundary/ImmersedBoundaries.cpp @@ -23,6 +23,7 @@ #include "bonded_interactions/bonded_interaction_data.hpp" #include "communication.hpp" #include "grid.hpp" +#include "interactions.hpp" #include diff --git a/src/core/immersed_boundary/ibm_tribend.cpp b/src/core/immersed_boundary/ibm_tribend.cpp index 4f5632d2009..bb5126500b1 100644 --- a/src/core/immersed_boundary/ibm_tribend.cpp +++ b/src/core/immersed_boundary/ibm_tribend.cpp @@ -20,8 +20,8 @@ #include "immersed_boundary/ibm_tribend.hpp" #include "bonded_interactions/bonded_interaction_data.hpp" -#include "communication.hpp" #include "grid.hpp" +#include "interactions.hpp" #include "particle_data.hpp" #include diff --git a/src/core/immersed_boundary/ibm_triel.cpp b/src/core/immersed_boundary/ibm_triel.cpp index a407c94f6cb..f3254e93300 100644 --- a/src/core/immersed_boundary/ibm_triel.cpp +++ b/src/core/immersed_boundary/ibm_triel.cpp @@ -20,8 +20,8 @@ #include "immersed_boundary/ibm_triel.hpp" #include "bonded_interactions/bonded_interaction_data.hpp" -#include "communication.hpp" #include "grid.hpp" +#include "interactions.hpp" #include "particle_data.hpp" #include diff --git a/src/core/integrate.cpp b/src/core/integrate.cpp index 1356fbda0c7..9bac61bfcbb 100644 --- a/src/core/integrate.cpp +++ b/src/core/integrate.cpp @@ -391,6 +391,29 @@ int python_integrate(int n_steps, bool recalc_forces, bool reuse_forces_par) { return ES_OK; } +static int mpi_steepest_descent_local(int steps, int) { + return integrate(steps, -1); +} +REGISTER_CALLBACK_MASTER_RANK(mpi_steepest_descent_local) + +int mpi_steepest_descent(int steps) { + return mpi_call(Communication::Result::master_rank, + mpi_steepest_descent_local, steps, 0); +} + +static int mpi_integrate_local(int n_steps, int reuse_forces) { + integrate(n_steps, reuse_forces); + + return check_runtime_errors_local(); +} + +REGISTER_CALLBACK_REDUCTION(mpi_integrate_local, std::plus()) + +int mpi_integrate(int n_steps, int reuse_forces) { + return mpi_call(Communication::Result::reduction, std::plus(), + mpi_integrate_local, n_steps, reuse_forces); +} + int integrate_set_steepest_descent(const double f_max, const double gamma, const double max_displacement) { if (f_max < 0.0) { @@ -509,3 +532,18 @@ double interaction_range() { auto const max_cut = maximal_cutoff(); return (max_cut > 0.) ? max_cut + skin : INACTIVE_CUTOFF; } + +void mpi_set_time_step_local(double dt) { + time_step = dt; + on_parameter_change(FIELD_TIMESTEP); +} + +REGISTER_CALLBACK(mpi_set_time_step_local) + +void mpi_set_time_step(double time_s) { + if (time_s <= 0.) + throw std::invalid_argument("time_step must be > 0."); + if (lb_lbfluid_get_lattice_switch() != ActiveLB::NONE) + check_tau_time_step_consistency(lb_lbfluid_get_tau(), time_s); + mpi_call_all(mpi_set_time_step_local, time_s); +} diff --git a/src/core/integrate.hpp b/src/core/integrate.hpp index 5c2a090fb9e..c6eb05f27df 100644 --- a/src/core/integrate.hpp +++ b/src/core/integrate.hpp @@ -113,6 +113,24 @@ int integrate(int n_steps, int reuse_forces); */ int python_integrate(int n_steps, bool recalc_forces, bool reuse_forces); +/** Start integrator. + * @param n_steps how many steps to do. + * @param reuse_forces whether to trust the old forces for the first half step + * @return nonzero on error + */ +int mpi_integrate(int n_steps, int reuse_forces); + +/** Steepest descent main integration loop + * + * Integration stops when the maximal force is lower than the user limit + * @ref SteepestDescentParameters::f_max "f_max" or when the maximal number + * of steps @p steps is reached. + * + * @param steps Maximal number of integration steps + * @return number of integrated steps + */ +int mpi_steepest_descent(int steps); + /** @brief Set the steepest descent integrator for energy minimization. * @retval ES_OK on success * @retval ES_ERROR on error @@ -148,4 +166,7 @@ int integrate_set_npt_isotropic(double ext_pressure, double piston, bool xdir_rescale, bool ydir_rescale, bool zdir_rescale, bool cubic_box); +/** Send new \ref time_step and rescale the velocities accordingly. */ +void mpi_set_time_step(double time_step); + #endif diff --git a/src/core/integrators/steepest_descent.cpp b/src/core/integrators/steepest_descent.cpp index 22ba9a72ba8..040f622c159 100644 --- a/src/core/integrators/steepest_descent.cpp +++ b/src/core/integrators/steepest_descent.cpp @@ -110,3 +110,9 @@ void steepest_descent_init(const double f_max, const double gamma, void mpi_bcast_steepest_descent_worker(int, int) { boost::mpi::broadcast(comm_cart, params, 0); } + +REGISTER_CALLBACK(mpi_bcast_steepest_descent_worker) + +void mpi_bcast_steepest_descent() { + mpi_call_all(mpi_bcast_steepest_descent_worker, -1, 0); +} diff --git a/src/core/integrators/steepest_descent.hpp b/src/core/integrators/steepest_descent.hpp index 917b850f955..dd0d1684e28 100644 --- a/src/core/integrators/steepest_descent.hpp +++ b/src/core/integrators/steepest_descent.hpp @@ -56,7 +56,8 @@ struct SteepestDescentParameters { */ void steepest_descent_init(double f_max, double gamma, double max_displacement); -void mpi_bcast_steepest_descent_worker(int, int); +/** Broadcast steepest descent parameters */ +void mpi_bcast_steepest_descent(); /** Steepest descent integrator * @return whether the maximum force/torque encountered is below the user diff --git a/src/core/interactions.cpp b/src/core/interactions.cpp new file mode 100644 index 00000000000..38ef0f7c51a --- /dev/null +++ b/src/core/interactions.cpp @@ -0,0 +1,81 @@ +/* + * Copyright (C) 2010-2020 The ESPResSo project + * Copyright (C) 2002,2003,2004,2005,2006,2007,2008,2009,2010 + * Max-Planck-Institute for Polymer Research, Theory Group + * + * This file is part of ESPResSo. + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ +#include "communication.hpp" + +#include "bonded_interactions/bonded_interaction_data.hpp" +#include "event.hpp" + +#include "serialization/IA_parameters.hpp" + +#include +#include + +#include + +inline bool is_tabulated_bond(BondedInteraction const type) { + return (type == BONDED_IA_TABULATED_DISTANCE or + type == BONDED_IA_TABULATED_ANGLE or + type == BONDED_IA_TABULATED_DIHEDRAL); +} + +void mpi_bcast_all_ia_params_slave() { + boost::mpi::broadcast(comm_cart, ia_params, 0); +} + +REGISTER_CALLBACK(mpi_bcast_all_ia_params_slave) + +void mpi_bcast_all_ia_params() { mpi_call_all(mpi_bcast_all_ia_params_slave); } + +void mpi_bcast_ia_params_slave(int i, int j) { + if (j >= 0) { + // non-bonded interaction parameters + boost::mpi::broadcast(comm_cart, *get_ia_param(i, j), 0); + } else { + // bonded interaction parameters + if (this_node) { + make_bond_type_exist(i); // realloc bonded_ia_params on slave nodes! + if (is_tabulated_bond(bonded_ia_params[i].type)) { + delete bonded_ia_params[i].p.tab.pot; + } + } + MPI_Bcast(&(bonded_ia_params[i]), sizeof(Bonded_ia_parameters), MPI_BYTE, 0, + comm_cart); + // for tabulated potentials we have to send the tables extra + if (is_tabulated_bond(bonded_ia_params[i].type)) { + if (this_node) { + bonded_ia_params[i].p.tab.pot = new TabulatedPotential(); + } + boost::mpi::broadcast(comm_cart, *bonded_ia_params[i].p.tab.pot, 0); + } + } + + on_short_range_ia_change(); +} + +REGISTER_CALLBACK(mpi_bcast_ia_params_slave) + +void mpi_bcast_ia_params(int i, int j) { + mpi_call_all(mpi_bcast_ia_params_slave, i, j); +} + +REGISTER_CALLBACK(realloc_ia_params) + +void mpi_realloc_ia_params(int ns) { mpi_call_all(realloc_ia_params, ns); } diff --git a/src/core/interactions.hpp b/src/core/interactions.hpp new file mode 100644 index 00000000000..549de87030d --- /dev/null +++ b/src/core/interactions.hpp @@ -0,0 +1,48 @@ +/* + * Copyright (C) 2010-2020 The ESPResSo project + * Copyright (C) 2002,2003,2004,2005,2006,2007,2008,2009,2010 + * Max-Planck-Institute for Polymer Research, Theory Group + * + * This file is part of ESPResSo. + * + * ESPResSo is free software: you can redistribute it and/or modify + * it under the terms of the GNU General Public License as published by + * the Free Software Foundation, either version 3 of the License, or + * (at your option) any later version. + * + * ESPResSo is distributed in the hope that it will be useful, + * but WITHOUT ANY WARRANTY; without even the implied warranty of + * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the + * GNU General Public License for more details. + * + * You should have received a copy of the GNU General Public License + * along with this program. If not, see . + */ +/** \file + * This file contains the asynchronous MPI communication for interactions. + */ +#ifndef _INTERACTIONS_HPP +#define _INTERACTIONS_HPP + +/** Broadcast \ref ia_params to all nodes. */ +void mpi_bcast_all_ia_params(); + +/** Send new IA params. + * Also calls \ref on_short_range_ia_change. + * + * Used for both bonded and non-bonded interaction parameters. Therefore + * @p i and @p j are used depending on their value: + * + * \param i particle type for non-bonded interaction parameters / + * bonded interaction type number. + * \param j if not negative: particle type for non-bonded interaction + * parameters / if negative: flag for bonded interaction + */ +void mpi_bcast_ia_params(int i, int j); + +/** Resize \ref ia_params. + * \param s the new size. + */ +void mpi_realloc_ia_params(int s); + +#endif diff --git a/src/core/nonbonded_interactions/bmhtf-nacl.cpp b/src/core/nonbonded_interactions/bmhtf-nacl.cpp index d59b612cbda..a0153881a5d 100644 --- a/src/core/nonbonded_interactions/bmhtf-nacl.cpp +++ b/src/core/nonbonded_interactions/bmhtf-nacl.cpp @@ -25,7 +25,7 @@ #include "bmhtf-nacl.hpp" #ifdef BMHTF_NACL -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/buckingham.cpp b/src/core/nonbonded_interactions/buckingham.cpp index 4b93d869bd4..10f0470a9db 100644 --- a/src/core/nonbonded_interactions/buckingham.cpp +++ b/src/core/nonbonded_interactions/buckingham.cpp @@ -25,7 +25,7 @@ #include "buckingham.hpp" #ifdef BUCKINGHAM -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/gaussian.cpp b/src/core/nonbonded_interactions/gaussian.cpp index 770b558cd72..87f83c48dcc 100644 --- a/src/core/nonbonded_interactions/gaussian.cpp +++ b/src/core/nonbonded_interactions/gaussian.cpp @@ -25,7 +25,7 @@ #include "gaussian.hpp" #ifdef GAUSSIAN -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/gay_berne.cpp b/src/core/nonbonded_interactions/gay_berne.cpp index b4508eef8ce..20049635951 100644 --- a/src/core/nonbonded_interactions/gay_berne.cpp +++ b/src/core/nonbonded_interactions/gay_berne.cpp @@ -25,7 +25,7 @@ #include "gay_berne.hpp" #ifdef GAY_BERNE -#include "communication.hpp" +#include "interactions.hpp" #include #include diff --git a/src/core/nonbonded_interactions/hat.cpp b/src/core/nonbonded_interactions/hat.cpp index 0dc9a0e474a..3559a2477dd 100644 --- a/src/core/nonbonded_interactions/hat.cpp +++ b/src/core/nonbonded_interactions/hat.cpp @@ -25,7 +25,7 @@ #include "hat.hpp" #ifdef HAT -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/hertzian.cpp b/src/core/nonbonded_interactions/hertzian.cpp index f51beea4610..faf0ba2be45 100644 --- a/src/core/nonbonded_interactions/hertzian.cpp +++ b/src/core/nonbonded_interactions/hertzian.cpp @@ -25,7 +25,7 @@ #include "hertzian.hpp" #ifdef HERTZIAN -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/lj.cpp b/src/core/nonbonded_interactions/lj.cpp index 296f9d9c078..8b028ccfa4b 100644 --- a/src/core/nonbonded_interactions/lj.cpp +++ b/src/core/nonbonded_interactions/lj.cpp @@ -25,7 +25,7 @@ #include "lj.hpp" #ifdef LENNARD_JONES -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/ljcos.cpp b/src/core/nonbonded_interactions/ljcos.cpp index b0d4e711dae..c96be9da29d 100644 --- a/src/core/nonbonded_interactions/ljcos.cpp +++ b/src/core/nonbonded_interactions/ljcos.cpp @@ -25,7 +25,7 @@ #include "ljcos.hpp" #ifdef LJCOS -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/ljcos2.cpp b/src/core/nonbonded_interactions/ljcos2.cpp index 68a4674a756..46f5fdc32a5 100644 --- a/src/core/nonbonded_interactions/ljcos2.cpp +++ b/src/core/nonbonded_interactions/ljcos2.cpp @@ -25,7 +25,7 @@ #include "ljcos2.hpp" #ifdef LJCOS2 -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/ljgen.cpp b/src/core/nonbonded_interactions/ljgen.cpp index 2e9e356db9e..4b5a95da4f6 100644 --- a/src/core/nonbonded_interactions/ljgen.cpp +++ b/src/core/nonbonded_interactions/ljgen.cpp @@ -25,7 +25,7 @@ #include "ljgen.hpp" #ifdef LENNARD_JONES_GENERIC -#include "communication.hpp" +#include "interactions.hpp" #include "lj.hpp" #include diff --git a/src/core/nonbonded_interactions/morse.cpp b/src/core/nonbonded_interactions/morse.cpp index b448dee51d5..75bc0488d7a 100644 --- a/src/core/nonbonded_interactions/morse.cpp +++ b/src/core/nonbonded_interactions/morse.cpp @@ -25,7 +25,7 @@ #include "morse.hpp" #ifdef MORSE -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/nonbonded_interaction_data.cpp b/src/core/nonbonded_interactions/nonbonded_interaction_data.cpp index ff7bb89c10e..c474ec66b43 100644 --- a/src/core/nonbonded_interactions/nonbonded_interaction_data.cpp +++ b/src/core/nonbonded_interactions/nonbonded_interaction_data.cpp @@ -24,8 +24,8 @@ #include "nonbonded_interactions/nonbonded_interaction_data.hpp" #include "bonded_interactions/bonded_interaction_data.hpp" #include "collision.hpp" -#include "communication.hpp" #include "grid.hpp" +#include "interactions.hpp" #include "serialization/IA_parameters.hpp" diff --git a/src/core/nonbonded_interactions/nonbonded_tab.cpp b/src/core/nonbonded_interactions/nonbonded_tab.cpp index b7f4e479281..7a76865e0ab 100644 --- a/src/core/nonbonded_interactions/nonbonded_tab.cpp +++ b/src/core/nonbonded_interactions/nonbonded_tab.cpp @@ -25,7 +25,7 @@ #include "nonbonded_interactions/nonbonded_tab.hpp" #ifdef TABULATED -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/smooth_step.cpp b/src/core/nonbonded_interactions/smooth_step.cpp index 7b5eb7427c2..e85dab70888 100644 --- a/src/core/nonbonded_interactions/smooth_step.cpp +++ b/src/core/nonbonded_interactions/smooth_step.cpp @@ -25,7 +25,7 @@ #include "smooth_step.hpp" #ifdef SMOOTH_STEP -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/soft_sphere.cpp b/src/core/nonbonded_interactions/soft_sphere.cpp index 4bcca6999f0..f15b29e97aa 100644 --- a/src/core/nonbonded_interactions/soft_sphere.cpp +++ b/src/core/nonbonded_interactions/soft_sphere.cpp @@ -26,7 +26,7 @@ #include "soft_sphere.hpp" #ifdef SOFT_SPHERE -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/thole.cpp b/src/core/nonbonded_interactions/thole.cpp index 4424782c655..641eb68bbea 100644 --- a/src/core/nonbonded_interactions/thole.cpp +++ b/src/core/nonbonded_interactions/thole.cpp @@ -25,7 +25,7 @@ #include "thole.hpp" #ifdef THOLE -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/nonbonded_interactions/wca.cpp b/src/core/nonbonded_interactions/wca.cpp index aebcce7656a..760d97b85fd 100644 --- a/src/core/nonbonded_interactions/wca.cpp +++ b/src/core/nonbonded_interactions/wca.cpp @@ -23,7 +23,7 @@ #include "wca.hpp" #ifdef WCA -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/npt.cpp b/src/core/npt.cpp index a736bdb2d24..d6d74596a12 100644 --- a/src/core/npt.cpp +++ b/src/core/npt.cpp @@ -39,6 +39,12 @@ void mpi_bcast_nptiso_geom_worker(int, int) { boost::mpi::broadcast(comm_cart, nptiso.non_const_dim, 0); } +REGISTER_CALLBACK(mpi_bcast_nptiso_geom_worker) + +void mpi_bcast_nptiso_geom() { + mpi_call_all(mpi_bcast_nptiso_geom_worker, -1, 0); +} + void npt_ensemble_init(const BoxGeometry &box) { if (integ_switch == INTEG_METHOD_NPT_ISO) { /* prepare NpT-integration */ diff --git a/src/core/npt.hpp b/src/core/npt.hpp index 5b607949846..8cb99e0f41a 100644 --- a/src/core/npt.hpp +++ b/src/core/npt.hpp @@ -86,10 +86,13 @@ extern nptiso_struct nptiso; /** @brief Synchronizes NpT state such as instantaneous and average pressure */ void synchronize_npt_state(); -void mpi_bcast_nptiso_geom_worker(int, int); void npt_ensemble_init(const BoxGeometry &box); void integrator_npt_sanity_checks(); void npt_reset_instantaneous_virials(); void npt_add_virial_contribution(const Utils::Vector3d &force, const Utils::Vector3d &d); + +/** Broadcast nptiso geometry parameters to all nodes. */ +void mpi_bcast_nptiso_geom(); + #endif diff --git a/src/core/object-in-fluid/oif_global_forces.cpp b/src/core/object-in-fluid/oif_global_forces.cpp index ff50fda95a5..e7c7451166e 100644 --- a/src/core/object-in-fluid/oif_global_forces.cpp +++ b/src/core/object-in-fluid/oif_global_forces.cpp @@ -21,9 +21,9 @@ #include "Particle.hpp" #include "bonded_interactions/bonded_interaction_data.hpp" -#include "communication.hpp" #include "grid.hpp" #include "grid_based_algorithms/lb_interface.hpp" +#include "interactions.hpp" #include using Utils::angle_btw_triangles; diff --git a/src/core/object-in-fluid/oif_local_forces.cpp b/src/core/object-in-fluid/oif_local_forces.cpp index 3ba22065a85..042b581d61d 100644 --- a/src/core/object-in-fluid/oif_local_forces.cpp +++ b/src/core/object-in-fluid/oif_local_forces.cpp @@ -19,7 +19,7 @@ #include "oif_local_forces.hpp" -#include "communication.hpp" +#include "interactions.hpp" #include diff --git a/src/core/particle_data.cpp b/src/core/particle_data.cpp index fe176369cfb..4cda890dad8 100644 --- a/src/core/particle_data.cpp +++ b/src/core/particle_data.cpp @@ -314,10 +314,7 @@ struct UpdateVisitor : public boost::static_visitor { }; } // namespace -/** - * @brief Callback for \ref mpi_send_update_message. - */ -void mpi_update_particle_slave(int node, int id) { +void mpi_send_update_message_local(int node, int id) { if (node == comm_cart.rank()) { UpdateMessage msg{}; comm_cart.recv(0, SOME_TAG, msg); @@ -327,6 +324,8 @@ void mpi_update_particle_slave(int node, int id) { on_particle_change(); } +REGISTER_CALLBACK(mpi_send_update_message_local) + /** * @brief Send a particle update message. * @@ -348,7 +347,7 @@ void mpi_update_particle_slave(int node, int id) { void mpi_send_update_message(int id, const UpdateMessage &msg) { auto const pnode = get_particle_node(id); - mpi_call(mpi_update_particle_slave, pnode, id); + mpi_call(mpi_send_update_message_local, pnode, id); /* If the particle is remote, send the * message to the target, otherwise we @@ -393,7 +392,7 @@ void add_exclusion(Particle *part, int part2); void auto_exclusion(int distance); -void mpi_who_has_slave(int, int) { +void mpi_who_has_local(int, int) { static std::vector sendbuf; auto local_particles = cell_structure.local_particles(); @@ -412,8 +411,10 @@ void mpi_who_has_slave(int, int) { MPI_Send(sendbuf.data(), n_part, MPI_INT, 0, SOME_TAG, comm_cart); } +REGISTER_CALLBACK(mpi_who_has_local) + void mpi_who_has() { - mpi_call(mpi_who_has_slave, -1, 0); + mpi_call(mpi_who_has_local, -1, 0); auto local_particles = cell_structure.local_particles(); @@ -509,7 +510,7 @@ const Particle &get_particle_data(int part) { return *cache_ptr; } -void mpi_get_particles_slave(int, int) { +void mpi_get_particles_local(int, int) { std::vector ids; boost::mpi::scatter(comm_cart, ids, 0); @@ -522,6 +523,8 @@ void mpi_get_particles_slave(int, int) { Utils::Mpi::gatherv(comm_cart, parts.data(), parts.size(), 0); } +REGISTER_CALLBACK(mpi_get_particles_local) + /** * @brief Get multiple particles at once. * @@ -529,10 +532,10 @@ void mpi_get_particles_slave(int, int) { * * @param ids The ids of the particles that should be returned. * - * @returns The particle data. + * @returns The particle list. */ std::vector mpi_get_particles(Utils::Span ids) { - mpi_call(mpi_get_particles_slave, 0, 0); + mpi_call(mpi_get_particles_local, 0, 0); /* Return value */ std::vector parts(ids.size()); @@ -597,15 +600,99 @@ void prefetch_particle_data(Utils::Span in_ids) { } } -int place_particle(int part, const double *pos) { +/** Move a particle to a new position. If it does not exist, it is created. + * The position must be on the local node! + * + * @param id the identity of the particle to move + * @param pos its new position + * @param _new if true, the particle is allocated, else has to exists already + * + * @return Pointer to the particle. + */ +Particle *local_place_particle(int id, const Utils::Vector3d &pos, int _new) { + auto pp = Utils::Vector3d{pos[0], pos[1], pos[2]}; + auto i = Utils::Vector3i{}; + fold_position(pp, i, box_geo); + + if (_new) { + Particle new_part; + new_part.p.identity = id; + new_part.r.p = pp; + new_part.l.i = i; + + return cell_structure.add_local_particle(std::move(new_part)); + } + + auto pt = cell_structure.get_local_particle(id); + pt->r.p = pp; + pt->l.i = i; + + return pt; +} + +boost::optional mpi_place_new_particle_local(int p_id, + Utils::Vector3d const &pos) { + auto p = local_place_particle(p_id, pos, 1); + on_particle_change(); + if (p) { + return comm_cart.rank(); + } + return {}; +} + +REGISTER_CALLBACK_ONE_RANK(mpi_place_new_particle_local) + +/** Create particle at a position on a node. + * Also calls \ref on_particle_change. + * \param p_id the particle to create. + * \param pos the particles position. + */ +int mpi_place_new_particle(int p_id, const Utils::Vector3d &pos) { + return mpi_call(Communication::Result::one_rank, mpi_place_new_particle_local, + p_id, pos); +} + +void mpi_place_particle_local(int pnode, int p_id) { + if (pnode == this_node) { + Utils::Vector3d pos; + comm_cart.recv(0, SOME_TAG, pos); + local_place_particle(p_id, pos, 0); + } + + cell_structure.set_resort_particles(Cells::RESORT_GLOBAL); + on_particle_change(); +} + +REGISTER_CALLBACK(mpi_place_particle_local) + +/** Move particle to a position on a node. + * Also calls \ref on_particle_change. + * \param node the node to attach it to. + * \param p_id the particle to move. + * \param pos the particles position. + */ +void mpi_place_particle(int node, int p_id, const Utils::Vector3d &pos) { + mpi_call(mpi_place_particle_local, node, p_id); + + if (node == this_node) + local_place_particle(p_id, pos, 0); + else { + comm_cart.send(node, SOME_TAG, pos); + } + + cell_structure.set_resort_particles(Cells::RESORT_GLOBAL); + on_particle_change(); +} + +int place_particle(int p_id, const double *pos) { Utils::Vector3d p{pos[0], pos[1], pos[2]}; - if (particle_exists(part)) { - mpi_place_particle(get_particle_node(part), part, p); + if (particle_exists(p_id)) { + mpi_place_particle(get_particle_node(p_id), p_id, p); return ES_PART_OK; } - particle_node[part] = mpi_place_new_particle(part, p); + particle_node[p_id] = mpi_place_new_particle(p_id, p); return ES_PART_CREATED; } @@ -848,6 +935,25 @@ const std::vector &get_particle_bonds(int part) { return ret; } +void mpi_remove_particle_local(int, int part) { + if (part != -1) { + cell_structure.remove_particle(part); + } else { + cell_structure.remove_all_particles(); + } + on_particle_change(); +} + +REGISTER_CALLBACK(mpi_remove_particle_local) + +/** Remove a particle. + * Also calls \ref on_particle_change. + * \param p_id the particle to remove, use -1 to remove all particles. + */ +void mpi_remove_particle(int, int p_id) { + mpi_call_all(mpi_remove_particle_local, -1, p_id); +} + void remove_all_particles() { mpi_remove_particle(-1, -1); clear_particle_node(); @@ -861,37 +967,18 @@ int remove_particle(int p_id) { remove_id_from_map(p_id, type); } - auto const pnode = get_particle_node(p_id); - particle_node[p_id] = -1; - mpi_remove_particle(pnode, p_id); + mpi_remove_particle(-1, p_id); particle_node.erase(p_id); return ES_OK; } -Particle *local_place_particle(int id, const Utils::Vector3d &pos, int _new) { - auto pp = Utils::Vector3d{pos[0], pos[1], pos[2]}; - auto i = Utils::Vector3i{}; - fold_position(pp, i, box_geo); - - if (_new) { - Particle new_part; - new_part.p.identity = id; - new_part.r.p = pp; - new_part.l.i = i; - - return cell_structure.add_local_particle(std::move(new_part)); - } - - auto pt = cell_structure.get_local_particle(id); - pt->r.p = pp; - pt->l.i = i; - - return pt; -} - +/** Locally rescale all particles on current node. + * @param dir direction to scale (0/1/2 = x/y/z, 3 = x+y+z isotropically) + * @param scale factor by which to rescale (>1: stretch, <1: contract) + */ void local_rescale_particles(int dir, double scale) { for (auto &p : cell_structure.local_particles()) { if (dir < 3) @@ -902,7 +989,33 @@ void local_rescale_particles(int dir, double scale) { } } +void mpi_rescale_particles_local(int, int dir) { + double scale = 0.0; + MPI_Recv(&scale, 1, MPI_DOUBLE, 0, SOME_TAG, comm_cart, MPI_STATUS_IGNORE); + local_rescale_particles(dir, scale); + on_particle_change(); +} + +REGISTER_CALLBACK(mpi_rescale_particles_local) + +void mpi_rescale_particles(int dir, double scale) { + mpi_call(mpi_rescale_particles_local, -1, dir); + for (int pnode = 0; pnode < n_nodes; pnode++) { + if (pnode == this_node) { + local_rescale_particles(dir, scale); + } else { + MPI_Send(&scale, 1, MPI_DOUBLE, pnode, SOME_TAG, comm_cart); + } + } + on_particle_change(); +} + #ifdef EXCLUSIONS +/** Locally add an exclusion to a particle. + * @param part1 the identity of the first exclusion partner + * @param part2 the identity of the second exclusion partner + * @param _delete if true, delete the exclusion instead of add + */ void local_change_exclusion(int part1, int part2, int _delete) { if (part1 == -1 && part2 == -1) { for (auto &p : cell_structure.local_particles()) { @@ -946,6 +1059,23 @@ void add_partner(std::vector &il, int i, int j, int distance) { } } // namespace +void mpi_send_exclusion_local(int part1, int part2, int _delete) { + local_change_exclusion(part1, part2, _delete); + on_particle_change(); +} + +REGISTER_CALLBACK(mpi_send_exclusion_local) + +/** Send exclusions. + * Also calls \ref on_particle_change. + * \param part1 identity of first particle of the exclusion. + * \param part2 identity of second particle of the exclusion. + * \param _delete if true, do not add the exclusion, rather delete it if found + */ +void mpi_send_exclusion(int part1, int part2, int _delete) { + mpi_call_all(mpi_send_exclusion_local, part1, part2, _delete); +} + int change_exclusion(int part1, int part2, int _delete) { if (particle_exists(part1) && particle_exists(part2)) { mpi_send_exclusion(part1, part2, _delete); diff --git a/src/core/particle_data.hpp b/src/core/particle_data.hpp index f2301b4b5fe..3991c7e0007 100644 --- a/src/core/particle_data.hpp +++ b/src/core/particle_data.hpp @@ -347,32 +347,8 @@ int remove_particle(int part); /** Remove all particles. */ void remove_all_particles(); -/** Used by \ref mpi_place_particle, should not be used elsewhere. - * Move a particle to a new position. If it does not exist, it is created. - * The position must be on the local node! - * - * @param id the identity of the particle to move - * @param pos its new position - * @param _new if true, the particle is allocated, else has to exists already - * - * @return Pointer to the particle. - */ -Particle *local_place_particle(int id, const Utils::Vector3d &pos, int _new); - -/** Used for example by \ref mpi_send_exclusion. - * Locally add an exclusion to a particle. - * @param part1 the identity of the first exclusion partner - * @param part2 the identity of the second exclusion partner - * @param _delete if true, delete the exclusion instead of add - */ -void local_change_exclusion(int part1, int part2, int _delete); - -/** Used by \ref mpi_rescale_particles, should not be used elsewhere. - * Locally rescale all particles on current node. - * @param dir direction to scale (0/1/2 = x/y/z, 3 = x+y+z isotropically) - * @param scale factor by which to rescale (>1: stretch, <1: contract) - */ -void local_rescale_particles(int dir, double scale); +/** Rescale all particle positions in direction @p dir by a factor @p scale. */ +void mpi_rescale_particles(int dir, double scale); /** Automatically add the next \ neighbors in each molecule to the * exclusion list. diff --git a/src/core/pressure.cpp b/src/core/pressure.cpp index 47710bc3d31..371b1696bed 100644 --- a/src/core/pressure.cpp +++ b/src/core/pressure.cpp @@ -124,7 +124,11 @@ void pressure_calc() { } } -void update_pressure() { mpi_gather_stats(GatherStats::pressure); } +void update_pressure_local(int, int) { pressure_calc(); } + +REGISTER_CALLBACK(update_pressure_local) + +void update_pressure() { mpi_call_all(update_pressure_local, -1, -1); } Utils::Vector9d observable_compute_pressure_tensor() { update_pressure(); diff --git a/src/core/rattle.cpp b/src/core/rattle.cpp index 31b615b914b..f5cafb08154 100644 --- a/src/core/rattle.cpp +++ b/src/core/rattle.cpp @@ -32,6 +32,7 @@ int n_rigidbonds = 0; #include "errorhandling.hpp" #include "global.hpp" #include "grid.hpp" +#include "interactions.hpp" #include diff --git a/src/core/rotate_system.cpp b/src/core/rotate_system.cpp index e4861c3236c..080e4398263 100644 --- a/src/core/rotate_system.cpp +++ b/src/core/rotate_system.cpp @@ -69,7 +69,7 @@ void local_rotate_system(double phi, double theta, double alpha, update_dependent_particles(); } -void mpi_rotate_system_slave(int, int) { +void mpi_rotate_system_local(int, int) { std::array params; mpi::broadcast(comm_cart, params, 0); @@ -77,8 +77,10 @@ void mpi_rotate_system_slave(int, int) { cell_structure.local_particles()); } +REGISTER_CALLBACK(mpi_rotate_system_local) + void mpi_rotate_system(double phi, double theta, double alpha) { - mpi_call(mpi_rotate_system_slave, 0, 0); + mpi_call(mpi_rotate_system_local, 0, 0); std::array params{{phi, theta, alpha}}; mpi::broadcast(comm_cart, params, 0); diff --git a/src/core/tuning.cpp b/src/core/tuning.cpp index d4481bbb72a..868faa43d2c 100644 --- a/src/core/tuning.cpp +++ b/src/core/tuning.cpp @@ -22,7 +22,6 @@ * Implementation of tuning.hpp. */ #include "cells.hpp" -#include "communication.hpp" #include "errorhandling.hpp" #include "global.hpp" #include "grid.hpp" diff --git a/src/python/espressomd/cellsystem.pxd b/src/python/espressomd/cellsystem.pxd index 87e8900bec0..798f777f89a 100644 --- a/src/python/espressomd/cellsystem.pxd +++ b/src/python/espressomd/cellsystem.pxd @@ -23,10 +23,7 @@ from libcpp.pair cimport pair from .utils cimport Vector3i cdef extern from "communication.hpp": - void mpi_bcast_cell_structure(int cs) - void mpi_set_use_verlet_lists(bool use_verlet_lists) int n_nodes - vector[int] mpi_resort_particles(int global_flag) cdef extern from "cells.hpp": int CELL_STRUCTURE_DOMDEC @@ -41,6 +38,9 @@ cdef extern from "cells.hpp": const DomainDecomposition * get_domain_decomposition() vector[pair[int, int]] mpi_get_pairs(double distance) + vector[int] mpi_resort_particles(int global_flag) + void mpi_bcast_cell_structure(int cs) + void mpi_set_use_verlet_lists(bool use_verlet_lists) cdef extern from "tuning.hpp": cdef void c_tune_skin "tune_skin" (double min_skin, double max_skin, double tol, int int_steps, bool adjust_max_skin) diff --git a/src/python/espressomd/electrostatic_extensions.pxd b/src/python/espressomd/electrostatic_extensions.pxd index fc6b8c74084..cc35493a54f 100644 --- a/src/python/espressomd/electrostatic_extensions.pxd +++ b/src/python/espressomd/electrostatic_extensions.pxd @@ -61,6 +61,4 @@ IF ELECTROSTATICS and P3M: iccp3m_struct iccp3m_cfg void iccp3m_alloc_lists() - - cdef extern from "communication.hpp": int mpi_iccp3m_init() diff --git a/src/python/espressomd/electrostatics.pxd b/src/python/espressomd/electrostatics.pxd index 0dab50df9c1..ca6bd85482b 100644 --- a/src/python/espressomd/electrostatics.pxd +++ b/src/python/espressomd/electrostatics.pxd @@ -34,7 +34,7 @@ cdef extern from "EspressoSystemInterface.hpp": IF ELECTROSTATICS: - cdef extern from "communication.hpp": + cdef extern from "electrostatics_magnetostatics/common.hpp": void mpi_bcast_coulomb_params() IF P3M: diff --git a/src/python/espressomd/galilei.pxd b/src/python/espressomd/galilei.pxd index d50ccdc9f8c..a6b37d0c7da 100644 --- a/src/python/espressomd/galilei.pxd +++ b/src/python/espressomd/galilei.pxd @@ -18,7 +18,7 @@ # from .utils cimport Vector3d -cdef extern from "communication.hpp": +cdef extern from "galilei.hpp": void mpi_kill_particle_motion(int rotation) void mpi_kill_particle_forces(int torque) Vector3d mpi_system_CMS() diff --git a/src/python/espressomd/globals.pxd b/src/python/espressomd/globals.pxd index 210b08819b5..214c034b3f1 100644 --- a/src/python/espressomd/globals.pxd +++ b/src/python/espressomd/globals.pxd @@ -43,15 +43,13 @@ cdef extern from "global.hpp": void mpi_bcast_parameter(int p) -cdef extern from "communication.hpp": - void mpi_set_time_step(double time_step) except + - cdef extern from "integrate.hpp": double time_step extern int integ_switch extern double sim_time extern double verlet_reuse extern double skin + void mpi_set_time_step(double time_step) except + cdef extern from "nonbonded_interactions/nonbonded_interaction_data.hpp": extern int max_seen_particle_type diff --git a/src/python/espressomd/integrate.pxd b/src/python/espressomd/integrate.pxd index 3d08d049a7b..041d81de518 100644 --- a/src/python/espressomd/integrate.pxd +++ b/src/python/espressomd/integrate.pxd @@ -27,6 +27,7 @@ cdef extern from "config.hpp": cdef extern from "integrate.hpp" nogil: cdef int python_integrate(int n_steps, cbool recalc_forces, int reuse_forces) + cdef int mpi_steepest_descent(int max_steps) cdef void integrate_set_sd() cdef void integrate_set_nvt() cdef int integrate_set_steepest_descent(const double f_max, const double gamma, @@ -66,6 +67,3 @@ IF(STOKESIAN_DYNAMICS or STOKESIAN_DYNAMICS_GPU): cdef inline int _integrate(int nSteps, cbool recalc_forces, int reuse_forces): with nogil: return python_integrate(nSteps, recalc_forces, reuse_forces) - -cdef extern from "communication.hpp": - int mpi_steepest_descent(int max_steps) diff --git a/src/python/espressomd/magnetostatics.pxd b/src/python/espressomd/magnetostatics.pxd index 73653277019..13e9409d927 100644 --- a/src/python/espressomd/magnetostatics.pxd +++ b/src/python/espressomd/magnetostatics.pxd @@ -17,7 +17,7 @@ include "myconfig.pxi" IF DIPOLES == 1: - cdef extern from "communication.hpp": + cdef extern from "electrostatics_magnetostatics/common.hpp": void mpi_bcast_coulomb_params() cdef extern from "electrostatics_magnetostatics/dipole.hpp": diff --git a/src/shapes/src/Slitpore.cpp b/src/shapes/src/Slitpore.cpp index e5e968e2c8a..a833f750b91 100644 --- a/src/shapes/src/Slitpore.cpp +++ b/src/shapes/src/Slitpore.cpp @@ -25,11 +25,11 @@ #include -using namespace std; - namespace Shapes { void Slitpore::calculate_dist(const Utils::Vector3d &pos, double &dist, Utils::Vector3d &vec) const { + using std::sqrt; + // the left circles Utils::Vector2d c11 = {dividing_plane() - m_pore_width / 2 - m_upper_smoothing_radius,