From d05263fec39b3bae983f5b27b73cfd0ea65ed593 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Thu, 26 Oct 2023 19:12:20 +0200 Subject: [PATCH] Encapsulate integrator --- src/core/integrate.cpp | 70 +++++++++---------- src/core/integrate.hpp | 39 ----------- src/core/system/System.hpp | 37 ++++++++++ src/core/tuning.cpp | 7 +- src/core/tuning.hpp | 4 +- .../EspressoSystemStandAlone_test.cpp | 4 +- src/core/unit_tests/Verlet_list_test.cpp | 6 +- src/python/espressomd/system.py | 2 +- .../integrators/BrownianDynamics.cpp | 2 +- .../integrators/BrownianDynamics.hpp | 7 +- .../integrators/Integrator.cpp | 5 +- .../integrators/Integrator.hpp | 11 ++- .../integrators/IntegratorHandle.cpp | 58 +++++++-------- .../integrators/IntegratorHandle.hpp | 15 ++-- .../integrators/SteepestDescent.cpp | 8 +-- .../integrators/SteepestDescent.hpp | 7 +- .../integrators/StokesianDynamics.cpp | 2 +- .../integrators/StokesianDynamics.hpp | 6 +- .../integrators/VelocityVerlet.cpp | 2 +- .../integrators/VelocityVerlet.hpp | 7 +- .../integrators/VelocityVerletIsoNPT.cpp | 5 +- .../integrators/VelocityVerletIsoNPT.hpp | 6 +- .../integrators/initialize.hpp | 5 +- src/script_interface/system/System.cpp | 2 + 24 files changed, 147 insertions(+), 170 deletions(-) diff --git a/src/core/integrate.cpp b/src/core/integrate.cpp index 5fa52ff2362..597f0dab6b6 100644 --- a/src/core/integrate.cpp +++ b/src/core/integrate.cpp @@ -58,7 +58,6 @@ #include "virtual_sites.hpp" #include -#include #ifdef CALIPER #include @@ -313,18 +312,17 @@ static void integrator_step_2(ParticleRange const &particles, double kT, } } -int integrate(System::System &system, int n_steps, int reuse_forces) { +namespace System { + +int System::integrate(int n_steps, int reuse_forces) { #ifdef CALIPER CALI_CXX_MARK_FUNCTION; #endif - auto &cell_structure = *system.cell_structure; - auto &box_geo = *system.box_geo; - auto &bond_breakage = *system.bond_breakage; - auto const time_step = system.get_time_step(); + auto const time_step = get_time_step(); // Prepare particle structure and run sanity checks of all active algorithms - system.on_integration_start(); + on_integration_start(); // If any method vetoes (e.g. P3M not initialized), immediately bail out if (check_runtime_errors(comm_cart)) @@ -343,13 +341,13 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { #endif // Communication step: distribute ghost positions - cells_update_ghosts(cell_structure, box_geo, global_ghost_flags()); + cells_update_ghosts(*cell_structure, *box_geo, global_ghost_flags()); - force_calc(system, time_step, temperature); + force_calc(*this, time_step, temperature); if (integ_switch != INTEG_METHOD_STEEPEST_DESCENT) { #ifdef ROTATION - convert_initial_torques(cell_structure.local_particles()); + convert_initial_torques(cell_structure->local_particles()); #endif } @@ -375,8 +373,8 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { auto lb_active = false; auto ek_active = false; if (integ_switch != INTEG_METHOD_STEEPEST_DESCENT) { - lb_active = system.lb.is_solver_set(); - ek_active = system.ek.is_ready_for_propagation(); + lb_active = lb.is_solver_set(); + ek_active = ek.is_ready_for_propagation(); } #ifdef VALGRIND @@ -392,19 +390,19 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { CALI_CXX_MARK_LOOP_ITERATION(integration_loop, step); #endif - auto particles = cell_structure.local_particles(); + auto particles = cell_structure->local_particles(); #ifdef BOND_CONSTRAINT if (n_rigidbonds) - save_old_position(particles, cell_structure.ghost_particles()); + save_old_position(particles, cell_structure->ghost_particles()); #endif - LeesEdwards::update_box_params(box_geo); + LeesEdwards::update_box_params(*box_geo); bool early_exit = integrator_step_1(particles, time_step); if (early_exit) break; - LeesEdwards::run_kernel(box_geo); + LeesEdwards::run_kernel(*box_geo); #ifdef NPT if (integ_switch != INTEG_METHOD_NPT_ISO) @@ -419,7 +417,7 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { #ifdef BOND_CONSTRAINT // Correct particle positions that participate in a rigid/constrained bond if (n_rigidbonds) { - correct_position_shake(cell_structure, box_geo); + correct_position_shake(*cell_structure, *box_geo); } #endif @@ -427,33 +425,31 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { virtual_sites()->update(); #endif - if (cell_structure.get_resort_particles() >= Cells::RESORT_LOCAL) + if (cell_structure->get_resort_particles() >= Cells::RESORT_LOCAL) n_verlet_updates++; // Communication step: distribute ghost positions - cells_update_ghosts(cell_structure, box_geo, global_ghost_flags()); + cells_update_ghosts(*cell_structure, *box_geo, global_ghost_flags()); - particles = cell_structure.local_particles(); + particles = cell_structure->local_particles(); - force_calc(system, time_step, temperature); + force_calc(*this, time_step, temperature); #ifdef VIRTUAL_SITES virtual_sites()->after_force_calc(time_step); #endif integrator_step_2(particles, temperature, time_step); - LeesEdwards::run_kernel(box_geo); + LeesEdwards::run_kernel(*box_geo); #ifdef BOND_CONSTRAINT // SHAKE velocity updates if (n_rigidbonds) { - correct_velocity_shake(cell_structure, box_geo); + correct_velocity_shake(*cell_structure, *box_geo); } #endif // propagate one-step functionalities if (integ_switch != INTEG_METHOD_STEEPEST_DESCENT) { if (lb_active and ek_active) { - auto &lb = system.lb; - auto &ek = system.ek; // assume that they are coupled, which is not necessarily true auto const md_steps_per_lb_step = calc_md_steps_per_tau(lb.get_tau()); auto const md_steps_per_ek_step = calc_md_steps_per_tau(ek.get_tau()); @@ -475,7 +471,6 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { } lb_lbcoupling_propagate(); } else if (lb_active) { - auto &lb = system.lb; auto const md_steps_per_lb_step = calc_md_steps_per_tau(lb.get_tau()); lb_skipped_md_steps += 1; if (lb_skipped_md_steps >= md_steps_per_lb_step) { @@ -484,7 +479,6 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { } lb_lbcoupling_propagate(); } else if (ek_active) { - auto &ek = system.ek; auto const md_steps_per_ek_step = calc_md_steps_per_tau(ek.get_tau()); ek_skipped_md_steps += 1; if (ek_skipped_md_steps >= md_steps_per_ek_step) { @@ -498,9 +492,9 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { #endif #ifdef COLLISION_DETECTION - handle_collisions(cell_structure); + handle_collisions(*cell_structure); #endif - bond_breakage.process_queue(system); + bond_breakage->process_queue(*this); } integrated_steps++; @@ -517,7 +511,7 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { } } // for-loop over integration steps - LeesEdwards::update_box_params(box_geo); + LeesEdwards::update_box_params(*box_geo); #ifdef CALIPER CALI_CXX_MARK_LOOP_END(integration_loop); #endif @@ -531,7 +525,7 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { #endif // Verlet list statistics - cell_structure.update_verlet_stats(n_steps, n_verlet_updates); + cell_structure->update_verlet_stats(n_steps, n_verlet_updates); #ifdef NPT if (integ_switch == INTEG_METHOD_NPT_ISO) { @@ -548,17 +542,17 @@ int integrate(System::System &system, int n_steps, int reuse_forces) { return integrated_steps; } -int integrate_with_signal_handler(System::System &system, int n_steps, - int reuse_forces, bool update_accumulators) { +int System::integrate_with_signal_handler(int n_steps, int reuse_forces, + bool update_accumulators) { assert(n_steps >= 0); // Override the signal handler so that the integrator obeys Ctrl+C SignalHandler sa(SIGINT, [](int) { ctrl_C = 1; }); /* if skin wasn't set, do an educated guess now */ - if (not system.cell_structure->is_verlet_skin_set()) { + if (not cell_structure->is_verlet_skin_set()) { try { - system.set_verlet_skin_heuristic(); + set_verlet_skin_heuristic(); } catch (...) { if (comm_cart.rank() == 0) { throw; @@ -568,7 +562,7 @@ int integrate_with_signal_handler(System::System &system, int n_steps, } if (not update_accumulators or n_steps == 0) { - return integrate(system, n_steps, reuse_forces); + return integrate(n_steps, reuse_forces); } using Accumulators::auto_update; @@ -579,7 +573,7 @@ int integrate_with_signal_handler(System::System &system, int n_steps, * end, depending on what comes first. */ auto const steps = std::min((n_steps - i), auto_update_next_update()); - auto const local_retval = integrate(system, steps, reuse_forces); + auto const local_retval = integrate(steps, reuse_forces); // make sure all ranks exit when one rank fails std::remove_const_t global_retval; @@ -599,6 +593,8 @@ int integrate_with_signal_handler(System::System &system, int n_steps, return 0; } +} // namespace System + double get_time_step() { return System::get_system().get_time_step(); } double get_sim_time() { return sim_time; } diff --git a/src/core/integrate.hpp b/src/core/integrate.hpp index bbe81d901a7..a2ebbabeae2 100644 --- a/src/core/integrate.hpp +++ b/src/core/integrate.hpp @@ -29,8 +29,6 @@ #include "config/config.hpp" -#include "system/System.hpp" - #ifdef WALBERLA #include @@ -83,43 +81,6 @@ void walberla_agrid_sanity_checks(std::string method, double agrid); #endif // WALBERLA -/** Integrate equations of motion - * @param system system to propagate - * @param n_steps Number of integration steps, can be zero - * @param reuse_forces Decide when to re-calculate forces - * - * @details This function calls two hooks for propagation kernels such as - * velocity verlet, velocity verlet + npt box changes, and steepest_descent. - * One hook is called before and one after the force calculation. - * It is up to the propagation kernels to increment the simulation time. - * - * This function propagates the system according to the choice of integrator - * stored in @ref integ_switch. The general structure is: - * - if reuse_forces is zero, recalculate the forces based on the current - * state of the system - * - Loop over the number of simulation steps: - * -# initialization (e.g., RATTLE) - * -# First hook for propagation kernels - * -# Update dependent particles and properties (RATTLE, virtual sites) - * -# Calculate forces for the current state of the system. This includes - * forces added by the Langevin thermostat and the - * Lattice-Boltzmann-particle coupling - * -# Second hook for propagation kernels - * -# Update dependent properties (Virtual sites, RATTLE) - * -# Run single step algorithms (Lattice-Boltzmann propagation, collision - * detection, NpT update) - * - Final update of dependent properties and statistics/counters - * - * High-level documentation of the integration and thermostatting schemes - * can be found in doc/sphinx/system_setup.rst and /doc/sphinx/running.rst - * - * @return number of steps that have been integrated, or a negative error code - */ -int integrate(System::System &system, int n_steps, int reuse_forces); - -int integrate_with_signal_handler(System::System &system, int n_steps, - int reuse_forces, bool update_accumulators); - /** Get time step */ double get_time_step(); diff --git a/src/core/system/System.hpp b/src/core/system/System.hpp index 42324fa677c..b69c31fc09d 100644 --- a/src/core/system/System.hpp +++ b/src/core/system/System.hpp @@ -133,6 +133,43 @@ class System { */ double particle_short_range_energy_contribution(int pid); + /** Integrate equations of motion + * @param n_steps Number of integration steps, can be zero + * @param reuse_forces Decide when to re-calculate forces + * + * @details This function calls two hooks for propagation kernels such as + * velocity verlet, velocity verlet + npt box changes, and steepest_descent. + * One hook is called before and one after the force calculation. + * It is up to the propagation kernels to increment the simulation time. + * + * This function propagates the system according to the choice of integrator + * stored in @ref integ_switch. The general structure is: + * - if reuse_forces is zero, recalculate the forces based on the current + * state of the system + * - Loop over the number of simulation steps: + * -# initialization (e.g., RATTLE) + * -# First hook for propagation kernels + * -# Update dependent particles and properties (RATTLE, virtual sites) + * -# Calculate forces for the current state of the system. This includes + * forces added by the Langevin thermostat and the + * Lattice-Boltzmann-particle coupling + * -# Second hook for propagation kernels + * -# Update dependent properties (Virtual sites, RATTLE) + * -# Run single step algorithms (Lattice-Boltzmann propagation, collision + * detection, NpT update) + * - Final update of dependent properties and statistics/counters + * + * High-level documentation of the integration and thermostatting schemes + * can be found in doc/sphinx/system_setup.rst and /doc/sphinx/running.rst + * + * @return number of steps that have been integrated, or a negative error + * code + */ + int integrate(int n_steps, int reuse_forces); + + int integrate_with_signal_handler(int n_steps, int reuse_forces, + bool update_accumulators); + /** \name Hook procedures * These procedures are called if several significant changes to * the system happen which may make a reinitialization of subsystems diff --git a/src/core/tuning.cpp b/src/core/tuning.cpp index a7a235b4234..acc79ee400b 100644 --- a/src/core/tuning.cpp +++ b/src/core/tuning.cpp @@ -66,7 +66,7 @@ static void check_statistics(Utils::Statistics::RunningAverage &acc) { } static void run_full_force_calc(System::System &system, int reuse_forces) { - auto const error_code = integrate(system, 0, reuse_forces); + auto const error_code = system.integrate(0, reuse_forces); if (error_code == INTEG_ERROR_RUNTIME) { throw TuningFailed{}; } @@ -107,15 +107,14 @@ double benchmark_integration_step(System::System &system, int int_steps) { */ static double time_calc(System::System &system, int int_steps) { auto const error_code_init = - integrate(system, 0, INTEG_REUSE_FORCES_CONDITIONALLY); + system.integrate(0, INTEG_REUSE_FORCES_CONDITIONALLY); if (error_code_init == INTEG_ERROR_RUNTIME) { return -1; } /* perform force calculation test */ auto const tick = MPI_Wtime(); - auto const error_code = - integrate(system, int_steps, INTEG_REUSE_FORCES_NEVER); + auto const error_code = system.integrate(int_steps, INTEG_REUSE_FORCES_NEVER); auto const tock = MPI_Wtime(); if (error_code == INTEG_ERROR_RUNTIME) { return -1; diff --git a/src/core/tuning.hpp b/src/core/tuning.hpp index a8f1df284ef..d56dca8bbe5 100644 --- a/src/core/tuning.hpp +++ b/src/core/tuning.hpp @@ -35,8 +35,8 @@ class TuningFailed : public std::runtime_error { /** * @brief Benchmark the integration loop. - * Call @ref integrate() several times and measure the elapsed time - * without propagating the system. It therefore doesn't include e.g. + * Call @ref System::System::integrate() several times and measure the elapsed + * time without propagating the system. It therefore doesn't include e.g. * Verlet list updates. * @param system The system to tune. * @param int_steps Number of integration steps. diff --git a/src/core/unit_tests/EspressoSystemStandAlone_test.cpp b/src/core/unit_tests/EspressoSystemStandAlone_test.cpp index a2591f11981..32f2a8042ec 100644 --- a/src/core/unit_tests/EspressoSystemStandAlone_test.cpp +++ b/src/core/unit_tests/EspressoSystemStandAlone_test.cpp @@ -347,7 +347,7 @@ BOOST_FIXTURE_TEST_CASE(espresso_system_stand_alone, ParticleFactory) { reset_particle_positions(); // recalculate forces without propagating the system - integrate(system, 0, INTEG_REUSE_FORCES_CONDITIONALLY); + system.integrate(0, INTEG_REUSE_FORCES_CONDITIONALLY); // particles are arranged in a right triangle auto const p1_opt = copy_particle_to_head_node(comm, pid1); @@ -390,7 +390,7 @@ BOOST_FIXTURE_TEST_CASE(espresso_system_stand_alone, ParticleFactory) { expected[pid] = p.pos(); } } - integrate(system, 1, INTEG_REUSE_FORCES_CONDITIONALLY); + system.integrate(1, INTEG_REUSE_FORCES_CONDITIONALLY); for (auto pid : pids) { auto const p_opt = copy_particle_to_head_node(comm, pid); if (rank == 0) { diff --git a/src/core/unit_tests/Verlet_list_test.cpp b/src/core/unit_tests/Verlet_list_test.cpp index 83bdf330d24..6d15c363a8d 100644 --- a/src/core/unit_tests/Verlet_list_test.cpp +++ b/src/core/unit_tests/Verlet_list_test.cpp @@ -205,7 +205,7 @@ BOOST_DATA_TEST_CASE_F(ParticleFactory, verlet_list_update, // integrate until both particles are closer than cutoff { - integrate(system, 11, INTEG_REUSE_FORCES_CONDITIONALLY); + system.integrate(11, INTEG_REUSE_FORCES_CONDITIONALLY); auto const p1_opt = copy_particle_to_head_node(comm, pid1); auto const p2_opt = copy_particle_to_head_node(comm, pid2); if (rank == 0) { @@ -217,7 +217,7 @@ BOOST_DATA_TEST_CASE_F(ParticleFactory, verlet_list_update, // check forces and Verlet update { - integrate(system, 1, INTEG_REUSE_FORCES_CONDITIONALLY); + system.integrate(1, INTEG_REUSE_FORCES_CONDITIONALLY); auto const p1_opt = copy_particle_to_head_node(comm, pid1); #ifdef EXTERNAL_FORCES auto const p2_opt = copy_particle_to_head_node(comm, pid2); @@ -254,7 +254,7 @@ BOOST_DATA_TEST_CASE_F(ParticleFactory, verlet_list_update, } } { - integrate(system, 0, INTEG_REUSE_FORCES_CONDITIONALLY); + system.integrate(0, INTEG_REUSE_FORCES_CONDITIONALLY); auto const p3_opt = copy_particle_to_head_node(comm, pid3); if (rank == 0) { auto const &p3 = *p3_opt; diff --git a/src/python/espressomd/system.py b/src/python/espressomd/system.py index 648d0670954..933bf5f5f2b 100644 --- a/src/python/espressomd/system.py +++ b/src/python/espressomd/system.py @@ -242,7 +242,7 @@ def _restore_object(cls, so_callback, so_callback_args, state): return so def __getstate__(self): - checkpointable_properties = ["integrator"] + checkpointable_properties = [] if has_features("VIRTUAL_SITES"): checkpointable_properties.append("_active_virtual_sites_handle") checkpointable_properties += [ diff --git a/src/script_interface/integrators/BrownianDynamics.cpp b/src/script_interface/integrators/BrownianDynamics.cpp index c61c8f1b3bd..ecc7e7699cd 100644 --- a/src/script_interface/integrators/BrownianDynamics.cpp +++ b/src/script_interface/integrators/BrownianDynamics.cpp @@ -26,7 +26,7 @@ namespace ScriptInterface { namespace Integrators { -void BrownianDynamics::activate() const { set_integ_switch(INTEG_METHOD_BD); } +void BrownianDynamics::activate() { set_integ_switch(INTEG_METHOD_BD); } } // namespace Integrators } // namespace ScriptInterface diff --git a/src/script_interface/integrators/BrownianDynamics.hpp b/src/script_interface/integrators/BrownianDynamics.hpp index f7379332328..eceb37c7b43 100644 --- a/src/script_interface/integrators/BrownianDynamics.hpp +++ b/src/script_interface/integrators/BrownianDynamics.hpp @@ -17,8 +17,7 @@ * along with this program. If not, see . */ -#ifndef ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_BROWNIAN_DYNAMICS_HPP -#define ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_BROWNIAN_DYNAMICS_HPP +#pragma once #include "Integrator.hpp" @@ -29,10 +28,8 @@ namespace ScriptInterface { namespace Integrators { class BrownianDynamics : public AutoParameters { - void activate() const override; + void activate() override; }; } // namespace Integrators } // namespace ScriptInterface - -#endif diff --git a/src/script_interface/integrators/Integrator.cpp b/src/script_interface/integrators/Integrator.cpp index 2aa263faf3b..8bc5b6d154d 100644 --- a/src/script_interface/integrators/Integrator.cpp +++ b/src/script_interface/integrators/Integrator.cpp @@ -49,9 +49,8 @@ Variant Integrator::integrate(VariantMap const ¶ms) { auto const reuse_forces = reuse_forces_flag ? INTEG_REUSE_FORCES_ALWAYS : INTEG_REUSE_FORCES_CONDITIONALLY; - auto &system = System::get_system(); - return ::integrate_with_signal_handler(system, steps, reuse_forces, - update_accumulators); + return get_system().integrate_with_signal_handler(steps, reuse_forces, + update_accumulators); } Variant Integrator::do_call_method(std::string const &name, diff --git a/src/script_interface/integrators/Integrator.hpp b/src/script_interface/integrators/Integrator.hpp index e2930901563..efd00ae173e 100644 --- a/src/script_interface/integrators/Integrator.hpp +++ b/src/script_interface/integrators/Integrator.hpp @@ -17,11 +17,11 @@ * along with this program. If not, see . */ -#ifndef ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_INTEGRATOR_HPP -#define ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_INTEGRATOR_HPP +#pragma once #include "script_interface/ObjectHandle.hpp" #include "script_interface/ScriptInterface.hpp" +#include "script_interface/system/Leaf.hpp" #include #include @@ -29,15 +29,14 @@ namespace ScriptInterface { namespace Integrators { -class Integrator : public ObjectHandle { +class Integrator : public System::Leaf { public: Variant do_call_method(std::string const &name, VariantMap const ¶ms) override; virtual Variant integrate(VariantMap const ¶ms); - virtual void activate() const = 0; + virtual void activate() = 0; + virtual void deactivate() { detach_system(); } }; } // namespace Integrators } // namespace ScriptInterface - -#endif diff --git a/src/script_interface/integrators/IntegratorHandle.cpp b/src/script_interface/integrators/IntegratorHandle.cpp index 789e3eec1da..59dacfe4830 100644 --- a/src/script_interface/integrators/IntegratorHandle.cpp +++ b/src/script_interface/integrators/IntegratorHandle.cpp @@ -38,9 +38,27 @@ namespace Integrators { IntegratorHandle::IntegratorHandle() { add_parameters({ + {"time_step", + [this](Variant const &v) { + context()->parallel_try_catch( + [&]() { get_system().set_time_step(get_value(v)); }); + }, + [this]() { return get_system().get_time_step(); }}, + {"time", [](Variant const &v) { set_time(get_value(v)); }, + []() { return get_sim_time(); }}, + {"force_cap", + [this](Variant const &v) { + get_system().set_force_cap(get_value(v)); + }, + [this]() { return get_system().get_force_cap(); }}, {"integrator", [this](Variant const &v) { + auto const old_instance = m_instance; m_instance = get_value>(v); + if (old_instance) { + old_instance->deactivate(); + } + m_instance->bind_system(m_system.lock()); m_instance->activate(); }, [this]() { @@ -68,40 +86,22 @@ IntegratorHandle::IntegratorHandle() { } } }}, - {"time_step", - [this](Variant const &v) { - context()->parallel_try_catch([&]() { - System::get_system().set_time_step(get_value(v)); - }); - }, - []() { return get_time_step(); }}, - {"time", [](Variant const &v) { set_time(get_value(v)); }, - []() { return get_sim_time(); }}, - {"force_cap", - [](Variant const &v) { - System::get_system().set_force_cap(get_value(v)); - }, - []() { return System::get_system().get_force_cap(); }}, }); } -static bool checkpoint_filter(typename VariantMap::value_type const &kv) { - /* When loading from a checkpoint file, defer the integrator object to last, - * and skip the time_step if it is -1 to avoid triggering sanity checks. - */ - return kv.first == "integrator" or - (kv.first == "time_step" and ::integ_switch == INTEG_METHOD_NVT and - get_time_step() == -1. and is_type(kv.second) and - get_value(kv.second) == -1.); -} - -void IntegratorHandle::do_construct(VariantMap const ¶ms) { - for (auto const &kv : params) { - if (not checkpoint_filter(kv)) { - do_set_parameter(kv.first, kv.second); +void IntegratorHandle::on_bind_system(::System::System &system) { + auto const ¶ms = *m_params; + for (auto const &key : get_parameter_insertion_order()) { + if (params.count(key) != 0ul) { + if (not(key == "time_step" and ::integ_switch == INTEG_METHOD_NVT and + system.get_time_step() == -1. and + is_type(params.at(key)) and + get_value(is_type(params.at(key))) == -1.)) { + do_set_parameter(key, params.at(key)); + } } } - do_set_parameter("integrator", params.at("integrator")); + m_params.reset(); } } // namespace Integrators diff --git a/src/script_interface/integrators/IntegratorHandle.hpp b/src/script_interface/integrators/IntegratorHandle.hpp index 732a49dc3b3..c21dc2b2289 100644 --- a/src/script_interface/integrators/IntegratorHandle.hpp +++ b/src/script_interface/integrators/IntegratorHandle.hpp @@ -17,11 +17,11 @@ * along with this program. If not, see . */ -#ifndef ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_INTEGRATOR_HANDLE_HPP -#define ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_INTEGRATOR_HANDLE_HPP +#pragma once #include "script_interface/ScriptInterface.hpp" #include "script_interface/auto_parameters/AutoParameters.hpp" +#include "script_interface/system/Leaf.hpp" #include "Integrator.hpp" @@ -31,16 +31,19 @@ namespace ScriptInterface { namespace Integrators { -class IntegratorHandle : public AutoParameters { +class IntegratorHandle : public AutoParameters { std::shared_ptr m_instance; + std::shared_ptr m_params; + + void on_bind_system(::System::System &system) override; public: IntegratorHandle(); - void do_construct(VariantMap const ¶ms) override; + void do_construct(VariantMap const ¶ms) override { + m_params = std::make_shared(params); + } }; } // namespace Integrators } // namespace ScriptInterface - -#endif diff --git a/src/script_interface/integrators/SteepestDescent.cpp b/src/script_interface/integrators/SteepestDescent.cpp index 90d95865677..c3a4b1bb02a 100644 --- a/src/script_interface/integrators/SteepestDescent.cpp +++ b/src/script_interface/integrators/SteepestDescent.cpp @@ -23,7 +23,6 @@ #include "core/integrate.hpp" #include "core/integrators/steepest_descent.hpp" -#include "core/system/System.hpp" #include @@ -42,9 +41,8 @@ Variant SteepestDescent::integrate(VariantMap const ¶ms) { throw std::domain_error("Parameter 'steps' must be positive"); } }); - auto &system = System::get_system(); - return ::integrate_with_signal_handler(system, steps, reuse_forces, - update_accumulators); + return get_system().integrate_with_signal_handler(steps, reuse_forces, + update_accumulators); } SteepestDescent::SteepestDescent() { @@ -69,7 +67,7 @@ void SteepestDescent::do_construct(VariantMap const ¶ms) { }); } -void SteepestDescent::activate() const { +void SteepestDescent::activate() { register_integrator(get_instance()); set_integ_switch(INTEG_METHOD_STEEPEST_DESCENT); } diff --git a/src/script_interface/integrators/SteepestDescent.hpp b/src/script_interface/integrators/SteepestDescent.hpp index 73d80da6dfc..ac78eba82b2 100644 --- a/src/script_interface/integrators/SteepestDescent.hpp +++ b/src/script_interface/integrators/SteepestDescent.hpp @@ -17,8 +17,7 @@ * along with this program. If not, see . */ -#ifndef ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_STEEPEST_DESCENT_HPP -#define ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_STEEPEST_DESCENT_HPP +#pragma once #include "Integrator.hpp" @@ -41,7 +40,7 @@ class SteepestDescent : public AutoParameters { void do_construct(VariantMap const ¶ms) override; Variant integrate(VariantMap const ¶ms) override; - void activate() const override; + void activate() override; ::SteepestDescentParameters const &get_instance() const { return *m_instance; @@ -50,5 +49,3 @@ class SteepestDescent : public AutoParameters { } // namespace Integrators } // namespace ScriptInterface - -#endif diff --git a/src/script_interface/integrators/StokesianDynamics.cpp b/src/script_interface/integrators/StokesianDynamics.cpp index c7dfa3569d0..54e7ab8c8f5 100644 --- a/src/script_interface/integrators/StokesianDynamics.cpp +++ b/src/script_interface/integrators/StokesianDynamics.cpp @@ -89,7 +89,7 @@ void StokesianDynamics::do_construct(VariantMap const ¶ms) { }); } -void StokesianDynamics::activate() const { +void StokesianDynamics::activate() { context()->parallel_try_catch([&]() { register_integrator(get_instance()); set_integ_switch(INTEG_METHOD_SD); diff --git a/src/script_interface/integrators/StokesianDynamics.hpp b/src/script_interface/integrators/StokesianDynamics.hpp index 27f20281ee1..7ce13d4a2f9 100644 --- a/src/script_interface/integrators/StokesianDynamics.hpp +++ b/src/script_interface/integrators/StokesianDynamics.hpp @@ -17,8 +17,7 @@ * along with this program. If not, see . */ -#ifndef ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_STOKESIAN_DYNAMICS_HPP -#define ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_STOKESIAN_DYNAMICS_HPP +#pragma once #include "config/config.hpp" @@ -44,7 +43,7 @@ class StokesianDynamics : public AutoParameters { StokesianDynamics(); void do_construct(VariantMap const ¶ms) override; - void activate() const override; + void activate() override; ::StokesianDynamicsParameters const &get_instance() const { return *m_instance; @@ -55,4 +54,3 @@ class StokesianDynamics : public AutoParameters { } // namespace ScriptInterface #endif // STOKESIAN_DYNAMICS -#endif diff --git a/src/script_interface/integrators/VelocityVerlet.cpp b/src/script_interface/integrators/VelocityVerlet.cpp index f32d0b11bd8..4d7c7e32ea7 100644 --- a/src/script_interface/integrators/VelocityVerlet.cpp +++ b/src/script_interface/integrators/VelocityVerlet.cpp @@ -26,7 +26,7 @@ namespace ScriptInterface { namespace Integrators { -void VelocityVerlet::activate() const { set_integ_switch(INTEG_METHOD_NVT); } +void VelocityVerlet::activate() { set_integ_switch(INTEG_METHOD_NVT); } } // namespace Integrators } // namespace ScriptInterface diff --git a/src/script_interface/integrators/VelocityVerlet.hpp b/src/script_interface/integrators/VelocityVerlet.hpp index 693b9039bdb..fcd9511e33c 100644 --- a/src/script_interface/integrators/VelocityVerlet.hpp +++ b/src/script_interface/integrators/VelocityVerlet.hpp @@ -17,8 +17,7 @@ * along with this program. If not, see . */ -#ifndef ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_VELOCITY_VERLET_HPP -#define ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_VELOCITY_VERLET_HPP +#pragma once #include "Integrator.hpp" @@ -29,10 +28,8 @@ namespace ScriptInterface { namespace Integrators { class VelocityVerlet : public AutoParameters { - void activate() const override; + void activate() override; }; } // namespace Integrators } // namespace ScriptInterface - -#endif diff --git a/src/script_interface/integrators/VelocityVerletIsoNPT.cpp b/src/script_interface/integrators/VelocityVerletIsoNPT.cpp index 98842c805f5..8464f75247c 100644 --- a/src/script_interface/integrators/VelocityVerletIsoNPT.cpp +++ b/src/script_interface/integrators/VelocityVerletIsoNPT.cpp @@ -27,7 +27,6 @@ #include "core/integrate.hpp" #include "core/npt.hpp" -#include "core/system/System.hpp" #include @@ -63,10 +62,10 @@ void VelocityVerletIsoNPT::do_construct(VariantMap const ¶ms) { }); } -void VelocityVerletIsoNPT::activate() const { +void VelocityVerletIsoNPT::activate() { ::nptiso = get_instance(); set_integ_switch(INTEG_METHOD_NPT_ISO); - System::get_system().on_thermostat_param_change(); + get_system().on_thermostat_param_change(); } } // namespace Integrators diff --git a/src/script_interface/integrators/VelocityVerletIsoNPT.hpp b/src/script_interface/integrators/VelocityVerletIsoNPT.hpp index 414b16ddb9f..f58aa21b201 100644 --- a/src/script_interface/integrators/VelocityVerletIsoNPT.hpp +++ b/src/script_interface/integrators/VelocityVerletIsoNPT.hpp @@ -17,8 +17,7 @@ * along with this program. If not, see . */ -#ifndef ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_NPT_HPP -#define ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_NPT_HPP +#pragma once #include "config/config.hpp" @@ -45,7 +44,7 @@ class VelocityVerletIsoNPT VelocityVerletIsoNPT(); void do_construct(VariantMap const ¶ms) override; - void activate() const override; + void activate() override; ::NptIsoParameters const &get_instance() const { return *m_instance; } }; @@ -54,4 +53,3 @@ class VelocityVerletIsoNPT } // namespace ScriptInterface #endif // NPT -#endif diff --git a/src/script_interface/integrators/initialize.hpp b/src/script_interface/integrators/initialize.hpp index 706621966a2..233e970b308 100644 --- a/src/script_interface/integrators/initialize.hpp +++ b/src/script_interface/integrators/initialize.hpp @@ -17,8 +17,7 @@ * along with this program. If not, see . */ -#ifndef ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_INITIALIZE_HPP -#define ESPRESSO_SRC_SCRIPT_INTERFACE_INTEGRATORS_INITIALIZE_HPP +#pragma once #include @@ -31,5 +30,3 @@ void initialize(Utils::Factory *om); } // namespace Integrators } // namespace ScriptInterface - -#endif diff --git a/src/script_interface/system/System.cpp b/src/script_interface/system/System.cpp index cd933ff8294..5d4d4e8ace2 100644 --- a/src/script_interface/system/System.cpp +++ b/src/script_interface/system/System.cpp @@ -66,6 +66,7 @@ static bool system_created = false; struct System::Leaves { Leaves() = default; std::shared_ptr cell_system; + std::shared_ptr integrator; std::shared_ptr analysis; std::shared_ptr comfixed; std::shared_ptr galilei; @@ -129,6 +130,7 @@ System::System() : m_instance{}, m_leaves{std::make_shared()} { {"max_oif_objects", ::max_oif_objects}, }); add_parameter("cell_system", &Leaves::cell_system); + add_parameter("integrator", &Leaves::integrator); add_parameter("analysis", &Leaves::analysis); add_parameter("comfixed", &Leaves::comfixed); add_parameter("galilei", &Leaves::galilei);