From 45dbc211eb1968e27906a7ec74e9fd4d4d1795ec Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Tue, 4 Oct 2022 18:43:59 +0200 Subject: [PATCH 1/4] tests: Check Monte Carlo displacement moves --- testsuite/python/CMakeLists.txt | 1 + testsuite/python/canonical_ensemble.py | 116 +++++++++++++++++++++++++ 2 files changed, 117 insertions(+) create mode 100644 testsuite/python/canonical_ensemble.py diff --git a/testsuite/python/CMakeLists.txt b/testsuite/python/CMakeLists.txt index 5adb9cf9581..bcba00d4747 100644 --- a/testsuite/python/CMakeLists.txt +++ b/testsuite/python/CMakeLists.txt @@ -177,6 +177,7 @@ python_test(FILE reaction_ensemble.py MAX_NUM_PROC 4) python_test(FILE widom_insertion.py MAX_NUM_PROC 1) python_test(FILE constant_pH.py MAX_NUM_PROC 1) python_test(FILE constant_pH_stats.py MAX_NUM_PROC 4 LABELS long) +python_test(FILE canonical_ensemble.py MAX_NUM_PROC 2) python_test(FILE writevtf.py MAX_NUM_PROC 4) python_test(FILE lb_stokes_sphere.py MAX_NUM_PROC 4 LABELS gpu long) python_test(FILE lb_pressure_tensor.py MAX_NUM_PROC 1 LABELS gpu long) diff --git a/testsuite/python/canonical_ensemble.py b/testsuite/python/canonical_ensemble.py new file mode 100644 index 00000000000..c062d712f4d --- /dev/null +++ b/testsuite/python/canonical_ensemble.py @@ -0,0 +1,116 @@ +# +# Copyright (C) 2022 The ESPResSo project +# +# This file is part of ESPResSo. +# +# ESPResSo is free software: you can redistribute it and/or modify +# it under the terms of the GNU General Public License as published by +# the Free Software Foundation, either version 3 of the License, or +# (at your option) any later version. +# +# ESPResSo is distributed in the hope that it will be useful, +# but WITHOUT ANY WARRANTY; without even the implied warranty of +# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the +# GNU General Public License for more details. +# +# You should have received a copy of the GNU General Public License +# along with this program. If not, see . +# + +import espressomd +import espressomd.constraints +import espressomd.observables +import espressomd.accumulators +import espressomd.reaction_methods + +import numpy as np +import scipy.optimize +import unittest as ut +import unittest_decorators as utx + + +class Test(ut.TestCase): + + system = espressomd.System(box_l=[1., 1., 1.]) + system.cell_system.skin = 0.4 + system.setup_type_map(type_list=[0]) + + def tearDown(self): + self.system.part.clear() + self.system.constraints.clear() + + @utx.skipIfMissingFeatures("ELECTROSTATICS") + def test_linear_potential(self): + """ + Set a particle in a box with a linear potential along the x-axis. + The particle distribution resulting from accepted Monte Carlo moves + should follow a Maxwell-Boltzmann distribution. + """ + + method = espressomd.reaction_methods.ReactionEnsemble( + kT=0.2, seed=42, exclusion_range=0., search_algorithm="order_n") + method.set_non_interacting_type(type=1) + + p = self.system.part.add(pos=[0., 0., 0.], q=1, type=0) + obs_pos = espressomd.observables.ParticlePositions(ids=(p.id,)) + obs_vel = espressomd.observables.ParticleVelocities(ids=(p.id,)) + acc_pos = espressomd.accumulators.TimeSeries(obs=obs_pos) + acc_vel = espressomd.accumulators.TimeSeries(obs=obs_vel) + + E = np.array([-1., 0., 0.]) + field = espressomd.constraints.LinearElectricPotential(E=E, phi0=0.) + self.system.constraints.add(field) + + for _ in range(5000): + accepted = method.displacement_mc_move_for_particles_of_type( + type_mc=0, particle_number_to_be_changed=1) + if accepted: + acc_pos.update() + acc_vel.update() + p.pos = [0., 0., 0.] + + # the x-position should follow an exponential distribution + # -> mean = kT, median = kT x ln(2), variance = kT^2 + series = acc_pos.time_series()[:, p.id, 0] + ydata, xbins = np.histogram(series, bins=15, range=[0., 1.]) + xdata = (xbins[1:] + xbins[:-1]) / 2. + ydata = ydata / float(ydata[0]) + (a, b, c), _ = scipy.optimize.curve_fit( + lambda x, a, b, c: a * np.exp(-b * x) + c, xdata, ydata) + # check histogram profile is roughly exponential + self.assertAlmostEqual(a, 1., delta=0.2) + self.assertAlmostEqual(b, 1. / method.kT, delta=0.3) + self.assertAlmostEqual(c, 0., delta=0.01) + # check distribution parameters with high accuracy + ln2 = np.log(2) + self.assertAlmostEqual(np.mean(series), method.kT, delta=0.02) + self.assertAlmostEqual(np.median(series) / ln2, method.kT, delta=0.02) + self.assertAlmostEqual(np.sqrt(np.var(series)), method.kT, delta=0.02) + + # the y- and z-position should follow a uniform distribution + for axis in (1, 2): + series = acc_pos.time_series()[:, p.id, axis] + ydata, _ = np.histogram(series, bins=10, range=[0., 1.]) + ydata = ydata / np.mean(ydata) + np.testing.assert_allclose(ydata, 1., atol=0.25) + + # the velocity vector should follow a normal distribution + # -> mean = 0, median = 0, variance = kT + for axis in (0, 1, 2): + series = acc_vel.time_series()[:, p.id, axis] + ydata, xbins = np.histogram(series, bins=25, range=[-1.5, 1.5]) + xdata = (xbins[1:] + xbins[:-1]) / 2. + ydata = ydata / len(series) + (_, b, c), _ = scipy.optimize.curve_fit( + lambda x, a, b, c: a * np.exp(-b * x**2) + c, xdata, ydata) + # check histogram profile is roughly gaussian + self.assertAlmostEqual(b, 0.5 / method.kT, delta=0.45) + self.assertAlmostEqual(c, 0., delta=0.002) + # check distribution parameters with high accuracy + self.assertAlmostEqual(np.mean(series), 0., delta=0.05) + self.assertAlmostEqual(np.median(series), 0., delta=0.025) + self.assertAlmostEqual(np.var(series), method.kT, delta=0.025) + + +if __name__ == "__main__": + ut.main() From a8f166ac637ae28c91e9e4e749b16857a0f3cf2c Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Thu, 6 Oct 2022 21:52:08 +0200 Subject: [PATCH 2/4] core: Fix Monte Carlo displacement bugs Restore original velocities when a displacement move is rejected. When one move in the sequence is rejected, skip the remaining move attempts. Fix infinite loop when not enough particles are available. --- .../reaction_methods/ReactionAlgorithm.cpp | 47 ++++++++++++------- .../reaction_methods/ReactionAlgorithm.hpp | 5 +- .../tests/ReactionAlgorithm_test.cpp | 16 ++++--- src/python/espressomd/reaction_methods.py | 15 ++++-- .../reaction_methods/ReactionAlgorithm.hpp | 2 +- testsuite/python/canonical_ensemble.py | 3 ++ 6 files changed, 55 insertions(+), 33 deletions(-) diff --git a/src/core/reaction_methods/ReactionAlgorithm.cpp b/src/core/reaction_methods/ReactionAlgorithm.cpp index 1f8330318be..304758c6b62 100644 --- a/src/core/reaction_methods/ReactionAlgorithm.cpp +++ b/src/core/reaction_methods/ReactionAlgorithm.cpp @@ -548,11 +548,11 @@ void ReactionAlgorithm::move_particle(int p_id, Utils::Vector3d const &new_pos, set_particle_v(p_id, vel); } -std::vector> +std::vector> ReactionAlgorithm::generate_new_particle_positions(int type, int n_particles) { - std::vector> old_positions; - old_positions.reserve(n_particles); + std::vector> old_state; + old_state.reserve(n_particles); // draw particle ids at random without replacement int p_id = -1; @@ -564,38 +564,49 @@ ReactionAlgorithm::generate_new_particle_positions(int type, int n_particles) { p_id = get_random_p_id(type, random_index); } drawn_pids.emplace_back(p_id); - // store original position + // store original state auto const &p = get_particle_data(p_id); - old_positions.emplace_back(std::pair{p_id, p.pos()}); - // write new position + old_state.emplace_back(std::tuple{ + p_id, p.pos(), p.v()}); + // write new position and new velocity auto const prefactor = std::sqrt(kT / p.mass()); auto const new_pos = get_random_position_in_box(); move_particle(p_id, new_pos, prefactor); check_exclusion_range(p_id); + if (particle_inside_exclusion_range_touched) { + break; + } } - return old_positions; + return old_state; } /** - * Performs a global MC move for a particle of the provided type. + * Performs a global MC move for particles of a given type. + * Particles are selected without replacement. + * @param type Type of particles to move. + * @param n_part Number of particles of that type to move. + * @returns true if all moves were accepted. */ -bool ReactionAlgorithm::do_global_mc_move_for_particles_of_type( - int type, int particle_number_of_type_to_be_changed) { +bool ReactionAlgorithm::displacement_move_for_particles_of_type(int type, + int n_part) { m_tried_configurational_MC_moves += 1; particle_inside_exclusion_range_touched = false; - auto const particle_number_of_type = number_of_particles_with_type(type); - if (particle_number_of_type == 0 or - particle_number_of_type_to_be_changed == 0) { + if (n_part == 0) { + // reject + return false; + } + + auto const n_particles_of_type = number_of_particles_with_type(type); + if (n_part > n_particles_of_type) { // reject return false; } auto const E_pot_old = mpi_calculate_potential_energy(); - auto const original_positions = generate_new_particle_positions( - type, particle_number_of_type_to_be_changed); + auto const original_state = generate_new_particle_positions(type, n_part); auto const E_pot_new = (particle_inside_exclusion_range_touched) ? std::numeric_limits::max() @@ -622,9 +633,11 @@ bool ReactionAlgorithm::do_global_mc_move_for_particles_of_type( m_accepted_configurational_MC_moves += 1; return true; } - // reject: restore original particle positions - for (auto const &item : original_positions) + // reject: restore original particle properties + for (auto const &item : original_state) { + set_particle_v(std::get<0>(item), std::get<2>(item)); place_particle(std::get<0>(item), std::get<1>(item)); + } return false; } diff --git a/src/core/reaction_methods/ReactionAlgorithm.hpp b/src/core/reaction_methods/ReactionAlgorithm.hpp index 8b559b6d3d4..93e8783abbd 100644 --- a/src/core/reaction_methods/ReactionAlgorithm.hpp +++ b/src/core/reaction_methods/ReactionAlgorithm.hpp @@ -133,8 +133,7 @@ class ReactionAlgorithm { reactions.erase(reactions.begin() + reaction_id); } - bool do_global_mc_move_for_particles_of_type(int type, - int particle_number_of_type); + bool displacement_move_for_particles_of_type(int type, int n_part); bool particle_inside_exclusion_range_touched = false; bool neighbor_search_order_n = true; @@ -165,7 +164,7 @@ class ReactionAlgorithm { std::tuple, std::vector, std::vector> make_reaction_attempt(SingleReaction const ¤t_reaction); - std::vector> + std::vector> generate_new_particle_positions(int type, int n_particles); void restore_properties(std::vector const &property_list); diff --git a/src/core/reaction_methods/tests/ReactionAlgorithm_test.cpp b/src/core/reaction_methods/tests/ReactionAlgorithm_test.cpp index adc96b10990..6d231db12f8 100644 --- a/src/core/reaction_methods/tests/ReactionAlgorithm_test.cpp +++ b/src/core/reaction_methods/tests/ReactionAlgorithm_test.cpp @@ -39,6 +39,7 @@ #include #include #include +#include #include #include @@ -164,14 +165,15 @@ BOOST_AUTO_TEST_CASE(ReactionAlgorithm_test) { BOOST_CHECK(r_algo.particle_inside_exclusion_range_touched); // check moves and bookkeeping for (auto const &item : bookkeeping) { - auto const pid = item.first; + auto const pid = std::get<0>(item); BOOST_REQUIRE(pid == 0 or pid == 1); auto const ref_old_pos = ref_positions[pid].first; auto const ref_old_vel = ref_positions[pid].second; auto const &p = get_particle_data(pid); auto const &new_pos = p.pos(); auto const &new_vel = p.v(); - BOOST_CHECK_EQUAL(item.second, ref_old_pos); + BOOST_CHECK_EQUAL(std::get<1>(item), ref_old_pos); + BOOST_CHECK_EQUAL(std::get<2>(item), ref_old_vel); BOOST_CHECK_GE(new_pos, Utils::Vector3d::broadcast(0.)); BOOST_CHECK_LE(new_pos, Utils::Vector3d::broadcast(box_l)); BOOST_CHECK_GE((new_pos - ref_old_pos).norm(), 0.1); @@ -193,14 +195,14 @@ BOOST_AUTO_TEST_CASE(ReactionAlgorithm_test) { set_particle_type(0, type_A); set_particle_type(1, type_A); // check early exit when a MC move cannot be performed - BOOST_REQUIRE(!r_algo.do_global_mc_move_for_particles_of_type(type_C, 1)); - BOOST_REQUIRE(!r_algo.do_global_mc_move_for_particles_of_type(type_B, 2)); - BOOST_REQUIRE(!r_algo.do_global_mc_move_for_particles_of_type(type_A, 0)); + BOOST_REQUIRE(!r_algo.displacement_move_for_particles_of_type(type_C, 1)); + BOOST_REQUIRE(!r_algo.displacement_move_for_particles_of_type(type_B, 2)); + BOOST_REQUIRE(!r_algo.displacement_move_for_particles_of_type(type_A, 0)); // force all MC moves to be rejected by picking particles inside // their exclusion radius r_algo.exclusion_range = box_l; r_algo.particle_inside_exclusion_range_touched = false; - BOOST_REQUIRE(!r_algo.do_global_mc_move_for_particles_of_type(type_A, 2)); + BOOST_REQUIRE(!r_algo.displacement_move_for_particles_of_type(type_A, 2)); // check none of the particles moved for (auto const pid : {0, 1}) { auto const ref_old_pos = ref_positions[pid]; @@ -211,7 +213,7 @@ BOOST_AUTO_TEST_CASE(ReactionAlgorithm_test) { // force a MC move to be accepted by using a constant Hamiltonian r_algo.exclusion_range = 0.; r_algo.particle_inside_exclusion_range_touched = false; - BOOST_REQUIRE(r_algo.do_global_mc_move_for_particles_of_type(type_A, 1)); + BOOST_REQUIRE(r_algo.displacement_move_for_particles_of_type(type_A, 1)); std::vector distances(2); // check that only one particle moved for (auto const pid : {0, 1}) { diff --git a/src/python/espressomd/reaction_methods.py b/src/python/espressomd/reaction_methods.py index 1cf361ebfc7..65bcb7db666 100644 --- a/src/python/espressomd/reaction_methods.py +++ b/src/python/espressomd/reaction_methods.py @@ -210,11 +210,15 @@ class ReactionAlgorithm(ScriptInterfaceHelper): The number of reactions to be performed at once, defaults to 1. displacement_mc_move_for_particles_of_type() - Performs a displacement Monte Carlo move for particles of given type. + Performs displacement Monte Carlo moves for particles of a given type. New positions of the displaced particles are chosen from the whole box - with a uniform probability distribution. If there are multiple types, - that are being moved in a simulation, they should be moved in a random - order to avoid artefacts. + with a uniform probability distribution and new velocities are + sampled from the Maxwell-Boltzmann distribution. + + The sequence of moves is only accepted if each individual move in + the sequence was accepted. Particles are sampled without replacement. + Therefore, calling this method once for 10 particles is not equivalent + to calling this method 10 times for 1 particle. Parameters ---------- @@ -222,11 +226,12 @@ class ReactionAlgorithm(ScriptInterfaceHelper): Particle type which should be moved particle_number_to_be_changed : :obj:`int` Number of particles to move, defaults to 1. + Particles are selected without replacement. Returns ------- :obj:`bool` - Whether the move was accepted. + Whether all moves were accepted. delete_particle() Deletes the particle of the given p_id and makes sure that the particle diff --git a/src/script_interface/reaction_methods/ReactionAlgorithm.hpp b/src/script_interface/reaction_methods/ReactionAlgorithm.hpp index 9f6b9939437..8029552df9c 100644 --- a/src/script_interface/reaction_methods/ReactionAlgorithm.hpp +++ b/src/script_interface/reaction_methods/ReactionAlgorithm.hpp @@ -128,7 +128,7 @@ class ReactionAlgorithm : public AutoParameters { } else if (name == "reaction") { RE()->do_reaction(get_value_or(parameters, "reaction_steps", 1)); } else if (name == "displacement_mc_move_for_particles_of_type") { - return RE()->do_global_mc_move_for_particles_of_type( + return RE()->displacement_move_for_particles_of_type( get_value(parameters, "type_mc"), get_value_or(parameters, "particle_number_to_be_changed", 1)); } else if (name == "check_reaction_method") { diff --git a/testsuite/python/canonical_ensemble.py b/testsuite/python/canonical_ensemble.py index c062d712f4d..391fe8add82 100644 --- a/testsuite/python/canonical_ensemble.py +++ b/testsuite/python/canonical_ensemble.py @@ -68,6 +68,9 @@ def test_linear_potential(self): acc_pos.update() acc_vel.update() p.pos = [0., 0., 0.] + p.v = [0., 0., 0.] + else: + self.assertAlmostEqual(np.linalg.norm(p.v), 0., delta=1e-12) # the x-position should follow an exponential distribution # -> mean = kT, median = kT x ln(2), variance = kT^2 From 292bf32a992cce80c838875ebbf21b814bd2214e Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Thu, 6 Oct 2022 22:19:55 +0200 Subject: [PATCH 3/4] core: Track default-constructed particles When no user-defined type is available, explicitly assign type 0. This is necessary for tracking particles ids in reaction methods. --- src/script_interface/particle_data/ParticleHandle.cpp | 3 +++ 1 file changed, 3 insertions(+) diff --git a/src/script_interface/particle_data/ParticleHandle.cpp b/src/script_interface/particle_data/ParticleHandle.cpp index 70b9330a93b..19e3f153230 100644 --- a/src/script_interface/particle_data/ParticleHandle.cpp +++ b/src/script_interface/particle_data/ParticleHandle.cpp @@ -483,6 +483,9 @@ void ParticleHandle::do_construct(VariantMap const ¶ms) { do_set_parameter(kv.first, kv.second); } } + if (params.count("type") == 0) { + do_set_parameter("type", 0); + } } } From 789323516335d15190617d7ae458b9ebaa8f4ff1 Mon Sep 17 00:00:00 2001 From: =?UTF-8?q?Jean-No=C3=ABl=20Grad?= Date: Thu, 13 Oct 2022 23:09:11 +0200 Subject: [PATCH 4/4] core: Add sanity checks in MC displacement moves --- src/core/reaction_methods/ReactionAlgorithm.cpp | 16 ++++++++++++---- src/core/reaction_methods/ReactionAlgorithm.hpp | 2 +- .../tests/ReactionAlgorithm_test.cpp | 14 ++++++++++---- testsuite/python/reaction_methods_interface.py | 8 ++++++++ 4 files changed, 31 insertions(+), 9 deletions(-) diff --git a/src/core/reaction_methods/ReactionAlgorithm.cpp b/src/core/reaction_methods/ReactionAlgorithm.cpp index 304758c6b62..902402ff5d5 100644 --- a/src/core/reaction_methods/ReactionAlgorithm.cpp +++ b/src/core/reaction_methods/ReactionAlgorithm.cpp @@ -50,7 +50,7 @@ namespace ReactionMethods { /** * Performs a randomly selected reaction in the reaction ensemble */ -int ReactionAlgorithm::do_reaction(int reaction_steps) { +void ReactionAlgorithm::do_reaction(int reaction_steps) { // calculate potential energy; only consider potential energy since we // assume that the kinetic part drops out in the process of calculating // ensemble averages (kinetic part may be separated and crossed out) @@ -59,7 +59,6 @@ int ReactionAlgorithm::do_reaction(int reaction_steps) { int reaction_id = i_random(static_cast(reactions.size())); generic_oneway_reaction(*reactions[reaction_id], current_E_pot); } - return 0; } /** @@ -590,14 +589,23 @@ ReactionAlgorithm::generate_new_particle_positions(int type, int n_particles) { */ bool ReactionAlgorithm::displacement_move_for_particles_of_type(int type, int n_part) { - m_tried_configurational_MC_moves += 1; - particle_inside_exclusion_range_touched = false; + + if (type < 0) { + throw std::domain_error("Parameter 'type_mc' must be >= 0"); + } + if (n_part < 0) { + throw std::domain_error( + "Parameter 'particle_number_to_be_changed' must be >= 0"); + } if (n_part == 0) { // reject return false; } + m_tried_configurational_MC_moves += 1; + particle_inside_exclusion_range_touched = false; + auto const n_particles_of_type = number_of_particles_with_type(type); if (n_part > n_particles_of_type) { // reject diff --git a/src/core/reaction_methods/ReactionAlgorithm.hpp b/src/core/reaction_methods/ReactionAlgorithm.hpp index 93e8783abbd..3ff739f0076 100644 --- a/src/core/reaction_methods/ReactionAlgorithm.hpp +++ b/src/core/reaction_methods/ReactionAlgorithm.hpp @@ -115,7 +115,7 @@ class ReactionAlgorithm { m_max_exclusion_range = max_exclusion_range; } - virtual int do_reaction(int reaction_steps); + void do_reaction(int reaction_steps); void check_reaction_method() const; void remove_constraint() { m_reaction_constraint = ReactionConstraint::NONE; } void set_cyl_constraint(double center_x, double center_y, double radius); diff --git a/src/core/reaction_methods/tests/ReactionAlgorithm_test.cpp b/src/core/reaction_methods/tests/ReactionAlgorithm_test.cpp index 6d231db12f8..f07ea69a343 100644 --- a/src/core/reaction_methods/tests/ReactionAlgorithm_test.cpp +++ b/src/core/reaction_methods/tests/ReactionAlgorithm_test.cpp @@ -36,6 +36,7 @@ #include #include +#include #include #include #include @@ -194,10 +195,15 @@ BOOST_AUTO_TEST_CASE(ReactionAlgorithm_test) { place_particle(1, ref_positions[1]); set_particle_type(0, type_A); set_particle_type(1, type_A); - // check early exit when a MC move cannot be performed - BOOST_REQUIRE(!r_algo.displacement_move_for_particles_of_type(type_C, 1)); - BOOST_REQUIRE(!r_algo.displacement_move_for_particles_of_type(type_B, 2)); - BOOST_REQUIRE(!r_algo.displacement_move_for_particles_of_type(type_A, 0)); + // check runtime errors when a MC move cannot be physically performed + auto displacement_move = + std::bind(&ReactionAlgorithm::displacement_move_for_particles_of_type, + &r_algo, std::placeholders::_1, std::placeholders::_2); + BOOST_REQUIRE(!displacement_move(type_C, 1)); + BOOST_REQUIRE(!displacement_move(type_B, 2)); + BOOST_REQUIRE(!displacement_move(type_A, 0)); + BOOST_CHECK_THROW(displacement_move(type_A, -2), std::domain_error); + BOOST_CHECK_THROW(displacement_move(-2, 1), std::domain_error); // force all MC moves to be rejected by picking particles inside // their exclusion radius r_algo.exclusion_range = box_l; diff --git a/testsuite/python/reaction_methods_interface.py b/testsuite/python/reaction_methods_interface.py index 8b35714dd8a..1551b9b3b35 100644 --- a/testsuite/python/reaction_methods_interface.py +++ b/testsuite/python/reaction_methods_interface.py @@ -103,6 +103,8 @@ def check_reaction_parameters(reactions, parameters): self.assertAlmostEqual(method.constant_pH, 8., delta=1e-10) self.assertFalse(method.displacement_mc_move_for_particles_of_type( type_mc=0, particle_number_to_be_changed=0)) + self.assertFalse(method.displacement_mc_move_for_particles_of_type( + type_mc=0, particle_number_to_be_changed=100000)) # check constraints method.set_wall_constraints_in_z_direction( @@ -291,6 +293,12 @@ def test_exceptions(self): with self.assertRaisesRegex(ValueError, "Invalid value for 'kT'"): espressomd.reaction_methods.ReactionEnsemble( kT=-1., exclusion_range=1., seed=12) + with self.assertRaisesRegex(ValueError, "Parameter 'particle_number_to_be_changed' must be >= 0"): + method.displacement_mc_move_for_particles_of_type( + type_mc=0, particle_number_to_be_changed=-1) + with self.assertRaisesRegex(ValueError, "Parameter 'type_mc' must be >= 0"): + method.displacement_mc_move_for_particles_of_type( + type_mc=-1, particle_number_to_be_changed=1) # check invalid exclusion ranges and radii with self.assertRaisesRegex(ValueError, "Invalid value for 'exclusion_range'"):