Generalized equilibria for color-gradient lattice Boltzmann model based on higher-order Hermite polynomials: A simplified implementation with central moments

We propose generalized equilibria of a three-dimensional color-gradient lattice Boltzmann model for two-component two-phase flows using higher-order Hermite polynomials. Although the resulting equilibrium distribution function, which includes a sixth-order term on the velocity, is computationally cumbersome, its equilibrium central moments (CMs) are velocity-independent and have a simplified form. Numerical experiments show that our approach, as in Wen et al. [Phys. Rev. E 100, 023301 (2019)] who consider terms up to third order, improves the Galilean invariance compared to that of the conventional approach. Dynamic problems can be solved with high accuracy at a density ratio of 10; however, the accuracy is still limited to a density ratio of 1000. For lower density ratios, the generalized equilibria benefit from the CM-based multiple-relaxation-time model, especially at very high Reynolds numbers, significantly improving the numerical stability.

Figure 1

Density field of droplet in a moving tube: (a) Model A, (b) Model B, (c) Model C, (d) Model D. Initially ($t=0$), a circular droplet is placed in a domain. Owing to the movement of the top and bottom walls, the droplet also begins to move. In Models A and B, the droplets transition to an elliptic shape over time owing to the lack of Galilean invariance. In Models C and D, the droplets remain nearly circular over time owing to improved Galilean invariance. To clarify whether ${\mathrm{\Phi }}_i$ or $\mathbf{Q}$ contributes to the improvement in Galilean invariance, (e) an additional simulation based on Model D but without $\mathbf{Q}$, is also conducted. In this case the droplet still deforms over time.

Figure 2

Comparison of the velocity field at $t=80000{\delta }_t$: (a) Model A and (b) Model D. The solid line in the figures represents the interface position. The correction of Galilean invariance significantly improves the discontinuity in the velocity distribution.

Figure 3

Velocity profiles of the simulated layered two-phase flow in a channel together with the analytical solution: (a) Case A, (b) Case B, and (c) Case C, as presented in Table 3. The four types of equilibria (Table 1) are examined. When a correction is applied (Models C and D), the velocity profiles obtained from the simulations agree well with the analytical solution.

Figure 4

Bubble rising simulation for Case 1 ($\text{Re}=35, \text{Eo}=10, {\rho }_r^0/{\rho }_b^0=10$, and ${\mu }_r/{\mu }_b=10$) in Table 4. Results from finite-element-method simulations with level-set-based approach (with FreeLIFE solver) [115] and phase-field-based approach [121] and the phase-field LB model [120] are also shown.

Figure 5

Bubble rising simulation for Case 2 ($\text{Re}=35, \text{Eo}=125, {\rho }_r^0/{\rho }_b^0=1000$, and ${\mu }_r/{\mu }_b=100$) in Table 4. Results from finite-element-method simulations with level-set-based approach (with FreeLIFE solver) [115] and phase-field-based approach [121] and the phase-field LB model [120] are also shown.

Figure 6

Interface evolution of three-dimensional Rayleigh-Taylor instability for $\text{Re}=1024$ and $\text{At}=0.5$. The original equilibria for the CG model are used for the top panels (Model A in Table 1), while the generalized equilibria proposed in this study are used for the bottom panels (Model D in Table 1).

Figure 7

Time histories of interface position for $\text{At}=0.5$ and $\text{Re}=1024$ with generalized (solid line) and original equilibria (broken line). Results from the improved CG RM-MRT model [62] and multiphase LB flux solver [124] are also shown.

Figure 8

Time history of spike position during Rayleigh-Taylor instability with images of interfaces at infinite Reynolds number with zero viscosity.