copy-edit - addressed comments by PE

This chapter dedicated to the lattice Boltzmann method (LBM), an attractive numerical approach for the simulation of fluid flow \cite{kruger:2017lattice}, chemical or thermal transport \cite{Sharma2020,Dapelo2021}, radiative transport \cite{Mink2020,Mink2022}, and other problems.

This chapter is dedicated to the lattice Boltzmann method (LBM), a popular numerical approach for the simulation of fluid flow \cite{kruger:2017lattice}, chemical or thermal transport \cite{Sharma2020,Dapelo2021}, radiative transport \cite{Mink2020,Mink2022}, and other problems.

The method can be easily and efficiently parallelized \cite{kruger:2017lattice} and the advent of general-purpose computing on graphics processing units (GPGPU) made large-scale and highly resolved numerical simulations of turbulent flows feasible \cite{kang2013,geier2015cumulant,kumar2018,peng2018,wittmann2013,zakirov2021}.

The content of the chapter focuses on the implementation of the method and performance optimizations for distributed parallel computing platforms, which is the author's contribution within the research group.

Hence, the method is described with the objective to formulate the computational algorithm and details such as rigorous derivation of the method are omitted.

Hence, LBM is described with the objective to formulate the computational algorithm and details such as rigorous derivation of the method are omitted.

The LBM implementation used in this work was developed, analyzed and validated in \cite{fucik2019,fucik2020,eichler2021a,eichler2021b,eichler2022,fucik:lbman,fucik:lbmat}.

The chapter is organized as follows.

\subsection{Lattice and velocity set}

\label{sec:lbm:lattice}

We consider the discretization of a cuboidal domain $\Omega \subset \mathbb R^D$ using an equidistant orthogonal grid represented by a finite set $\mathcal G_\Omega$ of grid nodes.

We consider the discretization of a cuboidal domain $\Omega \subset \mathbb R^3$ using an equidistant orthogonal grid represented by a finite set $\mathcal G_\Omega$ of grid nodes.

The extended grid $\mathcal G_{\overline \Omega}$ consists of the interior nodes from $\mathcal G_\Omega$ and an additional outer layer of nodes, which is added to facilitate the discretization of the domain boundary.

The Cartesian coordinates of grid nodes are defined as $\vec x_{i,j,k} = [i \delta_x + O_x, j \delta_x + O_y, k \delta_x + O_z]^T$, where $\delta_x$~[\si{\metre}] is the grid spacing parameter and $\vec O = [O_x, O_y, O_z]^T$~[\si{\metre}] is an offset vector that allows to fit the grid nodes to the physical coordinates of the domain.

The sets $\mathcal G_\Omega$ and $\mathcal G_{\overline \Omega}$ of interior and all grid nodes, respectively, are defined by specifying the indices $i,j,k$ for the points $\vec x_{i,j,k}$ as

\end{subequations}

where $\Delta x$ is the lattice spacing parameter that we set to $\Delta x = 1$ as usual in the LBM community~\cite{kruger:2017lattice}.

Note that a grid node $\vec x_{i,j,k}$ represents the same point in space as the corresponding lattice site $\Delta x [i, j, k]^T$, but in different unit systems.

Hence, we will use the notation $\vec x \in \mathcal G_\Omega$ to refer to a point using its physical coordinates and $\vec x \in \mathcal L_\Omega$ to refer to its lattice coordinates.

Hence, we will use the notation $\vec x \in \mathcal G_{\overline \Omega}$ to refer to a point using its physical coordinates and $\vec x \in \mathcal L_{\overline \Omega}$ to refer to its lattice coordinates.

Similarly to the space discretization, the time interval $\left[0, t_{\max}\right]$ is discretized by a finite set $\mathcal G_t = \left\{\ell  \delta_t \mid \ell \in \{0, \ldots, N_t\} \right\}$, where $\delta_t = t_{\max} / N_t$~[\si{\second}] is the time step in physical units and $N_t$ is the number of time steps.

The set of non-dimensional time levels is denoted as $\mathcal L_t = \left\{\ell  \Delta t \mid \ell \in \{0, \ldots, N_t\} \right\}$, where $\Delta t = 1$ is the non-dimensional time step.

\subsection{Collision operator}

The term $\mathcal C_q$ in \cref{eq:lbm:LBE} denotes the discrete collision operator; in this paper we use the cumulant operator proposed in~\cite{geier2015cumulant}.

The term $\mathcal C_q$ in \cref{eq:lbm:LBE} denotes the discrete collision operator; in this thesis we use the cumulant operator proposed in~\cite{geier2015cumulant}.

The operator has several parameters that drive the relaxation rates of different cumulants towards their equilibria.

The first relaxation rate, $\omega_1$, is linked to the kinematic viscosity $\nu$ in the Navier--Stokes equations via

The relaxation rate $\omega_1$ is linked to the kinematic viscosity $\nu$ in the Navier--Stokes equations via

\begin{equation}

    \nu = c_s^2 \left( \frac{1}{\omega_1} - \frac{1}{2} \right).

    \nu = c_s^2 \left( \frac{1}{\omega_1} - \frac{\Delta t}{2} \right) \frac{\delta_x^2}{\delta_t}.

\end{equation}

The other relaxation rates are set as suggested by the authors of the original paper.

The other relaxation rates are set as suggested by the authors of the original paper, including the limiters proposed in \cite{geier2017a,geier2017b}.

Furthermore, the approximations of the spatial velocity derivatives to reduce the artifacts due to the absence of higher-order cumulants are also incorporated in the collision operator according to \cite{geier2017a,geier2017b}.

\subsection{Equilibrium function}

In the kinetic theory of gases, the Maxwell--Boltzmann distribution $f^{\mathrm{eq}}(\vec x, \abs{\vec v}, t)$ describes the speeds of particles in an ideal gas moving with a mean velocity $\vec v$ in a thermodynamic equilibrium \cite{kruger:2017lattice}.

In the kinetic theory of gases, the Maxwell--Boltzmann distribution $f^{\mathrm{eq}}(\vec x, \abs{\vec v}, t)$ describes the speeds of particles in an ideal gas moving with a macroscopic velocity $\vec v$ in a thermodynamic equilibrium \cite{kruger:2017lattice}.

In the derivation of the discrete lattice Boltzmann equation, the equilibrium distribution function $f^{\mathrm{eq}}$ is discretized to obtain an approximation of the equilibrium state for the given discrete velocity set.

The discrete equilibrium distribution function $f_q^{\mathrm{eq}}$ for the cumulant collision operator can be derived in the space of cumulants following \cite{geier2015cumulant}.

The formulation of proper boundary conditions for LBM is a major complication in its applications, since the desired macroscopic conditions (e.g., given in terms of $\rho$ and $\vec v$) must be reproduced on the mesoscopic level using the discrete density distribution functions $f_q$ \cite{kruger:2017lattice}.

This problem does not have a unique solution and thus multiple boundary schemes for LBM exist that approximate the same macroscopic condition.

As of this writing, the LBM implementation used by the author contains simple boundary conditions based on the wet-node equilibrium approach for inflow \cite{kruger:2017lattice,latt2008,mohamad2009,haussmann2019}, extrapolation method for outflow, and full-way bounce-back scheme for solid walls \cite{kruger:2017lattice}.

As of this writing, the LBM implementation used by the author contains elementary boundary conditions based on the wet-node equilibrium approach for inflow \cite{kruger:2017lattice,latt2008,mohamad2009,haussmann2019}, extrapolation method for outflow, and full-way bounce-back scheme for solid walls \cite{kruger:2017lattice}.

Since the analysis of boundary conditions is not within the scope of this thesis, the details are omitted.

\subsection{Initialization}

Since the verification and validation of the LBM implementation was performed primarily by other research team members \cite{fucik2019,eichler2021a,eichler2021b,eichler2022}, this section presents only the performance evaluation, which is the author's original work.

The benchmark problem used in this section simulates a laminar flow in a rectangular duct for which an analytical solution can be derived \cite{white2005}.

The domain $\Omega = (0, L_x) \times (0, L_y) \times (0, L_z)$ has an inflow boundary on the left side, outflow boundary on the right side, and solid walls with the no-slip boundary condition on the remaining sides.

The domain $\Omega = (0, L_x) \times (0, L_y) \times (0, L_z)$ has an inflow boundary on the left hand side (i.e., $x = 0$), outflow boundary on the right hand side (i.e., $x = L_x$), and solid walls with the no-slip boundary condition on the remaining sides.

The inflow velocity profile is prescribed as $\vec v(0, y, z) = [v_{\mathrm{an},x}(y, z), 0, 0]^T$, where $v_{\mathrm{an},x}$ denotes the $x$-component of the velocity in the analytical solution.

The kinematic viscosity is set as $\nu = \SI{1.5e-5}{\metre\squared\per\second}$ and the final simulation time is $t_{\max} = \SI{100}{\second}$.