From 4e6cc22f43632663771a56eee3d0e7a1e455487a Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Fri, 21 Apr 2023 23:57:10 +0200 Subject: [PATCH] Rewrite integrators script interface as Python code The nogil attribute is required to make the integrator run smoothly with threads (e.g. from the visualizer). Functions marked with the nogil attribute can throw C++ exceptions since Cython 0.19. The Cython functions that handle signals can be replaced with the `signal.raise_signal()` function introduced in Python 3.8. --- src/core/EspressoSystemStandAlone.cpp | 2 +- src/core/integrate.cpp | 110 ++++++++---------- src/core/integrate.hpp | 42 ++----- src/core/tuning.cpp | 24 ++-- src/python/espressomd/integrate.pxd | 33 ------ .../{integrate.pyx => integrate.py} | 41 +++---- src/python/espressomd/lb.pxd | 3 + src/python/espressomd/lb.pyx | 2 +- src/python/espressomd/script_interface.pxd | 4 + src/python/espressomd/script_interface.pyx | 19 ++- src/python/espressomd/utils.pyx | 3 +- src/script_interface/ObjectHandle.hpp | 2 +- .../cell_system/CellSystem.cpp | 2 +- .../integrators/CMakeLists.txt | 1 + .../integrators/Integrator.cpp | 62 ++++++++++ .../integrators/Integrator.hpp | 3 + .../integrators/SteepestDescent.cpp | 13 +++ .../integrators/SteepestDescent.hpp | 1 + .../integrators/VelocityVerlet.hpp | 2 +- testsuite/python/integrator_exceptions.py | 11 +- .../python/integrator_steepest_descent.py | 4 + testsuite/python/interactions_bonded.py | 2 +- testsuite/python/scafacos_interface.py | 1 + testsuite/python/sigint.py | 2 +- testsuite/python/test_checkpoint.py | 3 + 25 files changed, 209 insertions(+), 183 deletions(-) delete mode 100644 src/python/espressomd/integrate.pxd rename src/python/espressomd/{integrate.pyx => integrate.py} (87%) create mode 100644 src/script_interface/integrators/Integrator.cpp diff --git a/src/core/EspressoSystemStandAlone.cpp b/src/core/EspressoSystemStandAlone.cpp index 3dee8e90b87..024fd403b33 100644 --- a/src/core/EspressoSystemStandAlone.cpp +++ b/src/core/EspressoSystemStandAlone.cpp @@ -62,5 +62,5 @@ void EspressoSystemStandAlone::set_time_step(double time_step) const { } void EspressoSystemStandAlone::set_skin(double new_skin) const { - mpi_set_skin_local(new_skin); + ::set_skin(new_skin); } diff --git a/src/core/integrate.cpp b/src/core/integrate.cpp index 575e30f53d1..6a5503ab1f9 100644 --- a/src/core/integrate.cpp +++ b/src/core/integrate.cpp @@ -93,14 +93,8 @@ static double verlet_reuse = 0.0; static int fluid_step = 0; -bool set_py_interrupt = false; namespace { volatile std::sig_atomic_t ctrl_C = 0; - -void notify_sig_int() { - ctrl_C = 0; // reset - set_py_interrupt = true; // global to notify Python -} } // namespace namespace LeesEdwards { @@ -259,7 +253,7 @@ int integrate(int n_steps, int reuse_forces) { // If any method vetoes (e.g. P3M not initialized), immediately bail out if (check_runtime_errors(comm_cart)) - return 0; + return INTEG_ERROR_RUNTIME; // Additional preparations for the first integration step if (reuse_forces == -1 || (recalc_forces && reuse_forces != 1)) { @@ -287,11 +281,17 @@ int integrate(int n_steps, int reuse_forces) { lb_lbcoupling_activate(); if (check_runtime_errors(comm_cart)) - return 0; + return INTEG_ERROR_RUNTIME; // Keep track of the number of Verlet updates (i.e. particle resorts) int n_verlet_updates = 0; + // Keep track of whether an interrupt signal was caught (only in singleton + // mode, since signal handlers are unreliable with more than 1 MPI rank) + auto const singleton_mode = comm_cart.size() == 1; + auto caught_sigint = false; + auto caught_error = false; + #ifdef VALGRIND_MARKERS CALLGRIND_START_INSTRUMENTATION; #endif @@ -384,12 +384,14 @@ int integrate(int n_steps, int reuse_forces) { integrated_steps++; - if (check_runtime_errors(comm_cart)) + if (check_runtime_errors(comm_cart)) { + caught_error = true; break; + } // Check if SIGINT has been caught. - if (ctrl_C == 1) { - notify_sig_int(); + if (singleton_mode and ctrl_C == 1) { + caught_sigint = true; break; } @@ -416,40 +418,50 @@ int integrate(int n_steps, int reuse_forces) { synchronize_npt_state(); } #endif + if (caught_sigint) { + ctrl_C = 0; + return INTEG_ERROR_SIGINT; + } + if (caught_error) { + return INTEG_ERROR_RUNTIME; + } return integrated_steps; } -int python_integrate(int n_steps, bool recalc_forces_par, - bool reuse_forces_par) { +int 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; }); - int reuse_forces = reuse_forces_par; - - if (recalc_forces_par) { - if (reuse_forces) { - runtimeErrorMsg() << "cannot reuse old forces and recalculate forces"; - } - reuse_forces = -1; + if (not update_accumulators or n_steps == 0) { + return integrate(n_steps, reuse_forces); } + auto const is_head_node = comm_cart.rank() == 0; + /* if skin wasn't set, do an educated guess now */ if (!skin_set) { auto const max_cut = maximal_cutoff(n_nodes); if (max_cut <= 0.0) { - runtimeErrorMsg() - << "cannot automatically determine skin, please set it manually"; - return ES_ERROR; + if (is_head_node) { + throw std::runtime_error( + "cannot automatically determine skin, please set it manually"); + } + return INTEG_ERROR_RUNTIME; } /* maximal skin that can be used without resorting is the maximal * range of the cell system minus what is needed for interactions. */ - auto const new_skin = - std::min(0.4 * max_cut, - *boost::min_element(cell_structure.max_cutoff()) - max_cut); - mpi_call_all(mpi_set_skin_local, new_skin); + auto const max_range = *boost::min_element(::cell_structure.max_cutoff()); + auto const new_skin = std::min(0.4 * max_cut, max_range - max_cut); + ::set_skin(new_skin); + } + + // re-acquire MpiCallbacks listener on worker nodes + if (not is_head_node) { + return 0; } using Accumulators::auto_update; @@ -459,8 +471,11 @@ int python_integrate(int n_steps, bool recalc_forces_par, /* Integrate to either the next accumulator update, or the * end, depending on what comes first. */ auto const steps = std::min((n_steps - i), auto_update_next_update()); - if (mpi_integrate(steps, reuse_forces)) - return ES_ERROR; + auto const retval = mpi_call(Communication::Result::main_rank, integrate, + steps, reuse_forces); + if (retval < 0) { + return retval; // propagate error code + } reuse_forces = 1; @@ -469,37 +484,10 @@ int python_integrate(int n_steps, bool recalc_forces_par, i += steps; } - if (n_steps == 0) { - if (mpi_integrate(0, reuse_forces)) - return ES_ERROR; - } - - return ES_OK; + return 0; } -static int mpi_steepest_descent_local(int steps) { - return integrate(steps, -1); -} - -REGISTER_CALLBACK_MAIN_RANK(mpi_steepest_descent_local) - -int mpi_steepest_descent(int steps) { - return mpi_call(Communication::Result::main_rank, mpi_steepest_descent_local, - steps); -} - -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); -} +REGISTER_CALLBACK_MAIN_RANK(integrate) double interaction_range() { /* Consider skin only if there are actually interactions */ @@ -524,14 +512,12 @@ void set_time_step(double value) { on_timestep_change(); } -void mpi_set_skin_local(double value) { +void set_skin(double value) { ::skin = value; - skin_set = true; + ::skin_set = true; on_skin_change(); } -REGISTER_CALLBACK(mpi_set_skin_local) - void set_time(double value) { ::sim_time = value; ::recalc_forces = true; diff --git a/src/core/integrate.hpp b/src/core/integrate.hpp index adfc53f0347..3c41e2f39b0 100644 --- a/src/core/integrate.hpp +++ b/src/core/integrate.hpp @@ -36,6 +36,12 @@ #define INTEG_METHOD_SD 7 /**@}*/ +/** \name Integrator error codes */ +/**@{*/ +#define INTEG_ERROR_RUNTIME -1 +#define INTEG_ERROR_SIGINT -2 +/**@}*/ + /** Switch determining which integrator to use. */ extern int integ_switch; @@ -45,9 +51,6 @@ extern double skin; /** If true, the forces will be recalculated before the next integration. */ extern bool recalc_forces; -/** Communicate signal handling to the Python interpreter */ -extern bool set_py_interrupt; - double interaction_range(); /** Check integrator parameters and incompatibilities between the integrator @@ -90,37 +93,12 @@ void integrator_sanity_checks(); * 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 + * @return number of steps that have been integrated, or a negative error code */ int integrate(int n_steps, int reuse_forces); -/** @brief Run the integration loop. Can be interrupted with Ctrl+C. - * - * @param n_steps Number of integration steps, can be zero - * @param recalc_forces Whether to recalculate forces - * @param reuse_forces Whether to re-use forces - * @retval ES_OK on success - * @retval ES_ERROR on error - */ -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); +int integrate_with_signal_handler(int n_steps, int reuse_forces, + bool update_accumulators); /** Get @c verlet_reuse */ double get_verlet_reuse(); @@ -138,7 +116,7 @@ void increment_sim_time(double amount); void set_time_step(double value); /** @brief Set new skin. */ -void mpi_set_skin_local(double value); +void set_skin(double value); /** @brief Set the simulation time. */ void set_time(double value); diff --git a/src/core/tuning.cpp b/src/core/tuning.cpp index bf9e198b681..255ef2d15c3 100644 --- a/src/core/tuning.cpp +++ b/src/core/tuning.cpp @@ -66,26 +66,22 @@ static void check_statistics(Utils::Statistics::RunningAverage &acc) { } } -static void throw_on_error() { - auto const n_errors = check_runtime_errors_local(); - if (boost::mpi::all_reduce(::comm_cart, n_errors, std::plus<>()) != 0) { - throw TuningFailed{}; - } -} - double benchmark_integration_step(int int_steps) { Utils::Statistics::RunningAverage running_average; - integrate(0, 0); - throw_on_error(); + if (integrate(0, 0) == INTEG_ERROR_RUNTIME) { + throw TuningFailed{}; + } /* perform force calculation test */ for (int i = 0; i < int_steps; i++) { auto const tick = MPI_Wtime(); - integrate(0, -1); + auto const error_code = integrate(0, -1); auto const tock = MPI_Wtime(); + if (error_code == INTEG_ERROR_RUNTIME) { + throw TuningFailed{}; + } running_average.add_sample((tock - tick)); - throw_on_error(); } if (this_node == 0) { @@ -147,10 +143,10 @@ void tune_skin(double min_skin, double max_skin, double tol, int int_steps, b = max_permissible_skin; while (fabs(a - b) > tol) { - mpi_set_skin_local(a); + ::set_skin(a); auto const time_a = time_calc(int_steps); - mpi_set_skin_local(b); + ::set_skin(b); auto const time_b = time_calc(int_steps); if (time_a > time_b) { @@ -160,5 +156,5 @@ void tune_skin(double min_skin, double max_skin, double tol, int int_steps, } } auto const new_skin = 0.5 * (a + b); - mpi_set_skin_local(new_skin); + ::set_skin(new_skin); } diff --git a/src/python/espressomd/integrate.pxd b/src/python/espressomd/integrate.pxd deleted file mode 100644 index 91c372f321d..00000000000 --- a/src/python/espressomd/integrate.pxd +++ /dev/null @@ -1,33 +0,0 @@ -# -# Copyright (C) 2013-2022 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 . -# -from libcpp cimport bool as cbool - -include "myconfig.pxi" - -cdef extern from "integrate.hpp": - double get_time_step() - -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 extern cbool set_py_interrupt - -cdef inline int _integrate(int nSteps, cbool recalc_forces, int reuse_forces): - with nogil: - return python_integrate(nSteps, recalc_forces, reuse_forces) diff --git a/src/python/espressomd/integrate.pyx b/src/python/espressomd/integrate.py similarity index 87% rename from src/python/espressomd/integrate.pyx rename to src/python/espressomd/integrate.py index 5c5b6e1ee7c..519dff9f702 100644 --- a/src/python/espressomd/integrate.pyx +++ b/src/python/espressomd/integrate.py @@ -16,12 +16,9 @@ # You should have received a copy of the GNU General Public License # along with this program. If not, see . # -from cpython.exc cimport PyErr_CheckSignals, PyErr_SetInterrupt -include "myconfig.pxi" -from . import utils -from . cimport integrate from .script_interface import ScriptInterfaceHelper, script_interface_register from .code_features import assert_features +import signal @script_interface_register @@ -115,22 +112,17 @@ def run(self, steps=1, recalc_forces=False, reuse_forces=False): Reuse the forces from previous time step. """ - utils.check_type_or_throw_except(steps, 1, int, "steps must be an int") - utils.check_type_or_throw_except( - recalc_forces, 1, bool, "recalc_forces has to be a bool") - utils.check_type_or_throw_except( - reuse_forces, 1, bool, "reuse_forces has to be a bool") - if steps < 0: - raise ValueError("steps must be positive") + status = self.call_method( + "integrate", + with_nogil=True, + steps=steps, + recalc_forces=recalc_forces, + reuse_forces=reuse_forces) + self.handle_sigint(status) - _integrate(steps, recalc_forces, reuse_forces) - - if integrate.set_py_interrupt: - PyErr_SetInterrupt() - integrate.set_py_interrupt = False - PyErr_CheckSignals() - - utils.handle_errors("Encountered errors during integrate") + def handle_sigint(self, error_code): + if error_code == -2: + signal.raise_signal(signal.Signals.SIGINT) @script_interface_register @@ -161,7 +153,7 @@ class SteepestDescent(Integrator): _so_name = "Integrators::SteepestDescent" _so_creation_policy = "GLOBAL" - def run(self, steps=1, **kwargs): + def run(self, steps=1): """ Run the steepest descent. @@ -176,12 +168,9 @@ def run(self, steps=1, **kwargs): Number of integrated steps. """ - utils.check_type_or_throw_except(steps, 1, int, "steps must be an int") - assert steps >= 0, "steps has to be positive" - - integrated = mpi_steepest_descent(steps) - - utils.handle_errors("Encountered errors during integrate") + integrated = self.call_method( + "integrate", with_nogil=True, steps=steps) + self.handle_sigint(integrated) return integrated diff --git a/src/python/espressomd/lb.pxd b/src/python/espressomd/lb.pxd index 98ebec7a146..3ac79be7c11 100644 --- a/src/python/espressomd/lb.pxd +++ b/src/python/espressomd/lb.pxd @@ -119,6 +119,9 @@ cdef extern from "grid_based_algorithms/lb_interpolation.hpp" namespace "Interpo cdef InterpolationOrder linear cdef InterpolationOrder quadratic +cdef extern from "integrate.hpp": + double get_time_step() + ############################################## # # Wrapper-functions to handle unit conversions diff --git a/src/python/espressomd/lb.pyx b/src/python/espressomd/lb.pyx index dceef548081..e2875002e33 100644 --- a/src/python/espressomd/lb.pyx +++ b/src/python/espressomd/lb.pyx @@ -28,7 +28,7 @@ from . import highlander from . import utils from . cimport utils from .utils cimport Vector3i, Vector3d, Vector6d, Vector19d -from .integrate cimport get_time_step +from .lb cimport get_time_step cdef class FluidActor: diff --git a/src/python/espressomd/script_interface.pxd b/src/python/espressomd/script_interface.pxd index 112a05de941..fc0104fff42 100644 --- a/src/python/espressomd/script_interface.pxd +++ b/src/python/espressomd/script_interface.pxd @@ -45,6 +45,10 @@ cdef extern from "script_interface/ScriptInterface.hpp" namespace "ScriptInterfa Variant get_parameter(const string & name) except + void set_parameter(const string & name, const Variant & value) except + Variant call_method(const string & name, const VariantMap & parameters) except + + # a function marked "nogil" can also be marked "except +": in the catch + # statement, the GIL will be automatically re-acquired by Cython + # https://github.com/cython/cython/commit/6795a2ba4761110c64c7a36401cb676614ba9f79 + Variant call_method_nogil "call_method"(const string & name, const VariantMap & parameters) nogil except + string_ref name() cdef extern from "script_interface/ContextManager.hpp" namespace "ScriptInterface::ContextManager": diff --git a/src/python/espressomd/script_interface.pyx b/src/python/espressomd/script_interface.pyx index b669b6e8abf..940c1f3e2cb 100644 --- a/src/python/espressomd/script_interface.pyx +++ b/src/python/espressomd/script_interface.pyx @@ -136,19 +136,23 @@ cdef class PScriptInterface: self.sip = sip - def call_method(self, method, handle_errors_message=None, **kwargs): + def call_method(self, method, handle_errors_message=None, + with_nogil=False, **kwargs): """ Call a method of the core class. Parameters ---------- - method : Creation policy. + method : :obj:`str` Name of the core method. handle_errors_message : :obj:`str`, optional Custom error message for runtime errors raised in a MPI context. + with_nogil : :obj:`bool`, optional + Run the core method without the GIL if ``True``. **kwargs Arguments for the method. """ + cdef ObjectHandle * handle = self.sip.get() cdef VariantMap parameters cdef Variant value @@ -156,7 +160,16 @@ cdef class PScriptInterface: parameters[utils.to_char_pointer(name)] = python_object_to_variant( kwargs[name]) - value = self.sip.get().call_method(utils.to_char_pointer(method), parameters) + # the internal buffer of a cython bytestring object can be accessed as + # a raw char pointer, but then the bytestring object must be kept alive + method_name_bytes_counted_reference = utils.to_char_pointer(method) + cdef char * method_name_char = method_name_bytes_counted_reference + + if with_nogil: + with nogil: + value = handle.call_method_nogil(method_name_char, parameters) + else: + value = handle.call_method(method_name_char, parameters) res = variant_to_python_object(value) if handle_errors_message is None: handle_errors_message = f"while calling method {method}()" diff --git a/src/python/espressomd/utils.pyx b/src/python/espressomd/utils.pyx index 987086c0baa..f1da1760064 100644 --- a/src/python/espressomd/utils.pyx +++ b/src/python/espressomd/utils.pyx @@ -70,7 +70,8 @@ cpdef check_type_or_throw_except(x, n, t, msg): def to_char_pointer(s): """ - Returns a char pointer which contains the information of the provided python string. + Returns a Cython bytes object which contains the information of the provided + Python string. Cython bytes objects implicitly cast to raw char pointers. Parameters ---------- diff --git a/src/script_interface/ObjectHandle.hpp b/src/script_interface/ObjectHandle.hpp index b28742d687c..a4325ee9129 100644 --- a/src/script_interface/ObjectHandle.hpp +++ b/src/script_interface/ObjectHandle.hpp @@ -133,7 +133,7 @@ class ObjectHandle { protected: /** - * @brief Local implementation of @c do_call_method. + * @brief Local implementation of @c call_method. * * If not overridden by the implementation, this does nothing. */ diff --git a/src/script_interface/cell_system/CellSystem.cpp b/src/script_interface/cell_system/CellSystem.cpp index 021b7479cb3..db39a9ad81b 100644 --- a/src/script_interface/cell_system/CellSystem.cpp +++ b/src/script_interface/cell_system/CellSystem.cpp @@ -96,7 +96,7 @@ CellSystem::CellSystem() { } throw Exception(""); } - mpi_set_skin_local(new_skin); + ::set_skin(new_skin); }, []() { return ::skin; }}, {"decomposition_type", AutoParameter::read_only, diff --git a/src/script_interface/integrators/CMakeLists.txt b/src/script_interface/integrators/CMakeLists.txt index a7d390ca37f..fe3c5359853 100644 --- a/src/script_interface/integrators/CMakeLists.txt +++ b/src/script_interface/integrators/CMakeLists.txt @@ -21,6 +21,7 @@ target_sources( espresso_script_interface PRIVATE ${CMAKE_CURRENT_SOURCE_DIR}/initialize.cpp ${CMAKE_CURRENT_SOURCE_DIR}/BrownianDynamics.cpp + ${CMAKE_CURRENT_SOURCE_DIR}/Integrator.cpp ${CMAKE_CURRENT_SOURCE_DIR}/IntegratorHandle.cpp ${CMAKE_CURRENT_SOURCE_DIR}/SteepestDescent.cpp ${CMAKE_CURRENT_SOURCE_DIR}/StokesianDynamics.cpp diff --git a/src/script_interface/integrators/Integrator.cpp b/src/script_interface/integrators/Integrator.cpp new file mode 100644 index 00000000000..b857b0f57af --- /dev/null +++ b/src/script_interface/integrators/Integrator.cpp @@ -0,0 +1,62 @@ +/* + * Copyright (C) 2022-2023 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 . + */ + +#include "Integrator.hpp" + +#include "core/integrate.hpp" + +#include "script_interface/ScriptInterface.hpp" + +#include +#include + +namespace ScriptInterface { +namespace Integrators { + +Variant Integrator::integrate(VariantMap const ¶ms) { + auto const steps = get_value(params, "steps"); + auto const recalc_forces_flag = get_value(params, "recalc_forces"); + auto const reuse_forces_flag = get_value(params, "reuse_forces"); + + context()->parallel_try_catch([&]() { + if (steps < 0) { + throw std::domain_error("Parameter 'steps' must be positive"); + } + if (recalc_forces_flag and reuse_forces_flag) { + throw std::invalid_argument( + "cannot reuse old forces and recalculate forces"); + } + }); + + auto constexpr update_accumulators = true; + auto const reuse_forces = static_cast(reuse_forces_flag); + return ::integrate_with_signal_handler(steps, reuse_forces, + update_accumulators); +} + +Variant Integrator::do_call_method(std::string const &name, + VariantMap const ¶ms) { + if (name == "integrate") { + return integrate(params); + } + return none; +} + +} // namespace Integrators +} // namespace ScriptInterface diff --git a/src/script_interface/integrators/Integrator.hpp b/src/script_interface/integrators/Integrator.hpp index 85f1f05af15..e2930901563 100644 --- a/src/script_interface/integrators/Integrator.hpp +++ b/src/script_interface/integrators/Integrator.hpp @@ -31,6 +31,9 @@ namespace Integrators { class Integrator : public ObjectHandle { public: + Variant do_call_method(std::string const &name, + VariantMap const ¶ms) override; + virtual Variant integrate(VariantMap const ¶ms); virtual void activate() const = 0; }; diff --git a/src/script_interface/integrators/SteepestDescent.cpp b/src/script_interface/integrators/SteepestDescent.cpp index 18f8bb180b3..9c395f54206 100644 --- a/src/script_interface/integrators/SteepestDescent.cpp +++ b/src/script_interface/integrators/SteepestDescent.cpp @@ -32,6 +32,19 @@ namespace ScriptInterface { namespace Integrators { +Variant SteepestDescent::integrate(VariantMap const ¶ms) { + auto constexpr reuse_forces = -1; + auto constexpr update_accumulators = false; + auto const steps = get_value(params, "steps"); + context()->parallel_try_catch([&]() { + if (steps < 0) { + throw std::domain_error("Parameter 'steps' must be positive"); + } + }); + return ::integrate_with_signal_handler(steps, reuse_forces, + update_accumulators); +} + SteepestDescent::SteepestDescent() { add_parameters({ {"f_max", AutoParameter::read_only, diff --git a/src/script_interface/integrators/SteepestDescent.hpp b/src/script_interface/integrators/SteepestDescent.hpp index 6e6aa6c80e3..73d80da6dfc 100644 --- a/src/script_interface/integrators/SteepestDescent.hpp +++ b/src/script_interface/integrators/SteepestDescent.hpp @@ -40,6 +40,7 @@ class SteepestDescent : public AutoParameters { SteepestDescent(); void do_construct(VariantMap const ¶ms) override; + Variant integrate(VariantMap const ¶ms) override; void activate() const override; ::SteepestDescentParameters const &get_instance() const { diff --git a/src/script_interface/integrators/VelocityVerlet.hpp b/src/script_interface/integrators/VelocityVerlet.hpp index 30ae3dc3d25..693b9039bdb 100644 --- a/src/script_interface/integrators/VelocityVerlet.hpp +++ b/src/script_interface/integrators/VelocityVerlet.hpp @@ -28,7 +28,7 @@ namespace ScriptInterface { namespace Integrators { -class VelocityVerlet : public Integrator { +class VelocityVerlet : public AutoParameters { void activate() const override; }; diff --git a/testsuite/python/integrator_exceptions.py b/testsuite/python/integrator_exceptions.py index 1e8a7c6b648..1461755e106 100644 --- a/testsuite/python/integrator_exceptions.py +++ b/testsuite/python/integrator_exceptions.py @@ -24,11 +24,11 @@ import unittest_decorators as utx -class Integrator_test(ut.TestCase): +class Test(ut.TestCase): system = espressomd.System(box_l=[1., 1., 1.]) system.time_step = 0.01 - msg = 'Encountered errors during integrate: ERROR: ' + msg = r'while calling method integrate\(\): ERROR: ' def setUp(self): self.system.part.add(pos=(0, 0, 0)) @@ -47,13 +47,14 @@ def test_00_common_interface(self): self.system.time_step = -2. with self.assertRaisesRegex(ValueError, 'time_step must be > 0.'): self.system.time_step = 0. - with self.assertRaisesRegex(ValueError, 'steps must be positive'): + with self.assertRaisesRegex(ValueError, "Parameter 'steps' must be positive"): self.system.integrator.run(steps=-1) - with self.assertRaisesRegex(Exception, self.msg + 'cannot automatically determine skin, please set it manually'): + with self.assertRaisesRegex(RuntimeError, 'cannot automatically determine skin, please set it manually'): self.system.integrator.run() self.system.cell_system.skin = 0.4 - with self.assertRaisesRegex(Exception, self.msg + 'cannot reuse old forces and recalculate forces'): + with self.assertRaisesRegex(ValueError, 'cannot reuse old forces and recalculate forces'): self.system.integrator.run(recalc_forces=True, reuse_forces=True) + self.assertIsNone(self.system.integrator.call_method("unknown")) if espressomd.has_features("WCA"): wca = self.system.non_bonded_inter[0, 0].wca wca.set_params(epsilon=1., sigma=0.01) diff --git a/testsuite/python/integrator_steepest_descent.py b/testsuite/python/integrator_steepest_descent.py index 1574329d11f..dfb5cf91965 100644 --- a/testsuite/python/integrator_steepest_descent.py +++ b/testsuite/python/integrator_steepest_descent.py @@ -152,6 +152,10 @@ def test_integrator_exceptions(self): with self.assertRaises(RuntimeError): self.system.integrator.set_steepest_descent( f_max=0, gamma=1, max_displacement=-1) + with self.assertRaisesRegex(ValueError, "Parameter 'steps' must be positive"): + self.system.integrator.set_steepest_descent( + f_max=0, gamma=1, max_displacement=0.01) + self.system.integrator.run(steps=-1) def test_integrator_recovery(self): # the system is still in a valid state after a failure diff --git a/testsuite/python/interactions_bonded.py b/testsuite/python/interactions_bonded.py index c4ff9afc7da..3584e2e587b 100644 --- a/testsuite/python/interactions_bonded.py +++ b/testsuite/python/interactions_bonded.py @@ -253,7 +253,7 @@ def run_test(self, bond_instance, force_func, energy_func, min_dist, p_tensor_sim, p_tensor_expected, atol=1E-12) if test_breakage: p2.pos = p1.pos + self.axis * cutoff * 1.01 - with self.assertRaisesRegex(Exception, "Encountered errors during integrate"): + with self.assertRaisesRegex(Exception, r"while calling method integrate\(\)"): self.system.integrator.run(recalc_forces=True, steps=0) if test_same_pos_exception: p2.pos = p1.pos diff --git a/testsuite/python/scafacos_interface.py b/testsuite/python/scafacos_interface.py index 5ad5715b717..8e26a66a773 100644 --- a/testsuite/python/scafacos_interface.py +++ b/testsuite/python/scafacos_interface.py @@ -383,6 +383,7 @@ def test_electrostatics_plus_magnetostatics(self): gamma=0.001, max_displacement=0.01) system.integrator.run(100) + system.integrator.set_vv() # compute forces and energies p3m_E_coulomb, p3m_E_dipoles, p3m_forces, p3m_torques = self.p3m_data() diff --git a/testsuite/python/sigint.py b/testsuite/python/sigint.py index 71553f51479..9467bfe7d7d 100644 --- a/testsuite/python/sigint.py +++ b/testsuite/python/sigint.py @@ -31,7 +31,7 @@ def setUp(self): self.process = subprocess.Popen([sys.executable, script]) def test_signal_handling(self): - self.process.send_signal(signal.SIGINT) + self.process.send_signal(signal.Signals.SIGINT) # Wait for the signal to arrive and one integration step to be finished time.sleep(1) self.assertIsNotNone(self.process.poll()) diff --git a/testsuite/python/test_checkpoint.py b/testsuite/python/test_checkpoint.py index 2174d3b5a53..4c99af5bd01 100644 --- a/testsuite/python/test_checkpoint.py +++ b/testsuite/python/test_checkpoint.py @@ -275,7 +275,10 @@ def test_shape_based_constraints_serialization(self): p3, p4 = system.part.by_ids([3, 4]) old_force = np.copy(p3.f) system.constraints.remove(system.constraints[0]) + old_integrator = system.integrator.integrator + system.integrator.set_vv() system.integrator.run(0, recalc_forces=True) + system.integrator.integrator = old_integrator np.testing.assert_allclose( np.copy(p3.f), -np.copy(p4.f), rtol=1e-4) self.assertGreater(np.linalg.norm(np.copy(p3.f) - old_force), 1e6)