2D LBM Simulation of Flow Behind Three Square Cylinders

A 2D LBM code was developed in this study for the flow behind three square cylinders with attached splitter plates at (g/D, l/D) = (0.5–4, 1–5) for a constant Re = 150. The 2D is used because of its low Reynolds number.
LBM is the latest numerical technique used to solve fluid dynamics problems. It was developed by Frisch et al. [33 (link)]. This method serves as a good alternative approach for modeling physical phenomena in fluid flow. The LBM is an approach to discrete kinetic theory that includes a mesoscale description of the fluid's microstructure compared to traditional numerical methods including the spectral, element method, FDM, FEM, and FVM, which are based on the discretization of the macroscopic continuum equations. Boltzmann's main idea is to bridge the gap between macro- and micro-scales by looking at the behaviour of a group of particles rather than individual particles. A distribution function represents the property of a collection of particles. This function serves as a representative of the collection of particles. The scale is known as the mesoscale.
Here, we will provide a brief summary of the method. In this study, the p₂q₉ (where p is space dimension and q is the particle number) model is applied. Each computing node in this model is composed of nine particles: one rest particle (i = 0) and eight moving particles (i = 1–8) make up each computing node (Fig. 2)).Fig. 2

LBM lattice velocities on square structure

The equation for density evolution is given by; $$f_i\left(\mathit{x}+e_i,t+1\right)=f_i\left(\mathit{x},t\right)-\left[f_i\left(\mathit{x},t\right)-f_i^{\left(\text{eq}\right)}\left(\mathit{x},t\right)\right]/\tau$$ where f_i and f_i^(eq) are particle distribution function and corresponding equilibrium distribution function, e_i are the velocity directions, t, x and τ are the dimensionless time, the particles position, and the single relaxation time, respectively.
The corresponding equilibrium distribution function is: $$f_i^{\left(\text{eq}\right)}=\rho {\omega }_i\left(1+3\left({\mathit{e}}_i·\mathit{u}\right)+4.5{\left({\mathit{e}}_i·\mathit{u}\right)}^2-1.5{\mathit{u}}^2\right)$$ Here ρ, u and ω_i are the fluid density, u is the instantaneous velocity, and ω_i is the corresponding weight functions (ω_i = 0 for i = 0, ω_i = 1/9 for i = 1–4, and ω_i = 1/36 for i = 5–8).
Equation (7) can be solved using two steps: collision, which makes use of a BGK operator [20 (link), 33 (link)], and propagation. These actions can be summed up as: $$\text{Collision}:f_i^{\left(\text{new}\right)}\left(\mathit{x},t\right)=f_i\left(\mathit{x},t\right)-\left[f_i\left(\mathit{x},t\right)-f_i^{\left(\text{eq}\right)}\left(\mathit{x},t\right)\right]/\tau$$ $$\text{Streaming}:f_i\left(\mathit{x}+{\mathit{e}}_i,t+1\right)=f_i^{\left(\text{new}\right)}\left(\mathit{x},t\right).$$
Following the streaming stage, boundary conditions are used, and the problem is then iteratively solved. At each computing node, the densities and velocities are determined using Eqs. (11) and (12). $$\rho =\sum _i{f}_i$$ $$\rho \mathit{u}=\sum _i{f}_i{e}_i$$
The equation of state (p = ρc_s², where c_s is the speed of sound and its value is 1/√3) is used to calculate pressure.