Poisseuille flow is the flow through a pipe or (in our case) a slit under a homogenous force density, e.g. gravity. In the limit of small Reynolds numbers, the flow can be described with the Stokes equation. We assume the slit being infinitely extended in $y$ and $z$ direction and a force density $f$ on the fluid in $y$ direction. No slip-boundary conditions (i.e. $\vec{u}=0$) are located at $z = \pm l/2$. Assuming invariance in $y$ and $z$ direction and a steady state the Stokes equation is simplified to: \begin{equation} \eta \partial_x^2 u_y = f \end{equation} where $f$ denotes the force density and $\eta$ the dynamic viscosity. This can be integrated twice and the integration constants are chosen so that $u_y=0$ at $z = \pm l/2$ and we obtain: \begin{equation} u_y = \frac{f}{2\eta} \left(l^2/4-x^2\right) \end{equation} With that knowledge investigate the script \texttt{poisseuille.py}. Note the use of the \texttt{lbboundaries} module. Two walls are created with normal vectors $\left(\pm 1, 0, 0 \right)$. An external force is applied to every node. After 5000 LB updates the steady state should be reached.

Task: Write a loop that prints the fluid velocity at the nodes (0,0,0) to (16,0,0) and the node position to a file. Use the `lbf.[node].quantity` method for that. Hint: to write
to a file, first open a file and then use the `write()` method to write into it. Do not forget to close the file afterwards. Example:

```
f = open("file.dat", "w")
f.write("Hello world!\n")
f.close()
```

Use the data to fit a parabolic function. Can you confirm the analytic solution?

In [1]:

```
## A script to simulate planar Poisseuille flow in Espresso
from espressomd import System, lb, shapes, lbboundaries
import numpy as np
# System setup
box_l = 16.0
system = System(box_l = [box_l, box_l, box_l])
system.set_random_state_PRNG()
np.random.seed(seed = system.seed)
system.time_step = 0.01
system.cell_system.skin = 0.2
lbf = lb.LBFluidGPU(agrid=1, dens=1, visc=1, tau=0.01, ext_force_density=[0, 0.001, 0])
system.actors.add(lbf)
system.thermostat.set_lb(kT=0)
# Setup boundaries
walls = [lbboundaries.LBBoundary() for k in range(2)]
walls[0].set_params(shape=shapes.Wall(normal=[1,0,0], dist = 1.5))
walls[1].set_params(shape=shapes.Wall(normal=[-1,0,0], dist = -14.5))
for wall in walls:
system.lbboundaries.add(wall)
## Perform enough iterations until the flow profile
## is static (5000 LB updates):
system.integrator.run(5000)
## Part of the solution
node_v_list = []
for i in range(int(box_l)):
node_v_list.append(lbf[i, 0, 0].velocity[1])
with open("lb_fluid_velocity.dat", "w") as f:
for line in node_v_list:
f.write(str(line)+"\n")
```

[1] S Succi. *The lattice Boltzmann equation for fluid dynamics and beyond. *Clarendon Press, Oxford, 2001.

[2] B. Dünweg and A. J. C. Ladd. *Advanced Computer Simulation Approaches for Soft Matter Sciences III*, chapter II, pages 89–166. Springer, 2009.

[3] B. Dünweg, U. Schiller, and A.J.C. Ladd. Statistical mechanics of the fluctuating lattice-boltzmann equation. *Phys. Rev. E*, 76:36704, 2007.

[4] P. G. de Gennes. *Scaling Concepts in Polymer Physics*. Cornell University Press, Ithaca, NY, 1979.

[5] M. Doi. *Introduction do Polymer Physics.* Clarendon Press, Oxford, 1996.

[6] Michael Rubinstein and Ralph H. Colby. *Polymer Physics.* Oxford University Press, Oxford, UK, 2003.

[7] Daan Frenkel and Berend Smit. *Understanding Molecular Simulation.* Academic Press, San Diego, second edition, 2002.