Multiple-relaxation-time color-gradient lattice Boltzmann model for simulating two-phase flows with high density ratio

In this paper we propose a color-gradient lattice Boltzmann (LB) model for simulating two-phase flows with high density ratio and high Reynolds number. The model applies a multirelaxation-time (MRT) collision operator to enhance the stability of the simulation. A source term, which is derived by the Chapman-Enskog analysis, is added into the MRT LB equation so that the Navier-Stokes equations can be exactly recovered. Also, a form of the equilibrium density distribution function is used to simplify the source term. To validate the proposed model, steady flows of a static droplet and the layered channel flow are first simulated with density ratios up to 1000. Small values of spurious velocities and interfacial tension errors are found in the static droplet test, and improved profiles of velocity are obtained by the present model in simulating channel flows. Then, two cases of unsteady flows, Rayleigh-Taylor instability and droplet splashing on a thin film, are simulated. In the former case, the density ratio of 3 and Reynolds numbers of 256 and 2048 are considered. The interface shapes and spike and bubble positions are in good agreement with the results of previous studies. In the latter case, the droplet spreading radius is found to obey the power law proposed in previous studies for the density ratio of 100 and Reynolds number up to 500.

Figure 1

Equilibrium droplet shape and velocity field for ${\rho }^R=10$ and ${\rho }^B=1$ at $\sigma =0.1$.

Figure 2

Schematic of the layered two-phase flow in a 2D channel. The red (wetting) fluid flows in the region of $\left|y\right|<a$ and the blue (wetting) fluid flows in the region of $a<\left|y\right|<b$.

Figure 3

Velocity profiles obtained by the present model and original color-gradient model for cases (a)–(f) in Table 2. The analytical velocity profile is also plotted for comparison.

Figure 4

Snapshots of the interface shape (${\rho }^N=0$) at different times for $\mathrm{At}=0.5$ and $\mathrm{Re}=256$, which are obtained by (a) the original color-gradient model and (b) the present model.

Figure 5

Interface shapes (${\rho }^N=0$) and contours of the $y$ component of the velocity at $t^{*}=2.5$, which are obtained by (a) the original color-gradient model and (b) the present model.

Figure 6

Snapshots of the interface shape (${\rho }^N=0$) at different times for $\mathrm{At}=0.5$ and $\mathrm{Re}=2048$, obtained by (a) the original color-gradient model and (b) the present model.

Figure 7

Positions of the spike bottom and bubble top obtained by the present model and the original color-gradient model at $\mathrm{At}=0.5$ for (a) $\mathrm{Re}=256$ and (b) $\mathrm{Re}=2048$. The previous numerical results of He et al. [44] and Chiappini et al. [48] are also presented for comparison.

Figure 8

Snapshots of the interface shape (${\rho }^N=0$) and the velocity field at (a) $t^{*}=0.5$, (b) $t^{*}=1.5$, and (c) $t^{*}=2.5$ for $\mathrm{Re}=40$.

Figure 9

Snapshots of the interface shape (${\rho }^N=0$) and the velocity field at (a) $t^{*}=0.5$, (b) $t^{*}=1.5$, and (c) $t^{*}=2.5$ for $\mathrm{Re}=100$.

Figure 10

Snapshots of the interface shape (${\rho }^N=0$) and the velocity field at (a) $t^{*}=0.25$, (b) $t^{*}=1$, (c) $t^{*}=1.75$, and (d) $t^{*}=2.5$ for $\mathrm{Re}=500$.

Figure 11

Dimensionless spreading radius as a function of the dimensionless time at $\mathrm{Re}=100$ and $500$. The spreading radius is normalized by $2R$.