Secondary instability of the spike-bubble structures induced by nonlinear Rayleigh-Taylor instability with a diffuse interface

Laminar-turbulent transition in Rayleigh-Taylor (RT) flows usually starts with infinitesimal perturbations, which evolve into the spike-bubble structures in the nonlinear saturation phase. It is well accepted that the emergence and rapid amplification of the small-scale perturbations are attributed to the Kelvin-Helmholtz-type secondary instability due to the high velocity shears induced by the stretch of the spike-bubble structures, however, there has been no quantitative description on such a secondary instability in literature. Moreover, the instability mechanism may not be that simple, because the acceleration or the “rising bubble” effect could also play a role. Therefore, based on the two-dimensional diffuse-interface RT nonlinear flows, the present paper employs the Arnoldi iteration and generalized Rayleigh quotient iteration methods to provide a quantitative study on the secondary instability. Both sinuous and varicose instability modes with high growth rates are observed, all of which are confirmed to be attributed to both the Rayleigh-Taylor and Kelvin-Helmholtz regimes. The former regime dominates the early-time instability due to the “rising bubble” effect, whereas the latter regime becomes more significant as time advances. Being similar to the primary RT instability [Yu et al., Phys. Rev. E 97, 013102 (2018), Dong et al., Phys. Rev. E 99, 013109 (2019), Fan and Dong, Phys. Rev. E 101, 063103 (2020)], the diffuse interface also leads to a multiplicity of the secondary instability modes and higher-order modes are found to exhibit more local extremes than the lower-order ones. Direct numerical simulations are carried out, which confirm the linear growth of the secondary instability modes with infinitesimal amplitudes and show their evolution to the turbulent-mixing state.

Figure 1

Contours of the density $\overline{\rho }$ (a) and velocity vector (b) for $t=7$, 8, and 9.

Figure 2

Contours of the vorticity $\overline{\mathrm{\Omega }}$ (a) and baroclinic torque (b) for $t=7$, 8, and 9.

Figure 3

Complex growth rate $\gamma$ of the secondary instability modes.

Figure 4

Eigenfunctions $\tilde{\rho }, \tilde{u}$, and $\tilde{v}$ of mode S1 for $t_0=7$. Left column: real part; right column: imaginary part.

Figure 5

Variation of the growth rate (a) and the frequency (b) with the acceleration $\overline{g}$.

Figure 6

Eigenfunctions $\tilde{\rho }, \tilde{u}$, and $\tilde{v}$ of mode S2 for $t_0=7$. Left column: real part; right column: imaginary part.

Figure 7

Eigenfunctions $\tilde{\rho }$ and $\tilde{u}$ of mode S3 for $t_0=7$. Left column: real part; right column: imaginary part.

Figure 8

Eigenfunctions $\tilde{\rho }$ and $\tilde{u}$ of mode S1 for $t_0=8$. Left column: real part; right column: imaginary part.

Figure 9

Eigenfunctions $\tilde{\rho }$ and $\tilde{u}$ of mode S1 for $t_0=9$. Left column: real part; right column: imaginary part.

Figure 10

Amplitude evolutions for case 1 and case 2.

Figure 11

The positions of the perturbation peak for case 1 and case 2.

Figure 12

Contours of the perturbation $\breve{u}$ at representative time instants for case 1.

Figure 13

Comparison of the energy spectra $E\left(k\right)$ among different time instants for case 1.

Figure 14

Contours of the perturbation $\breve{u}$ at representative time instants for case 2.

Figure 15

Contours of the density $\overline{\rho }$ for the multimode configuration.

Figure 16

Complex growth rate $\gamma$ of the secondary instability modes for a multimode configuration, where (b) is a zoom-in plot of (a).

Figure 17

Eigenfunctions of mode S1 for $t=6$ [(a) and (b)], 7 [(c) and (d)], and 8 [(e) and (f)], where the left and right columns are for $|\widehat{\rho }|$ and $|\widehat{u}|$, respectively.

Figure 18

Temporal evolution of the bubble front position $-h_b$.

Figure 19

Complex growth rate $\gamma$ of the secondary instability for filtered bubbles.

Figure 20

Variation of the growth rate with the acceleration $\overline{g}$ for modes BV1 (a) and BS1 (b).

Figure 21

Amplitude evolution under different mesh systems for case 1.