From 8dd2a640d3f4a2f58cc4fc5aa2347e3a66b67d3d Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Mon, 1 Mar 2021 13:47:44 +0100 Subject: [PATCH 1/7] core: Remove P3M member additional_mesh This member is a duplicate of the ELC member space_layer. --- src/core/electrostatics_magnetostatics/elc.cpp | 9 --------- src/core/electrostatics_magnetostatics/p3m-common.cpp | 7 +++++-- src/core/electrostatics_magnetostatics/p3m-common.hpp | 8 +++----- src/core/electrostatics_magnetostatics/p3m-dipolar.cpp | 2 +- src/core/electrostatics_magnetostatics/p3m.cpp | 8 +++++++- src/python/espressomd/magnetostatics.pyx | 2 +- src/python/espressomd/p3m_common.pxd | 1 - 7 files changed, 17 insertions(+), 20 deletions(-) diff --git a/src/core/electrostatics_magnetostatics/elc.cpp b/src/core/electrostatics_magnetostatics/elc.cpp index 7767f432fd3..41ded584289 100644 --- a/src/core/electrostatics_magnetostatics/elc.cpp +++ b/src/core/electrostatics_magnetostatics/elc.cpp @@ -1081,15 +1081,6 @@ void ELC_init() { runtimeErrorMsg() << "ELC auto-retuning failed"; } } - if (elc_params.dielectric_contrast_on) { - p3m.params.additional_mesh[0] = 0; - p3m.params.additional_mesh[1] = 0; - p3m.params.additional_mesh[2] = elc_params.space_layer; - } else { - p3m.params.additional_mesh[0] = 0; - p3m.params.additional_mesh[1] = 0; - p3m.params.additional_mesh[2] = 0; - } } int ELC_set_params(double maxPWerror, double gap_size, double far_cut, diff --git a/src/core/electrostatics_magnetostatics/p3m-common.cpp b/src/core/electrostatics_magnetostatics/p3m-common.cpp index 0bf10fc1a61..abbfffb0e58 100644 --- a/src/core/electrostatics_magnetostatics/p3m-common.cpp +++ b/src/core/electrostatics_magnetostatics/p3m-common.cpp @@ -115,14 +115,17 @@ double p3m_analytic_cotangent_sum(int n, double mesh_i, int cao) { void p3m_calc_local_ca_mesh(p3m_local_mesh &local_mesh, const P3MParameters ¶ms, - const LocalBox &local_geo, double skin) { + const LocalBox &local_geo, double skin, + double space_layer) { int i; int ind[3]; /* total skin size */ double full_skin[3]; for (i = 0; i < 3; i++) - full_skin[i] = params.cao_cut[i] + skin + params.additional_mesh[i]; + full_skin[i] = params.cao_cut[i] + skin; + + full_skin[2] += space_layer; /* inner left down grid point (global index) */ for (i = 0; i < 3; i++) diff --git a/src/core/electrostatics_magnetostatics/p3m-common.hpp b/src/core/electrostatics_magnetostatics/p3m-common.hpp index 7e002994c95..26e9fbaf3e5 100644 --- a/src/core/electrostatics_magnetostatics/p3m-common.hpp +++ b/src/core/electrostatics_magnetostatics/p3m-common.hpp @@ -147,14 +147,11 @@ typedef struct { /** number of points unto which a single charge is interpolated, i.e. * p3m.cao^3 */ int cao3 = 0; - /** additional points around the charge assignment mesh, for method like - * dielectric ELC creating virtual charges. */ - double additional_mesh[3] = {}; template void serialize(Archive &ar, long int) { ar &tuning &alpha_L &r_cut_iL &mesh; ar &mesh_off &cao &accuracy &epsilon &cao_cut; - ar &a &ai &alpha &r_cut &cao3 &additional_mesh; + ar &a &ai &alpha &r_cut &cao3; } } P3MParameters; @@ -186,7 +183,8 @@ double p3m_analytic_cotangent_sum(int n, double mesh_i, int cao); */ void p3m_calc_local_ca_mesh(p3m_local_mesh &local_mesh, const P3MParameters ¶ms, - const LocalBox &local_geo, double skin); + const LocalBox &local_geo, double skin, + double space_layer); /** Calculate the spatial position of the left down mesh * point of the local mesh, to be stored in diff --git a/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp b/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp index e1ac429d994..d2fe0cef479 100644 --- a/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp +++ b/src/core/electrostatics_magnetostatics/p3m-dipolar.cpp @@ -200,7 +200,7 @@ void dp3m_init() { * and the cutoff for charge assignment dp3m.params.cao_cut */ dp3m_init_a_ai_cao_cut(); - p3m_calc_local_ca_mesh(dp3m.local_mesh, dp3m.params, local_geo, skin); + p3m_calc_local_ca_mesh(dp3m.local_mesh, dp3m.params, local_geo, skin, 0.0); dp3m.sm.resize(comm_cart, dp3m.local_mesh); diff --git a/src/core/electrostatics_magnetostatics/p3m.cpp b/src/core/electrostatics_magnetostatics/p3m.cpp index 41ae1f8f3ad..58eb290e686 100644 --- a/src/core/electrostatics_magnetostatics/p3m.cpp +++ b/src/core/electrostatics_magnetostatics/p3m.cpp @@ -185,7 +185,13 @@ void p3m_init() { return; } - p3m_calc_local_ca_mesh(p3m.local_mesh, p3m.params, local_geo, skin); + double elc_layer = 0.0; + if (coulomb.method == COULOMB_ELC_P3M) { + elc_layer = elc_params.space_layer; + } + + p3m_calc_local_ca_mesh(p3m.local_mesh, p3m.params, local_geo, skin, + elc_layer); p3m.sm.resize(comm_cart, p3m.local_mesh); diff --git a/src/python/espressomd/magnetostatics.pyx b/src/python/espressomd/magnetostatics.pyx index 66a913d2828..5da552948aa 100644 --- a/src/python/espressomd/magnetostatics.pyx +++ b/src/python/espressomd/magnetostatics.pyx @@ -120,7 +120,7 @@ IF DP3M == 1: def valid_keys(self): return ["prefactor", "alpha_L", "r_cut_iL", "mesh", "mesh_off", "cao", "accuracy", "epsilon", "cao_cut", "a", "ai", - "alpha", "r_cut", "cao3", "additional_mesh", "tune", "verbose"] + "alpha", "r_cut", "cao3", "tune", "verbose"] def required_keys(self): return ["accuracy", ] diff --git a/src/python/espressomd/p3m_common.pxd b/src/python/espressomd/p3m_common.pxd index 3f7b685038c..aaf7f92b4d9 100644 --- a/src/python/espressomd/p3m_common.pxd +++ b/src/python/espressomd/p3m_common.pxd @@ -30,4 +30,3 @@ IF P3M == 1 or DP3M == 1: double a[3] double alpha double r_cut - double additional_mesh[3] From 0e527cebfe1acda7887a06be92c7a8331fa86d03 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Mon, 1 Mar 2021 14:29:59 +0100 Subject: [PATCH 2/7] python: Make ELC an actor via composition --- doc/sphinx/electrostatics.rst | 15 +- samples/visualization_elc.py | 8 +- .../espressomd/electrostatic_extensions.pxd | 18 --- .../espressomd/electrostatic_extensions.pyx | 118 --------------- src/python/espressomd/electrostatics.pxd | 18 +++ src/python/espressomd/electrostatics.pyx | 135 ++++++++++++++++++ testsuite/python/coulomb_mixed_periodicity.py | 11 +- testsuite/python/elc.py | 42 +++--- testsuite/python/elc_vs_analytic.py | 12 +- 9 files changed, 192 insertions(+), 185 deletions(-) diff --git a/doc/sphinx/electrostatics.rst b/doc/sphinx/electrostatics.rst index d39e7db170e..7c6061ab50c 100644 --- a/doc/sphinx/electrostatics.rst +++ b/doc/sphinx/electrostatics.rst @@ -234,7 +234,7 @@ using it. Electrostatic Layer Correction (ELC) ------------------------------------ -:class:`espressomd.electrostatic_extensions.ELC` +:class:`espressomd.electrostatics.ELC` *ELC* is an extension of the P3M electrostatics solver for explicit 2D periodic systems. It can account for different dielectric jumps on both sides of the @@ -260,8 +260,9 @@ Usage notes: *ELC* is an |es| actor and is used with:: - import espressomd.electrostatic_extensions - elc = electrostatic_extensions.ELC(gap_size=box_l * 0.2, maxPWerror=1e-3) + import espressomd.electrostatics + p3m = espressomd.electrostatics.P3M(prefactor=1, accuracy=1e-4) + elc = espressomd.electrostatics.ELC(p3m_actor=p3m, gap_size=box_l * 0.2, maxPWerror=1e-3) system.actors.add(elc) *ELC* can also be used to simulate 2D periodic systems with image charges, @@ -273,8 +274,8 @@ simulation region (*middle*) to *bottom* (at :math:`z=0`) and from *middle* to are :math:`\Delta_t=\frac{\varepsilon_m-\varepsilon_t}{\varepsilon_m+\varepsilon_t}` and :math:`\Delta_b=\frac{\varepsilon_m-\varepsilon_b}{\varepsilon_m+\varepsilon_b}`:: - elc = electrostatic_extensions.ELC(gap_size=box_l * 0.2, maxPWerror=1e-3, - delta_mid_top=0.9, delta_mid_bot=0.1) + elc = espressomd.electrostatics.ELC(p3m_actor=p3m, gap_size=box_l * 0.2, maxPWerror=1e-3, + delta_mid_top=0.9, delta_mid_bot=0.1) The fully metallic case :math:`\Delta_t=\Delta_b=-1` would lead to divergence of the forces/energies in *ELC* and is therefore only possible with the @@ -283,8 +284,8 @@ of the forces/energies in *ELC* and is therefore only possible with the Toggle ``const_pot`` on to maintain a constant electric potential difference ``pot_diff`` between the xy-planes at :math:`z=0` and :math:`z = L_z - h`:: - elc = electrostatic_extensions.ELC(gap_size=box_l * 0.2, maxPWerror=1e-3, - const_pot=True, delta_mid_bot=100.0) + elc = espressomd.electrostatics.ELC(p3m_actor=p3m, gap_size=box_l * 0.2, maxPWerror=1e-3, + const_pot=True, delta_mid_bot=100.0) This is done by countering the total dipole moment of the system with the electric field :math:`E_{\textrm{induced}}` and superposing a homogeneous diff --git a/samples/visualization_elc.py b/samples/visualization_elc.py index c8c67c73aea..950bc8b3708 100644 --- a/samples/visualization_elc.py +++ b/samples/visualization_elc.py @@ -27,7 +27,6 @@ import espressomd import espressomd.shapes from espressomd import electrostatics -from espressomd import electrostatic_extensions from espressomd import visualization required_features = ["P3M", "WCA"] @@ -75,11 +74,8 @@ system.thermostat.set_langevin(kT=0.1, gamma=1.0, seed=42) p3m = electrostatics.P3M(prefactor=1.0, accuracy=1e-2) - -system.actors.add(p3m) - -elc = electrostatic_extensions.ELC(maxPWerror=1.0, gap_size=elc_gap, - const_pot=True, pot_diff=potential_diff) +elc = electrostatics.ELC(p3m_actor=p3m, maxPWerror=1.0, gap_size=elc_gap, + const_pot=True, pot_diff=potential_diff) system.actors.add(elc) visualizer.run(1) diff --git a/src/python/espressomd/electrostatic_extensions.pxd b/src/python/espressomd/electrostatic_extensions.pxd index cc35493a54f..3eb19803613 100644 --- a/src/python/espressomd/electrostatic_extensions.pxd +++ b/src/python/espressomd/electrostatic_extensions.pxd @@ -20,28 +20,10 @@ include "myconfig.pxi" from .electrostatics cimport * from libcpp.vector cimport vector -from libcpp cimport bool from .utils cimport Vector3d IF ELECTROSTATICS and P3M: - cdef extern from "electrostatics_magnetostatics/elc.hpp": - ctypedef struct ELC_struct: - double maxPWerror - double gap_size - double far_cut - bool neutralize - double delta_mid_top, - double delta_mid_bot, - bool const_pot, - double pot_diff - - int ELC_set_params(double maxPWerror, double min_dist, double far_cut, - bool neutralize, double delta_mid_top, double delta_mid_bot, bool const_pot, double pot_diff) - - # links intern C-struct with python object - ELC_struct elc_params - cdef extern from "electrostatics_magnetostatics/icc.hpp": ctypedef struct iccp3m_struct: int n_ic diff --git a/src/python/espressomd/electrostatic_extensions.pyx b/src/python/espressomd/electrostatic_extensions.pyx index 2a3c14490fc..b1eb2c03cfd 100644 --- a/src/python/espressomd/electrostatic_extensions.pyx +++ b/src/python/espressomd/electrostatic_extensions.pyx @@ -31,124 +31,6 @@ IF ELECTROSTATICS and P3M: cdef class ElectrostaticExtensions(actors.Actor): pass - cdef class ELC(ElectrostaticExtensions): - """ - Electrostatics solver for systems with two periodic dimensions. - See :ref:`Electrostatic Layer Correction (ELC)` for more details. - - Parameters - ---------- - gap_size : :obj:`float`, required - The gap size gives the height :math:`h` of the empty region between - the system box and the neighboring artificial images. |es| does not - make sure that the gap is actually empty, this is the user's - responsibility. The method will run even if the condition is not - fulfilled, however, the error bound will not be reached. Therefore - you should really make sure that the gap region is empty (e.g. - with wall constraints). - maxPWerror : :obj:`float`, required - The maximal pairwise error sets the least upper bound (LUB) error - of the force between any two charges without prefactors (see the - papers). The algorithm tries to find parameters to meet this LUB - requirements or will throw an error if there are none. - delta_mid_top : :obj:`float`, optional - Dielectric contrast :math:`\\Delta_t` between the upper boundary - and the simulation box. - delta_mid_bottom : :obj:`float`, optional - Dielectric contrast :math:`\\Delta_b` between the lower boundary - and the simulation box. - const_pot : :obj:`bool`, optional - Activate a constant electric potential between the top and bottom - of the simulation box. - pot_diff : :obj:`float`, optional - If ``const_pot`` is enabled, this parameter controls the applied - voltage between the boundaries of the simulation box in the - *z*-direction (at :math:`z = 0` and :math:`z = L_z - h`). - neutralize : :obj:`bool`, optional - By default, *ELC* just as P3M adds a homogeneous neutralizing - background to the system in case of a net charge. However, unlike - in three dimensions, this background adds a parabolic potential - across the slab :cite:`ballenegger09a`. Therefore, under normal - circumstances, you will probably want to disable the neutralization - for non-neutral systems. This corresponds then to a formal - regularization of the forces and energies :cite:`ballenegger09a`. - Also, if you add neutralizing walls explicitly as constraints, you - have to disable the neutralization. When using a dielectric - contrast or full metallic walls (``delta_mid_top != 0`` or - ``delta_mid_bot != 0`` or ``const_pot=True``), ``neutralize`` is - overwritten and switched off internally. Note that the special - case of non-neutral systems with a *non-metallic* dielectric jump - (e.g. ``delta_mid_top`` or ``delta_mid_bot`` in ``]-1,1[``) is not - covered by the algorithm and will throw an error. - far_cut : :obj:`float`, optional - Cutoff radius, use with care, intended for testing purposes. When - setting the cutoff directly, the maximal pairwise error is ignored. - """ - - def validate_params(self): - default_params = self.default_params() - check_type_or_throw_except( - self._params["maxPWerror"], 1, float, "") - check_range_or_except( - self._params, "maxPWerror", 0, False, "inf", True) - check_type_or_throw_except(self._params["gap_size"], 1, float, "") - check_range_or_except( - self._params, "gap_size", 0, False, "inf", True) - check_type_or_throw_except(self._params["far_cut"], 1, float, "") - check_type_or_throw_except( - self._params["neutralize"], 1, type(True), "") - - def valid_keys(self): - return ["maxPWerror", "gap_size", "far_cut", "neutralize", - "delta_mid_top", "delta_mid_bot", "const_pot", "pot_diff", - "check_neutrality"] - - def required_keys(self): - return ["maxPWerror", "gap_size"] - - def default_params(self): - return {"maxPWerror": -1, - "gap_size": -1, - "far_cut": -1, - "delta_mid_top": 0, - "delta_mid_bot": 0, - "const_pot": False, - "pot_diff": 0.0, - "neutralize": True, - "check_neutrality": True} - - def _get_params_from_es_core(self): - params = {} - params.update(elc_params) - return params - - def _set_params_in_es_core(self): - if coulomb.method == COULOMB_P3M_GPU: - raise Exception("ELC is not set up to work with the GPU P3M") - - if self._params["const_pot"]: - self._params["delta_mid_top"] = -1 - self._params["delta_mid_bot"] = -1 - - if ELC_set_params( - self._params["maxPWerror"], - self._params["gap_size"], - self._params["far_cut"], - self._params["neutralize"], - self._params["delta_mid_top"], - self._params["delta_mid_bot"], - self._params["const_pot"], - self._params["pot_diff"]): - handle_errors("ELC tuning failed") - - def _activate_method(self): - check_neutrality(self._params) - self._set_params_in_es_core() - - def _deactivate_method(self): - raise Exception( - "Unable to remove ELC as the state of the underlying electrostatics method will remain unclear.") - cdef class ICC(ElectrostaticExtensions): """ Interface to the induced charge calculation scheme for dielectric diff --git a/src/python/espressomd/electrostatics.pxd b/src/python/espressomd/electrostatics.pxd index 0f60447daf4..7f9d5e90651 100644 --- a/src/python/espressomd/electrostatics.pxd +++ b/src/python/espressomd/electrostatics.pxd @@ -83,6 +83,24 @@ IF ELECTROSTATICS: cdef extern from "electrostatics_magnetostatics/p3m_gpu.hpp": void p3m_gpu_init(int cao, int * mesh, double alpha) except + + cdef extern from "electrostatics_magnetostatics/elc.hpp": + ctypedef struct ELC_struct: + double maxPWerror + double gap_size + double far_cut + bool neutralize + double delta_mid_top + double delta_mid_bot + bool const_pot + double pot_diff + + int ELC_set_params(double maxPWerror, double min_dist, double far_cut, + bool neutralize, double delta_mid_top, + double delta_mid_bot, bool const_pot, double pot_diff) + + # links intern C-struct with python object + ELC_struct elc_params + cdef extern from "electrostatics_magnetostatics/debye_hueckel.hpp": ctypedef struct Debye_hueckel_params: double r_cut diff --git a/src/python/espressomd/electrostatics.pyx b/src/python/espressomd/electrostatics.pyx index 860525ca503..65f4f44a0ac 100644 --- a/src/python/espressomd/electrostatics.pyx +++ b/src/python/espressomd/electrostatics.pyx @@ -25,9 +25,11 @@ IF SCAFACOS == 1: from .scafacos import ScafacosConnector from . cimport scafacos from .utils import is_valid_type, check_type_or_throw_except, to_str, handle_errors +from .utils cimport check_range_or_except from . cimport checks from .analyze cimport partCfg, PartCfg from .particle_data cimport particle +import sys IF ELECTROSTATICS == 1: @@ -395,6 +397,139 @@ IF P3M == 1: super()._set_params_in_es_core() handle_errors("P3M: initialization failed") + cdef class ELC(ElectrostaticInteraction): + """ + Electrostatics solver for systems with two periodic dimensions. + See :ref:`Electrostatic Layer Correction (ELC)` for more details. + + Parameters + ---------- + p3m_actor : :obj:`P3M`, required + Base P3M actor. + gap_size : :obj:`float`, required + The gap size gives the height :math:`h` of the empty region between + the system box and the neighboring artificial images. |es| does not + make sure that the gap is actually empty, this is the user's + responsibility. The method will run even if the condition is not + fulfilled, however, the error bound will not be reached. Therefore + you should really make sure that the gap region is empty (e.g. + with wall constraints). + maxPWerror : :obj:`float`, required + The maximal pairwise error sets the least upper bound (LUB) error + of the force between any two charges without prefactors (see the + papers). The algorithm tries to find parameters to meet this LUB + requirements or will throw an error if there are none. + delta_mid_top : :obj:`float`, optional + Dielectric contrast :math:`\\Delta_t` between the upper boundary + and the simulation box. + delta_mid_bottom : :obj:`float`, optional + Dielectric contrast :math:`\\Delta_b` between the lower boundary + and the simulation box. + const_pot : :obj:`bool`, optional + Activate a constant electric potential between the top and bottom + of the simulation box. + pot_diff : :obj:`float`, optional + If ``const_pot`` is enabled, this parameter controls the applied + voltage between the boundaries of the simulation box in the + *z*-direction (at :math:`z = 0` and :math:`z = L_z - h`). + neutralize : :obj:`bool`, optional + By default, *ELC* just as P3M adds a homogeneous neutralizing + background to the system in case of a net charge. However, unlike + in three dimensions, this background adds a parabolic potential + across the slab :cite:`ballenegger09a`. Therefore, under normal + circumstances, you will probably want to disable the neutralization + for non-neutral systems. This corresponds then to a formal + regularization of the forces and energies :cite:`ballenegger09a`. + Also, if you add neutralizing walls explicitly as constraints, you + have to disable the neutralization. When using a dielectric + contrast or full metallic walls (``delta_mid_top != 0`` or + ``delta_mid_bot != 0`` or ``const_pot=True``), ``neutralize`` is + overwritten and switched off internally. Note that the special + case of non-neutral systems with a *non-metallic* dielectric jump + (e.g. ``delta_mid_top`` or ``delta_mid_bot`` in ``]-1,1[``) is not + covered by the algorithm and will throw an error. + far_cut : :obj:`float`, optional + Cutoff radius, use with care, intended for testing purposes. When + setting the cutoff directly, the maximal pairwise error is ignored. + """ + + def validate_params(self): + # P3M + if CUDA: + if isinstance(self._params["p3m_actor"], P3MGPU): + raise ValueError( + "ELC is not set up to work with the GPU P3M") + check_type_or_throw_except( + self._params["p3m_actor"], 1, getattr( + sys.modules[__name__], "P3M"), + "p3m_actor has to be a P3M solver") + self._params["p3m_actor"]._params["epsilon"] = 0.0 + self._params["p3m_actor"].validate_params() + # ELC + check_type_or_throw_except( + self._params["maxPWerror"], 1, float, "") + check_range_or_except( + self._params, "maxPWerror", 0, False, "inf", True) + check_type_or_throw_except(self._params["gap_size"], 1, float, "") + check_range_or_except( + self._params, "gap_size", 0, False, "inf", True) + check_type_or_throw_except(self._params["far_cut"], 1, float, "") + check_type_or_throw_except( + self._params["neutralize"], 1, type(True), "") + + def valid_keys(self): + return ["p3m_actor", "maxPWerror", "gap_size", "far_cut", + "neutralize", "delta_mid_top", "delta_mid_bot", + "const_pot", "pot_diff", "check_neutrality"] + + def required_keys(self): + return ["p3m_actor", "maxPWerror", "gap_size"] + + def default_params(self): + return {"maxPWerror": -1, + "gap_size": -1, + "far_cut": -1, + "delta_mid_top": 0, + "delta_mid_bot": 0, + "const_pot": False, + "pot_diff": 0.0, + "neutralize": True, + "check_neutrality": True} + + def _get_params_from_es_core(self): + params = {} + params.update(elc_params) + params["p3m_actor"] = self._params["p3m_actor"] + return params + + def _set_params_in_es_core(self): + self._params["p3m_actor"]._set_params_in_es_core() + if coulomb.method == COULOMB_P3M_GPU: + raise Exception("ELC is not set up to work with the GPU P3M") + + if self._params["const_pot"]: + self._params["delta_mid_top"] = -1 + self._params["delta_mid_bot"] = -1 + + if ELC_set_params( + self._params["maxPWerror"], + self._params["gap_size"], + self._params["far_cut"], + self._params["neutralize"], + self._params["delta_mid_top"], + self._params["delta_mid_bot"], + self._params["const_pot"], + self._params["pot_diff"]): + handle_errors("ELC tuning failed") + + def tune(self, **tune_params_subset): + self._params["p3m_actor"].tune(**tune_params_subset) + + def _activate_method(self): + self._params["p3m_actor"]._activate_method() + check_neutrality(self._params) + self._set_params_in_es_core() + IF ELECTROSTATICS: cdef class MMM1D(ElectrostaticInteraction): """ diff --git a/testsuite/python/coulomb_mixed_periodicity.py b/testsuite/python/coulomb_mixed_periodicity.py index c275d3a9f7e..30f84005926 100644 --- a/testsuite/python/coulomb_mixed_periodicity.py +++ b/testsuite/python/coulomb_mixed_periodicity.py @@ -20,7 +20,7 @@ import unittest_decorators as utx import numpy as np import espressomd -from espressomd import electrostatics, electrostatic_extensions, scafacos +from espressomd import electrostatics, scafacos import tests_common @@ -44,7 +44,6 @@ def setUp(self): self.S.box_l = (10, 10, 10) self.S.time_step = 0.01 self.S.cell_system.skin = 0. - self.S.actors.clear() data = np.genfromtxt(tests_common.abspath( "data/coulomb_mixed_periodicity_system.data")) @@ -61,6 +60,7 @@ def setUp(self): def tearDown(self): self.S.part.clear() + self.S.actors.clear() def compare(self, method_name, energy=True): # Compare forces and energy now in the system to stored ones @@ -84,7 +84,7 @@ def compare(self, method_name, energy=True): # Tests for individual methods @utx.skipIfMissingFeatures(["P3M"]) - def test_zz_p3mElc(self): + def test_elc(self): # Make sure, the data satisfies the gap for p in self.S.part: if p.pos[2] < 0 or p.pos[2] > 9.: @@ -97,13 +97,11 @@ def test_zz_p3mElc(self): self.S.box_l = (10, 10, 10) p3m = electrostatics.P3M(prefactor=1, accuracy=1e-6, mesh=(64, 64, 64)) + elc = electrostatics.ELC(p3m_actor=p3m, maxPWerror=1E-6, gap_size=1) - self.S.actors.add(p3m) - elc = electrostatic_extensions.ELC(maxPWerror=1E-6, gap_size=1) self.S.actors.add(elc) self.S.integrator.run(0) self.compare("elc", energy=True) - self.S.actors.remove(p3m) @ut.skipIf(not espressomd.has_features("SCAFACOS") or 'p2nfft' not in scafacos.available_methods(), @@ -125,7 +123,6 @@ def test_scafacos_p2nfft(self): self.S.actors.add(scafacos) self.S.integrator.run(0) self.compare("scafacos_p2nfft", energy=True) - self.S.actors.remove(scafacos) if __name__ == "__main__": diff --git a/testsuite/python/elc.py b/testsuite/python/elc.py index ab058bccdce..32fe7754935 100644 --- a/testsuite/python/elc.py +++ b/testsuite/python/elc.py @@ -16,8 +16,7 @@ # along with this program. If not, see . import unittest as ut import unittest_decorators as utx -import espressomd -from espressomd import electrostatics, electrostatic_extensions +import espressomd.electrostatics import numpy as np @@ -38,24 +37,23 @@ def test_finite_potential_drop(self): p1 = system.part.add(pos=[0, 0, 1], q=+1) p2 = system.part.add(pos=[0, 0, 9], q=-1) - system.actors.add( - electrostatics.P3M( - # zero is not allowed - prefactor=1e-100, - mesh=32, - cao=5, - accuracy=1e-3, - )) - - system.actors.add( - electrostatic_extensions.ELC( - gap_size=GAP[2], - maxPWerror=1e-3, - delta_mid_top=-1, - delta_mid_bot=-1, - const_pot=1, - pot_diff=POTENTIAL_DIFFERENCE, - )) + p3m = espressomd.electrostatics.P3M( + # zero is not allowed + prefactor=1e-100, + mesh=32, + cao=5, + accuracy=1e-3, + ) + elc = espressomd.electrostatics.ELC( + p3m_actor=p3m, + gap_size=GAP[2], + maxPWerror=1e-3, + delta_mid_top=-1, + delta_mid_bot=-1, + const_pot=1, + pot_diff=POTENTIAL_DIFFERENCE, + ) + system.actors.add(elc) # Calculated energy U_elc = system.analysis.energy()['coulomb'] @@ -76,13 +74,13 @@ def test_finite_potential_drop(self): p1.pos = [BOX_L[0] / 2, BOX_L[1] / 2, BOX_L[2] - GAP[2] / 2] with self.assertRaises(Exception): self.system.analysis.energy() - with self.assertRaises(Exception, msg='entered ELC gap region'): + with self.assertRaisesRegex(Exception, 'entered ELC gap region'): self.system.integrator.run(2) # negative direction p1.pos = [BOX_L[0] / 2, BOX_L[1] / 2, -GAP[2] / 2] with self.assertRaises(Exception): self.system.analysis.energy() - with self.assertRaises(Exception, msg='entered ELC gap region'): + with self.assertRaisesRegex(Exception, 'entered ELC gap region'): self.system.integrator.run(2) diff --git a/testsuite/python/elc_vs_analytic.py b/testsuite/python/elc_vs_analytic.py index 8f186608061..d5bfbb28764 100644 --- a/testsuite/python/elc_vs_analytic.py +++ b/testsuite/python/elc_vs_analytic.py @@ -19,7 +19,6 @@ import espressomd import numpy as np import espressomd.electrostatics -from espressomd import electrostatic_extensions @utx.skipIfMissingFeatures(["P3M"]) @@ -62,12 +61,11 @@ def test_elc(self): accuracy=self.accuracy, mesh=[58, 58, 70], cao=4) - self.system.actors.add(p3m) - - elc = electrostatic_extensions.ELC(gap_size=self.elc_gap, - maxPWerror=self.accuracy, - delta_mid_bot=self.delta_mid_bot, - delta_mid_top=self.delta_mid_top) + elc = espressomd.electrostatics.ELC(p3m_actor=p3m, + gap_size=self.elc_gap, + maxPWerror=self.accuracy, + delta_mid_bot=self.delta_mid_bot, + delta_mid_top=self.delta_mid_top) self.system.actors.add(elc) elc_results = self.scan() From 7ba3b1374f8a13c1ac9dcc84d43f52a32506047b Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Mon, 1 Mar 2021 14:38:31 +0100 Subject: [PATCH 3/7] tests: Add DP3M/ELC/Scafacos checkpointing checks --- testsuite/python/CMakeLists.txt | 2 +- testsuite/python/save_checkpoint.py | 74 +++++++++++++++++---- testsuite/python/test_checkpoint.py | 99 +++++++++++++++++++++++++++-- 3 files changed, 155 insertions(+), 20 deletions(-) diff --git a/testsuite/python/CMakeLists.txt b/testsuite/python/CMakeLists.txt index 63bd46bb051..f13ae95562d 100644 --- a/testsuite/python/CMakeLists.txt +++ b/testsuite/python/CMakeLists.txt @@ -57,7 +57,7 @@ endfunction(PYTHON_TEST) # Separate features with hyphens, use a period to add an optional flag. foreach( TEST_COMBINATION - lb.cpu-p3m.cpu-lj-therm.lb;lb.gpu-p3m.cpu-lj-therm.lb;ek.gpu;lb.off-therm.npt-int.npt;lb.off-int.sd;lb.off-therm.langevin-int.nvt;lb.off-therm.dpd-int.nvt;lb.off-therm.bd-int.bd;lb.off-therm.sdm-int.sdm + lb.cpu-p3m.cpu-lj-therm.lb;lb.gpu-p3m.elc-lj-therm.lb;ek.gpu;lb.off-therm.npt-int.npt;lb.off-int.sd;lb.off-dp3m.cpu-therm.langevin-int.nvt;lb.off-therm.dpd-int.nvt;lb.off-scafacos-therm.bd-int.bd;lb.off-therm.sdm-int.sdm ) if(${TEST_COMBINATION} MATCHES "\\.gpu") set(TEST_LABELS "gpu") diff --git a/testsuite/python/save_checkpoint.py b/testsuite/python/save_checkpoint.py index f5726a3439a..ece4806e5d4 100644 --- a/testsuite/python/save_checkpoint.py +++ b/testsuite/python/save_checkpoint.py @@ -21,6 +21,7 @@ import espressomd import espressomd.checkpointing import espressomd.electrostatics +import espressomd.magnetostatics import espressomd.interactions import espressomd.virtual_sites import espressomd.accumulators @@ -36,8 +37,10 @@ modes = {x for mode in set("@TEST_COMBINATION@".upper().split('-')) for x in [mode, mode.split('.')[0]]} -# use a box with 3 different dimensions +# use a box with 3 different dimensions, unless DipolarP3M is used system = espressomd.System(box_l=[12.0, 14.0, 16.0]) +if 'DP3M' in modes: + system.box_l = 3 * [np.max(system.box_l)] system.cell_system.skin = 0.1 system.time_step = 0.01 system.min_global_cut = 2.0 @@ -89,15 +92,13 @@ ek.add_species(ek_species) system.actors.add(ek) -system.part.add(pos=[1.0] * 3) -system.part.add(pos=[1.0, 1.0, 2.0]) +p1 = system.part.add(pos=[1.0] * 3, q=1, dip=(1.3, 2.1, -6)) +p2 = system.part.add(pos=[1.0, 1.0, 2.0], q=-1, dip=(7.3, 6.1, -4)) if espressomd.has_features('EXCLUSIONS'): system.part.add(pos=[2.0] * 3, exclusions=[0, 1]) -if espressomd.has_features('P3M') and 'P3M.CPU' in modes: - system.part[0].q = 1 - system.part[1].q = -1 +if espressomd.has_features('P3M') and 'P3M' in modes: p3m = espressomd.electrostatics.P3M( prefactor=1.0, accuracy=0.1, @@ -106,7 +107,16 @@ alpha=1.0, r_cut=1.0, tune=False) - system.actors.add(p3m) + if 'P3M.CPU' in modes: + system.actors.add(p3m) + elif 'P3M.ELC' in modes: + elc = espressomd.electrostatics.ELC( + p3m_actor=p3m, + gap_size=6.0, + maxPWerror=0.1, + delta_mid_top=0.9, + delta_mid_bot=0.1) + system.actors.add(elc) obs = espressomd.observables.ParticlePositions(ids=[0, 1]) acc_mean_variance = espressomd.accumulators.MeanVarianceCalculator(obs=obs) @@ -117,7 +127,7 @@ acc_mean_variance.update() acc_time_series.update() acc_correlator.update() -system.part[0].pos = [1.0, 2.0, 3.0] +p1.pos = [1.0, 2.0, 3.0] acc_mean_variance.update() acc_time_series.update() acc_correlator.update() @@ -182,7 +192,7 @@ if espressomd.has_features(['VIRTUAL_SITES', 'VIRTUAL_SITES_RELATIVE']): system.virtual_sites = espressomd.virtual_sites.VirtualSitesRelative( have_quaternion=True) - system.part[1].vs_auto_relate_to(0) + p2.vs_auto_relate_to(p1) if espressomd.has_features(['LENNARD_JONES']) and 'LJ' in modes: system.non_bonded_inter[0, 0].lennard_jones.set_params( @@ -192,27 +202,65 @@ harmonic_bond = espressomd.interactions.HarmonicBond(r_0=0.0, k=1.0) system.bonded_inter.add(harmonic_bond) -system.part[1].add_bond((harmonic_bond, 0)) +p2.add_bond((harmonic_bond, p1)) if 'THERM.LB' not in modes: thermalized_bond = espressomd.interactions.ThermalizedBond( temp_com=0.0, gamma_com=0.0, temp_distance=0.2, gamma_distance=0.5, r_cut=2, seed=51) system.bonded_inter.add(thermalized_bond) - system.part[1].add_bond((thermalized_bond, 0)) + p2.add_bond((thermalized_bond, p1)) checkpoint.register("system") checkpoint.register("acc_mean_variance") checkpoint.register("acc_time_series") checkpoint.register("acc_correlator") # calculate forces system.integrator.run(0) -particle_force0 = np.copy(system.part[0].f) -particle_force1 = np.copy(system.part[1].f) +particle_force0 = np.copy(p1.f) +particle_force1 = np.copy(p2.f) checkpoint.register("particle_force0") checkpoint.register("particle_force1") if espressomd.has_features("COLLISION_DETECTION"): system.collision_detection.set_params( mode="bind_centers", distance=0.11, bond_centers=harmonic_bond) +if espressomd.has_features('DP3M') and 'DP3M' in modes: + dp3m = espressomd.magnetostatics.DipolarP3M( + prefactor=1., + epsilon=2., + mesh_off=[0.5, 0.5, 0.5], + r_cut=2.4, + cao=1, + mesh=[8, 8, 8], + alpha=12, + accuracy=0.01, + tune=False) + system.actors.add(dp3m) + +if espressomd.has_features('SCAFACOS') and 'SCAFACOS' in modes: + system.actors.add(espressomd.electrostatics.Scafacos( + prefactor=0.5, + method_name="p3m", + method_params={ + "p3m_r_cut": 1.0, + "p3m_grid": 64, + "p3m_cao": 7, + "p3m_alpha": 2.7})) + +if espressomd.has_features('SCAFACOS_DIPOLES') and 'SCAFACOS' in modes: + system.actors.add(espressomd.magnetostatics.Scafacos( + prefactor=1.2, + method_name='p2nfft', + method_params={ + "p2nfft_verbose_tuning": "0", + "pnfft_N": "32,32,32", + "pnfft_n": "32,32,32", + "pnfft_window_name": "bspline", + "pnfft_m": "4", + "p2nfft_ignore_tolerance": "1", + "pnfft_diff_ik": "0", + "p2nfft_r_cut": "11", + "p2nfft_alpha": "0.37"})) + if LB_implementation: m = np.pi / 12 nx = int(np.round(system.box_l[0] / lbf.get_params()["agrid"])) diff --git a/testsuite/python/test_checkpoint.py b/testsuite/python/test_checkpoint.py index e9626a081d4..119033e62e0 100644 --- a/testsuite/python/test_checkpoint.py +++ b/testsuite/python/test_checkpoint.py @@ -22,6 +22,8 @@ import espressomd import espressomd.checkpointing +import espressomd.electrostatics +import espressomd.magnetostatics import espressomd.virtual_sites import espressomd.integrate from espressomd.shapes import Sphere, Wall @@ -44,6 +46,16 @@ def setUpClass(cls): '.', '__'), checkpoint_path="@CMAKE_CURRENT_BINARY_DIR@") cls.checkpoint.load(0) + cls.ref_box_l = np.array([12.0, 14.0, 16.0]) + if 'DP3M' in modes: + cls.ref_box_l = np.array([16.0, 16.0, 16.0]) + + def get_active_actor_of_type(self, actor_type): + for actor in system.actors.active_actors: + if isinstance(actor, actor_type): + return actor + self.fail( + f"system doesn't have an actor of type {actor_type.__name__}") @ut.skipIf(not LB, "Skipping test due to missing mode.") def test_LB(self): @@ -120,6 +132,7 @@ def test_variables(self): self.assertEqual(system.cell_system.skin, 0.1) self.assertEqual(system.time_step, 0.01) self.assertEqual(system.min_global_cut, 2.0) + np.testing.assert_allclose(np.copy(system.box_l), self.ref_box_l) def test_part(self): np.testing.assert_allclose( @@ -306,18 +319,91 @@ def test_correlator(self): system.auto_update_accumulators[2].result(), expected) + @utx.skipIfMissingFeatures('DP3M') + @ut.skipIf('DP3M.CPU' not in modes, + "Skipping test due to missing combination.") + def test_dp3m(self): + actor = self.get_active_actor_of_type( + espressomd.magnetostatics.DipolarP3M) + state = actor.get_params() + reference = {'prefactor': 1.0, 'accuracy': 0.01, 'mesh': 3 * [8], + 'cao': 1, 'alpha': 12.0, 'r_cut': 2.4, 'tune': False, + 'mesh_off': [0.5, 0.5, 0.5], 'epsilon': 2.0} + for key in reference: + self.assertIn(key, state) + np.testing.assert_almost_equal(state[key], reference[key], + err_msg=f'for parameter {key}') + @utx.skipIfMissingFeatures('P3M') @ut.skipIf('P3M.CPU' not in modes, "Skipping test due to missing combination.") def test_p3m(self): - actor = system.actors.active_actors[-1] - self.assertTrue(isinstance(actor, espressomd.electrostatics.P3M)) + actor = self.get_active_actor_of_type(espressomd.electrostatics.P3M) state = actor.get_params() reference = {'prefactor': 1.0, 'accuracy': 0.1, 'mesh': 3 * [10], - 'cao': 1, 'alpha': 1.0, 'r_cut': 1.0} + 'cao': 1, 'alpha': 1.0, 'r_cut': 1.0, 'tune': False} + for key in reference: + self.assertIn(key, state) + np.testing.assert_almost_equal(state[key], reference[key], + err_msg=f'for parameter {key}') + + @utx.skipIfMissingFeatures('P3M') + @ut.skipIf('P3M.ELC' not in modes, + "Skipping test due to missing combination.") + def test_elc(self): + actor = self.get_active_actor_of_type(espressomd.electrostatics.ELC) + elc_state = actor.get_params() + p3m_state = elc_state['p3m_actor'].get_params() + p3m_reference = {'prefactor': 1.0, 'accuracy': 0.1, 'mesh': 3 * [10], + 'cao': 1, 'alpha': 1.0, 'r_cut': 1.0, 'tune': False} + elc_reference = {'gap_size': 6.0, 'maxPWerror': 0.1, + 'delta_mid_top': 0.9, 'delta_mid_bot': 0.1} + for key in elc_reference: + self.assertIn(key, elc_state) + np.testing.assert_almost_equal(elc_state[key], elc_reference[key], + err_msg=f'for parameter {key}') + for key in p3m_reference: + self.assertIn(key, p3m_state) + np.testing.assert_almost_equal(p3m_state[key], p3m_reference[key], + err_msg=f'for parameter {key}') + + @utx.skipIfMissingFeatures('SCAFACOS') + @ut.skipIf('SCAFACOS' not in modes, + "Skipping test due to missing combination.") + def test_scafacos(self): + actor = self.get_active_actor_of_type( + espressomd.electrostatics.Scafacos) + state = actor.get_params() + reference = {'prefactor': 0.5, 'method_name': 'p3m', + 'method_params': { + 'p3m_cao': '7', + 'p3m_r_cut': '1.0', + 'p3m_grid': '64', + 'p3m_alpha': '2.7'}} + for key in reference: + self.assertEqual(state[key], reference[key], msg=f'for {key}') + + @utx.skipIfMissingFeatures('SCAFACOS_DIPOLES') + @ut.skipIf('SCAFACOS' not in modes, + "Skipping test due to missing combination.") + def test_scafacos_dipoles(self): + actor = self.get_active_actor_of_type( + espressomd.magnetostatics.Scafacos) + state = actor.get_params() + reference = {'prefactor': 1.2, 'method_name': 'p2nfft', + 'method_params': { + "p2nfft_verbose_tuning": "0", + "pnfft_N": "32,32,32", + "pnfft_n": "32,32,32", + "pnfft_window_name": "bspline", + "pnfft_m": "4", + "p2nfft_ignore_tolerance": "1", + "pnfft_diff_ik": "0", + "p2nfft_r_cut": "11", + "p2nfft_alpha": "0.37"}} for key in reference: self.assertIn(key, state) - np.testing.assert_almost_equal(state[key], reference[key]) + self.assertEqual(state[key], reference[key], msg=f'for {key}') @utx.skipIfMissingFeatures('COLLISION_DETECTION') def test_collision_detection(self): @@ -345,6 +431,7 @@ def test_constraints(self): self.assertEqual(len(system.constraints), 8 - int(not espressomd.has_features("ELECTROSTATICS"))) c = system.constraints + ref_shape = self.ref_box_l.astype(int) + 2 self.assertIsInstance(c[0].shape, Sphere) self.assertAlmostEqual(c[0].shape.radius, 0.1, delta=1E-10) @@ -365,7 +452,7 @@ def test_constraints(self): self.assertAlmostEqual(c[4].gamma, 2.3, delta=1E-10) self.assertIsInstance(c[5], constraints.PotentialField) - self.assertEqual(c[5].field.shape, (14, 16, 18, 1)) + self.assertEqual(c[5].field.shape, tuple(list(ref_shape) + [1])) self.assertAlmostEqual(c[5].default_scale, 1.6, delta=1E-10) self.assertAlmostEqual(c[5].particle_scales[5], 6.0, delta=1E-10) np.testing.assert_allclose(np.copy(c[5].origin), [-0.5, -0.5, -0.5]) @@ -376,7 +463,7 @@ def test_constraints(self): atol=1e-10) self.assertIsInstance(c[6], constraints.ForceField) - self.assertEqual(c[6].field.shape, (14, 16, 18, 3)) + self.assertEqual(c[6].field.shape, tuple(list(ref_shape) + [3])) self.assertAlmostEqual(c[6].default_scale, 1.4, delta=1E-10) np.testing.assert_allclose(np.copy(c[6].origin), [-0.5, -0.5, -0.5]) np.testing.assert_allclose(np.copy(c[6].grid_spacing), np.ones(3)) From 944b9a7c571d7cc541cb2e5d74b6d83e9af53818 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Mon, 1 Mar 2021 14:57:25 +0100 Subject: [PATCH 4/7] tests: Check ELC exception mechanism --- testsuite/python/p3m_tuning_exceptions.py | 42 ++++++++++++++++++++++- 1 file changed, 41 insertions(+), 1 deletion(-) diff --git a/testsuite/python/p3m_tuning_exceptions.py b/testsuite/python/p3m_tuning_exceptions.py index 460376ad90b..e3fbc4c5d65 100644 --- a/testsuite/python/p3m_tuning_exceptions.py +++ b/testsuite/python/p3m_tuning_exceptions.py @@ -160,7 +160,7 @@ def test_03_non_cubic_box_dp3m_cpu(self): def check_invalid_params(self, solver_class, **custom_params): valid_params = { 'prefactor': 2, 'accuracy': .01, 'tune': False, 'cao': 1, - 'r_cut': 3.73e-01, 'alpha': 3.81e+00, 'mesh': (8, 8, 8), + 'r_cut': 0.373, 'alpha': 3.81, 'mesh': (8, 8, 8), 'mesh_off': [-1, -1, -1]} valid_params.update(custom_params) @@ -212,6 +212,46 @@ def test_04_invalid_params_dp3m_cpu(self): self.check_invalid_params(espressomd.magnetostatics.DipolarP3M) + @utx.skipIfMissingFeatures("P3M") + def test_04_invalid_params_p3m_elc_cpu(self): + import espressomd.electrostatics + + self.system.time_step = 0.01 + self.add_charged_particles() + + solver_p3m = espressomd.electrostatics.P3M( + prefactor=2, accuracy=0.01, tune=False, cao=1, + r_cut=0.373, alpha=3.81, mesh=(8, 8, 8)) + solver_elc = espressomd.electrostatics.ELC( + p3m_actor=solver_p3m, gap_size=1.2 * self.system.box_l[2], + maxPWerror=0.01) + with self.assertRaisesRegex(Exception, "gap size too large"): + self.system.actors.add(solver_elc) + + self.system.actors.clear() + solver_dh = espressomd.electrostatics.DH( + prefactor=1.2, kappa=0.8, r_cut=2.0) + solver_elc = espressomd.electrostatics.ELC( + p3m_actor=solver_dh, gap_size=1, maxPWerror=0.01) + with self.assertRaisesRegex(ValueError, "p3m_actor has to be a P3M solver"): + self.system.actors.add(solver_elc) + + @utx.skipIfMissingGPU() + @utx.skipIfMissingFeatures("P3M") + def test_04_invalid_params_p3m_elc_gpu(self): + import espressomd.electrostatics + + self.system.time_step = 0.01 + self.add_charged_particles() + + solver_p3m = espressomd.electrostatics.P3MGPU( + prefactor=2, accuracy=0.01, tune=False, cao=1, + r_cut=4.4434, alpha=0.3548, mesh=(28, 28, 28)) + solver_elc = espressomd.electrostatics.ELC( + p3m_actor=solver_p3m, gap_size=1, maxPWerror=0.01) + with self.assertRaisesRegex(ValueError, "ELC is not set up to work with the GPU P3M"): + self.system.actors.add(solver_elc) + ########################################################### # block of tests where tuning should not throw exceptions # ########################################################### From bf25f94c9c9cc6dfdaa445e76ee58a9833aee2f1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Mon, 1 Mar 2021 15:03:58 +0100 Subject: [PATCH 5/7] python: Improve ELC and MDLC docstrings Both methods check for the presence of particles in the gap region. Show variable names that failed validation. --- src/python/espressomd/electrostatics.pyx | 18 ++++++++++-------- .../espressomd/magnetostatic_extensions.pyx | 18 +++++++++++++----- 2 files changed, 23 insertions(+), 13 deletions(-) diff --git a/src/python/espressomd/electrostatics.pyx b/src/python/espressomd/electrostatics.pyx index 65f4f44a0ac..1fa9ae8a756 100644 --- a/src/python/espressomd/electrostatics.pyx +++ b/src/python/espressomd/electrostatics.pyx @@ -408,10 +408,8 @@ IF P3M == 1: Base P3M actor. gap_size : :obj:`float`, required The gap size gives the height :math:`h` of the empty region between - the system box and the neighboring artificial images. |es| does not - make sure that the gap is actually empty, this is the user's - responsibility. The method will run even if the condition is not - fulfilled, however, the error bound will not be reached. Therefore + the system box and the neighboring artificial images. |es| checks + the gap is empty and will throw an error if it isn't. Therefore you should really make sure that the gap region is empty (e.g. with wall constraints). maxPWerror : :obj:`float`, required @@ -467,15 +465,19 @@ IF P3M == 1: self._params["p3m_actor"].validate_params() # ELC check_type_or_throw_except( - self._params["maxPWerror"], 1, float, "") + self._params["maxPWerror"], 1, float, + "maxPWerror has to be a float") check_range_or_except( self._params, "maxPWerror", 0, False, "inf", True) - check_type_or_throw_except(self._params["gap_size"], 1, float, "") + check_type_or_throw_except(self._params["gap_size"], 1, float, + "gap_size has to be a float") check_range_or_except( self._params, "gap_size", 0, False, "inf", True) - check_type_or_throw_except(self._params["far_cut"], 1, float, "") + check_type_or_throw_except(self._params["far_cut"], 1, float, + "far_cut has to be a float") check_type_or_throw_except( - self._params["neutralize"], 1, type(True), "") + self._params["neutralize"], 1, type(True), + "neutralize has to be a bool") def valid_keys(self): return ["p3m_actor", "maxPWerror", "gap_size", "far_cut", diff --git a/src/python/espressomd/magnetostatic_extensions.pyx b/src/python/espressomd/magnetostatic_extensions.pyx index 69edd7ad57d..3f5cece537f 100644 --- a/src/python/espressomd/magnetostatic_extensions.pyx +++ b/src/python/espressomd/magnetostatic_extensions.pyx @@ -42,8 +42,11 @@ IF DIPOLES and DP3M: Parameters ---------- gap_size : :obj:`float` - Size of the empty gap. Note that DLC relies on the user to make - sure that this condition is fulfilled. + The gap size gives the height :math:`h` of the empty region between + the system box and the neighboring artificial images. |es| checks + the gap is empty and will throw an error if it isn't. Therefore + you should really make sure that the gap region is empty (e.g. + with wall constraints). maxPWerror : :obj:`float` Maximal pairwise error of the potential and force. far_cut : :obj:`float`, optional @@ -57,13 +60,18 @@ IF DIPOLES and DP3M: """ default_params = self.default_params() check_type_or_throw_except( - self._params["maxPWerror"], 1, float, "") + self._params["maxPWerror"], 1, float, + "maxPWerror has to be a float") check_range_or_except( self._params, "maxPWerror", 0, False, "inf", True) - check_type_or_throw_except(self._params["gap_size"], 1, float, "") + check_type_or_throw_except( + self._params["gap_size"], 1, float, + "gap_size has to be a float") check_range_or_except( self._params, "gap_size", 0, False, "inf", True) - check_type_or_throw_except(self._params["far_cut"], 1, float, "") + check_type_or_throw_except( + self._params["far_cut"], 1, float, + "far_cut has to be a float") def valid_keys(self): return ["maxPWerror", "gap_size", "far_cut"] From 5569e4b2af287de4b773dd187cd1d7c1d440dd8f Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Mon, 1 Mar 2021 15:52:37 +0100 Subject: [PATCH 6/7] tests: Fix regression --- testsuite/python/save_checkpoint.py | 12 ++++++++++-- 1 file changed, 10 insertions(+), 2 deletions(-) diff --git a/testsuite/python/save_checkpoint.py b/testsuite/python/save_checkpoint.py index ece4806e5d4..46a55ee7eca 100644 --- a/testsuite/python/save_checkpoint.py +++ b/testsuite/python/save_checkpoint.py @@ -92,8 +92,16 @@ ek.add_species(ek_species) system.actors.add(ek) -p1 = system.part.add(pos=[1.0] * 3, q=1, dip=(1.3, 2.1, -6)) -p2 = system.part.add(pos=[1.0, 1.0, 2.0], q=-1, dip=(7.3, 6.1, -4)) +p1 = system.part.add(pos=[1.0] * 3) +p2 = system.part.add(pos=[1.0, 1.0, 2.0]) + +if espressomd.has_features('ELECTROSTATICS'): + p1.q = 1 + p2.q = -1 + +if espressomd.has_features('DIPOLES'): + p1.dip = (1.3, 2.1, -6) + p2.dip = (7.3, 6.1, -4) if espressomd.has_features('EXCLUSIONS'): system.part.add(pos=[2.0] * 3, exclusions=[0, 1]) From b2dcdcee36430a0248c63ce82edada98f9edabbf Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Wed, 3 Mar 2021 11:48:00 +0100 Subject: [PATCH 7/7] docs: Apply suggested changes Co-authored-by: Alexander Reinauer --- src/python/espressomd/electrostatics.pyx | 2 +- src/python/espressomd/magnetostatic_extensions.pyx | 2 +- 2 files changed, 2 insertions(+), 2 deletions(-) diff --git a/src/python/espressomd/electrostatics.pyx b/src/python/espressomd/electrostatics.pyx index 1fa9ae8a756..f1af5f2dcbe 100644 --- a/src/python/espressomd/electrostatics.pyx +++ b/src/python/espressomd/electrostatics.pyx @@ -409,7 +409,7 @@ IF P3M == 1: gap_size : :obj:`float`, required The gap size gives the height :math:`h` of the empty region between the system box and the neighboring artificial images. |es| checks - the gap is empty and will throw an error if it isn't. Therefore + that the gap is empty and will throw an error if it isn't. Therefore you should really make sure that the gap region is empty (e.g. with wall constraints). maxPWerror : :obj:`float`, required diff --git a/src/python/espressomd/magnetostatic_extensions.pyx b/src/python/espressomd/magnetostatic_extensions.pyx index 3f5cece537f..e42eb695a05 100644 --- a/src/python/espressomd/magnetostatic_extensions.pyx +++ b/src/python/espressomd/magnetostatic_extensions.pyx @@ -44,7 +44,7 @@ IF DIPOLES and DP3M: gap_size : :obj:`float` The gap size gives the height :math:`h` of the empty region between the system box and the neighboring artificial images. |es| checks - the gap is empty and will throw an error if it isn't. Therefore + that the gap is empty and will throw an error if it isn't. Therefore you should really make sure that the gap region is empty (e.g. with wall constraints). maxPWerror : :obj:`float`