Skip to content

Commit

Permalink
Treat Off and None protocols of LEbc as identical
Browse files Browse the repository at this point in the history
  • Loading branch information
jngrad committed Sep 10, 2024
1 parent ed55184 commit d3a1057
Show file tree
Hide file tree
Showing 4 changed files with 58 additions and 20 deletions.
5 changes: 4 additions & 1 deletion src/core/lb/LBWalberla.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -38,6 +38,7 @@
#include <utils/math/int_pow.hpp>

#include <optional>
#include <variant>

namespace LB {

Expand Down Expand Up @@ -128,7 +129,9 @@ void LBWalberla::update_collision_model(LBWalberlaBase &lb,
LBWalberlaParams &params, double kT,
unsigned int seed) {
auto const &system = ::System::get_system();
if (auto le_protocol = system.lees_edwards->get_protocol()) {
auto le_protocol = system.lees_edwards->get_protocol();
if (le_protocol and
not std::holds_alternative<LeesEdwards::Off>(*le_protocol)) {
if (kT != 0.) {
throw std::runtime_error(
"Lees-Edwards LB doesn't support thermalization");
Expand Down
52 changes: 40 additions & 12 deletions src/script_interface/lees_edwards/LeesEdwards.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -45,15 +45,16 @@ class LeesEdwards : public AutoParameters<LeesEdwards, System::Leaf> {
add_parameters(
{{"protocol",
[this](Variant const &value) {
auto &system = get_system();
if (is_none(value)) {
auto const &system = get_system();
context()->parallel_try_catch([&system]() {
context()->parallel_try_catch([&]() {
auto constexpr invalid_dir = LeesEdwardsBC::invalid_dir;
system.lb.lebc_sanity_checks(invalid_dir, invalid_dir);
});
m_protocol = nullptr;
system.box_geo->set_lees_edwards_bc(LeesEdwardsBC{});
m_lees_edwards->unset_protocol();
system.on_lees_edwards_change();
return;
}
context()->parallel_try_catch([this]() {
Expand All @@ -67,8 +68,15 @@ class LeesEdwards : public AutoParameters<LeesEdwards, System::Leaf> {
"activation via set_boundary_conditions()");
}
});
m_protocol = get_value<std::shared_ptr<Protocol>>(value);
m_lees_edwards->set_protocol(m_protocol->protocol());
auto const protocol = get_value<std::shared_ptr<Protocol>>(value);
context()->parallel_try_catch([&]() {
try {
set_protocol(protocol);
} catch (...) {
set_protocol(m_protocol);
throw;
}
});
},
[this]() {
if (m_protocol)
Expand All @@ -89,13 +97,13 @@ class LeesEdwards : public AutoParameters<LeesEdwards, System::Leaf> {
VariantMap const &params) override {
if (name == "set_boundary_conditions") {
context()->parallel_try_catch([this, &params]() {
auto const protocol = params.at("protocol");
if (is_none(protocol)) {
do_set_parameter("protocol", protocol);
auto const variant = params.at("protocol");
if (is_none(variant)) {
do_set_parameter("protocol", variant);
return;
}
// check input arguments
m_protocol = get_value<std::shared_ptr<Protocol>>(protocol);
auto const protocol = get_value<std::shared_ptr<Protocol>>(variant);
auto const shear_direction = get_shear_axis(params, "shear_direction");
auto const shear_plane_normal =
get_shear_axis(params, "shear_plane_normal");
Expand All @@ -106,10 +114,18 @@ class LeesEdwards : public AutoParameters<LeesEdwards, System::Leaf> {
auto &system = get_system();
system.lb.lebc_sanity_checks(shear_direction, shear_plane_normal);
// update box geometry and cell structure
system.box_geo->set_lees_edwards_bc(
LeesEdwardsBC{0., 0., shear_direction, shear_plane_normal});
m_lees_edwards->set_protocol(m_protocol->protocol());
system.on_lees_edwards_change();
auto const old_shear_direction = get_lebc().shear_direction;
auto const old_shear_plane_normal = get_lebc().shear_plane_normal;
try {
system.box_geo->set_lees_edwards_bc(
LeesEdwardsBC{0., 0., shear_direction, shear_plane_normal});
set_protocol(protocol);
} catch (...) {
system.box_geo->set_lees_edwards_bc(LeesEdwardsBC{
0., 0., old_shear_direction, old_shear_plane_normal});
set_protocol(m_protocol);
throw;
}
});
}
return {};
Expand Down Expand Up @@ -159,6 +175,18 @@ class LeesEdwards : public AutoParameters<LeesEdwards, System::Leaf> {
}
m_params.reset();
}

void set_protocol(std::shared_ptr<Protocol> const &protocol) {
auto &system = get_system();
if (protocol) {
m_lees_edwards->set_protocol(protocol->protocol());
} else {
system.box_geo->set_lees_edwards_bc(LeesEdwardsBC{});
m_lees_edwards->unset_protocol();
}
system.on_lees_edwards_change();
m_protocol = protocol;
}
};

} // namespace LeesEdwards
Expand Down
3 changes: 2 additions & 1 deletion testsuite/python/integrator_exceptions.py
Original file line number Diff line number Diff line change
Expand Up @@ -170,7 +170,8 @@ def test_npt_integrator(self):
shear_direction="x", shear_plane_normal="y", protocol=protocol)
self.system.integrator.run(0)
with self.assertRaisesRegex(Exception, self.msg + 'The NpT integrator cannot use Lees-Edwards'):
self.system.lees_edwards.protocol = espressomd.lees_edwards.Off()
self.system.lees_edwards.protocol = espressomd.lees_edwards.LinearShear(
shear_velocity=0., initial_pos_offset=0., time_0=0.)
self.system.integrator.run(0)
self.system.lees_edwards.protocol = None
self.system.integrator.run(0)
Expand Down
18 changes: 12 additions & 6 deletions testsuite/python/lb_lees_edwards.py
Original file line number Diff line number Diff line change
Expand Up @@ -41,9 +41,13 @@ def __init__(self, **kwargs):

def __enter__(self):
self.lbf = espressomd.lb.LBFluidWalberla(
agrid=1., density=1., kinematic_viscosity=1., tau=system.time_step, **self.kwargs)
agrid=1., density=1., kinematic_viscosity=1.,
tau=system.time_step, **self.kwargs)
system.lb = self.lbf
system.thermostat.set_lb(LB_fluid=self.lbf, gamma=1.0)
thmst_kwargs = {}
if "seed" in self.kwargs:
thmst_kwargs["seed"] = self.kwargs["seed"]
system.thermostat.set_lb(LB_fluid=self.lbf, gamma=1.0, **thmst_kwargs)
return self.lbf

def __exit__(self, exc_type, exc_value, exc_traceback):
Expand Down Expand Up @@ -299,15 +303,17 @@ def test_lebc_mismatch(self):
"""
Check that MD LEbc and LB LEbc always agree.
"""
err_msg = "MD and LB Lees-Edwards boundary conditions disagree"
err_msg = "Lees-Edwards LB only supports shear_plane_normal=\"y\""
# MD and LB LEbc must match
with self.assertRaisesRegex(RuntimeError, err_msg):
with self.assertRaisesRegex(ValueError, err_msg):
with LBContextManager() as lbf:
LEContextManager('y', 'x', 1.).initialize()
system.lees_edwards.protocol = None
with LBContextManager() as lbf:
LEContextManager('x', 'y', 1.).initialize()
# when a LB actor with LEbc is active, the MD LEbc shear directions
# are immutable
err_msg = "MD and LB Lees-Edwards boundary conditions disagree"
with LEContextManager('x', 'y', 1.):
with LBContextManager() as lbf:
with self.assertRaisesRegex(RuntimeError, err_msg):
Expand All @@ -320,7 +326,7 @@ def test_lebc_mismatch(self):
self.assertEqual(system.lees_edwards.shear_plane_normal, "y")
# when de-activating and later re-activating a LB actor with LEbc,
# the MD LEbc must have the same shear directions
with self.assertRaisesRegex(Exception, err_msg):
with self.assertRaisesRegex(RuntimeError, err_msg):
with LEContextManager('z', 'y', 1.):
system.lb = lbf
self.assertIsNone(system.lb)
Expand Down Expand Up @@ -348,7 +354,7 @@ def test_lebc_mismatch(self):
self.assertIsNone(system.lb)
with self.assertRaisesRegex(RuntimeError, "Lees-Edwards LB doesn't support thermalization"):
with LBContextManager(kT=1., seed=42) as lbf:
LEContextManager('x', 'y', 1.)
LEContextManager('x', 'y', 1.).initialize()
self.assertIsNone(system.lb)

with self.assertRaisesRegex(ValueError, "Lees-Edwards sweep is implemented for a ghost layer of thickness 1"):
Expand Down

0 comments on commit d3a1057

Please sign in to comment.