Rewrite LB/EK interface

* remove the actor list * enforce SOLID principles * add fine-grained event hooks
espressomd · Aug 15, 2023 · 3c7a5d3 · 3c7a5d3
1 parent dc332cb
commit 3c7a5d3
Show file tree

Hide file tree

Showing 125 changed files with 2,527 additions and 1,842 deletions.
diff --git a/doc/sphinx/advanced_methods.rst b/doc/sphinx/advanced_methods.rst
@@ -460,7 +460,7 @@ Specification of fluid and movement
 
     lbf = espressomd.lb.LBFluidWalberla(agrid=1, density=1.0, kinematic_viscosity=1.5,
                                         tau=time_step, ext_force_density=[0.002, 0.0, 0.0])
-    system.actors.add(lbf)
+    self.system.lb = lbf
 
 This part of the script specifies the fluid that will get the system
 moving. Here ``agrid`` :math:`=\Delta x` is the spatial discretisation

diff --git a/doc/sphinx/ek.rst b/doc/sphinx/ek.rst
@@ -161,29 +161,34 @@ Here is a minimal working example::
     system.time_step = 0.01
     system.cell_system.skin = 1.0
 
-    ek_lattice = espressomd.electrokinetics.LatticeWalberla(agrid=0.5, n_ghost_layers=1)
-    ek_solver = espressomd.electrokinetics.EKNone(lattice=ek_lattice)
-    system.ekcontainer.solver = ek_solver
-    system.ekcontainer.tau = system.time_step
+    lattice = espressomd.electrokinetics.LatticeWalberla(agrid=0.5, n_ghost_layers=1)
+    ek_solver = espressomd.electrokinetics.EKNone(lattice=lattice)
+    system.ekcontainer = espressomd.electrokinetics.EKContainer(
+        solver=ek_solver, tau=system.time_step)
 
 where ``system.ekcontainer`` is the EK system, ``ek_solver`` is the Poisson
 solver (here ``EKNone`` doesn't actually solve the electrostatic field, but
-instead imposes a zero field), and ``ek_lattice`` contains the grid parameters.
+instead imposes a zero field), and ``lattice`` contains the grid parameters.
 In this setup, the EK system doesn't contain any species. The following
 sections will show how to add species that can diffuse, advect, react and/or
 electrostatically interact. An EK system can be set up at the same time as a
 LB system.
 
+To detach an EK system, use the following syntax::
+
+    system.ekcontainer = None
+
 .. _Diffusive species:
 
 Diffusive species
 ^^^^^^^^^^^^^^^^^
 ::
 
     ek_species = espressomd.electrokinetics.EKSpecies(
-        lattice=ek_lattice,
+        lattice=lattice,
         single_precision=False,
         kT=1.0,
+        tau=system.time_step,
         density=0.85,
         valency=0.0,
         diffusion=0.1,
@@ -234,8 +239,8 @@ Checkpointing
 
 ::
 
-    ek.save_checkpoint(path, binary)
-    ek.load_checkpoint(path, binary)
+    ek_species.save_checkpoint(path, binary)
+    ek_species.ekcontainer.load_checkpoint(path, binary)
 
 The first command saves all of the EK nodes' properties to an ASCII
 (``binary=False``) or binary (``binary=True``) format respectively.
@@ -260,7 +265,7 @@ one or multiple fluid field data into a single file using
     # create a VTK callback that automatically writes every 10 EK steps
     ek_vtk = espressomd.electrokinetics.VTKOutput(
         identifier="ek_vtk_automatic", observables=vtk_obs, delta_N=10)
-    ek.add_vtk_writer(vtk=ek_vtk)
+    ek_species.add_vtk_writer(vtk=ek_vtk)
     system.integrator.run(100)
     # can be deactivated
     ek_vtk.disable()
@@ -269,7 +274,7 @@ one or multiple fluid field data into a single file using
     # create a VTK callback that writes only when explicitly called
     ek_vtk_on_demand = espressomd.electrokinetics.VTKOutput(
         identifier="ek_vtk_now", observables=vtk_obs)
-    ek.add_vtk_writer(vtk=ek_vtk_on_demand)
+    ek_species.add_vtk_writer(vtk=ek_vtk_on_demand)
     ek_vtk_on_demand.write()
 
 Currently only supports the species density.
@@ -306,12 +311,14 @@ is available through :class:`~espressomd.io.vtk.VTKReader`::
     system.cell_system.skin = 0.4
     system.time_step = 0.1
 
-    lattice = espressomd.electrokinetics.LatticeWalberla(agrid=1.)
-    species = espressomd.electrokinetics.EKSpecies(
-            lattice=lattice, density=1., kT=1., diffusion=0.1, valency=0.,
-            advection=False, friction_coupling=False, tau=system.time_step)
-    system.ekcontainer.tau = species.tau
-    system.ekcontainer.add(species)
+    lattice = espressomd.electrokinetics.LatticeWalberla(agrid=1., n_ghost_layers=1)
+    ek_solver = espressomd.electrokinetics.EKNone(lattice=lattice)
+    ek_species = espressomd.electrokinetics.EKSpecies(
+        lattice=lattice, density=1., kT=1., diffusion=0.1, valency=0.,
+        advection=False, friction_coupling=False, tau=system.time_step)
+    system.ekcontainer = espressomd.electrokinetics.EKContainer(
+        solver=ek_solver, tau=ek_species.tau)
+    system.ekcontainer.add(ek_species)
     system.integrator.run(10)
 
     vtk_reader = espressomd.io.vtk.VTKReader()
@@ -327,15 +334,15 @@ is available through :class:`~espressomd.io.vtk.VTKReader`::
             identifier=label_vtk, delta_N=0,
             observables=["density"],
             base_folder=str(path_vtk_root))
-        species.add_vtk_writer(vtk=ek_vtk)
+        ek_species.add_vtk_writer(vtk=ek_vtk)
         ek_vtk.write()
 
         # read VTK file
         vtk_grids = vtk_reader.parse(path_vtk)
         vtk_density = vtk_grids[label_density]
 
         # check VTK values match node values
-        ek_density = np.copy(lbf[:, :, :].density)
+        ek_density = np.copy(ek_species[:, :, :].density)
         np.testing.assert_allclose(vtk_density, ek_density, rtol=1e-10, atol=0.)
 
 .. _Setting up EK boundary conditions:
@@ -361,17 +368,19 @@ One can set (or update) the boundary conditions of individual nodes::
     system.cell_system.skin = 0.1
     system.time_step = 0.01
     lattice = espressomd.electrokinetics.LatticeWalberla(agrid=0.5, n_ghost_layers=1)
+    ek_solver = espressomd.electrokinetics.EKNone(lattice=lattice)
     ek_species = espressomd.electrokinetics.EKSpecies(
-        kT=1.5, lattice=self.lattice, density=0.85, valency=0., diffusion=0.1,
+        kT=1.5, lattice=lattice, density=0.85, valency=0., diffusion=0.1,
         advection=False, friction_coupling=False, tau=system.time_step)
-    system.ekcontainer.tau = species.tau
+    system.ekcontainer = espressomd.electrokinetics.EKContainer(
+        solver=ek_solver, tau=ek_species.tau)
     system.ekcontainer.add(ek_species)
     # set node fixed density boundary conditions
-    lbf[0, 0, 0].boundary = espressomd.electrokinetics.DensityBoundary(1.)
+    ek_species[0, 0, 0].boundary = espressomd.electrokinetics.DensityBoundary(1.)
     # update node fixed density boundary conditions
-    lbf[0, 0, 0].boundary = espressomd.electrokinetics.DensityBoundary(2.)
+    ek_species[0, 0, 0].boundary = espressomd.electrokinetics.DensityBoundary(2.)
     # remove node boundary conditions
-    lbf[0, 0, 0].boundary = None
+    ek_species[0, 0, 0].boundary = None
 
 .. _Shape-based EK boundary conditions:
 
@@ -387,10 +396,12 @@ Adding a shape-based boundary is straightforward::
     system.cell_system.skin = 0.1
     system.time_step = 0.01
     lattice = espressomd.electrokinetics.LatticeWalberla(agrid=0.5, n_ghost_layers=1)
+    ek_solver = espressomd.electrokinetics.EKNone(lattice=lattice)
     ek_species = espressomd.electrokinetics.EKSpecies(
-        kT=1.5, lattice=self.lattice, density=0.85, valency=0.0, diffusion=0.1,
+        kT=1.5, lattice=lattice, density=0.85, valency=0.0, diffusion=0.1,
         advection=False, friction_coupling=False, tau=system.time_step)
-    system.ekcontainer.tau = species.tau
+    system.ekcontainer = espressomd.electrokinetics.EKContainer(
+        solver=ek_solver, tau=ek_species.tau)
     system.ekcontainer.add(ek_species)
     # set fixed density boundary conditions
     wall = espressomd.shapes.Wall(normal=[1., 0., 0.], dist=2.5)

diff --git a/doc/sphinx/integration.rst b/doc/sphinx/integration.rst
@@ -659,7 +659,7 @@ parameter ``gamma``. To enable the LB thermostat, use::
     import espressomd.lb
     system = espressomd.System(box_l=[1, 1, 1])
     lbf = espressomd.lb.LBFluidWalberla(agrid=1, density=1, kinematic_viscosity=1, tau=0.01)
-    system.actors.add(lbf)
+    self.system.lb = lbf
     system.thermostat.set_lb(LB_fluid=lbf, seed=123, gamma=1.5)
 
 No other thermostatting mechanism is necessary

diff --git a/doc/sphinx/lb.rst b/doc/sphinx/lb.rst
@@ -41,13 +41,13 @@ The following minimal example illustrates how to use the LBM in |es|::
     system = espressomd.System(box_l=[10, 20, 30])
     system.time_step = 0.01
     system.cell_system.skin = 0.4
-    lb = espressomd.lb.LBFluidWalberla(agrid=1.0, density=1.0, kinematic_viscosity=1.0, tau=0.01)
-    system.actors.add(lb)
+    lbf = espressomd.lb.LBFluidWalberla(agrid=1.0, density=1.0, kinematic_viscosity=1.0, tau=0.01)
+    system.lb = lbf
     system.integrator.run(100)
 
 To use the GPU-accelerated variant, replace line 6 in the example above by::
 
-    lb = espressomd.lb.LBFluidWalberlaGPU(agrid=1.0, density=1.0, kinematic_viscosity=1.0, tau=0.01)
+    lbf = espressomd.lb.LBFluidWalberlaGPU(agrid=1.0, density=1.0, kinematic_viscosity=1.0, tau=0.01)
 
 .. note:: Feature ``CUDA`` required for the GPU-accelerated variant
 
@@ -99,6 +99,11 @@ units to be applied to the fluid.
 Before running a simulation at least the following parameters must be
 set up: ``agrid``, ``tau``, ``kinematic_viscosity``, ``density``.
 
+To detach the LBM solver, use this syntax::
+
+    system.lb = None
+
+
 Performance considerations
 ^^^^^^^^^^^^^^^^^^^^^^^^^^
 
@@ -321,7 +326,7 @@ is available through :class:`~espressomd.io.vtk.VTKReader`::
 
     lbf = espressomd.lb.LBFluidWalberla(
         agrid=1., tau=0.1, density=1., kinematic_viscosity=1.)
-    system.actors.add(lbf)
+    system.lb = lbf
     system.integrator.run(10)
 
     vtk_reader = espressomd.io.vtk.VTKReader()
@@ -378,8 +383,8 @@ of the LBM in analogy to the example for the CPU given in section
     system = espressomd.System(box_l=[10, 20, 30])
     system.time_step = 0.01
     system.cell_system.skin = 0.4
-    lb = espressomd.lb.LBFluidWalberlaGPU(agrid=1.0, density=1.0, kinematic_viscosity=1.0, tau=0.01)
-    system.actors.add(lb)
+    lbf = espressomd.lb.LBFluidWalberlaGPU(agrid=1.0, density=1.0, kinematic_viscosity=1.0, tau=0.01)
+    system.lb = lbf
     system.integrator.run(100)
 
 .. _Electrohydrodynamics:
@@ -429,7 +434,7 @@ One can set (or update) the slip velocity of individual nodes::
     system.cell_system.skin = 0.1
     system.time_step = 0.01
     lbf = espressomd.lb.LBFluidWalberla(agrid=0.5, density=1.0, kinematic_viscosity=1.0, tau=0.01)
-    system.actors.add(lbf)
+    system.lb = lbf
     # make one node a boundary node with a slip velocity
     lbf[0, 0, 0].boundary = espressomd.lb.VelocityBounceBack([0, 0, 1])
     # update node for no-slip boundary conditions
@@ -450,7 +455,7 @@ Adding a shape-based boundary is straightforward::
     system.cell_system.skin = 0.1
     system.time_step = 0.01
     lbf = espressomd.lb.LBFluidWalberla(agrid=0.5, density=1.0, kinematic_viscosity=1.0, tau=0.01)
-    system.actors.add(lbf)
+    system.lb = lbf
     # set up shear flow between two sliding walls
     wall1 = espressomd.shapes.Wall(normal=[+1., 0., 0.], dist=2.5)
     lbf.add_boundary_from_shape(shape=wall1, velocity=[0., +0.05, 0.])

diff --git a/doc/tutorials/active_matter/active_matter.ipynb b/doc/tutorials/active_matter/active_matter.ipynb
@@ -903,7 +903,7 @@
     "                                    density=HYDRO_PARAMS['dens'],\n",
     "                                    kinematic_viscosity=HYDRO_PARAMS['visc'],\n",
     "                                    tau=HYDRO_PARAMS['time_step'])\n",
-    "system.actors.add(lbf)\n",
+    "system.lb = lbf\n",
     "system.thermostat.set_lb(LB_fluid=lbf, gamma=HYDRO_PARAMS['gamma'], seed=42)\n",
     "```"
    ]

diff --git a/doc/tutorials/electrokinetics/electrokinetics.ipynb b/doc/tutorials/electrokinetics/electrokinetics.ipynb
@@ -354,6 +354,7 @@
    "source": [
     "# Initializing espresso modules and the numpy package\n",
     "import espressomd\n",
+    "import espressomd.lb\n",
     "import espressomd.electrokinetics\n",
     "import espressomd.shapes\n",
     "\n",
@@ -451,8 +452,10 @@
    "outputs": [],
    "source": [
     "viscosity_kinematic = viscosity_dynamic / density_water\n",
-    "lbf = espressomd.lb.LBFluidWalberla(lattice=lattice, density=density_water, kinematic_viscosity=viscosity_kinematic, tau=dt, single_precision=single_precision)\n",
-    "system.actors.add(lbf)"
+    "lbf = espressomd.lb.LBFluidWalberla(lattice=lattice, density=density_water,\n",
+    "                                    kinematic_viscosity=viscosity_kinematic,\n",
+    "                                    tau=dt, single_precision=single_precision)\n",
+    "system.lb = lbf"
    ]
   },
   {
@@ -461,10 +464,11 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "eksolver = espressomd.electrokinetics.EKFFT(lattice=lattice, permittivity=permittivity, single_precision=single_precision)\n",
+    "eksolver = espressomd.electrokinetics.EKFFT(lattice=lattice,\n",
+    "                                            permittivity=permittivity,\n",
+    "                                            single_precision=single_precision)\n",
     "\n",
-    "system.ekcontainer.solver = eksolver\n",
-    "system.ekcontainer.tau = dt"
+    "system.ekcontainer = espressomd.electrokinetics.EKContainer(solver=eksolver, tau=dt)"
    ]
   },
   {
@@ -474,7 +478,10 @@
    "outputs": [],
    "source": [
     "density_counterions = -2.0 * sigma / width\n",
-    "ekspecies = espressomd.electrokinetics.EKSpecies(lattice=lattice, density=0.0, kT=kT, diffusion=D, valency=valency, advection=True, friction_coupling=True, ext_efield=ext_force_density, single_precision=single_precision, tau=dt)\n",
+    "ekspecies = espressomd.electrokinetics.EKSpecies(\n",
+    "    lattice=lattice, density=0.0, kT=kT, diffusion=D, valency=valency,\n",
+    "    advection=True, friction_coupling=True, ext_efield=ext_force_density,\n",
+    "    single_precision=single_precision, tau=dt)\n",
     "system.ekcontainer.add(ekspecies)"
    ]
   },
@@ -484,7 +491,10 @@
    "metadata": {},
    "outputs": [],
    "source": [
-    "ekwallcharge = espressomd.electrokinetics.EKSpecies(lattice=lattice, density=0.0, kT=kT, diffusion=0., valency=-valency, advection=False, friction_coupling=False, ext_efield=[0, 0, 0], single_precision=single_precision, tau=dt)\n",
+    "ekwallcharge = espressomd.electrokinetics.EKSpecies(\n",
+    "    lattice=lattice, density=0.0, kT=kT, diffusion=0., valency=-valency,\n",
+    "    advection=False, friction_coupling=False, ext_efield=[0, 0, 0],\n",
+    "    single_precision=single_precision, tau=dt)\n",
     "system.ekcontainer.add(ekwallcharge)"
    ]
   },
@@ -527,8 +537,10 @@
    "outputs": [],
    "source": [
     "for shape_obj in (wall_left, wall_right):\n",
-    "    ekspecies.add_boundary_from_shape(shape=shape_obj, value=[0., 0., 0.], boundary_type=espressomd.electrokinetics.FluxBoundary)\n",
-    "    ekspecies.add_boundary_from_shape(shape=shape_obj, value=0.0, boundary_type=espressomd.electrokinetics.DensityBoundary)\n",
+    "    ekspecies.add_boundary_from_shape(shape=shape_obj, value=[0., 0., 0.],\n",
+    "                                      boundary_type=espressomd.electrokinetics.FluxBoundary)\n",
+    "    ekspecies.add_boundary_from_shape(shape=shape_obj, value=0.0,\n",
+    "                                      boundary_type=espressomd.electrokinetics.DensityBoundary)\n",
     "    lbf.add_boundary_from_shape(shape=shape_obj, velocity=[0., 0., 0.])"
    ]
   },

diff --git a/doc/tutorials/lattice_boltzmann/lattice_boltzmann_poiseuille_flow.ipynb b/doc/tutorials/lattice_boltzmann/lattice_boltzmann_poiseuille_flow.ipynb
@@ -105,7 +105,7 @@
     "solution2_first": true
    },
    "source": [
-    "Create a lattice-Boltzmann actor and append it to the list of system actors. Use the GPU implementation of LB.\n",
+    "Create a lattice-Boltzmann actor and attach it to the system. Use the GPU implementation of LB.\n",
     "\n",
     "You can refer to section [setting up a LB fluid](https://espressomd.github.io/doc/lb.html#setting-up-a-lb-fluid)\n",
     "in the user guide."
@@ -123,7 +123,7 @@
     "                                    kinematic_viscosity=VISCOSITY,\n",
     "                                    tau=TIME_STEP,\n",
     "                                    ext_force_density=FORCE_DENSITY)\n",
-    "system.actors.add(lbf)\n",
+    "system.lb = lbf\n",
     "```"
    ]
   },

diff --git a/doc/tutorials/lattice_boltzmann/lattice_boltzmann_sedimentation.ipynb b/doc/tutorials/lattice_boltzmann/lattice_boltzmann_sedimentation.ipynb
@@ -329,7 +329,7 @@
     "\n",
     "**Exercise:**\n",
     "* create an unthermalized lattice-Boltzmann object with viscosity 1, density 1, grid size equal to\n",
-    "  ``spacing`` and LB time step equal to the MD time step and add this object to the system list of actors\n",
+    "  ``spacing`` and LB time step equal to the MD time step and add this object to the system\n",
     "  ([user guide](https://espressomd.github.io/doc/lb.html#setting-up-a-lb-fluid))\n",
     "* activate particle coupling to the fluid by setting the LB thermostat with $\\gamma=15$\n",
     "  ([user guide](https://espressomd.github.io/doc/lb.html#coupling-lb-to-a-md-simulation))"
@@ -346,7 +346,7 @@
     "                                    density=1.,\n",
     "                                    kinematic_viscosity=1.,\n",
     "                                    tau=system.time_step, kT=0.)\n",
-    "system.actors.add(lbf)\n",
+    "system.lb = lbf\n",
     "system.thermostat.set_lb(LB_fluid=lbf, gamma=15., seed=0)\n",
     "```"
    ]

diff --git a/doc/tutorials/lattice_boltzmann/lattice_boltzmann_theory.ipynb b/doc/tutorials/lattice_boltzmann/lattice_boltzmann_theory.ipynb
@@ -273,7 +273,7 @@
     "# Initialize an LBFluidWalberla with the minimum set of valid parameters.\n",
     "lbf = espressomd.lb.LBFluidWalberla(agrid=1, density=10, kinematic_viscosity=.1, tau=0.01)\n",
     "# Activate the LB by adding it to the System's actor list.\n",
-    "system.actors.add(lbf)\n",
+    "system.lb = lbf\n",
     "```"
    ]
   },

diff --git a/doc/tutorials/polymers/polymers.ipynb b/doc/tutorials/polymers/polymers.ipynb
@@ -300,7 +300,7 @@
     "    lbf = espressomd.lb.LBFluidWalberla(kT=kT, seed=42, agrid=1, density=1,\n",
     "                                        kinematic_viscosity=5, tau=system.time_step,\n",
     "                                        single_precision=True)\n",
-    "    system.actors.add(lbf)\n",
+    "    system.lb = lbf\n",
     "    system.thermostat.set_lb(LB_fluid=lbf, gamma=gamma, seed=42)"
    ]
   },
@@ -440,7 +440,7 @@
     "\n",
     "    # reset system\n",
     "    system.part.clear()\n",
-    "    system.actors.clear()\n",
+    "    system.lb = None\n",
     "    system.thermostat.turn_off()\n",
     "    system.auto_update_accumulators.clear()\n",
     "\n",

diff --git a/doc/tutorials/raspberry_electrophoresis/raspberry_electrophoresis.ipynb b/doc/tutorials/raspberry_electrophoresis/raspberry_electrophoresis.ipynb
@@ -675,7 +675,7 @@
    "outputs": [],
    "source": [
     "system.thermostat.turn_off()\n",
-    "system.actors.add(lb)\n",
+    "system.lb = lb\n",
     "system.thermostat.set_lb(LB_fluid=lb, seed=123, gamma=20.0)"
    ]
   },

diff --git a/maintainer/benchmarks/lb.py b/maintainer/benchmarks/lb.py
@@ -138,7 +138,7 @@
     lb_class = espressomd.lb.LBFluidWalberlaGPU
 lbf = lb_class(agrid=agrid, tau=system.time_step, kinematic_viscosity=1.,
                density=1., single_precision=args.single_precision)
-system.actors.add(lbf)
+system.lb = lbf
 
 
 # time integration loop