diff --git a/.gitlab-ci.yml b/.gitlab-ci.yml
index b4c3f70734..33d049e521 100644
--- a/.gitlab-ci.yml
+++ b/.gitlab-ci.yml
@@ -435,7 +435,8 @@ check_sphinx:
   <<: *global_job_definition
   stage: additional_checks
   needs:
-    - cuda12-maxset
+    - job: cuda12-maxset
+      artifacts: true
   when: on_success
   script:
     - cd ${CI_PROJECT_DIR}/build
@@ -458,7 +459,8 @@ run_tutorials:
   <<: *global_job_definition
   stage: additional_checks
   needs:
-    - cuda12-maxset
+    - job: cuda12-maxset
+      artifacts: true
   when: on_success
   script:
     - cd ${CI_PROJECT_DIR}/build
@@ -486,7 +488,8 @@ run_doxygen:
   <<: *global_job_definition
   stage: additional_checks
   needs:
-    - cuda12-maxset
+    - job: cuda12-maxset
+      artifacts: true
   when: on_success
   only:
     - python
@@ -508,7 +511,8 @@ maxset_no_gpu:
   stage: additional_checks
   when: on_success
   needs:
-    - cuda12-maxset
+    - job: cuda12-maxset
+      artifacts: true
   script:
     - export CUDA_VISIBLE_DEVICES=""
     - cd ${CI_PROJECT_DIR}/build
@@ -524,15 +528,32 @@ maxset_3_cores:
   stage: additional_checks
   when: on_success
   needs:
-    - cuda12-maxset
+    - job: cuda12-maxset
+      artifacts: false
+  variables:
+     CC: 'gcc-12'
+     CXX: 'g++-12'
+     GCOV: 'gcov-12'
+     myconfig: 'maxset'
+     with_cuda: 'true'
+     with_coverage: 'false'
+     with_coverage_python: 'false'
+     srcdir: '${CI_PROJECT_DIR}'
+     with_scafacos: 'true'
+     with_walberla: 'true'
+     with_walberla_avx: 'true'
+     with_stokesian_dynamics: 'true'
+     with_caliper: 'true'
+     check_odd_only: 'true'
+     cmake_params: '-D ESPRESSO_TEST_NP=3'
   script:
-    - cd ${CI_PROJECT_DIR}/build
-    - cmake -D ESPRESSO_TEST_NP=3 .
-    - make -t && make check_unit_tests && make check_python_parallel_odd
+    - bash maintainer/CI/build_cmake.sh
+  timeout: 90m
   tags:
     - espresso
     - cuda
     - numa
+    - avx2
     - reuse-artifacts-same-arch
 
 status_success:
diff --git a/src/core/lb/LBNone.hpp b/src/core/lb/LBNone.hpp
index 61248c5fe5..9ee33c07e8 100644
--- a/src/core/lb/LBNone.hpp
+++ b/src/core/lb/LBNone.hpp
@@ -67,6 +67,8 @@ struct LBNone {
   void on_node_grid_change() const { throw NoLBActive{}; }
   void on_timestep_change() const { throw NoLBActive{}; }
   void on_temperature_change() const { throw NoLBActive{}; }
+  void on_lees_edwards_change() const { throw NoLBActive{}; }
+  void update_collision_model() const { throw NoLBActive{}; }
 };
 
 } // namespace LB
diff --git a/src/core/lb/LBWalberla.cpp b/src/core/lb/LBWalberla.cpp
index d1dec13160..f31634f448 100644
--- a/src/core/lb/LBWalberla.cpp
+++ b/src/core/lb/LBWalberla.cpp
@@ -30,12 +30,15 @@
 #include "system/System.hpp"
 #include "thermostat.hpp"
 
+#include "lees_edwards/lees_edwards.hpp"
+
 #include <walberla_bridge/lattice_boltzmann/LBWalberlaBase.hpp>
 
 #include <utils/Vector.hpp>
 #include <utils/math/int_pow.hpp>
 
 #include <optional>
+#include <variant>
 
 namespace LB {
 
@@ -112,6 +115,45 @@ void LBWalberla::sanity_checks(System::System const &system) const {
                              system.get_time_step());
 }
 
+void LBWalberla::on_lees_edwards_change() { update_collision_model(); }
+
+void LBWalberla::update_collision_model() {
+  auto const energy_conversion =
+      Utils::int_pow<2>(lb_params->get_agrid() / lb_params->get_tau());
+  auto const kT = lb_fluid->get_kT() * energy_conversion;
+  auto const seed = lb_fluid->get_seed();
+  update_collision_model(*lb_fluid, *lb_params, kT, seed);
+}
+
+void LBWalberla::update_collision_model(LBWalberlaBase &lb,
+                                        LBWalberlaParams &params, double kT,
+                                        unsigned int seed) {
+  auto const &system = ::System::get_system();
+  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");
+    }
+    auto const &le_bc = system.box_geo->lees_edwards_bc();
+    auto lees_edwards_object = std::make_unique<LeesEdwardsPack>(
+        le_bc.shear_direction, le_bc.shear_plane_normal,
+        [&params, le_protocol, &system]() {
+          return get_pos_offset(system.get_sim_time(), *le_protocol) /
+                 params.get_agrid();
+        },
+        [&params, le_protocol, &system]() {
+          return get_shear_velocity(system.get_sim_time(), *le_protocol) *
+                 (params.get_tau() / params.get_agrid());
+        });
+    lb.set_collision_model(std::move(lees_edwards_object));
+  } else {
+    lb.set_collision_model(kT, seed);
+  }
+  lb.ghost_communication();
+}
+
 } // namespace LB
 
 #endif // WALBERLA
diff --git a/src/core/lb/LBWalberla.hpp b/src/core/lb/LBWalberla.hpp
index 14cf5d02f2..a5b6f77426 100644
--- a/src/core/lb/LBWalberla.hpp
+++ b/src/core/lb/LBWalberla.hpp
@@ -88,6 +88,11 @@ struct LBWalberla {
   }
   void on_timestep_change() const {}
   void on_temperature_change() const {}
+  void on_lees_edwards_change();
+  void update_collision_model();
+  static void update_collision_model(LBWalberlaBase &instance,
+                                     LBWalberlaParams &params, double kT,
+                                     unsigned int seed);
 };
 
 } // namespace LB
diff --git a/src/core/lb/Solver.cpp b/src/core/lb/Solver.cpp
index 009f953c57..69e733b4ed 100644
--- a/src/core/lb/Solver.cpp
+++ b/src/core/lb/Solver.cpp
@@ -140,6 +140,18 @@ void Solver::on_temperature_change() {
   }
 }
 
+void Solver::on_lees_edwards_change() {
+  if (impl->solver) {
+    std::visit([](auto &ptr) { ptr->on_lees_edwards_change(); }, *impl->solver);
+  }
+}
+
+void Solver::update_collision_model() {
+  if (impl->solver) {
+    std::visit([](auto &ptr) { ptr->update_collision_model(); }, *impl->solver);
+  }
+}
+
 bool Solver::is_gpu() const {
   check_solver(impl);
   return std::visit([](auto &ptr) { return ptr->is_gpu(); }, *impl->solver);
diff --git a/src/core/lb/Solver.hpp b/src/core/lb/Solver.hpp
index 6eefdc27cf..47d9d5e019 100644
--- a/src/core/lb/Solver.hpp
+++ b/src/core/lb/Solver.hpp
@@ -168,6 +168,8 @@ struct Solver : public System::Leaf<Solver> {
   void on_cell_structure_change();
   void on_timestep_change();
   void on_temperature_change();
+  void on_lees_edwards_change();
+  void update_collision_model();
   void veto_boxl_change() const;
 
 private:
diff --git a/src/core/system/System.cpp b/src/core/system/System.cpp
index 5243f82fc6..06cc3b9684 100644
--- a/src/core/system/System.cpp
+++ b/src/core/system/System.cpp
@@ -378,6 +378,8 @@ void System::on_observable_calc() {
   clear_particle_node();
 }
 
+void System::on_lees_edwards_change() { lb.on_lees_edwards_change(); }
+
 void System::update_local_geo() {
   *local_geo = LocalBox::make_regular_decomposition(
       box_geo->length(), ::communicator.calc_node_index(),
diff --git a/src/core/system/System.hpp b/src/core/system/System.hpp
index 55be810194..7c216a1633 100644
--- a/src/core/system/System.hpp
+++ b/src/core/system/System.hpp
@@ -253,6 +253,7 @@ class System : public std::enable_shared_from_this<System> {
    *  initialized immediately (P3M etc.).
    */
   void on_observable_calc();
+  void on_lees_edwards_change();
   void veto_boxl_change(bool skip_particle_checks = false) const;
   /**@}*/
 
diff --git a/src/core/unit_tests/lb_particle_coupling_test.cpp b/src/core/unit_tests/lb_particle_coupling_test.cpp
index 95c4244191..b42bea56eb 100644
--- a/src/core/unit_tests/lb_particle_coupling_test.cpp
+++ b/src/core/unit_tests/lb_particle_coupling_test.cpp
@@ -608,11 +608,13 @@ BOOST_AUTO_TEST_CASE(lb_exceptions) {
     BOOST_CHECK_THROW(lb.veto_kT(0.), NoLBActive);
     BOOST_CHECK_THROW(lb.lebc_sanity_checks(0u, 1u), NoLBActive);
     BOOST_CHECK_THROW(lb.propagate(), NoLBActive);
+    BOOST_CHECK_THROW(lb.update_collision_model(), NoLBActive);
     BOOST_CHECK_THROW(lb.on_cell_structure_change(), NoLBActive);
     BOOST_CHECK_THROW(lb.on_boxl_change(), NoLBActive);
     BOOST_CHECK_THROW(lb.on_node_grid_change(), NoLBActive);
     BOOST_CHECK_THROW(lb.on_timestep_change(), NoLBActive);
     BOOST_CHECK_THROW(lb.on_temperature_change(), NoLBActive);
+    BOOST_CHECK_THROW(lb.on_lees_edwards_change(), NoLBActive);
     BOOST_CHECK_THROW(lb_impl->get_density_at_pos(vec, true), NoLBActive);
     BOOST_CHECK_THROW(lb_impl->get_velocity_at_pos(vec, true), NoLBActive);
     BOOST_CHECK_THROW(lb_impl->add_force_at_pos(vec, vec), NoLBActive);
diff --git a/src/script_interface/lees_edwards/LeesEdwards.hpp b/src/script_interface/lees_edwards/LeesEdwards.hpp
index d4ddd81fb9..6efe545222 100644
--- a/src/script_interface/lees_edwards/LeesEdwards.hpp
+++ b/src/script_interface/lees_edwards/LeesEdwards.hpp
@@ -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]() {
@@ -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)
@@ -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");
@@ -103,12 +111,21 @@ class LeesEdwards : public AutoParameters<LeesEdwards, System::Leaf> {
           throw std::invalid_argument("Parameters 'shear_direction' and "
                                       "'shear_plane_normal' must differ");
         }
-        auto const &system = get_system();
+        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());
+        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 {};
@@ -158,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
diff --git a/src/script_interface/walberla/LBFluid.cpp b/src/script_interface/walberla/LBFluid.cpp
index 0cad46f999..ed4d08674c 100644
--- a/src/script_interface/walberla/LBFluid.cpp
+++ b/src/script_interface/walberla/LBFluid.cpp
@@ -158,7 +158,7 @@ void LBFluid::do_construct(VariantMap const &params) {
   auto const ext_f = get_value<Utils::Vector3d>(params, "ext_force_density");
   m_lb_params = std::make_shared<::LB::LBWalberlaParams>(agrid, tau);
   m_is_active = false;
-  m_seed = get_value<int>(params, "seed");
+  auto const seed = get_value<int>(params, "seed");
   context()->parallel_try_catch([&]() {
     if (tau <= 0.) {
       throw std::domain_error("Parameter 'tau' must be > 0");
@@ -175,7 +175,7 @@ void LBFluid::do_construct(VariantMap const &params) {
     auto const lb_dens = m_conv_dens * dens;
     auto const lb_kT = m_conv_energy * kT;
     auto const lb_ext_f = m_conv_force_dens * ext_f;
-    if (m_seed < 0) {
+    if (seed < 0) {
       throw std::domain_error("Parameter 'seed' must be >= 0");
     }
     if (lb_kT < 0.) {
@@ -188,28 +188,8 @@ void LBFluid::do_construct(VariantMap const &params) {
       throw std::domain_error("Parameter 'kinematic_viscosity' must be >= 0");
     }
     make_instance(params);
-    auto const &system = ::System::get_system();
-    if (auto le_protocol = system.lees_edwards->get_protocol()) {
-      if (lb_kT != 0.) {
-        throw std::runtime_error(
-            "Lees-Edwards LB doesn't support thermalization");
-      }
-      auto const &le_bc = system.box_geo->lees_edwards_bc();
-      auto lees_edwards_object = std::make_unique<LeesEdwardsPack>(
-          le_bc.shear_direction, le_bc.shear_plane_normal,
-          [this, le_protocol, &system]() {
-            return get_pos_offset(system.get_sim_time(), *le_protocol) /
-                   m_lb_params->get_agrid();
-          },
-          [this, le_protocol, &system]() {
-            return get_shear_velocity(system.get_sim_time(), *le_protocol) *
-                   (m_lb_params->get_tau() / m_lb_params->get_agrid());
-          });
-      m_instance->set_collision_model(std::move(lees_edwards_object));
-    } else {
-      m_instance->set_collision_model(lb_kT, m_seed);
-    }
-    m_instance->ghost_communication(); // synchronize ghost layers
+    ::LB::LBWalberla::update_collision_model(*m_instance, *m_lb_params, lb_kT,
+                                             static_cast<unsigned int>(seed));
     m_instance->set_external_force(lb_ext_f);
     for (auto &vtk : m_vtk_writers) {
       vtk->attach_to_lattice(m_instance, get_latice_to_md_units_conversion());
diff --git a/src/script_interface/walberla/LBFluid.hpp b/src/script_interface/walberla/LBFluid.hpp
index 036f408729..d8a9fab22f 100644
--- a/src/script_interface/walberla/LBFluid.hpp
+++ b/src/script_interface/walberla/LBFluid.hpp
@@ -51,7 +51,6 @@ class LBFluid : public LatticeModel<::LBWalberlaBase, LBVTKHandle> {
   using Base = LatticeModel<::LBWalberlaBase, LBVTKHandle>;
   std::shared_ptr<::LB::LBWalberlaParams> m_lb_params;
   bool m_is_active;
-  int m_seed;
   double m_conv_dist;
   double m_conv_visc;
   double m_conv_dens;
@@ -77,7 +76,8 @@ class LBFluid : public LatticeModel<::LBWalberlaBase, LBVTKHandle> {
          [this]() { return m_instance->get_lattice().get_grid_dimensions(); }},
         {"kT", AutoParameter::read_only,
          [this]() { return m_instance->get_kT() / m_conv_energy; }},
-        {"seed", AutoParameter::read_only, [this]() { return m_seed; }},
+        {"seed", AutoParameter::read_only,
+         [this]() { return static_cast<int>(m_instance->get_seed()); }},
         {"rng_state",
          [this](Variant const &v) {
            auto const rng_state = get_value<int>(v);
diff --git a/src/walberla_bridge/include/walberla_bridge/lattice_boltzmann/LBWalberlaBase.hpp b/src/walberla_bridge/include/walberla_bridge/lattice_boltzmann/LBWalberlaBase.hpp
index c0dd6406cc..bff71a129f 100644
--- a/src/walberla_bridge/include/walberla_bridge/lattice_boltzmann/LBWalberlaBase.hpp
+++ b/src/walberla_bridge/include/walberla_bridge/lattice_boltzmann/LBWalberlaBase.hpp
@@ -257,10 +257,13 @@ class LBWalberlaBase : public LatticeModel {
   /** @brief Get the fluid temperature (if thermalized). */
   virtual double get_kT() const noexcept = 0;
 
+  /** @brief Get the RNG seed (if thermalized). */
+  virtual unsigned int get_seed() const noexcept = 0;
+
   /** @brief Set the RNG counter (if thermalized). */
   [[nodiscard]] virtual std::optional<uint64_t> get_rng_state() const = 0;
 
-  /** @brief Set the rng state of thermalized LBs */
+  /** @brief Set the RNG counter (if thermalized). */
   virtual void set_rng_state(uint64_t counter) = 0;
 
   /** @brief get the velocity field id */
diff --git a/src/walberla_bridge/src/lattice_boltzmann/LBWalberlaImpl.hpp b/src/walberla_bridge/src/lattice_boltzmann/LBWalberlaImpl.hpp
index 1f1a8bb24e..6105b996c6 100644
--- a/src/walberla_bridge/src/lattice_boltzmann/LBWalberlaImpl.hpp
+++ b/src/walberla_bridge/src/lattice_boltzmann/LBWalberlaImpl.hpp
@@ -257,6 +257,7 @@ class LBWalberlaImpl : public LBWalberlaBase {
   FloatType m_viscosity; /// kinematic viscosity
   FloatType m_density;
   FloatType m_kT;
+  unsigned int m_seed;
 
   // Block data access handles
   BlockDataID m_pdf_field_id;
@@ -401,7 +402,7 @@ class LBWalberlaImpl : public LBWalberlaBase {
   LBWalberlaImpl(std::shared_ptr<LatticeWalberla> lattice, double viscosity,
                  double density)
       : m_viscosity(FloatType_c(viscosity)), m_density(FloatType_c(density)),
-        m_kT(FloatType{0}), m_lattice(std::move(lattice)) {
+        m_kT(FloatType{0}), m_seed(0u), m_lattice(std::move(lattice)) {
 
     auto const &blocks = m_lattice->get_blocks();
     auto const n_ghost_layers = m_lattice->get_ghost_layers();
@@ -592,6 +593,7 @@ class LBWalberlaImpl : public LBWalberlaBase {
     auto const omega_odd = odd_mode_relaxation_rate(omega);
     auto const blocks = get_lattice().get_blocks();
     m_kT = FloatType_c(kT);
+    m_seed = seed;
     auto obj = CollisionModelThermalized(m_last_applied_force_field_id,
                                          m_pdf_field_id, m_kT, omega, omega,
                                          omega_odd, omega, seed, uint32_t{0u});
@@ -1366,6 +1368,10 @@ class LBWalberlaImpl : public LBWalberlaBase {
     return numeric_cast<double>(m_kT);
   }
 
+  [[nodiscard]] unsigned int get_seed() const noexcept override {
+    return m_seed;
+  }
+
   [[nodiscard]] std::optional<uint64_t> get_rng_state() const override {
     auto const cm = std::get_if<CollisionModelThermalized>(&*m_collision_model);
     if (!cm or m_kT == 0.) {
diff --git a/testsuite/python/integrator_exceptions.py b/testsuite/python/integrator_exceptions.py
index 6e11c04c36..77953f8487 100644
--- a/testsuite/python/integrator_exceptions.py
+++ b/testsuite/python/integrator_exceptions.py
@@ -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)
diff --git a/testsuite/python/lb_lees_edwards.py b/testsuite/python/lb_lees_edwards.py
index 9fd79a6366..077586f716 100644
--- a/testsuite/python/lb_lees_edwards.py
+++ b/testsuite/python/lb_lees_edwards.py
@@ -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):
@@ -191,8 +195,6 @@ def test_velocity_shift_from_particles(self):
             with LBContextManager() as lbf:
                 for profile in self.sample_lb_velocities(lbf):
                     self.check_profile(profile, stencil_D2Q8, '', 'SNWE', tol)
-
-        # order of instantiation doesn't matter
         with LBContextManager() as lbf:
             with LEContextManager('x', 'y', 0):
                 for profile in self.sample_lb_velocities(lbf):
@@ -208,6 +210,13 @@ def test_velocity_shift_from_particles(self):
             with LBContextManager() as lbf:
                 for profile in self.sample_lb_velocities(lbf):
                     self.check_profile(profile, stencil, 'SN', 'WE', tol)
+        with LBContextManager() as lbf:
+            with LEContextManager('x', 'y', le_offset):
+                stencil = {'N~': (8 - le_offset, 0),
+                           'S~': (8 + le_offset, 16),
+                           **stencil_D2Q8}
+                for profile in self.sample_lb_velocities(lbf):
+                    self.check_profile(profile, stencil, 'SN', 'WE', tol)
 
         # TODO: re-enable this check once LB can be sheared in any direction
 #        # East and West are sheared vertically
@@ -253,8 +262,6 @@ def create_impulse(lbf, stencil):
                 create_impulse(lbf, stencil_D2Q8)
                 for profile in self.sample_lb_velocities(lbf):
                     self.check_profile(profile, stencil_D2Q8, '', 'SNWE', tol)
-
-        # order of instantiation doesn't matter
         with LBContextManager() as lbf:
             with LEContextManager('x', 'y', 0):
                 create_impulse(lbf, stencil_D2Q8)
@@ -272,6 +279,14 @@ def create_impulse(lbf, stencil):
                 create_impulse(lbf, stencil_D2Q8)
                 for profile in self.sample_lb_velocities(lbf):
                     self.check_profile(profile, stencil, 'SN', 'WE', tol)
+        with LBContextManager() as lbf:
+            with LEContextManager('x', 'y', le_offset):
+                stencil = {'N~': (8 - le_offset, 1),
+                           'S~': (8 + le_offset, 15),
+                           **stencil_D2Q8}
+                create_impulse(lbf, stencil_D2Q8)
+                for profile in self.sample_lb_velocities(lbf):
+                    self.check_profile(profile, stencil, 'SN', 'WE', tol)
 
         # TODO: re-enable this check once LB can be sheared in any direction
 #        # East and West are sheared vertically
@@ -288,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):
@@ -309,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)
@@ -329,12 +346,27 @@ def test_lebc_mismatch(self):
             system.lb = lbf
             system.lb = None
         # no thermalization
+        system.lees_edwards.protocol = None
         with self.assertRaisesRegex(RuntimeError, "Lees-Edwards LB doesn't support thermalization"):
             with LEContextManager('x', 'y', 1.):
                 system.lb = espressomd.lb.LBFluidWalberla(
                     agrid=1., density=1., kinematic_viscosity=1., kT=1., seed=42,
                     tau=system.time_step)
         self.assertIsNone(system.lb)
+        system.lees_edwards.protocol = None
+        with self.assertRaisesRegex(RuntimeError, "Lees-Edwards LB doesn't support thermalization"):
+            with LBContextManager(kT=1., seed=42) as lbf:
+                LEContextManager('x', 'y', 1.).initialize()
+        self.assertIsNone(system.lb)
+        system.lees_edwards.protocol = None
+        with self.assertRaisesRegex(RuntimeError, "Lees-Edwards LB doesn't support thermalization"):
+            with LBContextManager(kT=1., seed=42) as lbf:
+                system.lees_edwards.set_boundary_conditions(
+                    shear_direction='x', shear_plane_normal='y',
+                    protocol=espressomd.lees_edwards.Off())
+                system.lees_edwards.protocol = espressomd.lees_edwards.LinearShear(
+                    shear_velocity=0., initial_pos_offset=0., time_0=0.)
+        self.assertIsNone(system.lb)
 
         with self.assertRaisesRegex(ValueError, "Lees-Edwards sweep is implemented for a ghost layer of thickness 1"):
             lattice = espressomd.lb.LatticeWalberla(agrid=1., n_ghost_layers=2)