Kelvin-Helmholtz and Holmboe instabilities of a diffusive interface between miscible phases

The stability of a shear flow imposed along a diffusive interface that separates two miscible liquids (a heavier liquid lies underneath) is studied using direct numerical simulations. The phase-field approach is employed for description of a thermo- and hydrodynamic evolution of a heterogeneous binary mixture. The approach takes into account the dynamic interfacial stresses at a miscible interface and uses the extended Fick's law for setting the diffusion transport (the diffusion flux is proportional to the gradient of chemical potential). The shear flow is unstable to two kinds of instabilities: (1) the Kelvin-Helmholtz instability, with an immovable vortex formed in the middle of an interface (in the vertical direction) and (2) the Holmboe instability, with traveling waves along the interfacial boundary. The development of the Holmboe instability results in a stronger enhancement of molecular mixing between the mixture components. Earlier, the boundaries of these instabilities were determined using the linear stability analysis and employing the concept of a “frozen interface.” In the current work, through the solution of full equations, we obtain the stability boundaries for several sets of governing parameters, showing a greater variety of the possible shapes of the stability diagrams. The Kelvin-Helmholtz instability always occurs at lower gravity effects (lower density contrasts), while the Holmboe instability occurs when gravity is stronger. We show that for some parameters these two instabilities are separated by a zone where the shear flow is stable, and this zone disappears for the other sets of parameters.

Figure 1

The classical part of a free energy function (a) and the shapes of the phase diagrams (b) given by formulas (8), solid lines, and (7), dashed lines. The dots in panel (b) indicate the initial state of the mixture, and the arrows indicate the directions of the thermodynamic transformations experienced by the mixture in the current work.

Figure 2

The time dependencies of the total kinetic energy (a), (c) and of the length of interface (b), (d) are plotted for the numerical runs performed with the use of different grids, with $250\times 250$ (dash-dotted lines), $500\times 500$ (dashed lines), and $750\times 750$ (solid lines) nodes. The other parameters for these runs were $k=6.28, \text{Pe}={10}^6, \text{Cn}=4\times {10}^{-4}, {\delta }_0={\delta }_U=0.028$, and $U_0=0.2$ (a, b) and $U_0=0,4$ (c), (d).

Figure 3

Development of the Kelvin-Helmholtz instability. The fields of concentration (isolines) and velocity (vectors) at different time moments. The data are obtained for $\delta ={\delta }_U=0.2, U_0=1, \text{Re}=100, \text{Cn}={10}^{-3}, \text{Pe}={10}^6, \text{Gr}=0.5, k=4.19$, resolution $750\times 500$ nodes.

Figure 4

Development of the Holmboe instability. The fields of concentration (isolines) and velocity (vectors) at different time moments. The data are obtained for $\delta ={\delta }_U=0.2, U_0=1, \text{Re}=100, \text{Cn}={10}^{-3}, \text{Pe}={10}^6, \text{Gr}=10, k=4.19$, resolution $750\times 500$ nodes.

Figure 5

(a) The total kinetic energy, (b) the thickness of the interface, (c) the length of the interface, (d) the surface tension coefficient, and (e), (f) the average concentrations in each phase vs time. The data are obtained for $\delta ={\delta }_U=0.2, U_0=1, \text{Re}=100, \text{Cn}=0.001, \text{Pe}={10}^6, k=4.19, \text{Gr}=0.5$ (solid line) and $\text{Gr}=10$ (dashed line).

Figure 6

The neutral curves defining the zones of the Kelvin-Helmoltz (marked by “KHI”) and Holmboe (“HI”) instabilities. (a) The results of the linear stability analysis (see Ref. [26] for more details) obtained for an immiscible interface ($\text{Pe}=\infty$) with no surface tension effects ($\text{Cn}=0$) that separates two semi-infinite inviscid liquid domains ($\text{Re}=\infty$); other parameters are $U_0=1, {\delta }_0=0.2$. The thickness of the velocity profile, ${\delta }_U$, was used as a length scale. (b)–(d) The results of the direct numerical simulations. The data are obtained for $\text{Pe}={10}^6, \text{Re}=100, U_0=1$, and (b) ${\delta }_0=0.04, {\delta }_U=0.2, \text{Cn}=4\times {10}^{-5}$; (c) ${\delta }_0={\delta }_U=0.2, \text{Cn}=0.001$; (d) ${\delta }_0=0.04, {\delta }_U=0.2, \text{Cn}=0.001$. The unit of length is the height of the plane layer. The axes in (b–d) are rescaled to simplify the comparisons with (a). Square symbols correspond to the phase-field simulations, and circles are the points obtained using the classical approach.