Collision Operators#

The collision operators adopted in the current solver are all based on BGK collision operator. While the original BGK collision operator lacks stability at values of \(\tau\) close to 0.5, it is possible to compensate this weakness by following some strategies that result in the recursive regularized-BGK (RR-BGK) collision operator. First it is important to highlight that populations can be written in terms of its equilibrium and non-equilibrium portions as \(f_{i}=f_{i}^{\mathrm{eq}}+f_{i}^{\mathrm{neq}}\).

Note

The following development brings certain complexity and might require knowledge of some specific topics to its fully understanding.

Recursive Regularized-BGK (RR-BGK)#

The functional form of the equilibrium distribution function \(f_i^\mathrm{eq}\) is taken as a mesoscopic velocity expansion in Hermite polynomials up to a \(N\) order, as suggested by Malaspinas[1] and Mattila et al.[2]:

(1)#\[f_{i}^{\mathrm{eq}} = w_{i}\sum_{n=0}^{N}\frac{1}{n!}\mathbf{H}_{i}^{(n)}:\mathbf{a}^{(n),\mathrm{eq}}\]

in which “:” stands for a full index contraction, and \(\mathbf{H}_{i}^{(n)}\) is the 3 dimensional Hermite polynomial associated to the mesoscopic velocity vector \(\xi_{i,\alpha}=c_{i,\alpha}/c_{s}\). The Hermite polynomials up to a order \(N\) are given by:

(2)#\[H^{(0)}\left(\mathbf{\xi}_{i}\right) = 1,\]
(3)#\[H_{i,\alpha}^{(1)}\left(\mathbf{\xi}_{i}\right) = \xi_{i,\alpha},\]
(4)#\[H_{i,\alpha_{1}\alpha_{2}}^{(2)}\left(\mathbf{\xi}_{i}\right) = \xi_{i,\alpha_{1}}\xi_{i,\alpha_{2}}-\delta_{\alpha_{1}\alpha_{2}},\]
(5)#\[H_{i,\alpha_{1}\alpha_{2}\alpha_{3}}^{(3)}\left(\mathbf{\xi}_{i}\right) = \xi_{i,\alpha_{1}}\xi_{i,\alpha_{2}}\xi_{i,\alpha_{3}} -\left(\delta_{\alpha_{1}\alpha_{2}}\xi_{i,\alpha_{3}} + \delta_{\alpha_{1}\alpha_{3}}\xi_{i,\alpha_{2}} + \delta_{\alpha_{2}\alpha_{3}}\xi_{i,\alpha_{1}}\right),\]
(6)#\[H_{i,\alpha_{0}...\alpha_{n}}^{(n+1)}\left(\mathbf{ \xi_{i}}\right)=\xi_{i,\alpha_{0}}H_{i,\alpha_{1}...\alpha_{n}}^{(n)}\left(\mathbf{ \xi_{i}}\right)-\sum_{i=1}^{n}\delta_{\alpha_{0}\alpha_{i}}H_{i,\alpha_{1}...\alpha_{i-1}\alpha_{i+1}...\alpha_{n}}^{(n-1)}\left(\mathbf{\xi}_{i}\right)\]

Due to orthogonality of Hermite polynomials, the value of the coefficients \(\mathbf{a}^{(n),\mathrm{eq}}\) can be found with:

(7)#\[\mathbf{a}^{(n),\mathrm{eq}}=\sum_{i=0}^{q-1}\mathbf{H}_{i}^{(n)}f_{i}^{\mathrm{eq}},\]

where \(q\) is the number of velocity directions of the velocity set adopted. Knowing that the moments of equilibrium populations are :

(8)#\[\sum_{i=0}^{q-1}f_{i}^{\mathrm{eq}}=\rho,\]
(9)#\[\sum_{i=0}^{q-1}f_{i}^{\mathrm{eq}}\xi_{i,\alpha}=\frac{1}{c_{s}}\rho u_{\alpha},\]
(10)#\[\sum_{i=0}^{q-1}f_{i}^{\mathrm{eq}}\xi_{i,\alpha}\xi_{i,\beta}=\frac{1}{c_{s}^{2}}\rho u_{\alpha}u_{\beta} + \frac{1}{c_{s}^{2}} \rho c_{s}^{2} \delta_{\alpha\beta}.\]

The values of \(\mathbf{a}^{(n),\mathrm{eq}}\) can be written recursively as:

(11)#\[a^{(0),\mathrm{eq}}=\rho,\]
(12)#\[a^{(n),\mathrm{eq}}_{\alpha_{1}...\alpha_{n}}=\frac{1}{c_{s}}u_{\alpha_{n}}a^{(n-1),\mathrm{eq}}_{\alpha_{1}...\alpha_{n-1}}.\]

We adopt a expansion up to third-order to represent the equilibrium populations, which gives:

(13)#\[f_{i}^{\mathrm{eq}}=w_{i}\left[\rho + \frac{\mathbf{c}_{i}\cdot \left(\rho\mathbf{u}\right)}{c_{s}^{2}} + \frac{1}{2}\mathbf{H}_{i}^{(2)}:\mathbf{a}^{(2),\mathrm{eq}} + \frac{1}{6}\mathbf{H}_{i}^{(3)}:\mathbf{a}^{(3),\mathrm{eq}}\right]\]

or:

(14)#\[f_{i}^{\mathrm{eq}}=w_{i}\left\{\rho + \frac{c_{i,\alpha} \left(\rho u_{\alpha}\right)}{c_{s}^{2}} + \frac{1}{2}\frac{\left(c_{i,\alpha}c_{i,\beta}-c_{s}^{2}\delta_{\alpha\beta}\right)\left(\rho u_{\alpha}u_{\beta}\right)}{c_{s}^{4}} + \frac{1}{6}\frac{\left[c_{i,\alpha}c_{i,\beta}c_{i,\gamma}-c_{s}^{2}\left(\delta_{\alpha\beta}c_{i,\gamma} + \delta_{\alpha\gamma}c_{i,\beta} + \delta_{\beta\gamma}c_{i,\alpha}\right)\right]\left[\rho u_{\alpha}u_{\beta}u_{\gamma}\right]}{c_{s}^{6}}\right\}\]

Following the same approach, the non-equilibrium populations are written as:

(15)#\[f_{i}^{\mathrm{neq}} = w_{i}\sum_{n=0}^{N}\frac{1}{n!}\mathbf{H}_{i}^{(n)}:\mathbf{a}^{(n),\mathrm{neq}}\]

Being the first moments of non-equilibrium populations:

(16)#\[\sum_{i=0}^{q-1}f_{i}^{\mathrm{neq}}=0,\]
(17)#\[\sum_{i=0}^{q-1}f_{i}^{\mathrm{neq}}\xi_{i,\alpha}=-\frac{\Delta t}{2c_{s}}F_{\alpha},\]
(18)#\[\sum_{i=0}^{q-1}f_{i}^{\mathrm{neq}}\xi_{i,\alpha}\xi_{i,\beta}=\frac{1}{c_{s}^{2}}\Pi_{\alpha\beta}^{\mathrm{neq}},\]

where \(\Pi_{\alpha\beta}^{\mathrm{neq}}=\sum_{i=0}^{q-1}f_{i}^{\mathrm{neq}}c_{i,\alpha}c_{i,\beta}\). From those moments, the coefficients \(\mathbf{a}^{(n),\mathrm{neq}}\) can be written recursively as:

(19)#\[a^{(0),\mathrm{neq}}=0,\]
(20)#\[a_{\alpha_{1}}^{(1),\mathrm{neq}}=-\frac{\Delta t}{2c_{s}}F_{\alpha_{1}},\]
(21)#\[a_{\alpha_{1}\alpha_{2}}^{(2),\mathrm{neq}}=\frac{1}{c_{s}^{2}}\sum_{i}f_{i}^{\mathrm{neq}}e_{i,\alpha_{1}}e_{i,\alpha_{2}},\]
(22)#\[a_{\alpha_{1}\alpha_{2}\alpha_{3}}^{(3),\mathrm{neq}} = \frac{1}{c_{s}}u_{\alpha_{3}}a_{\alpha_{1}\alpha_{2}}^{(2),\mathrm{neq}}+ \frac{1}{c_{s}}\left[ u_{\alpha_{1}}a_{\alpha_{2}\alpha_{3}}^{(2),\mathrm{neq}}+ u_{\alpha_{2}}a_{\alpha_{1}\alpha_{3}}^{(2),\mathrm{neq}} \right]\]
(23)#\[a_{\alpha_{1}...\alpha_{n}}^{(n),\mathrm{neq}}=\frac{1}{c_{s}}u_{\alpha_{n}}a_{\alpha_{1}...\alpha_{n-1}}^{(n-1),\mathrm{neq}}+ \frac{1}{c_{s}^{(n-2)}}\left[ u_{\alpha_{1}}...u_{\alpha_{n-2}}a_{\alpha_{n-1}\alpha_{n}}^{(2),\mathrm{neq}}+\mathrm{perm}({\alpha_{n}}) \right]\]

where the operator “perm” stands for all the cyclic index permutations with \(u_{\alpha_{1}}...u_{\alpha_{n-2}}a_{\alpha_{n-1}\alpha_{n}}^{(2),\mathrm{neq}}\) of indexes from \(\alpha_{1}\) to \(\alpha_{n-1}\) (\(\alpha_{n}\) is never permuted).

The non-equilibrium populations truncated up to a third order Hermite polynomial expansion will be given by:

(24)#\[f_{i}^{\mathrm{neq}}=w_{i}\left[\mathbf{H}_{i}^{(1)}:\mathbf{a}^{(1),\mathrm{neq}} + \frac{1}{2}\mathbf{H}_{i}^{(2)}:\mathbf{a}^{(2),\mathrm{neq}} + \frac{1}{6}\mathbf{H}_{i}^{(3)}:\mathbf{a}^{(3),\mathrm{neq}} \right]\]

or:

(25)#\[f_{i}^{\mathrm{neq}}=w_{i}\left\{\frac{c_{i,\alpha}}{c_{s}}\left[-\frac{\Delta t}{2 c_{s}}F_{\alpha}\right] + \frac{1}{2}\left[\frac{c_{i,\alpha}c_{i,\beta}-c_{s}^{2}\delta_{\alpha\beta}}{c_{s}^{2}}\right]\left[\frac{1}{c_{s}^{2}}\Pi_{\alpha\beta}^{\mathrm{neq}}\right] + \frac{1}{6}\left[\frac{c_{i,\alpha}c_{i,\beta}c_{i,\gamma} - c_{s}^{2}\left(\delta_{\alpha\beta}c_{i,\gamma} + \delta_{\alpha\gamma}c_{i,\beta} + \delta_{\beta\gamma}c_{i,\alpha}\right)}{c_{s}^{3}}\right]\left[\frac{u_{\gamma}\Pi_{\alpha\beta}^{\mathrm{neq}} + \left(u_{\alpha}\Pi_{\beta\gamma}^{\mathrm{neq}} + u_{\beta}\Pi_{\alpha\gamma}^{\mathrm{neq}}\right)}{c_s^{3}}\right] \right\}\]

Important

Through this projection in an Hermite polynomial basis, the stability of LBM is greatly increased.

The discrete form of lattice Boltzmann equation is then written as:

(26)#\[f_i(\boldsymbol{x} + \boldsymbol{c}_i \Delta t, t+\Delta t) = f_i(\boldsymbol{x}, t)^{\mathrm{eq}} + \left(1-\omega\right)f_i(\boldsymbol{x}, t)^{\mathrm{neq}} + \Delta t\left(1-\frac{\omega}{2}\right)F_i(\boldsymbol{x}, t)\]

which describes the flow evolution.

Note

A simpler formulation of this development is the regularized-BGK (R-BGK) collision operator, in which the truncation of hermite polynomials is performed only up to second order. This operator is also available in the present solver.

Hybrid Recursive Regularized-BGK (HRR-BGK)#

The RR-BGK presents instabilities in LES at very low relaxation times. To overcome this issue, it is reccomendable to adopt the Hybrid Recursive Regularized-BGK collision operator (HRR-BGK) with LES (Jacob et al.[3]).

This modification consists on hybridizing the rate-of-strain tensor, considering it as a composition of the second-order centered finnite difference calculated rate-of-strain \(S_{\alpha\beta}^{\mathrm{FD}}\) and the value given by the mesoscopic terms \(S_{\alpha\beta}\), hence:

(27)#\[S_{\alpha\beta}^{\mathrm{HRR}}=\sigma S_{\alpha\beta} + \left(1-\sigma\right)S_{\alpha\beta}^{\mathrm{FD}}\]

where the relaxation parameter \(\sigma\) is adjusted dynamically with:

(28)#\[\sigma = \frac{1}{6\nu_{\mathrm{SGS}}\frac{\left(|\nabla \mathbf{u}|/|\nabla^{2} \mathbf{u}|\right)^{2}}{\Delta x^{2}c_{s}^{2}\tau}+1}\]

This adjustment assures that the dissipation of kinetic energy associated with the error of finnite difference scheme formulation recovers a perfectly controlled amount of local dissipation. The HRR-BGK presents a high stability and it allows the use of D3Q19 velocity set for high Reynolds LES simulations.

Note

The LES associated dissipation is given by \(\nu_{\mathrm{SGS}}|\nabla \mathbf{u}|^{2}\).