34#include "communication.hpp"
41#include "thermostat.hpp"
47#include <boost/mpi/collectives/all_reduce.hpp>
61 auto handle = std::make_shared<System>(Private());
67 box_geo = std::make_shared<BoxGeometry>();
68 local_geo = std::make_shared<LocalBox>();
69 cell_structure = std::make_shared<CellStructure>(*box_geo);
70 propagation = std::make_shared<Propagation>();
71 thermostat = std::make_shared<Thermostat::Thermostat>();
72 nonbonded_ias = std::make_shared<InteractionsNonBonded>();
73 comfixed = std::make_shared<ComFixed>();
74 galilei = std::make_shared<Galilei>();
75 bond_breakage = std::make_shared<BondBreakage::BondBreakage>();
76 lees_edwards = std::make_shared<LeesEdwards::LeesEdwards>();
84void System::initialize() {
85 auto handle = shared_from_this();
86 cell_structure->bind_system(handle);
87 lees_edwards->bind_system(handle);
88 thermostat->bind_system(handle);
90 gpu.bind_system(handle);
93 lb.bind_system(handle);
94 ek.bind_system(handle);
109 throw std::domain_error(
"time_step must be > 0.");
110 if (lb.is_solver_set()) {
111 lb.veto_time_step(value);
113 if (ek.is_solver_set()) {
114 ek.veto_time_step(value);
117 on_timestep_change();
121 if (lb.is_solver_set()) {
124 if (ek.is_solver_set()) {
131 propagation->recalc_forces =
true;
135 min_global_cut = value;
136 on_verlet_skin_change();
141 cell_structure->set_regular_decomposition(get_interaction_range());
143 cell_structure->set_atom_decomposition();
148 std::as_const(*cell_structure).decomposition());
149 cell_structure->set_hybrid_decomposition(
150 old_hybrid_decomposition.get_cutoff_regular(),
151 old_hybrid_decomposition.get_n_square_types());
156 set_cell_structure_topology(cell_structure->decomposition_type());
161 rebuild_cell_structure();
164 if (not skip_method_adaption) {
168 coulomb.on_boxl_change();
171 dipoles.on_boxl_change();
178 if (not skip_particle_checks) {
179 auto const n_part = boost::mpi::all_reduce(
180 ::comm_cart, cell_structure->local_particles().size(), std::plus<>());
182 throw std::runtime_error(
183 "Cannot reset the box length when particles are present");
187 lb.veto_boxl_change();
188 ek.veto_boxl_change();
193 lb.on_node_grid_change();
194 ek.on_node_grid_change();
196 coulomb.on_node_grid_change();
199 dipoles.on_node_grid_change();
201 rebuild_cell_structure();
206 coulomb.on_periodicity_change();
210 dipoles.on_periodicity_change();
213#ifdef STOKESIAN_DYNAMICS
215 if (box_geo->periodic(0
u) or box_geo->periodic(1u) or box_geo->periodic(2u))
217 <<
"(False, False, False)\n";
220 on_verlet_skin_change();
225 lb.on_cell_structure_change();
226 ek.on_cell_structure_change();
228 coulomb.on_cell_structure_change();
231 dipoles.on_cell_structure_change();
238 rebuild_cell_structure();
240 coulomb.on_coulomb_change();
243 dipoles.on_dipoles_change();
245 on_short_range_ia_change();
249 lb.on_temperature_change();
250 ek.on_temperature_change();
254 lb.on_timestep_change();
255 ek.on_timestep_change();
256 on_thermostat_param_change();
260 rebuild_cell_structure();
261 propagation->recalc_forces =
true;
265 nonbonded_ias->recalc_maximal_cutoffs();
266 rebuild_cell_structure();
267 propagation->recalc_forces =
true;
272 coulomb.on_coulomb_change();
274 on_short_range_ia_change();
279 dipoles.on_dipoles_change();
281 on_short_range_ia_change();
287 propagation->recalc_forces =
true;
291 cell_structure->update_ghosts_and_resort_particle(get_global_ghost_flags());
292 propagation->recalc_forces =
true;
302 coulomb.on_particle_change();
305 dipoles.on_particle_change();
307 propagation->recalc_forces =
true;
315 coulomb.on_particle_change();
321#ifdef VIRTUAL_SITES_RELATIVE
324 cell_structure->update_ghosts_and_resort_particle(get_global_ghost_flags());
328 update_icc_particles();
340 cell_structure->update_ghosts_and_resort_particle(get_global_ghost_flags());
341 update_dependent_particles();
344 coulomb.on_observable_calc();
348 dipoles.on_observable_calc();
362 max_cut = std::max(max_cut, get_min_global_cut());
364 max_cut = std::max(max_cut, coulomb.cutoff());
367 max_cut = std::max(max_cut, dipoles.cutoff());
374 max_cut = std::max(max_cut, nonbonded_ias->maximal_cutoff());
383 coulomb.sanity_checks();
386 dipoles.sanity_checks();
388 }
catch (std::runtime_error
const &err) {
396 auto const max_cut = maximal_cutoff();
397 auto const verlet_skin = cell_structure->get_verlet_skin();
403 box_geo->set_length(box_l);
409 integrator_sanity_checks();
413 long_range_interactions_sanity_checks();
419 thermostat->recalc_prefactors(time_step);
420 reinit_thermo =
false;
421 propagation->recalc_forces =
true;
432#ifdef ADDITIONAL_CHECKS
438 auto const &actor = coulomb.impl->solver;
446 auto const &actor = dipoles.impl->solver;
454 on_observable_calc();
466 if (lb.is_solver_set())
477#ifdef COLLISION_DETECTION
CellStructureType
Cell structure topology.
@ NSQUARE
Atom decomposition (N-square).
@ HYBRID
Hybrid decomposition.
@ REGULAR
Regular decomposition.
Vector implementation and trait types for boost qvm interoperability.
double maximal_cutoff_bonded()
Calculate the maximal cutoff of bonded interactions, required to determine the cell size for communic...
Data structures for bonded interactions.
Hybrid decomposition cell system.
void init_volume_conservation(CellStructure &cs)
Initialize volume conservation.
static LocalBox make_regular_decomposition(Utils::Vector3d const &box_l, Utils::Vector3i const &node_index, Utils::Vector3i const &node_grid)
void on_periodicity_change()
double get_interaction_range() const
Get the interaction range.
double maximal_cutoff() const
Calculate the maximal cutoff of all interactions.
void on_constraint_change()
Called every time a constraint is changed.
void on_cell_structure_change()
unsigned get_global_ghost_flags() const
Returns the ghost flags required for running pair kernels for the global state, e....
void on_boxl_change(bool skip_method_adaption=false)
Called when the box length has changed.
void set_box_l(Utils::Vector3d const &box_l)
Change the box dimensions.
void on_temperature_change()
void on_thermostat_param_change()
void on_integration_start()
void set_force_cap(double value)
Set force_cap.
void on_non_bonded_ia_change()
void on_verlet_skin_change()
void update_dependent_particles()
Update particles with properties depending on other particles, namely virtual sites and ICC charges.
void on_lb_boundary_conditions_change()
Called when the LB boundary conditions change (geometry, slip velocity, or both).
void veto_boxl_change(bool skip_particle_checks=false) const
void set_time_step(double value)
Set time_step.
void on_particle_change()
Called every time a particle property changes.
bool long_range_interactions_sanity_checks() const
Check electrostatic and magnetostatic methods are properly initialized.
void on_node_grid_change()
void on_particle_charge_change()
Called every time a particle charge changes.
void on_particle_local_change()
Called every time a particle local property changes.
void on_observable_calc()
called before calculating observables, i.e.
void on_timestep_change()
void check_kT(double value) const
Veto temperature change.
static std::shared_ptr< System > create()
void set_min_global_cut(double value)
Set min_global_cut.
void set_cell_structure_topology(CellStructureType topology)
Change cell structure topology.
void on_short_range_ia_change()
void rebuild_cell_structure()
Rebuild cell lists.
Collision_parameters collision_params
Parameters for collision detection.
@ OFF
Deactivate collision detection.
double collision_detection_cutoff()
Communicator communicator
boost::mpi::communicator comm_cart
The communicator.
This file contains the defaults for ESPResSo.
This file contains the errorhandling code for severe errors, like a broken bond or illegal parameter ...
#define runtimeErrorMsg()
ICC is a method that allows to take into account the influence of arbitrarily shaped dielectric inter...
ImmersedBoundaries immersed_boundaries
@ DATA_PART_MOMENTUM
Particle::m.
@ DATA_PART_PROPERTIES
Particle::p.
@ DATA_PART_BONDS
Particle::bonds.
@ DATA_PART_POSITION
Particle::r.
Constraints< ParticleRange, Constraint > constraints
static std::shared_ptr< System > instance
void set_system(std::shared_ptr< System > new_instance)
bool all_compare(boost::mpi::communicator const &comm, T const &value)
Compare values on all nodes.
constexpr double INACTIVE_CUTOFF
Cutoff for deactivated interactions.
void npt_ensemble_init(Utils::Vector3d const &box_l, bool recalc_forces)
void integrator_npt_sanity_checks()
Exports for the NpT code.
void invalidate_fetch_cache()
Invalidate the fetch cache for get_particle_data.
void clear_particle_node()
Invalidate particle_node.
Particles creation and deletion.
void vs_relative_update_particles(CellStructure &cell_structure, BoxGeometry const &box_geo)
Utils::Vector3i calc_node_index() const
Calculate the node index in the Cartesian topology.
Utils::Vector3i node_grid
Routines to thermalize the center of mass and distance of a particle pair.