Phase-field method based on discrete unified gas-kinetic scheme for large-density-ratio two-phase flows

In this paper, a phase-field method under the framework of discrete unified gas-kinetic scheme (DUGKS) for incompressible multiphase fluid flows is proposed. Two kinetic models are constructed to solve the conservative Allen-Cahn equation that accounts for the interface behavior and the incompressible hydrodynamic equations that govern the flow field, respectively. With a truncated equilibrium distribution function as well as a temporal derivative added to the source term, the macroscopic governing equations can be exactly recovered from the kinetic models through the Chapman-Enskog analysis. Calculation of source terms involving high-order derivatives existed in the quasi-incompressible model is simplified. A series of benchmark cases including four interface-capturing tests and four binary flow tests are carried out. Results compared to that of the lattice Boltzmann method (LBM) have been obtained. A convergence rate of second order can be guaranteed in the test of interface diagonal translation. The capability of the present method to track the interface that undergoes a severe deformation has been verified. Stationary bubble and spinodal decomposition problems, both with a density ratio as high as 1000, are conducted and reliable solutions have been provided. The layered Poiseuille flow with a large viscosity ratio is simulated and numerical results agree well with the analytical solutions. Variation of positions of the bubble front and spike tip during the evolution of Rayleigh-Taylor instability has been predicted precisely. However, the detailed depiction of complicated interface patterns appearing during the evolution process is failed, which is mainly caused by the relatively large numerical dissipation of DUGKS compared to that of LBM. A high-order DUGKS is needed to overcome this problem.

Figure 1

Convergence rate of DUGKS and LBM, Pe $=$ 128, Cn $=$ 1/32, ${\mathrm{M}}_{\varphi }=0.02, \mathrm{\Delta }t=0.5$.

Figure 2

Initial state of Zalesak's disk, $\mathrm{Cn}=4/256, {\mathrm{M}}_{\varphi }=0.02$.

Figure 3

Results of Zalesak's disk after one period at $\mathrm{Pe}=128$.

Figure 4

Results of Zalesak's disk after one period at $\mathrm{Pe}=256$.

Figure 5

Results of Zalesak's disk after one period at $\mathrm{Pe}=512$.

Figure 6

Results of Zalesak's disk after one period at $\mathrm{Pe}=1024$.

Figure 7

Results of interface elongation in a shear flow at $t=0.25T, {\mathrm{M}}_{\varphi }=0.01, \mathrm{Pe}=256$.

Figure 8

Results of interface elongation in a shear flow at $t=0.5T, {\mathrm{M}}_{\varphi }=0.01, \mathrm{Pe}=256$.

Figure 9

Results of interface elongation in a shear flow at $t=0.75T, {\mathrm{M}}_{\varphi }=0.01, \mathrm{Pe}=256$.

Figure 10

Results of interface elongation in a shear flow at $t=T, {\mathrm{M}}_{\varphi }=0.01, \mathrm{Pe}=256$.

Figure 11

Results of interface elongation in a shear flow at $t=0.5T, {\mathrm{M}}_{\varphi }=0.1, \mathrm{Pe}=256$.

Figure 12

Results of interface deformation in a smoothed flow at $t=0.25T, {\mathrm{M}}_{\varphi }=0.02, \mathrm{Pe}=2048$.

Figure 13

Results of interface deformation in a smoothed flow at $t=0.5T, {\mathrm{M}}_{\varphi }=0.02, \mathrm{Pe}=2048$.

Figure 14

Results of interface deformation in a smoothed flow at $t=0.75T, {\mathrm{M}}_{\varphi }=0.02, \mathrm{Pe}=2048$.

Figure 15

Results of interface deformation in a smoothed flow at $t=T, {\mathrm{M}}_{\varphi }=0.02, \mathrm{Pe}=2048$.

Figure 16

Validation of Laplace's law with a density ratio of 1000 and ${\mathrm{M}}_{\varphi }=0.1$.

Figure 17

Density profile along the vertical center line with a density ratio of 1000.

Figure 18

Maximum magnitude of spurious velocity with various La number and density ratios.

Figure 19

Velocity profile of layered Poiseuille flow with various viscosity ratios.

Figure 20

Relative errors of layered Poiseuille flow with ${\mu }^{*}=10$.

Figure 21

Velocity profile (a) and relaxation time (b) with different interpolation schemes (${\mu }^{*}=10$).

Figure 22

Contours of density distribution at various moments in the process of phase separation.

Figure 23

Mass variation along with the process of phase separation.

Figure 24

Time evolution of interface patterns of Rayleigh-Taylor instability at $\mathrm{At}=0.5({\rho }^{*}=3), \mathrm{Re}=3000$.

Figure 25

Time evolution of bubble front and spike tip positions. Comparison with the results of Zhang et al. [41], Li et al. [47], Ren et al. [46], and Ding et al. [58].

Figure 26

Time evolution of interface patterns of Rayleigh-Taylor instability at $\mathrm{At}=0.1({\rho }^{*}=1.1/0.9), \mathrm{Re}=150$.

Figure 27

Time evolution of interface patterns of Rayleigh-Taylor instability at $\mathrm{At}=0.1({\rho }^{*}=1.1/0.9), \mathrm{Re}=3000$.

Figure 28

Time evolution of interface patterns of Rayleigh-Taylor instability by the HCZ model at $\mathrm{At}=0.1({\rho }^{*}=1.1/0.9), \mathrm{Re}=3000$.

Figure 29

Time evolution of interface patterns of Rayleigh-Taylor instability at (a) $\mathrm{At}=0.98({\rho }^{*}=99), \mathrm{Re}=150$; (b) $\mathrm{At}=0.998({\rho }^{*}=999), \mathrm{Re}=50$.