Effects of the initial perturbations on the Rayleigh–Taylor–Kelvin–Helmholtz instability system

Abstract

The effects of initial perturbations on the Rayleigh–Taylor instability (RTI), Kelvin–Helmholtz instability (KHI), and the coupled Rayleigh–Taylor–Kelvin–Helmholtz instability (RTKHI) systems are investigated using a multiple-relaxation-time discrete Boltzmann model. Six different perturbation interfaces are designed to study the effects of the initial perturbations on the instability systems. It is found that the initial perturbation has a significant influence on the evolution of RTI. The sharper the interface, the faster the growth of bubble or spike. While the influence of initial interface shape on KHI evolution can be ignored. Based on the mean heat flux strength D_3,1, the effects of initial interfaces on the coupled RTKHI are examined in detail. The research is focused on two aspects: (i) the main mechanism in the early stage of the RTKHI, (ii) the transition point from KHI-like to RTI-like for the case where the KHI dominates at earlier time and the RTI dominates at later time. It is found that the early main mechanism is related to the shape of the initial interface, which is represented by both the bilateral contact angle θ₁ and the middle contact angle θ₂. The increase of θ₁ and the decrease of θ₂ have opposite effects on the critical velocity. When θ₂ remains roughly unchanged at 90 degrees, if θ₁ is greater than 90 degrees (such as the parabolic interface), the critical shear velocity increases with the increase of θ₁, and the ellipse perturbation is its limiting case; If θ₁ is less than 90 degrees (such as the inverted parabolic and the inverted ellipse disturbances), the critical shear velocities are basically the same, which is less than that of the sinusoidal and sawtooth disturbances. The influence of inverted parabolic and inverted ellipse perturbations on the transition point of the RTKHI system is greater than that of other interfaces: (i) For the same amplitude, the smaller the contact angle θ₁, the later the transition point appears; (ii) For the same interface morphology, the disturbance amplitude increases, resulting in a shorter duration of the linear growth stage, so the transition point is greatly advanced.

Graphical abstract

Keywords

discrete Boltzmann method / hydrodynamic instability / non-equilibrium characteristic / initial perturbation

1 Introduction

Hydrodynamic instability-driven mixing processes are ubiquitous in nature, and of great significance both for fundamental research and for engineering application in a number of research fields, including Inertial Confinement Fusion (ICF), astrophysical phenomena, and supersonic combustion. In the last several decades, much attention has been paid to the challenging problem [1–20].

In view of its direct consequence in various applications, the effects of initial perturbations are being explored extensively [21–38]. Miles et al. [21] investigated the effect of initial conditions on two-dimensional (2D) Rayleigh–Taylor instability (RTI) and transition to turbulence in planar blast-wave-driven systems, and found that the initial conditions have a strong effect on the time to transition to the quasi-self-similar regime. Dimonte [22] studied the dependency of the self-similar Rayleigh–Taylor bubble acceleration constant on the initial perturbation amplitude. Ramaprabhu et al. [23] investigated the dependence of the RTI growth coefficient and the self-similar parameter on the amplitude, the spectral shape, the longest wavelength imposed, and mode-coupling effects. Gowardhan et al. [25] demonstrated that the initial material interface morphology controls the evolution characteristics of Richtmyer–Meshkov instability (RMI). Wei et al. [27] investigated the role of the initial perturbation shape on the RTI development. It is found that the time when the instability reaches a mean quadratic growth depends on the initial perturbation shape, and the subsequent vortical interactions are also very sensitive to details of the initial perturbation shape. Liu et al. [29] found that the temporal evolution of the bubble tip velocity is sensitively dependent on the Atwood numbers, the initial perturbation amplitude and the initial perturbation velocity in classical RTI. McFarland et al. [30] investigated the effects of inclination angle and incident shock Mach number on the inclined interface RMI. Zhai et al. studied the RMI problems of three light gas interfaces [31] and six heavy gas interfaces [32], including the different shock wave refraction on the interface and the influence of initial interface shape on the interface features. Dell et al. [33] systematically studied the effect of the initial perturbation on RMI, and observed that the initial growth rate of RMI is a non-monotone function of the initial perturbation amplitude. Xiao et al. [34] studied the effects of an initial perturbation on RMI, and found that the evolution of the interface with large initial amplitude in a low-density nonuniform area is fastest, while that with a small initial amplitude in a high-density nonuniform area is slower. Xie et al. [35] showed that the growth of the RTI mixing zone is substantially retarded by superimposing an optimized additional mode on its random initial perturbations. Liang et al. [8] analyzed the effects of the initial conditions in terms of the perturbation wavelength and amplitude. It is found that the instability undergoes a faster growth at the intermediate stage for a larger wavelength, while the late-time bubble and spike growth rates are insensitive to the changes of the initially perturbed wavelength and amplitude. Kord et al. [36] investigated the sensitivity of objective functions to the initial perturbation amplitudes at various stages of the RTI, such as the mole fraction, kinetic energy norm. Ding et al. [37] demonstrated that controlling the layer thickness is an effective way to modulate the late-stage instability growth, which may be useful for the target design. Si et al. [38] investigated the mode-composition effect on the three-dimensional (3D) interface RMI development, and found that the mode-coupling has an evident influence on the bubble evolution.

It is generally believed that in the ICF, in the late stage of RTI, the materials on both sides penetrate each other to form shear flow, and Kelvin–Helmholtz instability (KHI) occurs at the interface. However, in nature or practical engineering applications, such as in the atmosphere and oceans, air/fuel mixing in combustion chambers, the outer region of supernovae, and the compression of the fuel capsule in ICF, there is rarely a single instability. Due to the complexity of the interface, KHI often occurs simultaneously when RTI develops. It would be interesting to study the effect of the combined occurrence of two instability mechanisms. A complete understanding of the relation between the RTI and the KHI is important in understanding the hydrodynamic instability-driven mixing processes. In addition to the growth rate and final state, some studies have been made on the early linear and early nonlinear evolution of the coupled instability [39–49]. Both the instability modes and devolution processes that may be manifested when two such mechanisms occur simultaneously are of current interest. Ye et al. [39, 40] investigated the competitions between RTI and KHI in 2D incompressible fluids within a linear growth regime. It is found that the competition between the RTI and the KHI is dependent on the Froude number, the density ratio of the two fluids, and the thicknesses of the density transition layer and the velocity shear layer. Mandal et al. [41] investigated the nonlinear evolution of bubble and spike due to the combined action of RTI and KHI. Olson et al. [42] studied the coupled RTKHI in the early nonlinear regime, and indicated a complex and non-monotonic behavior where small amounts of shear in fact decrease the growth rate. Akula et al. [43] showed that the superposition of shear on RTI at small Atwood number increases the mixing width and growth rate at early times. Vadivukkarasan et al. [44–46] described the 3D destabilization characteristics of cylindrical and annular interfaces under the combined RTKHI, and the effects of various parameters on the most unstable wavenumbers were studied. Sarychev et al. [47] found that an undulating topography on the interface coating/base material is resulted from the coupled RTKHI. Brizzolara et al. [48] observed that the cross-over time can correctly predict the transition from shear- to buoyancy-driven turbulence, in terms of turbulent kinetic energy production, energy spectra scaling and mixing layer thickness. In 2020, we investigated the coupled RTKHI system with a Multiple-Relaxation Time (MRT) Discrete Boltzmann Model (DBM) [49]. To quantitatively analyse the coupled RTKHI process, we resort to morphological and non-equilibrium analysis techniques. We found that both the total boundary length L of the condensed temperature field and the mean heat flux strength D_3,1 can be used to quantitatively judge the main mechanism in the early stage of the RTKHI system. For the case where the KHI dominates at earlier time and the RTI dominates at later time, the ending point of linear increasing L or D_3,1 can work as a geometric or physical criterion for discriminating the two stages.

In this paper, we investigate the effects of initial perturbations on the RTI, KHI and the coupled RTKHI systems with the MRT-DBM. The paper is organized as follows: Section 2 briefly introduces the methodology of Discrete Boltzmann Modeling (DBM) method [50]. Systematic numerical simulations on the effects of initial perturbations are shown and analyzed in Section 3. A brief conclusion is given in Section 4.

2 Discrete Boltzmann modeling method

The DBM [51–57] is developed from a hybrid of Lattice Boltzmann method (LBM) [58–76] and the phase space description method of statistical physics. It is one of the concrete applications of statistical physics coarse-grained modeling method in the field of fluid mechanics, and is a further development of the description method of statistical physical phase space in the form of discrete Boltzmann equation. Its idea originated from a research review published by Xu et al. in 2012 [51]. In the process of development, it was inspired by the morphological phase space description method [55, 57, 77–79]. The methodology of DBM is as follows: It selects a perspective to study a set of kinetic properties of the system according to research requirements. Therefore, the kinetic moments describing this set of properties are required to maintain their values in the model simplification. Based on the independent components of the non-conserved kinetic moments of (f−f ^eq ), construct the phase space, and the phase space and its subspaces are used to describe the non-equilibrium behaviors of the system. The research perspective and modeling accuracy should be adjusted as the research progresses, where f ^eq is the corresponding equilibrium distribution function of f .

The establishment of a DBM model needs to go through three steps. The first step is the linearization of the collision term, and the second step is the discretization of the particle velocity space. The principle for coarse-grained modeling is that the physical quantities we choose to measure the system must keep the same values after simplification. The third step is the purpose and core of DBM modeling, and the specific scheme of non-equilibrium state and behavior description should be given. DBM is mainly aimed at the “mesoscale” and “dilemma” situations where continuum modeling fails or physical functions are insufficient and molecular dynamics method is unable to do due to the limited applicable scale. The physical information provided by DBM is between macroscopic continuous description and microscopic molecular dynamics. Compared with macroscopic description, DBM observes the system from a wider perspective. The necessity and benefit of non-conserved moment description in DBM increase with the increase of non-equilibrium degree.

Generally, different from the traditional LBM, whose main function is to recover macroscopic fluid equations, the DBM does not aim to solve the macro fluid equations, but mainly to describe and measure both hydrodynamic and thermodynamic nonequilibrium effects beyond the macro fluid equations. DBM only puts forward a physical requirement for the used discrete velocities and does not include specific discrete format. The physical requirement is that the kinetic moments concerned must keep the value unchanged after being converted into sums for calculation. DBM does not use the “lattice gas” image of “virtual particle propagation + collision”. This is also one reason why it should no longer be called the “lattice” Boltzmann method. A more detailed introduction can be seen in our recent review articles [55–57].

The DBM equation with external force term can be written as

\begin{array}{}\text{(1)}& \begin{array}{ll}\frac{\partial f_i}{\partial t}+v_{i\alpha }\frac{\partial f_i}{\partial x_{\alpha }}=\hfill & -M_{il}^{-1}[{\widehat{S}}_{lk}({\widehat{f}}_k-{\widehat{f}}_k^{eq})+{\widehat{A}}_l]\hfill \\ \hfill & -g_{\alpha }\frac{(v_{i\alpha }-u_{\alpha })}{RT}{f}_i^{eq}\mathrm{,}\hfill \end{array}\end{array}

where the variable t and x_α are the time and spatial coordinates, T is the temperature, g_α and u_α denote the acceleration and velocity in the x_α direction, v_iα is the discrete particle velocity, i = 1, . . ., N , and the subscript α indicates the x, y, or z component. f_i and {\widehat{f}}_i ( f_i^{eq} and {\widehat{f}}_i^{eq}) are the particle (equilibrium) distribution functions in velocity space and kinetic moment space, respectively; the mapping between moment space and velocity space is defined by the linear transformation M_ij, i.e., {\widehat{f}}_i=M_{ij}{f}_j and f_i=M_{ij}^{-1}{\widehat{f}}_j . The matrix \widehat{S}=\text{diag}(s_1\mathrm{,}{s}_2\mathrm{,}\cdots \mathrm{,}{s}_N) is the diagonal relaxation matrix. {\widehat{A}}_l is the correction term used for recovering reasonable hydrodynamic behaviors.

In the paper, the following two-dimensional discrete velocity model is used:

\begin{array}{}\text{(2)}& (v_{ix\mathrm{,}}{v}_{iy})=\{\begin{array}{ll}\mathbf{\text{c}}\mathbf{\text{y}}\mathbf{\text{c}}:c\left(\pm \mathrm{1,0}\right)\mathrm{,}\hfill & \text{for}1\le i\le \mathrm{4,}\hfill \\ c\left(\pm \mathrm{1,}\pm 1\right)\mathrm{,}\hfill & \text{for}5\le i\le \mathrm{8,}\hfill \\ \mathbf{\text{c}}\mathbf{\text{y}}\mathbf{\text{c}}:2c\left(\pm \mathrm{1,0}\right)\mathrm{,}\hfill & \text{for}9\le i\le \mathrm{12,}\hfill \\ 2c\left(\pm \mathrm{1,}\pm 1\right)\mathrm{,}\hfill & \text{for}13\le i\le \mathrm{16,}\hfill \end{array}\end{array}

and η_i = η₀ for i = 1, · · · , 4, and η_i = 0 for i = 5, · · · , 16, which is introduced to control the specific-heat-ratio γ. “cyc” indicates the cyclic permutation, and c and η₀ are two free parameters, which are adjusted to optimize the properties of the model. The specific forms of the correction term, transformation matrix and the corresponding equilibrium distribution functions in kinetic moment space can be seen in the appendix for details.

The kinetic moments of f_i-f_i^{eq} can be used to check and measure the Thermodynamic Non-Equilibrium (TNE) state and effects. {\text{Δ}}_i^{*}=M_{ij}^{*}(f_j-f_j^{eq}), where M_{ij}^{*}(f_i) represents kinetic central moments of f_i, and its specific expression can be obtained by replacing the variable v_iα by v_iα−u_α in M_ij . At present, the commonly used nonequilibrium characteristic quantities are {\text{Δ}}_{2\alpha \beta }^{*} , {\text{Δ}}_{(\mathrm{3,1})\alpha }^{*}, {\text{Δ}}_{3\alpha \beta \gamma }^{*} and {\text{Δ}}_{(\mathrm{4,2})\alpha \beta }^{*}\mathrm{.}{\text{Δ}}_{2\alpha \beta }^{*} indicates the viscous stress tensor and {\text{Δ}}_{(\mathrm{3,1})\alpha }^{*} indicates the heat flux. {\text{Δ}}_{3\alpha \beta \gamma }^{*} and {\text{Δ}}_{(\mathrm{4,2})\alpha \beta }^{*} indicate the flux of viscous stress and of heat flux, respectively. Where, the subscripts “2” and “3” indicate tensors of second order and third order, the subscript “3,1” means a first-order tensor contracted from a third-order tensor, and the subscript “4,2” means a second-order tensor contracted from a fourth-order tensor.

In the phase space opened by the non-conserved kinetic moments and their subspaces, the origin of coordinates corresponds to a thermodynamic equilibrium state, and any other point corresponds to a specific thermodynamic non-equilibrium state. The corresponding nonequilibrium strength can be defined by means of the distance from the origin. For example, the mean heat flux strength can be defined as D_{(\mathrm{3,1})}=\overline{\sqrt{{\text{Δ}}_{(\mathrm{3,1})\alpha }^{*2}}} . These more abundant but previously poorly understood characteristics of non-equilibrium behavior contain a large number of physical functions to be explored. At the same time, it should be pointed out that in addition to the idea of DBM modeling, a series of data processing techniques, characteristic scale extraction technology and structure analysis techniques developed in DBM application process can be directly used for data processing and feature analysis of various complex configuration and dynamic physical field.

3 Numerical simulations

In recent decades, numerous investigations have shown that the initial conditions have important effects on the evolution of interface instability. In this section, we investigate the effects of initial perturbations on the RTI, KHI and the coupled RTKHI systems with the multiple-relaxation-time DBM model, which has been validated by some well-known benchmark tests [49, 80, 81].

Six different initial conditions, shown in Fig. 1, are designed to study the effects of the initial perturbation shape on the instability system, including (a) sinusoidal perturbation, (b) sawtooth perturbation, (c) parabolic y = x², (d) elliptic disturbance, (e) inverted parabolic y =−x², and (f) inverted elliptic disturbance interface. The angle between the tangent at the end of the disturbance interface and the boundary is called the bilateral contact angle θ₁, and the angle between the centerline and the tangent at the center point is called the middle contact angle θ₂, shown in Fig. 2(a). The two contact angles together describe the basic morphology of the disturbance interface. Except for the sawtooth perturbation, the middle contact angle θ₂ of the other five interfaces is approximately equal to 90◦. These six interfaces basically include the continuous curve shape from the boundary to the center line determined by θ₁ and θ₂. Their related understanding and recognition can be helpful to analyze the multimode disturbance in the future.

In the simulation, a two-dimensional computational domain with height H = 80 and width W = 20 is adopted, and divided into N X × N Y = 101 × 400 mesh-cells. The computational domain can be divided into upper and lower parts by the disturbed interface. In the two parts, two different homogeneous temperatures are fixed, and the corresponding hydrostatic density profiles ρ(y) follow the regularity, ∂_y p(y) =−g_y ρ(y), and p(y) = T ρ(y) in each part. With fixed T , the solution has an exponentially decaying behavior in the two half volumes. The initial hydrostatic unstable configuration is therefore given by

where ρ_b and T_b (ρ_u and T_u) are the density and temperature of the light (heavy) fluid, u₀ is the shear velocity, g_y is the gravitational acceleration, y_s is the initial disturbance interface. For these six interfaces, the specific expressions of initial disturbances are as follows:

\begin{array}{}\text{(4a)}& y_s^a=H/2+a_0\cos (2\pi x/\lambda )\mathrm{,}\end{array}

\begin{array}{}\text{(4b)}& y_s^b=H/2-a_0+2a_0\frac{|x-\lambda /2|}{\lambda /2}\mathrm{,}\end{array}

\begin{array}{}\text{(4c)}& y_s^c=H/2-a_0+2a_0\frac{{(x-\lambda /2)}^2}{{(\lambda /2)}^2}\mathrm{,}\end{array}

\begin{array}{}\text{(4d)}& y_s^d=H/2+a_0-2a_0{(1-\frac{{(x-\lambda /2)}^2}{{(\lambda /2)}^2})}^{0.5}\mathrm{,}\end{array}

\begin{array}{}\text{(4e)}& y_s^e=H/2+a_0-2a_0\frac{{(x-\lambda /2)}^2}{{(\lambda /2)}^2}\mathrm{,}\end{array}

\begin{array}{}\text{(4f)}& y_s^f=H/2-a_0+2a_0{(1-\frac{{(x-\lambda /2)}^2}{{(\lambda /2)}^2})}^{0.5}\mathrm{,}\end{array}

where λ = W is the wavelength of the perturbation, a₀ is the initial amplitude. To have a finite width of the initial interface, a smooth interpolation between the two half-volumes is plused on the initial configuration. When the relative velocity u₀ is set to zero, this is a pure RTI system; when the acceleration g_y is set to zero, it is a pure KHI; with various tangential velocities and accelerations, the coupled RTKHI systems are gotten. The fifth-order weighted essentially non-oscillatory (WENO) scheme is adopted for space discretization, and the third-order Runge–Kutta scheme is adopted for time discretization. The left and right boundaries are periodic boundary conditions. In the bottom and the top boundaries, we apply the solid wall boundary conditions, i.e., (ρ, U_x, T )|_ix,−2 = (ρ, U_x, T )|_ix,−1 = (ρ, U_x, T )|_ix,0 = (ρ, U_x, T )|_ix,1, U_y |_ix,−2 = U_y |_ix,−1 = U_y |_ix,0 = 0, where (ix,−2), (ix,−1), and (ix, 0) are the indexes of ghost nodes out of the bottom boundary (ix, 1) when the fifth-order WENO scheme is used. The distribution functions f_i on the ghost nodes are approximated by theirs equilibrium distribution functions f_i^{eq}. On the top side, we can operate in a similar way.

We first investigate the influence of the initial interface shape on the evolution of RTI. Figure 3 shows the temperature patterns of the above six different initial disturbances evolving to t = 300 under the action of RTI. Figure 4 shows the time evolutions of bubbles, spikes, disturbance amplitudes and morphological boundary lengths for different initial interface systems. The definitions of bubble and spike amplitudes are shown in Fig. 2(b). The disturbance amplitude of RTI is defined as half of the sum of the bubble amplitude and the spike amplitude. The morphological boundary length L is the total length of the interface between high and low temperature (light and heavy) fluids. The methodology of morphological analysis technique is as follows: By defining a temperature threshold T_th, the value higher than the threshold is defined as the high temperature area (white area), and the value below the threshold is defined as the low temperature area (black area). In this way, the original continuous temperature field is transformed into a temperature pattern of black and white pixels, the total length of the dividing lines between white and black regions is defined as the morphological boundary length L. The specific calculation method of L can refer to our previous work [49, 56, 77, 78]. Both disturbance amplitude and morphological boundary length are important parameters that describe the degree of instability development and material mixing. They have different perspectives, but they are related and cannot replace each other. In the simulation, the dimensionless parameters are set as follows: ρ_b = 1, T_b = 1.4, ρ_u = 2.33333, T_u = 0.6, g_x = 0, g_y = 0.005, u₀ = 0, c = 1, η₀ = 3, temperature threshold T_th = 1.0, the specific-heat-ratio γ = 1.4, the initial amplitude a₀ = 2, and spatial step dx = dy = 0.2, temporal step dt = 10⁻³.

It can be seen from Figs. 3 and 4 that the initial disturbance interface shape has a great influence on the evolution of RTI. In the evolution of RTI, the heavy and light fluids gradually penetrate into each other as time progresses, with the light fluid rising to form a bubble and the heavy fluid falling to generate a spike. The sharper the interface (the interface has a higher curvature), the faster the growth of bubble or spike. For example, in the system with an inverted elliptic interface, the heavy fluid falls along both sides of the computing domain to form spikes. Compared with other interfaces where spikes are formed, the interface in this place is the sharpest, and the spike falls at the fastest rate; The light fluid rises in the middle part of the interface to form bubbles, compared with other interfaces to form bubbles, the interface here is the bluntest and the bubble rises at the slowest speed. The difference of amplitudes of RTI systems with different initial interfaces is small [Fig. 4(c)], but the morphological boundary length differs greatly in the later stage [Fig. 4(d)], indicating that the initial interface shape has substantial influence on the mixing degree of different media in the RTI system.

Figure 5 shows the temperature patterns of the six different initial disturbances evolving to t = 300 under the action of KHI. Figure 6 shows the time evolutions of amplitudes and morphological boundary lengths for different KHI systems. The amplitude of KHI systems is defined as half of the mixing width. The initial parameters are g_x = 0, g_y = 0, u₀ = 0.1, and the other parameters are the same as Fig. 3. It can be seen from Figs. 5 and 6 that the initial interface shape has little influence on the evolution of KHI. There is no significant difference in the final degree of material mixing in KHI systems with different initial interface shapes.

Figure 7 shows the evolutions of the coupled RTKHI systems with different initial disturbances at times t = 0, 50, 100, 150, and 200. The initial conditions are g_x = 0, g_y = 0.005, u₀ = 0.15, and the other parameters remain the same. As shown in the figure, under the action of initial perturbation and tangential velocity, the perturbations gradually grow to rolled-up vortexes, which are the main characteristic of free shear flows, so KHI plays a major role at the initial stage. For the inverted parabolic and elliptic interfaces, the vortexes are rolled up on the right side of the computational domain, which is different from the other four interfaces. As time progresses, more fluid is entrained in the vortical structure, the secondary RTI develops along the vortex arms, and plays a major role. It should be pointed out that in the systems with different disturbance interfaces, the secondary RTI develops at different times and with different speeds. Among them, the first three systems develop relatively fast, while the last three develop relatively slowly. In other words, the shape of the perturbed interface has a certain degree of influence on the evolution of the coupled RTKHI system. The focal points of the following research are: (i) the influence of perturbed interface on the transition point of the coupled system; (ii) the influence of perturbed interface on the main mechanism in the early stage of the RTKHI system.

In the previous work, we have shown that the main mechanism of the coupled RTKHI system in the early stage depends on the comparison of buoyancy and shear strength, i.e., gravity acceleration g and shear velocity u₀, and both the morphological boundary length L and the mean heat flux strength D_3,1 can be used to quantitatively judge the main mechanism in the early stage. The judgment method is as follows: Comparing the morphological boundary length L or the non-equilibrium strength D_3,1 of pure RTI with acceleration g and KHI systems with shear velocity u₀, if L^{KHI}=L^{RTI}(D_{\mathrm{3,1}}^{KHI}=D_{\mathrm{3,1}}^{RTI}), the buoyancy and shear effects are balanced in the early stages of the corresponding coupled system; The shear velocity u₀ is defined as the critical velocity u_C of the coupled system. It stands to reason that, keeping g constant, if the shear velocity increases, the shear effect is stronger than the buoyancy effect; when the shear velocity decreases, the buoyancy effect is stronger than the shear effect.

Figure 8 shows the influence of different initial disturbance interface shapes on the early mechanism judgment of the coupled system. In the previous work, we have proved that, the two quantities, L and D_3,1, always show a high correlation, especially in the early stage. Therefore, the simulation results of D_3,1 are only shown. In Figs. 8(a) and (b), the black curve of the pure RTI system (g = 0.005) and the green curve of the pure KHI system with shear velocity u₀ = 0.1 intersect at a point, and remain approximately coincidence until this point, i.e., these two systems have the highest similar degree of medium mixing in the early stages. The critical shear velocity is u_C = 0.1. In figures 8(c) and (d), as the interface shape changes, the bilateral contact angle θ₁ increases, the black curve of the pure RTI system deviates from the green curve, which indicates that the critical shear velocity u_C increases. For the elliptic structure [Fig. 8(d)], the critical shear velocity increases to u_C = 0.12. Figures 8(e) and (f) show similar trends. The black curve of the pure RTI system with g = 0.005 is approximately coincident with the red curve of the pure KHI system with shear velocity u₀ = 0.08 in the early stage, and then gradually deviates from it, and intersects with the green and blue lines (shear velocity u₀ = 0.1, 0.12). In general, for a given gravitational acceleration g, the critical shear velocity u_C is related to the shape of the initial interface.

For the case where the KHI dominates at earlier time and the RTI dominates at later time, the evolution process can be roughly divided into two stages. In the previous work, we have pointed out that, before the transition point of the two stages, both L and D_3,1 initially increase exponentially, and then increase linearly. Hence, the ending point of linear increasing L or D_3,1 can work as a geometric or physical criterion for discriminating the two stages.

Figure 9 shows the effects of initial disturbance shapes on the transition points of the coupled system. The green box denotes the mean heat flux strength D_3,1, the red circle represents the dD_3,1/dt, and the blue vertical line marks the position of the transition point. As can be seen from the figure, for the initial disturbance interface whose bilateral contact angle θ₁ is greater than 90◦[Figs. 9(b), (c), and (d)], the delay of transition point is relatively slow and the magnitude is almost negligible, as the bilateral contact angle increases; For the inverted parabolic and elliptic interfaces [Figs. 9(e) and (f)], the transition point has a relatively rapid and substantial delay, as the bilateral contact angle decreases; The inverted parabolic and elliptic interfaces [Figs. 9(e) and (f)] have a longer linear growth stage, and the transition point is later than that of the parabolic and elliptic interfaces [Figs. 9(c) and (d)].

Figure 10 shows the effect of disturbance amplitude on the early main mechanism of the coupled system. From top to bottom, the corresponding initial disturbance amplitudes are a₀ = 2, 2.8, and 3.6, respectively. Figures 10(a)–(d) correspond to sinusoidal, parabolic, inverted parabolic and sawtooth disturbances, respectively. From the curves we can see that:

(i) For the sinusoidal disturbance, the contact angles (θ₁ and θ₂) of the interface are basically unchanged with the increase of amplitude, so the change of amplitude has no significant impact on the judgment of the early mechanism of the RTKHI system [Figs. 10(a1)–(a3)]; For ellipse and inverted ellipse perturbation interfaces, there are consistent conclusions.

(ii) For the parabolic disturbance, with the increase of the disturbance amplitude, the middle contact angle θ₂ almost remains unchanged, but the bilateral contact angle θ₁ further increases, the black curve of pure RTI system with g = 0.005 gradually approaches the blue curve of pure KHI system with u₀ = 0.12, that is, the critical velocity tends to increase with the increase of bilateral contact angle θ₁ [Figs. 10(b1)–(b3)]. This is consistent with the conclusion of Figs. 8(c) and (d).

(iii) For the inverted parabolic disturbance (θ₂ = 90◦ and θ₁ < 90◦), the amplitude increases and the bilateral contact angle θ₁ decreases, but it does not affect the judgment of the early mechanism of the RTKHI system [Figs. 10(c1)–(c3)], and the images are similar to the system with inverted ellipse interface. This is also consistent with the conclusion of Figs. 8(e) and (f).

(iv) For the sawtooth disturbance, the amplitude increases, the bilateral contact angle θ₁ increases, but the middle contact angle θ₂ decreases, and the critical velocity is basically unchanged due to the combined action of the two factors [Figs. 10 (d1)–(d3)].

Figure 11 shows the effect of the initial disturbance amplitude on the transition point of the coupled system. The green and red curves correspond to the cases with a₀ = 2 and a₀ = 3.6, respectively. The virtual vertical lines indicate the positions of transition points. It can be seen from Fig. 11 that, (i) For the RTKHI systems with different initial interfaces, the effects of amplitude on the transition point are not the same; (ii) For the sinusoidal and sawtooth [Figs. 11(a) and (b)], the increase of amplitude has no significant effect on the transition point (marked by blue vertical lines). (iii) For the parabolic and ellipse interfaces [Figs. 11(c) and (d)], as the amplitude increases, the transition points will be brought forward by a smaller amplitude. The results of a₀ = 2 and a₀ = 3.6 are marked with green and red vertical lines, respectively. (iv) For the inverted parabolic and inverted ellipse perturbation interfaces [Figs. 11(e) and (f)], the increase of amplitude has a significant effect on the evolution process of the system and makes the transition point of the coupled system advance significantly. The effect of amplitude on the inverted ellipse interface is the greatest.

4 Conclusions

In this paper, we design six typical initial disturbance interfaces, and compare their effects on the RTI, KHI, and the coupled RTKHI systems. It is found that the initial perturbation has a great influence on the evolution of RTI. The sharper the interface, the faster the growth of bubbles or spikes. The influence of initial interface shape on KHI evolution can be ignored. Based on the mean heat flux strength D_3,1, the effects of initial interfaces on the coupled RTKHI system are investigated, and the researches focus on two aspects: (i) how the perturbed interface affects the transition point of the coupled system; (ii) how the perturbed interface affects the main mechanism in the early stage of the RTKHI system.

It is found that, in the coupled RTKHI system, for a given gravitational acceleration g, the magnitude of the critical shear velocity u_C is closely related to the shape of the initial disturbance interface. The interface shape can be described by two contact angles, the bilateral contact angle θ₁ and the middle contact angle θ₂. For the interface that θ₁ and θ₂ are basically constant, the critical velocity is basically constant. The increase of θ₁ and the decrease of θ₂ have opposite effects on the critical velocity. For the parabolic perturbation interface (θ₂ = 90◦, θ₁ > 90◦), the critical shear velocity increases as the increase of bilateral contact angle θ₁. The ellipse perturbation interface is its limiting case. For the inverted parabolic and the inverted ellipse disturbance (θ₂ = 90◦, θ₁ < 90◦), the critical shear velocity is basically the same, which is less than that of the sinusoidal and sawtooth disturbances. The influence of inverted parabolic and inverted ellipse perturbations on the transition point of the RTKHI system is greater than that of other interfaces. When the amplitude is constant, the bilateral contact angle θ₁ of the interface decreases, the transition point will be greatly delayed, and the transition point of the inverted ellipse structure system appears at the latest. For the same interface morphology, the disturbance amplitude increases, resulting in a shorter duration of the linear growth stage, and the transition point is greatly advanced.

References

[1] Y. Zhou, Rayleigh–Taylor and Richtmyer–Meshkov instability induced flow, turbulence, and mixing (I), Phys. Rep. 720, 1 (2017)

[2] Y. Zhou, Rayleigh–Taylor and Richtmyer–Meshkov instability induced flow, turbulence, and mixing (II), Phys. Rep. 723-725, 1 (2017)

[3] Y. Zhou, T. T. Clark, D. S. Clark, G. S. Gail, S. M. Aaron, C. M. Huntington, O. A. Hurricane, A. M. Dimits, and B. A. Remington, Turbulent mixing and transition criteria of flows induced by hydrodynamic instabilities, Phys. Plasmas 26(8), 080901 (2019)

[4] H. Li, B. Tian, Z. He, and Y. Zhang, Growth mechanism of interfacial fluid mixing width induced by successive nonlinear wave interactions, Phys. Rev. E 103(5), 053109 (2021)

[5] L. F. Wang, C. Xue, W. H. Ye, and Y. J. Li, Destabilizing effect of density gradient on the Kelvin–Helmholtz instability, Phys. Plasmas 16(11), 112104 (2009)

[6] F. Chen, A. G. Xu, G. C. Zhang, Y. J. Li, and S. Succi, Multiple-relaxation-time lattice Boltzmann approach to compressible flows with flexible specific-heat ratio and Prandtl number, Europhys. Lett. 90(5), 54003 (2010)

[7] H. Liang, Q. X. Li, B. C. Shi, and Z. H. Chai, Lattice Boltzmann simulation of three-dimensional Rayleigh–Taylor instability, Phys. Rev. E 93(3), 033113 (2016)

[8] H. Liang, X. L. Hu, X. F. Huang, and J. R. Xu, Direct numerical simulations of multi-mode immiscible Rayleigh–Taylor instability with high Reynolds numbers, Phys. Fluids 31(11), 112104 (2019)

[9] H. Liang, Z. H. Xia, and H. W. Huang, Late-time description of immiscible Rayleigh–Taylor instability: A lattice Boltzmann study, Phys. Fluids 33(8), 082103 (2021)

[10] Z. Zhai, L. Zou, Q. Wu, and X. Luo, Review of experimental Richtmyer–Meshkov instability in shock tube: From simple to complex, J. Mech. Eng. Sci. 232(16), 2830 (2018)

[11] L. Zou, J. Liu, S. Liao, X. Zheng, Z. Zhai, and X. Luo, Richtmyer–Meshkov instability of a flat interface subjected to a rippled shock wave, Phys. Rev. E 95(1), 013107 (2017)

[12] L. Zou, M. Al-Marouf, W. Cheng, R. Samtaney, J. Ding, and X. Luo, Richtmyer–Meshkov instability of an unperturbed interface subjected to a diffracted convergent shock, J. Fluid Mech. 879, 448 (2019)

[13] A. Ravid, R. I. Citron, and R. Jeanloz, Hydrodynamic instability at impact interfaces and planetary implications, Nat. Commun. 12(1), 2104 (2021)

[14] Y. W. Bin, M. J. Xiao, Y. P. Shi, Y. S. Zhang, and S. Y. Chen, A new idea to predict reshocked Richtmyer–Meshkov mixing: Constrained large-eddy simulation, J. Fluid Mech. 918, R1 (2021)

[15] H. Y. Ye, H. L. Lai, D. M. Li, Y. B. Gan, C. D. Lin, L. Chen, and A. G. Xu, Knudsen number effects on two-dimensional Rayleigh–Taylor instability in compressible fluid: Based on a discrete Boltzmann method, Entropy(Basel) 22(5), 500 (2020)

[16] L. Chen, H. L. Lai, C. D. Lin, and D. M. Li, Specific heat ratio effects of compressible Rayleigh–Taylor instability studied by discrete Boltzmann method, Front. Phys. 16(5), 52500 (2021)

[17] J. G. Tang, F. Zhang, X. S. Luo, and Z. G. Zhai, Effect of Atwood number on convergent Richtmyer–Meshkov instability, Acta Mech. Sin. 37(3), 434 (2021)

[18] C. D. Lin, K. H. Luo, Y. B. Gan, and Z. P. Liu, Kinetic simulation of nonequilibrium Kelvin–Helmholtz instability, Commun. Theor. Phys. 71(1), 132 (2019)

[19] R. H. Zeng, J. J. Tao, and Y. B. Sun, Three-dimensional viscous Rayleigh–Taylor instability at the cylindrical interface, Phys. Rev. E 102(2), 023112 (2020)

[20] Y. B. Sun, R. H. Zeng, and J. J. Tao, Effects of viscosity and elasticity on Rayleigh–Taylor instability in a cylindrical geometry, Phys. Plasmas 28(6), 062701 (2021)

[21] G. Dimonte, Dependence of turbulent Rayleigh–Taylor (RT) instability on initial perturbations, Phys. Rev. E 69(5), 056305 (2004)

[22] A. R. Miles, M. J. Edwards, and J. A. Greenough, Effect of initial conditions on two-dimensional Rayleigh–Taylor instability and transition to turbulence in planar blast-wave-driven systems, Phys. Plasmas 11(11), 5278 (2004)

[23] P. Ramaprabhu, G. Dimonte, and M. J. Andrews, A numerical study of the influence of initial perturbations on the turbulent Rayleigh–Taylor instability, J. Fluid Mech. 536, 285 (2005)

[24] D. H. Olson and J. W. Jacobs, Experimental study of Rayleigh–Taylor instability with a complex initial perturbation, Phys. Fluids 21(3), 034103 (2009)

[25] A. A. Gowardhan, J. R. Ristorcelli, and F. F. Grinstein, The bipolar behavior of the Richtmyer–Meshkov instability, Phys. Fluids 23(7), 071701 (2011)

[26] Y. Doron and A. Duggleby, Optical density measurements and analysis for single-mode initial-condition buoyancy-driven mixing, J. Fluids Eng. 133(10), 101204 (2011)

[27] T. Wei and D. Livescu, Late-time quadratic growth in single-mode Rayleigh–Taylor instability, Phys. Rev. E 86(4), 046405 (2012)

[28] S. Kuchibhatla and D. Ranjan, Effect of initial conditions on Rayleigh–Taylor mixing: Modal interaction, Phys. Scr. T155, 014057 (2013)

[29] W. H. Liu, L. F. Wang, W. H. Ye, and X. T. He, Temporal evolution of bubble tip velocity in classical Rayleigh–Taylor instability at arbitrary Atwood numbers, Phys. Plasmas 20(6), 062101 (2013)

[30] J. A. Mc Farland, J. A. Greenough, and D. Ranjan, Investigation of the initial perturbation amplitude for the inclined interface Richtmyer–Meshkov instability, Phys. Scr. T155, 014014 (2013)

[31] Z. G. Zhai, M. H. Wang, T. Si, and X. S. Luo, On the interaction of a planar shock with a light polygonal interface, J. Fluid Mech. 757, 800 (2014)

[32] X. S. Luo, M. H. Wang, T. Si, and Z. G. Zhai, On the interaction of a planar shock with an SF₆ polygon, J. Fluid Mech. 773, 366 (2015)

[33] Z. Dell, R. F. Stellingwerf, and S. I. Abarzhi, Effect of initial perturbation amplitude on Richtmyer–Meshkov flows induced by strong shocks, Phys. Plasmas 22(9), 092711 (2015)

[34] J. X. Xiao, J. S. Bai, and T. Wang, Numerical study of initial perturbation effects on Richtmyer–Meshkov instability in non-uniform flows, Phys. Rev. E 94(1), 013112 (2016)

[35] C. Y. Xie, J. J. Tao, Z. L. Sun, and J. Li, Retarding viscous Rayleigh–Taylor mixing by an optimized additional mode, Phys. Rev. E 95(2), 023109 (2017)

[36] A. Kord and J. Capecelatro, Optimal perturbations for controlling the growth of a Rayleigh–Taylor instability, J. Fluid Mech. 876, 150 (2019)

[37] R. Sun, J. C. Ding, Z. G. Zhai, T. Si, and X. S. Luo, Convergent Richtmyer–Meshkov instability of heavy gas layer with perturbed inner surface, J. Fluid Mech. 902, A3 (2020)

[38] Y. Liang, L. L. Liu, Z. G. Zhai, T. Si, and X. S. Luo, Universal perturbation growth of Richtmyer–Meshkov instability for minimum-surface featured interface induced by weak shock waves, Phys. Fluids 33(3), 032110 (2021)

[39] L. F. Wang, W. H. Ye, and Y. J. Li, Combined effect of the density and velocity gradients in the combination of Kelvin–Helmholtz and Rayleigh–Taylor instabilities, Phys. Plasmas 17(4), 042103 (2010)

[40] W. H. Ye, L. F. Wang, C. Xue, Z. F. Fan, and X. T. He, Competitions between Rayleigh–Taylor instability and Kelvin–Helmholtz instability with continuous density and velocity profiles, Phys. Plasmas 18(2), 022704 (2011)

[41] L. Mandal, S. Roy, R. Banerjee, M. Khan, and M. R. Gupta, Evolution of nonlinear interfacial structure induced by combined effect of Rayleigh–Taylor and Kelvin–Helmholtz instability, Nucl. Instr. Meth. Phys. Res. A 653(1), 103 (2011)

[42] B. J. Olson, J. Larsson, S. K. Lele, and A. W. Cook, Nonlinear effects in the combined Rayleigh–Taylor/Kelvin–Helmholtz instability, Phys. Fluids 23(11), 114107 (2011)

[43] B. Akula, M. J. Andrews, and D. Ranjan, Effect of shear on Rayleigh–Taylor mixing at small Atwood number, Phys. Rev. E 87(3), 033013 (2013)

[44] M. Vadivukkarasan and M. V. Panchagnula, Helical modes in combined Rayleigh–Taylor and Kelvin–Helmholtz instability of a cylindrical interface, Int. J. Spray Combust. 8(4), 219 (2016)

[45] M. Vadivukkarasan and M. V. Panchagnula, Combined Rayleigh–Taylor and Kelvin–Helmholtz instabilities on an annular liquid sheet, J. Fluid Mech. 812, 152 (2017)

[46] M. Vadivukkarasan, Temporal instability characteristics of Rayleigh–Taylor and Kelvin–Helmholtz mechanisms of an inviscid cylindrical interface, Meccanica 56(1), 117 (2021)

[47] V. D. Sarychev, S. A. Nevskii, A. Y. Granovskii, S. V. Konovalov, and V. E. Gromov, Combined Rayleigh–Taylor-Kelvin–Helmholtz instability and its role in the formation of the surface relief of the coating/substrate, AIP Conf. Proc. 2167, 020307 (2019)

[48] S. Brizzolara, J. Mollicone, M. Van Reeuwijk, A. Mazzino, and M. Holzner, Transition from shear-dominated to Rayleigh–Taylor turbulence, J. Fluid Mech. 924, A10 (2021)

[49] F. Chen, A. G. Xu, Y. D. Zhang, and Q. K. Zeng, Morphological and nonequilibrium analysis of coupled Rayleigh–Taylor-Kelvin–Helmholtz instability, Phys. Fluids 32(10), 104111 (2020)

[50] Without causing misunderstanding, DBM is used as an abbreviation of discrete Boltzmann Model/Modeling/Method.

[51] A. G. Xu, G. C. Zhang, Y. B. Gan, F. Chen, and X. Yu, Lattice Boltzmann modeling and simulation of compressible flows, Front. Phys. 7(5), 582 (2012)

[52] A. G. Xu, G. C. Zhang, and Y. J. Ying, Progress of discrete Boltzmann modeling and simulation of combustion system, Acta Phys. Sin. 64(18), 184701 (2015)

[53] A. Xu, G. Zhang, and Y. Gan, Progress in studies on discrete Boltzmann modeling of phase separation process, Mech. Eng. 38, 361 (2016)

[54] A. G. Xu, G. C. Zhang, and Y. D. Zhang, Discrete Boltzmann Modeling of Compressible Flows, Chapter 2 in Kinetic Theory, edited by G. Z. Kyzas and A. C. Mitropoulos, Rijeka: In Tech, 2018

[55] A. G. Xu, J. Chen, J. H. Song, D. W. Chen, and Z. H. Chen, Progress of discrete Boltzmann study on multiphase complex flows, Acta Aerodyn. Sin. 39(3), 138 (2021)

[56] A. G. Xu, J. H. Song, F. Chen, K. Xie, and Y. J. Ying, Modeling and analysis methods for complex fields based on phase space, Chinese J. Comput. Phys. 38(6), 631 (2021) (in Chinese)

[57] A. G. Xu, Y. M. Shan, F. Chen, Y. B. Gan, and C. D. Lin, Progress of mesoscale modeling and investigation of combustion multiphase flow, Acta Aero. Astro. Sin. 42(12), 625842 (2021)

[58] S. Succi, The Lattice Boltzmann Equation for Fluid Dynamics and Beyond, New York: Oxford University Press, 2001

[59] R. Benzi, S. Succi, and M. Vergassola, The lattice Boltzmann equation: Theory and applications, Phys. Rep. 222(3), 145 (1992)

[60] X. Shan and H. Chen, Lattice Boltzmann model for simulating flows with multiple phases and components, Phys. Rev. E 47(3), 1815 (1993)

[61] Y. Zhang, R. Qin, and D. Emerson, Lattice Boltzmann simulation of rarefied gas flows in microchannels, Phys. Rev. E 71(4), 047702 (2005)

[62] V. E. Ambruç and V. Sofonea, Quadrature-Based Lattice Boltzmann Models for Rarefied Gas Flow, edited by F. Toschi and M. Sega, Springer, 2019

[63] Y. B. Li and X. W. Shan, Lattice Boltzmann method for adiabatic acoustics, Phil. Trans. R. Soc. A 369(1944), 2371 (2011)

[64] Q. Li, K. H. Luo, Y. J. Gao, and Y. L. He, Additional interfacial force in lattice Boltzmann models for incompressible multiphase flows, Phys. Rev. E 85(2), 026704 (2012)

[65] Z. Wang, Y. Wei, and Y. Qian, A simple direct heating thermal immersed boundary-lattice Boltzmann method for its application in incompressible flow, Comput. Math. Appl. 80(6), 1633 (2020)

[66] Z. Chen, C. Shu, and D. Tan, Highly accurate simplified lattice Boltzmann method, Phys. Fluids 30(10), 103605 (2018)

[67] F. B. Tian, H. Luo, L. Zhu, J. C. Liao, and X. Y. Lu, An efficient immersed boundary-lattice Boltzmann method for the hydrodynamic interaction of elastic filaments, J. Comput. Phys. 230(19), 7266 (2011)

[68] F. B. Tian, Y. Wang, H. Liu, and Y. Zhang, The lattice Boltzmann method and its applications in complex flows and fluid-structure interactions, Inst. Mech. Eng. C 232(3), 403 (2018)

[69] H. Liang, B. C. Shi, Z. L. Guo, and Z. H. Chai, Phase-field-based multiple relaxation-time lattice Boltzmann model for incompressible multiphase flows, Phys. Rev. E 89(5), 053320 (2014)

[70] Y. Wang, C. Zhong, C. Zhuo, and S. Liu, A simplified finite volume lattice Boltzmann method for simulations of fluid flows from laminar to turbulent regime, Part I: Numerical framework and its application to laminar flow simulation, Comput. Math. Appl. 79(5), 1590 (2020)

[71] K. Pasieczynski and B. X. Chen, Multipseudopotential interaction models for thermal lattice Boltzmann method simulations, Phys. Rev. E 102(1), 013311 (2020)

[72] R. Qiu, Y. Bao, T. Zhou, H. Che, R. Chen, and Y. You, Study of regular reflection shock waves using a mesoscopic kinetic approach: Curvature pattern and effects of viscosity, Phys. Fluids 32(10), 106106 (2020)

[73] R. Qiu, T. Zhou, Y. Bao, K. Zhou, H. Che, and Y. You, Mesoscopic kinetic approach for studying nonequilibrium hydrodynamic and thermodynamic effects of shock wave, contact discontinuity, and rarefaction wave in the unsteady shock tube, Phys. Rev. E 103(5), 053113 (2021)

[74] D. K. Sun, A discrete kinetic scheme to model anisotropic liquid-solid phase transitions, Appl. Math. Lett. 103, 106222 (2020)

[75] D. K. Sun, H. Xing, X. L. Dong, and Y. S. Han, An anisotropic lattice Boltzmann-phase field scheme for numerical simulations of dendritic growth with melt convection, Int. J. Heat Mass Tran. 133, 1240 (2019)

[76] C. J. Zhan, Z. H. Chai, and B. C. Shi, A lattice Boltzmann model for the coupled cross-diffusion-fluid system, Appl. Math. Comput. 400, 126105 (2021)

[77] A. Xu, G. Zhang, X. Pan, P. Zhang, and J. Zhu, Morphological characterization of shocked porous material, J. Phys. D 42(7), 075409 (2009)

[78] A. G. Xu, G. C. Zhang, H. Li, Y. J. Ying, X. J. Yu, and J. S. Zhu, Temperature pattern dynamics in shocked porous materials, Sci. China Phys. Mech. Astron. 53(8), 1466 (2010)

[79] A. G. Xu, G. C. Zhang, Y. J. Ying, and C. Wang, Complex fields in heterogeneous materials under shock: Modeling, simulation and analysis, Sci. China Phys. Mech. Astron. 59(5), 650501 (2016)

[80] F. Chen, A. Xu, and G. Zhang, Viscosity, heat conductivity, and Prandtl number effects in the Rayleigh–Taylor Instability, Front. Phys. 11(6), 114703 (2016)

[81] F. Chen, A. Xu, and G. Zhang, Collaboration and competition between Richtmyer–Meshkov instability and Rayleigh–Taylor instability, Phys. Fluids 30(10), 102105 (2018)

Rights & permissions

Higher Education Press 2022