Stabilizing the thermal lattice Boltzmann method by spatial filtering

We propose to stabilize the thermal lattice Boltzmann method by filtering the second- and third-order moments of the collision operator. By means of the Chapman-Enskog expansion, we show that the additional numerical diffusivity diminishes in the low-wavnumber limit. To demonstrate the enhanced stability, we consider a three-dimensional thermal lattice Boltzmann system involving 33 discrete velocities. Filtering extends the linear stability of this thermal lattice Boltzmann method to 10-fold smaller transport coefficients. We further demonstrate that the filtering does not compromise the accuracy of the hydrodynamics by comparing simulation results to reference solutions for a number of standardized test cases, including natural convection in two dimensions.

Figure 1

Viscosity of conserved modes ${\nu }_{\alpha }$ as a function of absolute wave vector $k$ along $\mathit{k}=k{\mathbf{\delta }}_x$ (a, d), $\mathit{k}=k\frac{1}{\sqrt{2}}\left({\mathbf{\delta }}_x+{\mathbf{\delta }}_y\right)$ (b, e), and $\mathit{k}=k\frac{1}{\sqrt{3}}\left({\mathbf{\delta }}_x+{\mathbf{\delta }}_y+{\mathbf{\delta }}_z\right)$ (c, f). Data are shown for the unfiltered model (solid gray lines) and the filtered model (dashed black lines). In all cases the background temperature equals $T_0=\frac{2}{5}$. In (a)–(c) ${\mathit{u}}_0=\mathit{0}$ and $\tau -\frac{1}{2}=0.01$. In (d)–(f) ${\mathit{u}}_0=\frac{c_s\mathit{k}}{4k}$ and $\tau -\frac{1}{2}=0.02$.

Figure 2

(a) Relative viscosity variation $\left({\nu }_{\alpha }-{\langle {\nu }_{\alpha }\rangle }_{\theta }\right)/{\langle {\nu }_{\alpha }\rangle }_{\theta }$ of the stress flux modes as a function of the orientation angle $\theta$ of the wave vector $\mathit{k}=k\left[cos\left(\theta \right){\mathbf{\delta }}_x+sin\left(\theta \right){\mathbf{\delta }}_y\right]$ at a fixed absolute value $k=0.04$ of the wave vector. Here ${\langle ·\rangle }_{\theta }$ is the average of $·$ over $\theta$. (b) Root-mean-square of the orientational viscosity variation $\sqrt{{\langle {\nu }_{\alpha }^2\rangle }_{\theta }-{\langle {\nu }_{\alpha }\rangle }_{\theta }^2}/{\langle {\nu }_{\alpha }\rangle }_{\theta }$ of the stress flux modes as a function of the absolute value $k$ of the wave vector $\mathit{k}=k\left[cos\left(\theta \right){\mathbf{\delta }}_x+sin\left(\theta \right){\mathbf{\delta }}_y\right]$. The orientational variations diminish as $k^2$.

Figure 3

Discrete Fourier transformed filter $1+{\nu }^F\widehat{\psi }$ as a function of absolute wavenumber $k$ using ${\nu }^F=\frac{1}{15}$ for $\mathit{k}=k{\mathbf{\delta }}_x$ (dashed line), $\mathit{k}=k\frac{1}{\sqrt{2}}\left({\mathbf{\delta }}_x+{\mathbf{\delta }}_y\right)$ (dotted line), and $\mathit{k}=k\frac{1}{\sqrt{3}}\left({\mathbf{\delta }}_x+{\mathbf{\delta }}_y+{\mathbf{\delta }}_z\right)$ (solid line).

Figure 4

The relative difference between the conserved mode viscosity ${\nu }_{\alpha }$ and the viscosity ${\nu }_0$ at zero wave number: $|{\nu }_{\alpha }-{\nu }_0|/{\nu }_0$ as a function of the absolute wave vector $k$ along $\mathit{k}=k{\mathbf{\delta }}_x$ (a, d), $\mathit{k}=k\frac{1}{\sqrt{2}}\left({\mathbf{\delta }}_x+{\mathbf{\delta }}_y\right)$ (b, e), and $\mathit{k}=k\frac{1}{\sqrt{3}}\left({\mathbf{\delta }}_x+{\mathbf{\delta }}_y+{\mathbf{\delta }}_z\right)$ (c, f). Data are shown for the unfiltered model (solid gray lines) and the filtered model (dashed black lines). In all cases the background temperature equals $T_0=\frac{2}{5}$. In (a)–(c) ${\mathit{u}}_0=\mathit{0}$ and $\tau -\frac{1}{2}=0.01$. In (d)–(f) ${\mathit{u}}_0=\frac{c_s\mathit{k}}{4k}$ and $\tau -\frac{1}{2}=0.02$.

Figure 5

Density $\rho$ and temperature $T$ as functions of the distance $x$ to the center of Sod's shock tube at time $t=200$ after removing the diaphragm, computed with the filtered thermal lattice Boltzmann model (markers) and reference solution (lines) presented in Ref. [40].

Figure 6

Fully developed velocity $u_x$ and temperature $T$ as functions of the distance $y$ from the center line in thermal Poiseuille flow, computed with the filtered thermal lattice Boltzmann model (markers) and reference solution (lines) presented in Appendix pp4.

Figure 7

Fully developed density $\rho$ and temperature $T$ as functions of the distance $x$ from the hot wall in the stationary Rayleigh Bénard setup, computed with the filtered thermal lattice Boltzmann model (markers) and reference solution (lines) presented in Appendix pp5.

Figure 8

Temperature distribution in Rayleigh Bénard convection for three Rayleigh numbers of $\mathrm{Ra}=5000$ (top), $\mathrm{Ra}=10000$ (middle), and $\mathrm{Ra}=20000$ (bottom), visualized using 10 temperature contours that are equally spaced between the bottom temperature of 0.402 and the top temperature of 0.398. The domain dimensions and other parameters are reported in Table 1.

Figure 9

The Nusselt number Nu as a function of the Rayleigh number Ra in Rayleigh Bénard convection, obtained from the present simulations (crosses), reference values [41] (circles), and an empirical formula $\mathrm{Nu}=1.56{(\mathrm{Ra}/{\mathrm{Ra}}_c)}^{0.296}$, where the critical Rayleigh number ${\mathrm{Ra}}_c=1708$ [8].