From 83a4c9820cb37c9c092ccf1aa4b38c5ad6b08e4b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Fri, 9 Sep 2022 16:35:54 +0200 Subject: [PATCH] Decythonize non-bonded interactions --- src/core/nonbonded_interactions/thole.hpp | 2 - src/python/espressomd/interactions.pxd | 6 - src/python/espressomd/interactions.pyx | 1499 ++++++++--------- .../interactions/NonBondedInteraction.hpp | 73 +- .../interactions/NonBondedInteractions.hpp | 13 +- testsuite/python/dpd.py | 2 +- .../python/interactions_bonded_interface.py | 13 + .../interactions_non-bonded_interface.py | 70 +- testsuite/python/save_checkpoint.py | 9 + testsuite/python/tabulated.py | 4 + testsuite/python/test_checkpoint.py | 11 +- 11 files changed, 795 insertions(+), 907 deletions(-) diff --git a/src/core/nonbonded_interactions/thole.hpp b/src/core/nonbonded_interactions/thole.hpp index 9c9f4acdebb..041840b8bb3 100644 --- a/src/core/nonbonded_interactions/thole.hpp +++ b/src/core/nonbonded_interactions/thole.hpp @@ -23,8 +23,6 @@ /** \file * Routines to calculate the Thole damping potential between particle pairs. * See @cite thole81a. - * - * Implementation in \ref thole.cpp. */ #include "config.hpp" diff --git a/src/python/espressomd/interactions.pxd b/src/python/espressomd/interactions.pxd index 6d5ca5f26fc..e8dd0dd7013 100644 --- a/src/python/espressomd/interactions.pxd +++ b/src/python/espressomd/interactions.pxd @@ -23,12 +23,6 @@ from libc cimport stdint from .thermostat cimport thermalized_bond -include "myconfig.pxi" - -# force include of config.hpp -cdef extern from "config.hpp": - pass - cdef extern from "script_interface/interactions/bonded.hpp": int bonded_ia_params_zero_based_type(int bond_id) except + diff --git a/src/python/espressomd/interactions.pyx b/src/python/espressomd/interactions.pyx index 36c2b25d1df..c1510cc96be 100644 --- a/src/python/espressomd/interactions.pyx +++ b/src/python/espressomd/interactions.pyx @@ -18,8 +18,8 @@ # cimport cpython.object -include "myconfig.pxi" from . import utils +from . import code_features from .script_interface import ScriptObjectMap, ScriptInterfaceHelper, script_interface_register @@ -36,45 +36,34 @@ class NonBondedInteraction(ScriptInterfaceHelper): _so_bind_methods = ("deactivate",) def __init__(self, **kwargs): + code_features.assert_features(self.__class__.__dict__["_so_feature"]) if "sip" in kwargs: super().__init__(**kwargs) else: - self._validate(kwargs) params = self.default_params() params.update(kwargs) super().__init__(**params) - self._validate_post(params) def __str__(self): return f'{self.__class__.__name__}({self.get_params()})' - def _validate(self, params): - for key in self.required_keys(): - if key not in params: - raise RuntimeError( - f"{self.__class__.__name__} parameter '{key}' is missing") - - def _validate_post(self, params): - valid_parameters = self.valid_keys() - for key in params: - if key != "_types" and key not in valid_parameters: - raise RuntimeError( - f"Parameter '{key}' is not a valid {self.__class__.__name__} parameter") - def set_params(self, **kwargs): """Set new parameters. """ - self._validate(kwargs) params = self.default_params() params.update(kwargs) - self._validate_post(params) err_msg = f"setting {self.__class__.__name__} raised an error" self.call_method("set_params", handle_errors_message=err_msg, **params) - def _serialize(self): - return (self.__class__.__name__, self.get_params()) + def __reduce__(self): + return (NonBondedInteraction._restore_object, + (self.__class__, self.get_params())) + + @classmethod + def _restore_object(cls, derived_class, kwargs): + return derived_class(**kwargs) def default_params(self): """Virtual method. @@ -83,709 +72,607 @@ class NonBondedInteraction(ScriptInterfaceHelper): raise Exception( "Subclasses of NonBondedInteraction must define the default_params() method.") - def valid_keys(self): - """All parameters that can be set. - """ - return set(self._valid_parameters()) +@script_interface_register +class LennardJonesInteraction(NonBondedInteraction): + """ + Standard 6-12 Lennard-Jones potential. - def required_keys(self): - """Virtual method. + Methods + ------- + set_params() + Set new parameters for the interaction. - """ - raise Exception( - "Subclasses of NonBondedInteraction must define the required_keys() method.") + Parameters + ---------- + epsilon : :obj:`float` + Magnitude of the interaction. + sigma : :obj:`float` + Interaction length scale. + cutoff : :obj:`float` + Cutoff distance of the interaction. + shift : :obj:`float` or :obj:`str` \{'auto'\} + Constant shift of the potential. If ``'auto'``, a default value + is computed from ``sigma`` and ``cutoff``. The LJ potential + will be shifted by :math:`4\\epsilon\\cdot\\text{shift}`. + offset : :obj:`float`, optional + Offset distance of the interaction. + min : :obj:`float`, optional + Restricts the interaction to a minimal distance. + + """ + _so_name = "Interactions::InteractionLJ" + _so_feature = "LENNARD_JONES" -IF LENNARD_JONES == 1: + def default_params(self): + """Python dictionary of default parameters. - @script_interface_register - class LennardJonesInteraction(NonBondedInteraction): """ - Standard 6-12 Lennard-Jones potential. + return {"offset": 0., "min": 0.} - Methods - ------- - set_params() - Set new parameters for the interaction. - Parameters - ---------- - epsilon : :obj:`float` - Magnitude of the interaction. - sigma : :obj:`float` - Interaction length scale. - cutoff : :obj:`float` - Cutoff distance of the interaction. - shift : :obj:`float` or :obj:`str` \{'auto'\} - Constant shift of the potential. If ``'auto'``, a default value - is computed from ``sigma`` and ``cutoff``. The LJ potential - will be shifted by :math:`4\\epsilon\\cdot\\text{shift}`. - offset : :obj:`float`, optional - Offset distance of the interaction. - min : :obj:`float`, optional - Restricts the interaction to a minimal distance. +@script_interface_register +class WCAInteraction(NonBondedInteraction): + """ + Standard 6-12 Weeks-Chandler-Andersen potential. - """ + Methods + ------- + set_params() + Set new parameters for the interaction. - _so_name = "Interactions::InteractionLJ" + Parameters + ---------- + epsilon : :obj:`float` + Magnitude of the interaction. + sigma : :obj:`float` + Interaction length scale. - def default_params(self): - """Python dictionary of default parameters. + """ - """ - return {"offset": 0., "min": 0.} + _so_name = "Interactions::InteractionWCA" + _so_feature = "WCA" - def required_keys(self): - """Parameters that have to be set. + def default_params(self): + """Python dictionary of default parameters. - """ - return {"epsilon", "sigma", "cutoff", "shift"} + """ + return {} -IF WCA == 1: + @property + def cutoff(self): + return self.call_method("get_cutoff") - @script_interface_register - class WCAInteraction(NonBondedInteraction): - """ - Standard 6-12 Weeks-Chandler-Andersen potential. - Methods - ------- - set_params() - Set new parameters for the interaction. +@script_interface_register +class GenericLennardJonesInteraction(NonBondedInteraction): + """ + Generalized Lennard-Jones potential. - Parameters - ---------- - epsilon : :obj:`float` - Magnitude of the interaction. - sigma : :obj:`float` - Interaction length scale. + Methods + ------- + set_params() + Set new parameters for the interaction. - """ + Parameters + ---------- + epsilon : :obj:`float` + Magnitude of the interaction. + sigma : :obj:`float` + Interaction length scale. + cutoff : :obj:`float` + Cutoff distance of the interaction. + shift : :obj:`float` or :obj:`str` \{'auto'\} + Constant shift of the potential. If ``'auto'``, a default value + is computed from the other parameters. The LJ potential + will be shifted by :math:`\\epsilon\\cdot\\text{shift}`. + offset : :obj:`float` + Offset distance of the interaction. + e1 : :obj:`int` + Exponent of the repulsion term. + e2 : :obj:`int` + Exponent of the attraction term. + b1 : :obj:`float` + Prefactor of the repulsion term. + b2 : :obj:`float` + Prefactor of the attraction term. + delta : :obj:`float`, optional + ``LJGEN_SOFTCORE`` parameter delta. Allows control over how + smoothly the potential drops to zero as lambda approaches zero. + lam : :obj:`float`, optional + ``LJGEN_SOFTCORE`` parameter lambda. Tune the strength of the + interaction. - _so_name = "Interactions::InteractionWCA" + """ - def default_params(self): - """Python dictionary of default parameters. + _so_name = "Interactions::InteractionLJGen" + _so_feature = "LENNARD_JONES_GENERIC" - """ - return {} + def default_params(self): + """Python dictionary of default parameters. - def required_keys(self): - """Parameters that have to be set. + """ + if code_features.has_features("LJGEN_SOFTCORE"): + return {"delta": 0., "lam": 1.} + return {} - """ - return {"epsilon", "sigma"} - @property - def cutoff(self): - return self.call_method("get_cutoff") +@script_interface_register +class LennardJonesCosInteraction(NonBondedInteraction): + """Lennard-Jones cosine interaction. -IF LENNARD_JONES_GENERIC == 1: + Methods + ------- + set_params() + Set or update parameters for the interaction. + Parameters marked as required become optional once the + interaction has been activated for the first time; + subsequent calls to this method update the existing values. - @script_interface_register - class GenericLennardJonesInteraction(NonBondedInteraction): - """ - Generalized Lennard-Jones potential. + Parameters + ---------- + epsilon : :obj:`float` + Magnitude of the interaction. + sigma : :obj:`float` + Interaction length scale. + cutoff : :obj:`float` + Cutoff distance of the interaction. + offset : :obj:`float`, optional + Offset distance of the interaction. + """ - Methods - ------- - set_params() - Set new parameters for the interaction. + _so_name = "Interactions::InteractionLJcos" + _so_feature = "LJCOS" - Parameters - ---------- - epsilon : :obj:`float` - Magnitude of the interaction. - sigma : :obj:`float` - Interaction length scale. - cutoff : :obj:`float` - Cutoff distance of the interaction. - shift : :obj:`float` or :obj:`str` \{'auto'\} - Constant shift of the potential. If ``'auto'``, a default value - is computed from the other parameters. The LJ potential - will be shifted by :math:`\\epsilon\\cdot\\text{shift}`. - offset : :obj:`float` - Offset distance of the interaction. - e1 : :obj:`int` - Exponent of the repulsion term. - e2 : :obj:`int` - Exponent of the attraction term. - b1 : :obj:`float` - Prefactor of the repulsion term. - b2 : :obj:`float` - Prefactor of the attraction term. - delta : :obj:`float`, optional - ``LJGEN_SOFTCORE`` parameter delta. Allows control over how - smoothly the potential drops to zero as lambda approaches zero. - lam : :obj:`float`, optional - ``LJGEN_SOFTCORE`` parameter lambda. Tune the strength of the - interaction. + def default_params(self): + """Python dictionary of default parameters. """ + return {"offset": 0.} - _so_name = "Interactions::InteractionLJGen" - def default_params(self): - """Python dictionary of default parameters. +@script_interface_register +class LennardJonesCos2Interaction(NonBondedInteraction): + """Second variant of the Lennard-Jones cosine interaction. - """ - IF LJGEN_SOFTCORE: - return {"delta": 0., "lam": 1.} - ELSE: - return {} + Methods + ------- + set_params() + Set new parameters for the interaction. - def required_keys(self): - """Parameters that have to be set. + Parameters + ---------- + epsilon : :obj:`float` + Magnitude of the interaction. + sigma : :obj:`float` + Interaction length scale. + offset : :obj:`float`, optional + Offset distance of the interaction. + width : :obj:`float` + Width of interaction. + """ - """ - return {"epsilon", "sigma", "cutoff", - "shift", "offset", "e1", "e2", "b1", "b2"} + _so_name = "Interactions::InteractionLJcos2" + _so_feature = "LJCOS2" -IF LJCOS == 1: + def default_params(self): + """Python dictionary of default parameters. - @script_interface_register - class LennardJonesCosInteraction(NonBondedInteraction): - """Lennard-Jones cosine interaction. + """ + return {"offset": 0.} - Methods - ------- - set_params() - Set or update parameters for the interaction. - Parameters marked as required become optional once the - interaction has been activated for the first time; - subsequent calls to this method update the existing values. + @property + def cutoff(self): + return self.call_method("get_cutoff") - Parameters - ---------- - epsilon : :obj:`float` - Magnitude of the interaction. - sigma : :obj:`float` - Interaction length scale. - cutoff : :obj:`float` - Cutoff distance of the interaction. - offset : :obj:`float`, optional - Offset distance of the interaction. - """ - _so_name = "Interactions::InteractionLJcos" - - def default_params(self): - """Python dictionary of default parameters. - - """ - return {"offset": 0.} - - def required_keys(self): - """Parameters that have to be set. - - """ - return {"epsilon", "sigma", "cutoff"} - -IF LJCOS2 == 1: - - @script_interface_register - class LennardJonesCos2Interaction(NonBondedInteraction): - """Second variant of the Lennard-Jones cosine interaction. - - Methods - ------- - set_params() - Set new parameters for the interaction. - - Parameters - ---------- - epsilon : :obj:`float` - Magnitude of the interaction. - sigma : :obj:`float` - Interaction length scale. - offset : :obj:`float`, optional - Offset distance of the interaction. - width : :obj:`float` - Width of interaction. - """ +@script_interface_register +class HatInteraction(NonBondedInteraction): + """Hat interaction. - _so_name = "Interactions::InteractionLJcos2" + Methods + ------- + set_params() + Set new parameters for the interaction. - def default_params(self): - """Python dictionary of default parameters. + Parameters + ---------- + F_max : :obj:`float` + Magnitude of the interaction. + cutoff : :obj:`float` + Cutoff distance of the interaction. - """ - return {"offset": 0.} + """ - def required_keys(self): - """Parameters that have to be set. + _so_name = "Interactions::InteractionHat" + _so_feature = "HAT" - """ - return {"epsilon", "sigma", "width"} + def default_params(self): + return {} + + +@script_interface_register +class GayBerneInteraction(NonBondedInteraction): + """Gay--Berne interaction. - @property - def cutoff(self): - return self.call_method("get_cutoff") + Methods + ------- + set_params() + Set new parameters for the interaction. -IF HAT == 1: + Parameters + ---------- + eps : :obj:`float` + Potential well depth. + sig : :obj:`float` + Interaction range. + cut : :obj:`float` + Cutoff distance of the interaction. + k1 : :obj:`float` or :obj:`str` + Molecular elongation. + k2 : :obj:`float`, optional + Ratio of the potential well depths for the side-by-side + and end-to-end configurations. + mu : :obj:`float`, optional + Adjustable exponent. + nu : :obj:`float`, optional + Adjustable exponent. - @script_interface_register - class HatInteraction(NonBondedInteraction): - """Hat interaction. + """ - Methods - ------- - set_params() - Set new parameters for the interaction. + _so_name = "Interactions::InteractionGayBerne" + _so_feature = "GAY_BERNE" - Parameters - ---------- - F_max : :obj:`float` - Magnitude of the interaction. - cutoff : :obj:`float` - Cutoff distance of the interaction. + def default_params(self): + """Python dictionary of default parameters. """ + return {} - _so_name = "Interactions::InteractionHat" - def default_params(self): - return {} +@script_interface_register +class TabulatedNonBonded(NonBondedInteraction): + """Tabulated interaction. - def required_keys(self): - return {"F_max", "cutoff"} + Methods + ------- + set_params() + Set new parameters for the interaction. -IF GAY_BERNE: + Parameters + ---------- + min : :obj:`float`, + The minimal interaction distance. + max : :obj:`float`, + The maximal interaction distance. + energy: array_like of :obj:`float` + The energy table. + force: array_like of :obj:`float` + The force table. - @script_interface_register - class GayBerneInteraction(NonBondedInteraction): - """Gay--Berne interaction. + """ - Methods - ------- - set_params() - Set new parameters for the interaction. + _so_name = "Interactions::InteractionTabulated" + _so_feature = "TABULATED" - Parameters - ---------- - eps : :obj:`float` - Potential well depth. - sig : :obj:`float` - Interaction range. - cut : :obj:`float` - Cutoff distance of the interaction. - k1 : :obj:`float` or :obj:`str` - Molecular elongation. - k2 : :obj:`float`, optional - Ratio of the potential well depths for the side-by-side - and end-to-end configurations. - mu : :obj:`float`, optional - Adjustable exponent. - nu : :obj:`float`, optional - Adjustable exponent. + def default_params(self): + """Python dictionary of default parameters. """ + return {} - _so_name = "Interactions::InteractionGayBerne" + @property + def cutoff(self): + return self.call_method("get_cutoff") - def default_params(self): - """Python dictionary of default parameters. - """ - return {} +@script_interface_register +class DPDInteraction(NonBondedInteraction): + """DPD interaction. - def required_keys(self): - """Parameters that have to be set. + Methods + ------- + set_params() + Set new parameters for the interaction. - """ - return {"eps", "sig", "cut", "k1", "k2", "mu", "nu"} + Parameters + ---------- + weight_function : :obj:`int`, \{0, 1\} + The distance dependence of the parallel part, + either 0 (constant) or 1 (linear) + gamma : :obj:`float` + Friction coefficient of the parallel part + k : :obj:`float` + Exponent in the modified weight function + r_cut : :obj:`float` + Cutoff of the parallel part + trans_weight_function : :obj:`int`, \{0, 1\} + The distance dependence of the orthogonal part, + either 0 (constant) or 1 (linear) + trans_gamma : :obj:`float` + Friction coefficient of the orthogonal part + trans_r_cut : :obj:`float` + Cutoff of the orthogonal part -IF TABULATED: + """ - @script_interface_register - class TabulatedNonBonded(NonBondedInteraction): - """Tabulated interaction. + _so_name = "Interactions::InteractionDPD" + _so_feature = "DPD" - Methods - ------- - set_params() - Set new parameters for the interaction. + def default_params(self): + return { + "weight_function": 0, + "gamma": 0.0, + "k": 1.0, + "r_cut": -1.0, + "trans_weight_function": 0, + "trans_gamma": 0.0, + "trans_r_cut": -1.0} - Parameters - ---------- - min : :obj:`float`, - The minimal interaction distance. - max : :obj:`float`, - The maximal interaction distance. - energy: array_like of :obj:`float` - The energy table. - force: array_like of :obj:`float` - The force table. - """ +@script_interface_register +class SmoothStepInteraction(NonBondedInteraction): + """Smooth step interaction. - _so_name = "Interactions::InteractionTabulated" + Methods + ------- + set_params() + Set new parameters for the interaction. - def default_params(self): - """Python dictionary of default parameters. + Parameters + ---------- + d : :obj:`float` + Short range repulsion parameter. + n : :obj:`int`, optional + Exponent of short range repulsion. + eps : :obj:`float` + Magnitude of the second (soft) repulsion. + k0 : :obj:`float`, optional + Exponential factor in second (soft) repulsion. + sig : :obj:`float`, optional + Length scale of second (soft) repulsion. + cutoff : :obj:`float` + Cutoff distance of the interaction. - """ - return {} + """ - def required_keys(self): - """Parameters that have to be set. + _so_name = "Interactions::InteractionSmoothStep" + _so_feature = "SMOOTH_STEP" - """ - return {"min", "max", "energy", "force"} + def default_params(self): + """Python dictionary of default parameters. - @property - def cutoff(self): - return self.call_method("get_cutoff") + """ + return {"n": 10, "k0": 0., "sig": 0.} -IF DPD: - @script_interface_register - class DPDInteraction(NonBondedInteraction): - """DPD interaction. - - Methods - ------- - set_params() - Set new parameters for the interaction. - - Parameters - ---------- - weight_function : :obj:`int`, \{0, 1\} - The distance dependence of the parallel part, - either 0 (constant) or 1 (linear) - gamma : :obj:`float` - Friction coefficient of the parallel part - k : :obj:`float` - Exponent in the modified weight function - r_cut : :obj:`float` - Cutoff of the parallel part - trans_weight_function : :obj:`int`, \{0, 1\} - The distance dependence of the orthogonal part, - either 0 (constant) or 1 (linear) - trans_gamma : :obj:`float` - Friction coefficient of the orthogonal part - trans_r_cut : :obj:`float` - Cutoff of the orthogonal part +@script_interface_register +class BMHTFInteraction(NonBondedInteraction): + """BMHTF interaction. - """ + Methods + ------- + set_params() + Set new parameters for the interaction. - _so_name = "Interactions::InteractionDPD" + Parameters + ---------- + a : :obj:`float` + Magnitude of exponential part of the interaction. + b : :obj:`float` + Exponential factor of the interaction. + c : :obj:`float` + Magnitude of the term decaying with the sixth power of r. + d : :obj:`float` + Magnitude of the term decaying with the eighth power of r. + sig : :obj:`float` + Shift in the exponent. + cutoff : :obj:`float` + Cutoff distance of the interaction. - def default_params(self): - return { - "weight_function": 0, - "gamma": 0.0, - "k": 1.0, - "r_cut": -1.0, - "trans_weight_function": 0, - "trans_gamma": 0.0, - "trans_r_cut": -1.0} + """ - def required_keys(self): - return set() - -IF SMOOTH_STEP == 1: - - @script_interface_register - class SmoothStepInteraction(NonBondedInteraction): - """Smooth step interaction. - - Methods - ------- - set_params() - Set new parameters for the interaction. + _so_name = "Interactions::InteractionBMHTF" + _so_feature = "BMHTF_NACL" - Parameters - ---------- - d : :obj:`float` - Short range repulsion parameter. - n : :obj:`int`, optional - Exponent of short range repulsion. - eps : :obj:`float` - Magnitude of the second (soft) repulsion. - k0 : :obj:`float`, optional - Exponential factor in second (soft) repulsion. - sig : :obj:`float`, optional - Length scale of second (soft) repulsion. - cutoff : :obj:`float` - Cutoff distance of the interaction. + def default_params(self): + """Python dictionary of default parameters. """ + return {} - _so_name = "Interactions::InteractionSmoothStep" - def default_params(self): - """Python dictionary of default parameters. +@script_interface_register +class MorseInteraction(NonBondedInteraction): + """Morse interaction. - """ - return {"n": 10, "k0": 0., "sig": 0.} - - def required_keys(self): - """Parameters that have to be set. - - """ - return {"d", "eps", "cutoff"} - -IF BMHTF_NACL == 1: - - @script_interface_register - class BMHTFInteraction(NonBondedInteraction): - """BMHTF interaction. - - Methods - ------- - set_params() - Set new parameters for the interaction. - - Parameters - ---------- - a : :obj:`float` - Magnitude of exponential part of the interaction. - b : :obj:`float` - Exponential factor of the interaction. - c : :obj:`float` - Magnitude of the term decaying with the sixth power of r. - d : :obj:`float` - Magnitude of the term decaying with the eighth power of r. - sig : :obj:`float` - Shift in the exponent. - cutoff : :obj:`float` - Cutoff distance of the interaction. - """ - - _so_name = "Interactions::InteractionBMHTF" - - def default_params(self): - """Python dictionary of default parameters. + Methods + ------- + set_params() + Set new parameters for the interaction. - """ - return {} + Parameters + ---------- + eps : :obj:`float` + The magnitude of the interaction. + alpha : :obj:`float` + Stiffness of the Morse interaction. + rmin : :obj:`float` + Distance of potential minimum + cutoff : :obj:`float`, optional + Cutoff distance of the interaction. - def required_keys(self): - """Parameters that have to be set. + """ - """ - return {"a", "b", "c", "d", "sig", "cutoff"} + _so_name = "Interactions::InteractionMorse" + _so_feature = "MORSE" -IF MORSE == 1: - - @script_interface_register - class MorseInteraction(NonBondedInteraction): - """Morse interaction. - - Methods - ------- - set_params() - Set new parameters for the interaction. + def default_params(self): + """Python dictionary of default parameters. - Parameters - ---------- - eps : :obj:`float` - The magnitude of the interaction. - alpha : :obj:`float` - Stiffness of the Morse interaction. - rmin : :obj:`float` - Distance of potential minimum - cutoff : :obj:`float`, optional - Cutoff distance of the interaction. """ + return {"cutoff": 0.} - _so_name = "Interactions::InteractionMorse" - - def default_params(self): - """Python dictionary of default parameters. - """ - return {"cutoff": 0.} +@script_interface_register +class BuckinghamInteraction(NonBondedInteraction): + """Buckingham interaction. - def required_keys(self): - """Parameters that have to be set. + Methods + ------- + set_params() + Set new parameters for the interaction. - """ - return {"eps", "alpha", "rmin"} + Parameters + ---------- + a : :obj:`float` + Magnitude of the exponential part of the interaction. + b : :obj:`float`, optional + Exponent of the exponential part of the interaction. + c : :obj:`float` + Prefactor of term decaying with the sixth power of distance. + d : :obj:`float` + Prefactor of term decaying with the fourth power of distance. + discont : :obj:`float` + Distance below which the potential is linearly continued. + cutoff : :obj:`float` + Cutoff distance of the interaction. + shift: :obj:`float`, optional + Constant potential shift. -IF BUCKINGHAM == 1: + """ - @script_interface_register - class BuckinghamInteraction(NonBondedInteraction): - """Buckingham interaction. + _so_name = "Interactions::InteractionBuckingham" + _so_feature = "BUCKINGHAM" - Methods - ------- - set_params() - Set new parameters for the interaction. + def default_params(self): + """Python dictionary of default parameters. - Parameters - ---------- - a : :obj:`float` - Magnitude of the exponential part of the interaction. - b : :obj:`float`, optional - Exponent of the exponential part of the interaction. - c : :obj:`float` - Prefactor of term decaying with the sixth power of distance. - d : :obj:`float` - Prefactor of term decaying with the fourth power of distance. - discont : :obj:`float` - Distance below which the potential is linearly continued. - cutoff : :obj:`float` - Cutoff distance of the interaction. - shift: :obj:`float`, optional - Constant potential shift. """ + return {"b": 0., "shift": 0.} - _so_name = "Interactions::InteractionBuckingham" - - def default_params(self): - """Python dictionary of default parameters. - """ - return {"b": 0., "shift": 0.} +@script_interface_register +class SoftSphereInteraction(NonBondedInteraction): + """Soft sphere interaction. - def required_keys(self): - """Parameters that have to be set. + Methods + ------- + set_params() + Set new parameters for the interaction. - """ - return {"a", "c", "d", "discont", "cutoff"} + Parameters + ---------- + a : :obj:`float` + Magnitude of the interaction. + n : :obj:`float` + Exponent of the power law. + cutoff : :obj:`float` + Cutoff distance of the interaction. + offset : :obj:`float`, optional + Offset distance of the interaction. -IF SOFT_SPHERE == 1: + """ - @script_interface_register - class SoftSphereInteraction(NonBondedInteraction): - """Soft sphere interaction. + _so_name = "Interactions::InteractionSoftSphere" + _so_feature = "SOFT_SPHERE" - Methods - ------- - set_params() - Set new parameters for the interaction. + def default_params(self): + """Python dictionary of default parameters. - Parameters - ---------- - a : :obj:`float` - Magnitude of the interaction. - n : :obj:`float` - Exponent of the power law. - cutoff : :obj:`float` - Cutoff distance of the interaction. - offset : :obj:`float`, optional - Offset distance of the interaction. """ + return {"offset": 0.} - _so_name = "Interactions::InteractionSoftSphere" - - def default_params(self): - """Python dictionary of default parameters. - """ - return {"offset": 0.} +@script_interface_register +class HertzianInteraction(NonBondedInteraction): + """Hertzian interaction. - def required_keys(self): - """Parameters that have to be set. + Methods + ------- + set_params() + Set new parameters for the interaction. - """ - return {"a", "n", "cutoff"} + Parameters + ---------- + eps : :obj:`float` + Magnitude of the interaction. + sig : :obj:`float` + Interaction length scale. -IF HERTZIAN == 1: + """ - @script_interface_register - class HertzianInteraction(NonBondedInteraction): - """Hertzian interaction. + _so_name = "Interactions::InteractionHertzian" + _so_feature = "HERTZIAN" - Methods - ------- - set_params() - Set new parameters for the interaction. + def default_params(self): + """Python dictionary of default parameters. - Parameters - ---------- - eps : :obj:`float` - Magnitude of the interaction. - sig : :obj:`float` - Interaction length scale. """ + return {} - _so_name = "Interactions::InteractionHertzian" - - def default_params(self): - """Python dictionary of default parameters. - """ - return {} +@script_interface_register +class GaussianInteraction(NonBondedInteraction): + """Gaussian interaction. - def required_keys(self): - """Parameters that have to be set. + Methods + ------- + set_params() + Set new parameters for the interaction. - """ - return {"eps", "sig"} + Parameters + ---------- + eps : :obj:`float` + Overlap energy epsilon. + sig : :obj:`float` + Variance sigma of the Gaussian interaction. + cutoff : :obj:`float` + Cutoff distance of the interaction. -IF GAUSSIAN == 1: + """ - @script_interface_register - class GaussianInteraction(NonBondedInteraction): - """Gaussian interaction. + _so_name = "Interactions::InteractionGaussian" + _so_feature = "GAUSSIAN" - Methods - ------- - set_params() - Set new parameters for the interaction. + def default_params(self): + """Python dictionary of default parameters. - Parameters - ---------- - eps : :obj:`float` - Overlap energy epsilon. - sig : :obj:`float` - Variance sigma of the Gaussian interaction. - cutoff : :obj:`float` - Cutoff distance of the interaction. """ + return {} - _so_name = "Interactions::InteractionGaussian" - - def default_params(self): - """Python dictionary of default parameters. - - """ - return {} - - def required_keys(self): - """Parameters that have to be set. - - """ - return {"eps", "sig", "cutoff"} - -IF THOLE: - @script_interface_register - class TholeInteraction(NonBondedInteraction): - """Thole interaction. +@script_interface_register +class TholeInteraction(NonBondedInteraction): + """Thole interaction. - Methods - ------- - set_params() - Set new parameters for the interaction. + Methods + ------- + set_params() + Set new parameters for the interaction. - Parameters - ---------- - scaling_coeff : :obj:`float` - The factor used in the Thole damping function between - polarizable particles i and j. Usually calculated by - the polarizabilities :math:`\\alpha_i`, :math:`\\alpha_j` - and damping parameters :math:`a_i`, :math:`a_j` via - :math:`s_{ij} = \\frac{(a_i+a_j)/2}{((\\alpha_i\\cdot\\alpha_j)^{1/2})^{1/3}}` - q1q2: :obj:`float` - Charge factor of the involved charges. Has to be set because - it acts only on the portion of the Drude core charge that is - associated to the dipole of the atom. For charged, polarizable - atoms that charge is not equal to the particle charge property. - """ + Parameters + ---------- + scaling_coeff : :obj:`float` + The factor used in the Thole damping function between + polarizable particles i and j. Usually calculated by + the polarizabilities :math:`\\alpha_i`, :math:`\\alpha_j` + and damping parameters :math:`a_i`, :math:`a_j` via + :math:`s_{ij} = \\frac{(a_i+a_j)/2}{((\\alpha_i\\cdot\\alpha_j)^{1/2})^{1/3}}` + q1q2: :obj:`float` + Charge factor of the involved charges. Has to be set because + it acts only on the portion of the Drude core charge that is + associated to the dipole of the atom. For charged, polarizable + atoms that charge is not equal to the particle charge property. - _so_name = "Interactions::InteractionThole" + """ - def default_params(self): - return {} + _so_name = "Interactions::InteractionThole" + _so_feature = "THOLE" - def required_keys(self): - return {"scaling_coeff", "q1q2"} + def default_params(self): + return {} @script_interface_register @@ -802,29 +689,27 @@ class NonBondedInteractionHandle(ScriptInterfaceHelper): return globals()[obj.__class__.__name__]( _types=self.call_method("get_types"), **obj.get_params()) - def __init__(self, *args, **kwargs): - if "sip" in kwargs: - super().__init__(**kwargs) - return - if len(args): - _type1, _type2 = args - else: - _type1, _type2 = kwargs.pop("_types") - if not (utils.is_valid_type(_type1, int) - and utils.is_valid_type(_type2, int)): - raise TypeError("The particle types have to be of type integer.") - super().__init__(_types=[_type1, _type2], **kwargs) - def _serialize(self): serialized = [] for name, obj in self.get_params().items(): - serialized.append((name, *obj._serialize())) + serialized.append((name, obj.__reduce__()[1])) return serialized def reset(self): for key in self._valid_parameters(): getattr(self, key).deactivate() + @classmethod + def _restore_object(cls, types, kwargs): + objects = {} + for name, (obj_class, obj_params) in kwargs: + objects[name] = obj_class(**obj_params) + return NonBondedInteractionHandle(_types=types, **objects) + + def __reduce__(self): + return (NonBondedInteractionHandle._restore_object, + (self.call_method("get_types"), self._serialize())) + @script_interface_register class NonBondedInteractions(ScriptInterfaceHelper): @@ -870,10 +755,7 @@ class NonBondedInteractions(ScriptInterfaceHelper): def __setstate__(self, params): for types, kwargs in params["state"]: - objects = {} - for name, class_name, obj_params in kwargs: - objects[name] = globals()[class_name](**obj_params) - obj = NonBondedInteractionHandle(_types=types, **objects) + obj = NonBondedInteractionHandle._restore_object(types, kwargs) self.call_method("insert", key=types, object=obj) @classmethod @@ -908,6 +790,10 @@ class BondedInteraction(ScriptInterfaceHelper): kwargs = args[0] args = [] + feature = self.__class__.__dict__.get("_so_feature") + if feature is not None: + code_features.assert_features(feature) + if not 'sip' in kwargs: if len(args) == 1 and utils.is_valid_type(args[0], int): # create a new script interface object for a bond that already @@ -1019,46 +905,6 @@ class BondedInteraction(ScriptInterfaceHelper): raise NotImplementedError("only equality comparison is supported") -class BondedInteractionNotDefined: - - def __init__(self, *args, **kwargs): - raise Exception( - self.__class__.__name__ + " not compiled into ESPResSo core") - - def type_number(self): - raise Exception(f"{self.name} has to be defined in myconfig.hpp.") - - def type_name(self): - """Name of interaction type. - - """ - raise Exception(f"{self.name} has to be defined in myconfig.hpp.") - - def valid_keys(self): - """All parameters that can be set. - - """ - raise Exception(f"{self.name} has to be defined in myconfig.hpp.") - - def required_keys(self): - """Parameters that have to be set. - - """ - raise Exception(f"{self.name} has to be defined in myconfig.hpp.") - - def get_default_params(self): - """Gets default values of optional parameters. - - """ - raise Exception(f"{self.name} has to be defined in myconfig.hpp.") - - def _get_params_from_es_core(self): - raise Exception(f"{self.name} has to be defined in myconfig.hpp.") - - def _set_params_in_es_core(self): - raise Exception(f"{self.name} has to be defined in myconfig.hpp.") - - @script_interface_register class FeneBond(BondedInteraction): @@ -1129,67 +975,65 @@ class HarmonicBond(BondedInteraction): return {"r_cut": 0.} -if ELECTROSTATICS: - - @script_interface_register - class BondedCoulomb(BondedInteraction): - - """ - Bonded Coulomb bond. +@script_interface_register +class BondedCoulomb(BondedInteraction): - Parameters - ---------- + """ + Bonded Coulomb bond. - prefactor : :obj:`float` - Coulomb prefactor of the bonded Coulomb interaction. - """ + Parameters + ---------- - _so_name = "Interactions::BondedCoulomb" + prefactor : :obj:`float` + Coulomb prefactor of the bonded Coulomb interaction. + """ - def __init__(self, *args, **kwargs): - super().__init__(*args, **kwargs) + _so_name = "Interactions::BondedCoulomb" + _so_feature = "ELECTROSTATICS" - def type_number(self): - return BONDED_IA_BONDED_COULOMB + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) - def type_name(self): - return "BONDED_COULOMB" + def type_number(self): + return BONDED_IA_BONDED_COULOMB - def get_default_params(self): - """Gets default values of optional parameters. + def type_name(self): + return "BONDED_COULOMB" - """ - return {} + def get_default_params(self): + """Gets default values of optional parameters. + """ + return {} -if ELECTROSTATICS: - @script_interface_register - class BondedCoulombSRBond(BondedInteraction): +@script_interface_register +class BondedCoulombSRBond(BondedInteraction): - """ - Bonded Coulomb short range bond. Calculates the short range part of - Coulomb interactions. + """ + Bonded Coulomb short range bond. Calculates the short range part of + Coulomb interactions. - Parameters - ---------- + Parameters + ---------- - q1q2 : :obj:`float` - Charge factor of the involved particle pair. Note the - particle charges are used to allow e.g. only partial subtraction - of the involved charges. - """ + q1q2 : :obj:`float` + Charge factor of the involved particle pair. Note the + particle charges are used to allow e.g. only partial subtraction + of the involved charges. + """ - _so_name = "Interactions::BondedCoulombSR" + _so_name = "Interactions::BondedCoulombSR" + _so_feature = "ELECTROSTATICS" - def type_number(self): - return BONDED_IA_BONDED_COULOMB_SR + def type_number(self): + return BONDED_IA_BONDED_COULOMB_SR - def type_name(self): - return "BONDED_COULOMB_SR" + def type_name(self): + return "BONDED_COULOMB_SR" - def get_default_params(self): - return {} + def get_default_params(self): + return {} @script_interface_register @@ -1262,49 +1106,44 @@ class ThermalizedBond(BondedInteraction): return {"r_cut": 0., "seed": None} -IF BOND_CONSTRAINT == 1: - - @script_interface_register - class RigidBond(BondedInteraction): - - """ - Rigid bond. +@script_interface_register +class RigidBond(BondedInteraction): - Parameters - ---------- - r : :obj:`float` - Bond length. - ptol : :obj:`float`, optional - Tolerance for positional deviations. - vtop : :obj:`float`, optional - Tolerance for velocity deviations. + """ + Rigid bond. - """ + Parameters + ---------- + r : :obj:`float` + Bond length. + ptol : :obj:`float`, optional + Tolerance for positional deviations. + vtop : :obj:`float`, optional + Tolerance for velocity deviations. - _so_name = "Interactions::RigidBond" + """ - def __init__(self, *args, **kwargs): - super().__init__(*args, **kwargs) + _so_name = "Interactions::RigidBond" + _so_feature = "BOND_CONSTRAINT" - def type_number(self): - return BONDED_IA_RIGID_BOND + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) - def type_name(self): - """Name of interaction type. + def type_number(self): + return BONDED_IA_RIGID_BOND - """ - return "RIGID" + def type_name(self): + """Name of interaction type. - def get_default_params(self): - """Gets default values of optional parameters. + """ + return "RIGID" - """ - # TODO rationality of Default Parameters has to be checked - return {"ptol": 0.001, "vtol": 0.001} + def get_default_params(self): + """Gets default values of optional parameters. -ELSE: - class RigidBond(BondedInteractionNotDefined): - name = "RIGID" + """ + # TODO rationality of Default Parameters has to be checked + return {"ptol": 0.001, "vtol": 0.001} @script_interface_register @@ -1350,164 +1189,148 @@ class Dihedral(BondedInteraction): return {} -IF TABULATED: - class _TabulatedBase(BondedInteraction): - - """ - Parent class for tabulated bonds. +@script_interface_register +class TabulatedDistance(BondedInteraction): - Parameters - ---------- + """ + Tabulated bond length. - min : :obj:`float` - The minimal interaction distance. Has to be 0 for angles and dihedrals. - max : :obj:`float` - The maximal interaction distance. Has to be pi for angles and 2pi for - dihedrals. - energy: array_like of :obj:`float` - The energy table. - force: array_like of :obj:`float` - The force table. + Parameters + ---------- - """ + min : :obj:`float` + The minimal interaction distance. + max : :obj:`float` + The maximal interaction distance. + energy: array_like of :obj:`float` + The energy table. + force: array_like of :obj:`float` + The force table. - def __init__(self, *args, **kwargs): - super().__init__(*args, **kwargs) + """ - def type_number(self): - return "BONDED_IA_TABULATED" + _so_name = "Interactions::TabulatedDistanceBond" + _so_feature = "TABULATED" - def type_name(self): - """Name of interaction type. + def __init__(self, *args, **kwargs): + super().__init__(*args, **kwargs) - """ - return "TABULATED_BOND" + def type_number(self): + return BONDED_IA_TABULATED_DISTANCE - def get_default_params(self): - """Gets default values of optional parameters. + def type_name(self): + """Name of interaction type. - """ - return {} + """ + return "TABULATED_DISTANCE" - @script_interface_register - class TabulatedDistance(_TabulatedBase): + def get_default_params(self): + """Gets default values of optional parameters. """ - Tabulated bond length. - - Parameters - ---------- + return {} - min : :obj:`float` - The minimal interaction distance. - max : :obj:`float` - The maximal interaction distance. - energy: array_like of :obj:`float` - The energy table. - force: array_like of :obj:`float` - The force table. - """ +@script_interface_register +class TabulatedAngle(BondedInteraction): - _so_name = "Interactions::TabulatedDistanceBond" + """ + Tabulated bond angle. - def __init__(self, *args, **kwargs): - super().__init__(*args, **kwargs) + Parameters + ---------- - def type_number(self): - return BONDED_IA_TABULATED_DISTANCE + energy: array_like of :obj:`float` + The energy table for the range :math:`0-\\pi`. + force: array_like of :obj:`float` + The force table for the range :math:`0-\\pi`. - def type_name(self): - """Name of interaction type. + """ - """ - return "TABULATED_DISTANCE" + _so_name = "Interactions::TabulatedAngleBond" + _so_feature = "TABULATED" - @script_interface_register - class TabulatedAngle(_TabulatedBase): + pi = 3.14159265358979 - """ - Tabulated bond angle. + def __init__(self, *args, **kwargs): + if len(args) == 0: + kwargs.update({"min": 0., "max": self.pi}) + super().__init__(*args, **kwargs) - Parameters - ---------- + def type_number(self): + return BONDED_IA_TABULATED_ANGLE - energy: array_like of :obj:`float` - The energy table for the range :math:`0-\\pi`. - force: array_like of :obj:`float` - The force table for the range :math:`0-\\pi`. + def type_name(self): + """Name of interaction type. """ + return "TABULATED_ANGLE" - _so_name = "Interactions::TabulatedAngleBond" - - pi = 3.14159265358979 + def validate_params(self, params): + """Check that parameters are valid. - def __init__(self, *args, **kwargs): - if len(args) == 0: - kwargs.update({"min": 0., "max": self.pi}) - super().__init__(*args, **kwargs) + """ + phi = [params["min"], params["max"]] + if abs(phi[0] - 0.) > 1e-5 or abs(phi[1] - self.pi) > 1e-5: + raise ValueError(f"Tabulated angle expects forces/energies " + f"within the range [0, pi], got {phi}") - def type_number(self): - return BONDED_IA_TABULATED_ANGLE + def get_default_params(self): + """Gets default values of optional parameters. - def type_name(self): - """Name of interaction type. + """ + return {} - """ - return "TABULATED_ANGLE" - def validate_params(self, params): - """Check that parameters are valid. +@script_interface_register +class TabulatedDihedral(BondedInteraction): - """ - phi = [params["min"], params["max"]] - if abs(phi[0] - 0.) > 1e-5 or abs(phi[1] - self.pi) > 1e-5: - raise ValueError(f"Tabulated angle expects forces/energies " - f"within the range [0, pi], got {phi}") + """ + Tabulated bond dihedral. - @script_interface_register - class TabulatedDihedral(_TabulatedBase): + Parameters + ---------- - """ - Tabulated bond dihedral. + energy: array_like of :obj:`float` + The energy table for the range :math:`0-2\\pi`. + force: array_like of :obj:`float` + The force table for the range :math:`0-2\\pi`. - Parameters - ---------- + """ - energy: array_like of :obj:`float` - The energy table for the range :math:`0-2\\pi`. - force: array_like of :obj:`float` - The force table for the range :math:`0-2\\pi`. + _so_name = "Interactions::TabulatedDihedralBond" + _so_feature = "TABULATED" - """ + pi = 3.14159265358979 - _so_name = "Interactions::TabulatedDihedralBond" + def __init__(self, *args, **kwargs): + if len(args) == 0: + kwargs.update({"min": 0., "max": 2. * self.pi}) + super().__init__(*args, **kwargs) - pi = 3.14159265358979 + def type_number(self): + return BONDED_IA_TABULATED_DIHEDRAL - def __init__(self, *args, **kwargs): - if len(args) == 0: - kwargs.update({"min": 0., "max": 2. * self.pi}) - super().__init__(*args, **kwargs) + def type_name(self): + """Name of interaction type. - def type_number(self): - return BONDED_IA_TABULATED_DIHEDRAL + """ + return "TABULATED_DIHEDRAL" - def type_name(self): - """Name of interaction type. + def validate_params(self, params): + """Check that parameters are valid. - """ - return "TABULATED_DIHEDRAL" + """ + phi = [params["min"], params["max"]] + if abs(phi[0] - 0.) > 1e-5 or abs(phi[1] - 2 * self.pi) > 1e-5: + raise ValueError(f"Tabulated dihedral expects forces/energies " + f"within the range [0, 2*pi], got {phi}") - def validate_params(self, params): - """Check that parameters are valid. + def get_default_params(self): + """Gets default values of optional parameters. - """ - phi = [params["min"], params["max"]] - if abs(phi[0] - 0.) > 1e-5 or abs(phi[1] - 2 * self.pi) > 1e-5: - raise ValueError(f"Tabulated dihedral expects forces/energies " - f"within the range [0, 2*pi], got {phi}") + """ + return {} @script_interface_register @@ -1924,8 +1747,13 @@ class QuarticBond(BondedInteraction): bonded_interaction_classes = { int(BONDED_IA_FENE): FeneBond, int(BONDED_IA_HARMONIC): HarmonicBond, + int(BONDED_IA_BONDED_COULOMB): BondedCoulomb, + int(BONDED_IA_BONDED_COULOMB_SR): BondedCoulombSRBond, int(BONDED_IA_RIGID_BOND): RigidBond, int(BONDED_IA_DIHEDRAL): Dihedral, + int(BONDED_IA_TABULATED_DISTANCE): TabulatedDistance, + int(BONDED_IA_TABULATED_ANGLE): TabulatedAngle, + int(BONDED_IA_TABULATED_DIHEDRAL): TabulatedDihedral, int(BONDED_IA_VIRTUAL_BOND): Virtual, int(BONDED_IA_ANGLE_HARMONIC): AngleHarmonic, int(BONDED_IA_ANGLE_COSINE): AngleCosine, @@ -1938,17 +1766,6 @@ bonded_interaction_classes = { int(BONDED_IA_THERMALIZED_DIST): ThermalizedBond, int(BONDED_IA_QUARTIC): QuarticBond, } -IF ELECTROSTATICS: - bonded_interaction_classes[int(BONDED_IA_BONDED_COULOMB)] = BondedCoulomb - bonded_interaction_classes[ - int(BONDED_IA_BONDED_COULOMB_SR)] = BondedCoulombSRBond -IF TABULATED: - bonded_interaction_classes[ - int(BONDED_IA_TABULATED_DISTANCE)] = TabulatedDistance - bonded_interaction_classes[ - int(BONDED_IA_TABULATED_ANGLE)] = TabulatedAngle - bonded_interaction_classes[ - int(BONDED_IA_TABULATED_DIHEDRAL)] = TabulatedDihedral def get_bonded_interaction_type_from_es_core(bond_id): diff --git a/src/script_interface/interactions/NonBondedInteraction.hpp b/src/script_interface/interactions/NonBondedInteraction.hpp index effd77a3d26..6ed08ce4906 100644 --- a/src/script_interface/interactions/NonBondedInteraction.hpp +++ b/src/script_interface/interactions/NonBondedInteraction.hpp @@ -32,10 +32,13 @@ #include "core/event.hpp" #include "core/nonbonded_interactions/nonbonded_interaction_data.hpp" +#include #include #include #include +#include #include +#include #include #include #include @@ -55,6 +58,7 @@ class InteractionPotentialInterface protected: using AutoParameters>::context; + using AutoParameters>::valid_parameters; /** @brief Managed object. */ std::shared_ptr m_ia_si; /** @brief Pointer to the corresponding member in a handle. */ @@ -72,12 +76,39 @@ class InteractionPotentialInterface [this, ptr]() { return m_ia_si.get()->*ptr; }}; } +private: + std::set get_valid_parameters() const { + auto const vec = valid_parameters(); + auto valid_keys = std::set(); + std::transform(vec.begin(), vec.end(), + std::inserter(valid_keys, valid_keys.begin()), + [](auto const &key) { return std::string{key}; }); + return valid_keys; + } + + void check_valid_parameters(VariantMap const ¶ms) const { + auto const valid_keys = get_valid_parameters(); + for (auto const &key : valid_keys) { + if (params.count(std::string(key)) == 0) { + throw std::runtime_error("Parameter '" + key + "' is missing"); + } + } + for (auto const &kv : params) { + if (valid_keys.count(kv.first) == 0) { + throw std::runtime_error("Parameter '" + kv.first + + "' is not recognized"); + } + } + } + public: Variant do_call_method(std::string const &name, VariantMap const ¶ms) override { if (name == "set_params") { - context()->parallel_try_catch( - [this, ¶ms]() { make_new_instance(params); }); + context()->parallel_try_catch([this, ¶ms]() { + check_valid_parameters(params); + make_new_instance(params); + }); if (m_types[0] != -1) { copy_si_to_core(); on_non_bonded_ia_change(); @@ -92,6 +123,9 @@ class InteractionPotentialInterface } return {}; } + if (name == "is_registered") { + return m_types[0] != -1; + } if (name == "bind_types") { auto types = get_value>(params, "_types"); if (types[0] > types[1]) { @@ -124,8 +158,10 @@ class InteractionPotentialInterface inactive_cutoff()) < 1e-9) { m_ia_si = std::make_shared(); } else { - context()->parallel_try_catch( - [this, ¶ms]() { make_new_instance(params); }); + context()->parallel_try_catch([this, ¶ms]() { + check_valid_parameters(params); + make_new_instance(params); + }); } } } @@ -160,6 +196,7 @@ class InteractionWCA : public InteractionPotentialInterface<::WCA_Parameters> { }); } +private: std::string inactive_parameter() const override { return "sigma"; } double inactive_cutoff() const override { return 0.; } @@ -168,6 +205,7 @@ class InteractionWCA : public InteractionPotentialInterface<::WCA_Parameters> { params, "epsilon", "sigma"); } +public: Variant do_call_method(std::string const &name, VariantMap const ¶ms) override { if (name == "get_cutoff") { @@ -198,6 +236,7 @@ class InteractionLJ : public InteractionPotentialInterface<::LJ_Parameters> { }); } +private: void make_new_instance(VariantMap const ¶ms) override { auto new_params = params; auto const *shift_string = boost::get(¶ms.at("shift")); @@ -245,6 +284,7 @@ class InteractionLJGen }); } +private: void make_new_instance(VariantMap const ¶ms) override { auto new_params = params; auto const *shift_string = boost::get(¶ms.at("shift")); @@ -291,6 +331,7 @@ class InteractionLJcos }); } +private: void make_new_instance(VariantMap const ¶ms) override { m_ia_si = make_shared_from_args( @@ -317,6 +358,7 @@ class InteractionLJcos2 }); } +private: std::string inactive_parameter() const override { return "sigma"; } double inactive_cutoff() const override { return 0.; } @@ -326,6 +368,7 @@ class InteractionLJcos2 params, "epsilon", "sigma", "offset", "width"); } +public: Variant do_call_method(std::string const &name, VariantMap const ¶ms) override { if (name == "get_cutoff") { @@ -353,6 +396,7 @@ class InteractionHertzian }); } +private: std::string inactive_parameter() const override { return "sig"; } void make_new_instance(VariantMap const ¶ms) override { @@ -379,6 +423,7 @@ class InteractionGaussian }); } +private: void make_new_instance(VariantMap const ¶ms) override { m_ia_si = make_shared_from_args( params, "eps", "sig", "cutoff"); @@ -406,6 +451,7 @@ class InteractionBMHTF }); } +private: void make_new_instance(VariantMap const ¶ms) override { m_ia_si = make_shared_from_args( @@ -432,6 +478,7 @@ class InteractionMorse }); } +private: void make_new_instance(VariantMap const ¶ms) override { m_ia_si = make_shared_from_args( @@ -461,6 +508,7 @@ class InteractionBuckingham }); } +private: void make_new_instance(VariantMap const ¶ms) override { m_ia_si = make_shared_from_args( @@ -487,6 +535,7 @@ class InteractionSoftSphere }); } +private: void make_new_instance(VariantMap const ¶ms) override { m_ia_si = make_shared_from_args( @@ -510,6 +559,7 @@ class InteractionHat : public InteractionPotentialInterface<::Hat_Parameters> { }); } +private: void make_new_instance(VariantMap const ¶ms) override { m_ia_si = make_shared_from_args( params, "F_max", "cutoff"); @@ -538,6 +588,7 @@ class InteractionGayBerne }); } +private: std::string inactive_parameter() const override { return "cut"; } void make_new_instance(VariantMap const ¶ms) override { @@ -566,6 +617,7 @@ class InteractionTabulated }); } +private: std::string inactive_parameter() const override { return "max"; } void make_new_instance(VariantMap const ¶ms) override { @@ -574,6 +626,7 @@ class InteractionTabulated params, "min", "max", "force", "energy"); } +public: Variant do_call_method(std::string const &name, VariantMap const ¶ms) override { if (name == "get_cutoff") { @@ -613,6 +666,7 @@ class InteractionDPD : public InteractionPotentialInterface<::DPD_Parameters> { std::ignore = get_ptr_offset(); // for code coverage } +private: std::string inactive_parameter() const override { return "r_cut"; } void make_new_instance(VariantMap const ¶ms) override { @@ -620,6 +674,14 @@ class InteractionDPD : public InteractionPotentialInterface<::DPD_Parameters> { double, double, double, double>( params, "gamma", "k", "r_cut", "weight_function", "trans_gamma", "trans_r_cut", "trans_weight_function"); + if (m_ia_si->radial.wf != 0 and m_ia_si->radial.wf != 1) { + throw std::domain_error( + "DPDInteraction parameter 'weight_function' must be 0 or 1"); + } + if (m_ia_si->trans.wf != 0 and m_ia_si->trans.wf != 1) { + throw std::domain_error( + "DPDInteraction parameter 'trans_weight_function' must be 0 or 1"); + } } }; #endif // DPD @@ -640,6 +702,7 @@ class InteractionThole }); } +private: std::string inactive_parameter() const override { return "scaling_coeff"; } double inactive_cutoff() const override { return 0.; } @@ -669,6 +732,8 @@ class InteractionSmoothStep make_autoparameter(&CoreInteraction::k0, "k0"), }); } + +private: void make_new_instance(VariantMap const ¶ms) override { m_ia_si = make_shared_from_args( diff --git a/src/script_interface/interactions/NonBondedInteractions.hpp b/src/script_interface/interactions/NonBondedInteractions.hpp index 9680d37419d..945ac1be2f1 100644 --- a/src/script_interface/interactions/NonBondedInteractions.hpp +++ b/src/script_interface/interactions/NonBondedInteractions.hpp @@ -69,17 +69,7 @@ class NonBondedInteractions : public ObjectHandle { void do_construct(VariantMap const ¶ms) override { auto const size = ::max_seen_particle_type; - { - // when reloading from a checkpoint file, need to resize IA lists - auto const new_size = ::max_seen_particle_type; - auto const n_pairs = new_size * (new_size + 1) / 2; - ::nonbonded_ia_params.resize(n_pairs); - for (auto &ia_params : ::nonbonded_ia_params) { - if (ia_params == nullptr) { - ia_params = std::make_shared<::IA_parameters>(); - } - } - } + make_particle_type_exist_local(size); for (int i = 0; i < size; i++) { for (int j = i; j < size; j++) { auto const key = Utils::upper_triangular(i, j, size); @@ -99,6 +89,7 @@ class NonBondedInteractions : public ObjectHandle { } if (name == "insert") { auto const types = get_value>(params.at("key")); + make_particle_type_exist_local(std::max(types[0], types[1])); auto const key = get_ia_param_key(std::min(types[0], types[1]), std::max(types[0], types[1])); auto obj_ptr = get_value>( diff --git a/testsuite/python/dpd.py b/testsuite/python/dpd.py index 7d8ce3b6641..64f9dfc003a 100644 --- a/testsuite/python/dpd.py +++ b/testsuite/python/dpd.py @@ -59,7 +59,7 @@ def reset_particles(): gamma = 1.5 # No seed should throw exception - with self.assertRaises(ValueError): + with self.assertRaisesRegex(ValueError, "A seed has to be given as keyword argument on first activation of the thermostat"): system.thermostat.set_dpd(kT=kT) system.thermostat.set_dpd(kT=kT, seed=41) diff --git a/testsuite/python/interactions_bonded_interface.py b/testsuite/python/interactions_bonded_interface.py index c5b0c256232..240b9d551cc 100644 --- a/testsuite/python/interactions_bonded_interface.py +++ b/testsuite/python/interactions_bonded_interface.py @@ -22,6 +22,7 @@ import espressomd import espressomd.interactions +import espressomd.code_features import numpy as np @@ -336,6 +337,18 @@ def test_exceptions(self): self.assertEqual(len(self.system.bonded_inter), 2) self.system.bonded_inter.clear() + def test_feature_checks(self): + base_class = espressomd.interactions.BondedInteraction + FeaturesError = espressomd.code_features.FeaturesError + for class_bond in espressomd.interactions.__dict__.values(): + if isinstance(class_bond, type) and issubclass( + class_bond, base_class) and class_bond != base_class: + feature = getattr(class_bond, "_so_feature", None) + if feature is not None and not espressomd.code_features.has_features( + feature): + with self.assertRaisesRegex(FeaturesError, f"Missing features {feature}"): + class_bond() + if __name__ == "__main__": ut.main() diff --git a/testsuite/python/interactions_non-bonded_interface.py b/testsuite/python/interactions_non-bonded_interface.py index 793dac07572..1bdaf1f9a6c 100644 --- a/testsuite/python/interactions_non-bonded_interface.py +++ b/testsuite/python/interactions_non-bonded_interface.py @@ -22,6 +22,7 @@ import espressomd import espressomd.interactions +import espressomd.code_features class Test(ut.TestCase): @@ -38,40 +39,6 @@ def intersMatch(self, inType, outInter, inParams, outParams, msg_long): self.assertIsInstance(outInter, inType) tests_common.assert_params_match(self, inParams, outParams, msg_long) - def parameterKeys(self, interObject): - """ - Check :meth:`~espressomd.interactions.NonBondedInteraction.valid_keys` - and :meth:`~espressomd.interactions.NonBondedInteraction.required_keys` - return sets, and that - :meth:`~espressomd.interactions.NonBondedInteraction.default_params` - returns a dictionary with the correct keys. - - Parameters - ---------- - interObject: instance of a class derived from :class:`espressomd.interactions.NonBondedInteraction` - Object of the interaction to test, e.g. - :class:`~espressomd.interactions.LennardJonesInteraction` - """ - classname = interObject.__class__.__name__ - valid_keys = interObject.valid_keys() - required_keys = interObject.required_keys() - default_keys = set(interObject.default_params().keys()) - self.assertIsInstance(valid_keys, set, - "{}.valid_keys() must return a set".format( - classname)) - self.assertIsInstance(required_keys, set, - "{}.required_keys() must return a set".format( - classname)) - self.assertTrue(default_keys.issubset(valid_keys), - "{}.default_params() has unknown parameters: {}".format( - classname, default_keys.difference(valid_keys))) - self.assertTrue(default_keys.isdisjoint(required_keys), - "{}.default_params() has extra parameters: {}".format( - classname, default_keys.intersection(required_keys))) - self.assertSetEqual(default_keys, valid_keys - required_keys, - "{}.default_params() should have keys: {}, got: {}".format( - classname, valid_keys - required_keys, default_keys)) - def generateTestForNon_bonded_interaction( _partType1, _partType2, _interClass, _params, _interName): """Generates test cases for checking interaction parameters set and @@ -115,7 +82,6 @@ def func(self): f"{interClass.__name__}: value set and value gotten back " f"differ for particle types {partType1} and {partType2}: " f"{params} vs. {outParams}") - self.parameterKeys(outInter) return func @@ -204,16 +170,16 @@ def test_set_params(self): @utx.skipIfMissingFeatures("LENNARD_JONES") def test_exceptions(self): - with self.assertRaisesRegex(RuntimeError, "LennardJonesInteraction parameter 'shift' is missing"): + with self.assertRaisesRegex(RuntimeError, "Parameter 'shift' is missing"): espressomd.interactions.LennardJonesInteraction( epsilon=1., sigma=2., cutoff=2.) - with self.assertRaisesRegex(RuntimeError, "LennardJonesInteraction parameter 'shift' is missing"): + with self.assertRaisesRegex(RuntimeError, "Parameter 'shift' is missing"): self.system.non_bonded_inter[0, 0].lennard_jones.set_params( epsilon=1., sigma=2., cutoff=2.) - with self.assertRaisesRegex(RuntimeError, "Parameter 'unknown' is not a valid LennardJonesInteraction parameter"): + with self.assertRaisesRegex(RuntimeError, "Parameter 'unknown' is not recognized"): espressomd.interactions.LennardJonesInteraction( epsilon=1., sigma=2., cutoff=3., shift=4., unknown=5.) - with self.assertRaisesRegex(RuntimeError, "Parameter 'unknown' is not a valid LennardJonesInteraction parameter"): + with self.assertRaisesRegex(RuntimeError, "Parameter 'unknown' is not recognized"): self.system.non_bonded_inter[0, 0].lennard_jones.set_params( epsilon=1., sigma=2., cutoff=3., shift=4., unknown=5.) with self.assertRaisesRegex(ValueError, "LJ parameter 'shift' has to be 'auto' or a float"): @@ -236,11 +202,23 @@ def test_exceptions(self): self.assertEqual(inter_00_obj.call_method("get_types"), [0, 0]) self.assertIsNone(inter_00_obj.call_method("unknown")) - def check_potential_exceptions(self, ia_class, params, check_keys): + def test_feature_checks(self): + base_class = espressomd.interactions.NonBondedInteraction + FeaturesError = espressomd.code_features.FeaturesError + for class_ia in espressomd.interactions.__dict__.values(): + if isinstance(class_ia, type) and issubclass( + class_ia, base_class) and class_ia != base_class: + feature = class_ia._so_feature + if not espressomd.code_features.has_features(feature): + with self.assertRaisesRegex(FeaturesError, f"Missing features {feature}"): + class_ia() + + def check_potential_exceptions( + self, ia_class, params, check_keys, invalid_value=-0.1): for key in check_keys: with self.assertRaisesRegex(ValueError, f"parameter '{key}'"): invalid_params = params.copy() - invalid_params[key] = -0.1 + invalid_params[key] = invalid_value ia_class(**invalid_params) @utx.skipIfMissingFeatures("WCA") @@ -358,6 +336,16 @@ def test_smooth_step_exceptions(self): ("eps", "sig", "cutoff") ) + @utx.skipIfMissingFeatures("DPD") + def test_dpd_exceptions(self): + self.check_potential_exceptions( + espressomd.interactions.DPDInteraction, + {"weight_function": 1, "gamma": 2., "r_cut": 2., + "trans_weight_function": 1, "trans_gamma": 1., "trans_r_cut": 2.}, + ("weight_function", "trans_weight_function"), + invalid_value=-1 + ) + if __name__ == "__main__": ut.main() diff --git a/testsuite/python/save_checkpoint.py b/testsuite/python/save_checkpoint.py index 1dd3047feef..c6710cac1c1 100644 --- a/testsuite/python/save_checkpoint.py +++ b/testsuite/python/save_checkpoint.py @@ -233,6 +233,15 @@ epsilon=1.2, sigma=1.7, cutoff=2.0, shift=0.1) system.non_bonded_inter[1, 7].lennard_jones.set_params( epsilon=1.2e5, sigma=1.7, cutoff=2.0, shift=0.1) +if espressomd.has_features(['DPD']): + dpd_params = {"weight_function": 1, "gamma": 2., "trans_r_cut": 2., "k": 2., + "trans_weight_function": 0, "trans_gamma": 1., "r_cut": 2.} + dpd_ia = espressomd.interactions.DPDInteraction(**dpd_params) + handle_ia = espressomd.interactions.NonBondedInteractionHandle( + _types=(0, 0)) + checkpoint.register("dpd_ia") + checkpoint.register("dpd_params") + checkpoint.register("handle_ia") # bonded interactions harmonic_bond = espressomd.interactions.HarmonicBond(r_0=0.0, k=1.0) diff --git a/testsuite/python/tabulated.py b/testsuite/python/tabulated.py index 198fdd34cc1..73304f865af 100644 --- a/testsuite/python/tabulated.py +++ b/testsuite/python/tabulated.py @@ -59,6 +59,8 @@ def check(self): def test_non_bonded(self): self.system.non_bonded_inter[0, 0].tabulated.set_params( min=self.min_, max=self.max_, energy=self.energy, force=self.force) + self.assertEqual( + self.system.non_bonded_inter[0, 0].tabulated.cutoff, self.max_) params = self.system.non_bonded_inter[0, 0].tabulated.get_params() np.testing.assert_allclose(params['force'], self.force) @@ -70,6 +72,8 @@ def test_non_bonded(self): self.system.non_bonded_inter[0, 0].tabulated.set_params( min=-1, max=-1, energy=[], force=[]) + self.assertEqual( + self.system.non_bonded_inter[0, 0].tabulated.cutoff, -1.) with self.assertRaisesRegex(ValueError, "TabulatedPotential parameter 'max' must be larger than or equal to parameter 'min'"): espressomd.interactions.TabulatedNonBonded( diff --git a/testsuite/python/test_checkpoint.py b/testsuite/python/test_checkpoint.py index e9182c9bc7d..2d5b874722b 100644 --- a/testsuite/python/test_checkpoint.py +++ b/testsuite/python/test_checkpoint.py @@ -330,7 +330,9 @@ def test_integrator_SDM(self): @utx.skipIfMissingFeatures('LENNARD_JONES') @ut.skipIf('LJ' not in modes, "Skipping test due to missing mode.") - def test_non_bonded_inter(self): + def test_non_bonded_inter_lj(self): + self.assertTrue( + system.non_bonded_inter[0, 0].lennard_jones.call_method("is_registered")) params1 = system.non_bonded_inter[0, 0].lennard_jones.get_params() params2 = system.non_bonded_inter[3, 0].lennard_jones.get_params() reference1 = {'shift': 0.1, 'sigma': 1.3, 'epsilon': 1.2, @@ -339,6 +341,13 @@ def test_non_bonded_inter(self): 'cutoff': 2.0, 'offset': 0.0, 'min': 0.0} self.assertEqual(params1, reference1) self.assertEqual(params2, reference2) + self.assertTrue(handle_ia.lennard_jones.call_method("is_registered")) + self.assertEqual(handle_ia.lennard_jones.get_params(), reference1) + + @utx.skipIfMissingFeatures('DPD') + def test_non_bonded_inter_dpd(self): + self.assertEqual(dpd_ia.get_params(), dpd_params) + self.assertFalse(dpd_ia.call_method("is_registered")) def test_bonded_inter(self): # check the ObjectHandle was correctly initialized (including MPI)