Kelvin-Helmholtz instability under stabilizing parallel magnetic field in nonhomogeneous compressible MHD flows

G. Peyrichon,¹ H. El-Rabii,² J.-F. Ripoll,¹,³ A. Casner,¹ A. Michael,⁴ V. G. Merkin,⁴ and O. A. Hurricane⁵

¹) CEA, DAM, DIF, F-91297 Arpajon, France

²) Institut Pprime, UPR 3346 CNRS, 86073 Poitiers, France

³) Université Paris-Saclay, CEA, LMCE, F-91680 Bruyères-le-Châtel, France

⁴) The Johns Hopkins University Applied Physics Laboratory, Laurel, Maryland 20723, United States

⁵) Lawrence Livermore National Laboratory, 7000 East Avenue, Livermore, California 94550, USA

(*Electronic mail: guillaume.peyrichon@cea.fr, hazem.el-rabii@cnrs.fr)

We study the Kelvin-Helmholtz instability (KHI) for the general case of a compressible, nonhomogeneous, magnetized plasma flow. The study is limited to a vortex sheet interface with an imposed parallel magnetic field. We introduce a new formalism based on a convective Mach number $M_c$, a convective Alfvénic Mach number $M_{Ac}$, and a total convective Mach number that combines the two. We derive an analytic expression of the KHI growth rate for a homogeneous flow (i.e. zero Atwood number, $A = 0$) that converges toward both the expression for unmagnetized compressible flow and Chandrasekhar's expression for magnetized incompressible flow. Otherwise, the dispersion relation is solved numerically and allows deriving general stability diagrams of magnetized KHI for the triplet ($A, M_c, \beta$ - plasma) parameters. We show these parameters uniquely define all configurations for a parallel magnetic field. We also construct diagrams with respect to the convective Alfvénic Mach number, the $\beta$-plasma parameter or the magnetic field showing which magnetic field strength is required for stabilizing a given shear flow. The theoretical growth rates are compared with 18 simulations made with the GAMERA code, currently used for 3D magnetospheric simulations. Finally, we apply our results to the analysis of a past KHI experiment performed at the OMEGA laser facility, showing linear theory succeeds to provide accurate estimates of the growth rate at early times. We further discuss how our results can inform future experiments in the high-Mach magnetized regime at the National Ignition Facility. Possible limitations of the study due to resistive, mixing, or turbulence effects are discussed.

# I. INTRODUCTION

The Kelvin-Helmholtz instability (KHI) is a well-known hydrodynamic and magnetohydrodynamics (MHD) instability that develops at the interface between two fluids flowing in opposite directions¹. This instability is fundamental for the dynamics of vortex and shear flows and, more generally, in the transition to turbulence. Interactions of turbulence and KHI development generate edge modes in Tokamak plasmas²,³. KHI has been extensively studied in the context of astrophysics flows, including the outer flanks of the Earth's magnetosphere⁴⁻⁷, solar jets⁸, comet tails⁹, the solar atmosphere¹⁰, geophysical flows¹¹ and distant astrophysical jets¹². An extended review on magnetized KHI in the Earth's magnetosphere is offered in Faganello and Califano¹³.

Linear theory has been carried out on the temporal stability of magnetospheric sheared flows for the vortex sheet approximation, with an infinitely thin interface between the two sheared flows,¹⁴⁻²⁰ and for shear layers with a finite thickness²¹⁻²³. In the linear case of an unmagnetized compressible flow and assuming the vortex sheet approximation, Miles²⁴ showed that compressibility stabilizes KHI for a homogeneous plasma where the Atwood number, $A$, is zero, with $A = (\rho_h - \rho_l)/(\rho_h + \rho_l)$ where $\rho_h$ is the heavier fluid and $\rho_l$ is the lighter fluid.

For a homogeneous, unmagnetized flow, Choudhury²⁵ rewrote the dispersion relation by introducing a reduced Mach number in the convective frame to obtain a useful expression

of the KHI growth rate as a function of this reduced Mach number. This approach will be used in our article. More generally, for nonhomogeneous, unmagnetized fluids, Shimony et al.²⁶ show numerically that the dispersion relation, a polynomial of 6th order, only varies depending on two parameters: the Atwood number and the convective Mach number $M_c$ defined as $\Delta U/(a_{01} + a_{02})$ where $\Delta U$ is the shear velocity and $a_{0i}$ is the speed of sound in region i.

Ong and Roderick²⁷ derived stability diagrams from linear theory for a homogeneous, magnetized, incompressible flow with a finite shear layer thickness. They showed the stabilization of the layer due to compressibility effects when a perpendicular magnetic field to the flow is applied. In that case and for a parallel magnetic field, defining the Alfvénic Mach number, $M_A = \Delta U/v_A$, with $v_A$ the Alfvénic speed, Chandrasekhar¹ showed the flow is fully stabilized for $M_A \le 2$. For $2 < M_A < 5$, the flow is linearly unstable and generates KHI with the resulting structures remaining only wavelike due to the tension of magnetic field lines²⁸. For $M_A > 5$, vortices fold the interface between the two flowing plasmas and roll-up magnetic field lines, creating magnetic inversion layer in the 2D plane, where magnetic reconnection occurs as “Type II Vortex Induced Reconnection” (VIR) if some non-ideal MHD term is active²⁹,³⁰.

Miura and Pritchett²² treated the more general case of a magnetized compressible flow for arbitrary magnetic field orientations with a finite shear layer thickness. They established numerically seminal stability diagrams for a large range of flow and magnetic field parameters. The analysis,
however, was performed for a homogeneous plasma, leading to $A = 0$ and to a Mach number $M = 2M_c$. For finite shear layers, the growth rate follows a bell-shaped profile with respect to the wavenumber, $k$, of the perturbation, both in the incompressible and compressible regimes$^{22,27}$. However, the growth rate is linear with $k$ for a vortex sheet$^1$. In all these studies, the linear theory approach for magnetized KHI has demonstrated its capacity for understanding how the magnetic field stabilizes the flow at the first order.

González and Gratton$^{16,17}$ presented a complete development of KHI for nonhomogeneous, compressible plasmas in the presence of a magnetic field of any orientation. They solved an 8th order polynomial dispersion relation using a graphical method (rather than a 6th order polynomial for nonhomogeneous, unmagnetized fluids). They show that the KHI can be overstable ($\omega_r \neq 0$ and $\gamma > 0$) or monotonically unstable (zero real angular frequency $\omega_r$ and positive growth rate $\gamma$) for homogeneous plasmas while it is always overstable for nonzero Atwood numbers. They introduced the notion of main and secondary KH waves (MKH and SKH, respectively) that develop depending on the shear velocity range. MKH waves are defined by higher growth rates than SKH waves. SKH waves only appear in compressible plasmas for a lower range of velocity, which is stable in the incompressible limit. Turkakin, Rankin, and Mann$^{20}$ followed González and Gratton$^{16,17}$'s work and applied it to the Earth's magnetosphere, assuming a vortex sheet and a bounded domain. They revealed that SKH are predominant in the near-Earth magnetosphere flanks while MKH waves are especially present on the nightside magnetopause. Among the two, MKH is the most important wave generation mechanism for particle and energy transport in the magnetosphere. However, the configuration with a parallel magnetic field to the plasma flow is barely mentioned in their work and not investigated, leaving gaps our study will fill with a theoretical and numerical approach.

In the nonlinear, magnetized regime, KHI has been investigated theoretically by Baek and Sohn$^{31}$ for a parallel magnetic field and numerically for a magnetic field of any orientation$^{28-30,32}$. The nonlinear regime is not investigated in our study.

In our study, as well as in a vast majority of the literature, hydromagnetic equations in a collision-dominated regime are used for studying the Kelvin-Helmoltz instability. In that case, the gas pressure is approximated by a scalar. However, there exist a few studies which address the question of the effect of an anisotropic pressure, mostly for collisionless plasmas. In that case, the pressure can then be modeled by a vector, with a state law closing each component of the pressure vector$^{33}$. A more rigorous but numerically costly alternative is to solve for more hydrodynamic moments Kuznetsov and Dzhalilov$^{34}$.

Our theoretical and numerical investigations, despite being related to future laser experiments, are not specially devoted to Inertial Confinement Fusion (ICF) experiments. In current ignition experiments, such as the ones performed with NIF, no external B-fields are applied$^{35}$. The presence of an external magnetic field could be envisioned to boost

the yield in the future$^{36}$. General questions remain therefore to be addressed, such as the effect of the magnetic field on the Rayleigh-Taylor instability (RTI)$^{37,38}$. Regarding RTI at the ablation front, an extensive experimental work has already been performed in indirect or direct drive$^{39,40}$. RTI and Richmyer-Meshkov instability (RMI) are the principal instabilities in ICF since the density gradients in ICF capsules are primarily aligned with the direction of acceleration/deceleration. Shear flows and KH instability appear as secondary instability$^{41}$. KHI has been observed at the interface in an indirect drive and deteriorating the gold/gas mixing in ICF experiments$^{42}$. KHI vortex suppression in ICF experiments through imposing a constant magnetic field is of great interest. However, it has not been deeply studied in the literature, at least in the fully compressible regime. KHI has, however, been studied in laser experiments performed at the OMEGA laser facility$^{43-48}$ at low convective Mach number ($M_c < 1$). In these experiments, the velocity shear is created by the ablation of foam and plastic of different densities after the passage of a laser-driven blast-wave. The laser-driven blast-wave, namely shock, creates a baroclinic effect, which is the production term of the hydrodynamic vorticity$^{49}$. The laser experiment in Hurricane et al.$^{43}$ is particularly interesting because of the possibility to reach the compressible regime. Our article aims to create a link between linear theory and the suppression of KH vortices by applying an external magnetic field to a nonhomogeneous, compressible, magnetized plasma flow. To this end, we limit the study to a vortex sheet interface with an imposed parallel magnetic field. This will be especially useful for laser experiments, or any similar shear flows. However, we focus upon fundamental understanding of KHI and the potential applications are just motivation. The questions we target are: what physical quantities drive KHI? Which physical mechanisms are involved in the development and stabilization of KHI for these plasmas? What magnetic field intensity is required such that the shock induced shear flow is partially or fully stabilized during a KH ablative laser experiment?

The article is organized as follows. In section II, we show that the dispersion relation found by González and Gratton$^{16,17}$ can be rewritten for such a plasmas as a function of four dimensionless variables that drives the stability of the flow: the Atwood number, the convective Mach number, the plasma-β parameter, and the angle $\phi$ formed by the wave vector $k$ with the direction of the flow. The strength of this new formalism lies in the use of a convective Mach number, a convective Alfvénic Mach number, and a total convective Mach number, which combines both of them, which provides a generalization of the unmagnetized, compressible case$^{24-26}$. We derive an expression of the KHI growth rate which is analytic for a homogeneous flow with $A = 0$. This expression converges both toward Choudhury$^{25}$ expression for an unmagnetized, compressible flow and Chandrasekhar$^{1}$ expression for a magnetized, incompressible flow. Otherwise, the dispersion relation (of 8th order) has to be solved numerically, which is done in section III. There, we compute general stability diagrams of magnetized KHI for the entire range of the $(A, M_c, \beta)$ parameters.
ters and discuss the different stability regimes of the flow. This provides new theoretical and general insights of the KHI growth for magnetized/unmagnetized flows. In section IV, numerical simulations are performed using the 3D MHD GAMERA⁵⁰ code, which is currently used for global simulations of the magnetosphere⁵¹⁻⁵⁵. GAMERA is also specifically used for the simulation of magnetized KHI⁵⁶,⁵⁷. Recently, Briard et al.⁵⁸ showed that GAMERA is able to reproduce the diagrams of Miura and Pritchett²² for A ≤ 0.3. Comparisons of theoretical and simulated growth rates are performed in this study to verify the theoretical linear growth rates, assess the numerical accuracy of the code, and for future assessment of differences with nonlinear regimes. Finally, in section V, we apply our results to the analysis of the KHI experiment performed by Hurricane et al.⁴³ at the OMEGA laser facility, showing linear theory succeeds to provide accurate estimates. We also discuss how our results can be applied to the design of future experiments, particularly in the high-Mach magnetized regime, which new experiments planned at the National Ignition Facility (NIF) reside. The extended laser capabilities of NIF enable sustaining high energy density plasma flows for longer periods of time, offering the capabilities to explore unchartered, turbulent-like plasma regimes³⁹,⁴⁰.

## II. MAGNETIZED KHI IN THE VORTEX SHEET APPROXIMATION FOR NONHOMOGENEOUS AND COMPRESSIBLE PLASMA FLOW

### A. Ideal MHD equation for compressible plasma

KHI is derived from the linearization of the ideal MHD equations:

$$ \frac{\partial \rho}{\partial t} = -\nabla \cdot (\rho U), \qquad (1) $$

$$ \frac{\partial \rho U}{\partial t} = -\nabla \cdot (\rho U \otimes U + \bar{I}P) - \nabla \cdot (\bar{I} \frac{B^2}{2} - B \otimes B), \qquad (2) $$

$$ \frac{\partial B}{\partial t} + (U \cdot \nabla)B = (B \cdot \nabla)U - B(\nabla \cdot U). \qquad (3) $$

By making the assumption that the flow is homoentropic, we can write

$$ \frac{d(P\rho^{-\Gamma})}{dt} = 0, \qquad (4) $$

where $\Gamma = 5/3$ is the specific heat coefficient for monoatomic gas.

We define the real magnetic field as $\tilde{B}$ and introduce the reduced magnetic field, $B = \tilde{B}/\sqrt{\mu_0}$, which we use going forward to ease the calculus.

<figure><img src="image_1.png" /><figcaption>FIG. 1. Scheme of the vortex sheet configuration with $U_1$ and $U_2$ the velocities on each side of the discontinuity interface, $\rho_1$ and $\rho_2$ the densities on each side of the discontinuity interface, $k = (k_x, 0, k_z)$ the wavevector in the plane $(x, z)$ and $\phi$ the angle which $k$ makes with the $x$ direction. $B_0$ is the magnetic field and $P_0$ is the constant pressure, both imposed over the whole domain.</figcaption></figure>

In perturbation theory, each variable is written as the sum of a background (zero order) term and a perturbative term of order 1, so that each physical quantity is developed in the following equations

$$ U(x,t) = U_0(y)\vec{x} + u(x,t), \qquad (5) $$

$$ \rho(x,t) = \rho_0(y) + \rho'(x,t), \qquad (6) $$

$$ P(x,t) = P_0 + p(x,t), \qquad (7) $$

$$ B(x,t) = B_0\vec{x} + b(x,t), \qquad (8) $$

where $u = (u, v, w)$, $b = (b_x, b_y, b_z)$, and $x = (x, y, z)$. In this article, $P_0$ and $B_0$ are chosen constants. More generally, we call $F_0$ the background flow variables and $F'$ (or $f$) the perturbed flow variables.

Fourier analysis allows us to understand the wave dynamics of the flow, so that for each perturbative term, we can write $f = \hat{f} \exp i(k \cdot x - \omega t)$, with $k$ the wavevector, $\hat{f}$ the amplitude of the perturbed flow variables $f$, and $\omega$ the complex angular frequency. KHI is an initial value problem, for which each initial amplitude $\hat{f}$ has a very low value (a seed such as $\hat{f}/\lambda \ll F_0/L$ where $\lambda$ is the wavelength and $L$ is the reference length). The seed then grows exponentially over time. This translates mathematically into the fact that the angular frequency $\omega$ is complex, $\omega = \omega_r + i\gamma$, in which $\omega_r$ is the real part of the angular frequency, $\gamma$ is the imaginary part representing the growth rate, and the wavenumber $k$ is real.

### B. The incompressible KHI case in the vortex sheet limit

The development of KHI in an incompressible plasma flow is a simple configuration to work on and corresponds
to real configurations in space physics and magnetospheric plasma flows$^{27,59}$. For a parallel magnetic field, Fig. 1 depicts schematically the vortex sheet approximation, i.e. an infinitely thin interface between the two sheared flows, that is investigated here.

In the incompressible regime, KHI is studied by linearizing equations (1), (2) and (3) and setting $\nabla \cdot U = 0$. The dispersion relation of this system was first derived by Chandrasekhar$^1$ in the vortex sheet approximation for a nonhomogeneous incompressible flow. He gave the theoretical relation that exists between the incompressible growth rate $\gamma_*$ with the shear velocity, the density ratio, and the magnetic field in a perpendicular or a parallel configuration. For a parallel magnetic field, assuming the flow is only perturbed by propagating waves in the flow direction $k = k_x$, $\gamma_*$ is

$$ \gamma_* = k_x \sqrt{\frac{(1 - A^2)\Delta U^2}{4} - \frac{2B^2}{(\rho_h + \rho_l)}}, \quad (9) $$

where $k_x$ the wavenumber in the direction parallel to the flow, $\Delta U$ is the shear velocity, $B$ is the reduced magnetic field, $\rho_h$ and $\rho_l$ are the highest and lowest density throughout the flow interface and A is the Atwood number defined by

$$ A = \frac{\rho_h - \rho_l}{\rho_h + \rho_l}. \quad (10) $$

For the next paragraphs, we set $\rho_1 = \rho_h$ and $\rho_2 = \rho_l$ such that $\rho_1 \ge \rho_2$. We use the *, either as a lower or an upper script, to designate that variable in the incompressible regime.

We now express the growth rate, $\gamma_*$, in terms of dimensionless quantities, such as the Mach number $M$ and the Alfvénic Mach number $M_A$, in order to address magnetized, incompressible flows. Because of the nonuniformity of the flow due to nonzero Atwood numbers, we choose to use the convective Mach number, $M_c$, defined from the convective velocity, $U$, by $U = \Delta U/2 = (U_1 - U_2)/2$, such that $M_c = 2U/(a_{01} + a_{02})$ where $a_{0i} = \sqrt{\Gamma P_0 / \rho_i}$ is the sound speed in each region $i$. By extension of the convective Mach number, we introduce another dimensionless parameter, the convective Alfvénic Mach number, $M_{Ac}$ defined by $M_{Ac} = 2U/(v_{A1} + v_{A2})$, where $v_{Ai} = B_0 / \sqrt{\rho_i}$ is the Alfvénic speed in region $i$. Using the set of parameters ($A$, $M_c$, $M_{Ac}$), the incompressible growth rate $\gamma_*$ can be conveniently written as

$$ \frac{\gamma_*^2}{\gamma_U^2} = (1 - A^2) \left( - \frac{4}{M_{Ac}^2 (\sqrt{1 - A} + \sqrt{1 + A})^2} + 1 \right), \quad (11) $$

where the growth rate $\gamma_*$ is given by the principal square root of (11) normalized by $\gamma_U = k_x U$. The meaning of the normalization by $\gamma_U$ can be understood for an homogeneous, unmagnetized, incompressible flow, for which we have $\gamma_*(A = 0, M_{Ac}^{-1} = 0) = \gamma_U(A = 0, M_c \ll 1)$. We can also define $\gamma_* = k_x c_i^*$, where $c_i$ is the imaginary phase velocity. $\gamma_*/\gamma_U$ then simplifies to the velocity ratio $c_i^*/U$. Note that

$\gamma_*/\gamma_U$ in Eq. (11) depends only on the A and $M_{Ac}$ numbers and not the convective Mach number. The dependence in velocity is thus only carried by $\gamma_U$, see section V F for further discussion. By setting $\gamma_*/\gamma_U = \varepsilon$ in the expression (11), we can easily find the limiting convective Alfvénic Mach number $M_{Ac,lim}^*$ in the incompressible regime below which the flow is stabilized, given by

$$ M_{Ac,lim}^* = \frac{2}{(1 - \varepsilon) (\sqrt{1 - A} + \sqrt{1 + A})}. \quad (12) $$

For $A = 0$, this gives $M_{Ac,lim}^* = 1/(1 - \varepsilon)$. $M_{Ac}$ above 1 (resp. $M_{Ac}^{-1} < 1$) stabilizes the flow until full stabilization for $M_{Ac} \le 1$ (resp. $M_{Ac}^{-1} \ge 1$). This is equivalent to the well-known criteria$^1$ $B^2 \ge \rho \Delta U^2/4$ or $2M_{Ac,lim}^* = M_A = 2$.

### C. Extension to unmagnetized compressible KHI

The hydrodynamic compressible flow satisfies equations (1), (2), and (3), the compressibility Eq. (4) and $B = 0$. Miles$^{24}$, equation (5.9) derived the dispersion relation for a homogeneous, compressible flow in the vortex sheet configuration. In addition, Shimony et al.$^{26}$ have shown that it is possible to generalize the formulation of Choudhury$^{25}$, equation (28) based on a reduced Mach number, C, for homogeneous plasma flows ($C = M_c(A = 0) = M/2$) to the same expression but with a convective Mach number $M_c$ to account for a nonzero Atwood number. Therefore, rewriting Miles$^{24}$, equation (5.9) using the convective Mach number $M_c$ gives

$$ \frac{\gamma_*^2}{\gamma_U^2} = - \left( 1 + \frac{1}{M_c^2} \right) + \frac{2}{M_c} \sqrt{1 + \frac{1}{4M_c^2}}. \quad (13) $$

$\gamma/\gamma_U$ is given by the principal square root of (13), a 6th order polynomial function of $\gamma$. Eq. (13) turns out to be a simple expression of $\gamma$, which helps to understand the dynamics of compressible KHI and becomes a reference limit when A starts to vary. In the nonhomogeneous plasma flow ($A \ne 0$) configuration, the dimensionless growth rate $\gamma/\gamma_U$ depends only on two dimensionless parameters: the A and $M_c$ numbers. For this case, the dispersion relation is a 6th order polynomial function in $\gamma$. A simple theoretical relation has yet to be found, therefore, solutions will be computed numerically, though using a new formulation of the general dispersion relation derived below.

### D. Nonhomogeneous magnetized compressible KHI with a parallel magnetic field

By linearizing the ideal MHD equations for compressible plasmas (cf. subsection II A), we find a system of 8 equations and 8 unknowns in the Fourier domain. After a series
of calculations (not presented here), we reach the differential equation on the transverse velocity perturbation $\hat{v}$ given by

$$ \left[-\rho_0\Omega + \frac{k_x^2 B_0^2}{\Omega}\right] \hat{v} + D\{DU_0 k_x \rho_0 f \hat{v}\} + D\{\rho_0 \Omega f D \hat{v}\} = 0, \quad (14) $$

with

$$ f(M, M_A, \omega/\omega_{DS}) = \frac{H_1(M_A, \omega/\omega_{DS})}{k^2 - k_x^2 H_2(M, M_A, \omega/\omega_{DS})}, \quad (15) $$

$$ H_1(M_A, \omega/\omega_{DS}) = 1 - \frac{1}{M_A^2 (1 - \omega/\omega_{DS})^2}, \quad (16) $$

$$ H_2(M, M_A, \omega/\omega_{DS}) = \frac{M^2(1 - \omega/\omega_{DS})^2}{(M/M_A)^2 + H_1(M_A, \omega/\omega_{DS})}, \quad (17) $$

and

$$ \Omega = \omega - \omega_{DS}, \quad (18) $$

where $D = d/dy$, $k^2 = k_x^2 + k_z^2$ is the square of the wavenumber, $M = U_0/a_0$ is the Mach number, $M_A = U_0/v_A$ is the Alfvénic Mach number and $\Omega$ is the Doppler shifted angular frequency for which $\omega_{DS} = k_x U_0$ is the Doppler shift that refers to the difference between the frequency measured by an observer at rest and the frequency measured by an observer moving with the mean flow $U_0$. We can also interpret the Doppler shifted angular frequency as a growth rate normalization term since $\gamma_U = k_x U$, as done in the previous section (II B). In this derivation, gravity is not taken into account.

It is possible to numerically solve the partial differential Eq. (14), but this is not attempted here. We rather stay in the limit of the vortex sheet approximation for which the shear layer's thickness, $\sigma$, is infinitely small compared with the wavenumber $k$, i.e. $k\sigma \ll 1$. For the sake of clarity, this configuration is presented in Fig. 1. This leads to a reduced set of differential equations for regions 1 ($y > 0$) and 2 ($y < 0$) given by

$$ D^2 \hat{v} - \kappa^2 \hat{v} = 0. \quad (19) $$

The dispersion relation is then found directly and written as

$$ \left(1 - \frac{k_x^2 v_{A1}^2}{\Omega_1^2}\right) \rho_1 \frac{\Omega_1^2}{\kappa_1} = - \left(1 - \frac{k_x^2 v_{A2}^2}{\Omega_2^2}\right) \rho_2 \frac{\Omega_2^2}{\kappa_2}, \quad (20) $$

where

$$ \kappa_{1,2}^2 = \frac{k^2 \left( 1 - \left( \frac{\Omega_{1,2}}{ka_{01,2}} \right)^2 + \left( \frac{v_{A1,2}}{a_{01,2}} \right)^2 - \left( \frac{v_{A1,2}k_x}{\Omega_{1,2}} \right)^2 \right)}{\left( 1 + \left( \frac{v_{A1,2}}{a_{01,2}} \right)^2 - \left( \frac{v_{A1,2}k_x}{\Omega_{1,2}} \right)^2 \right)}. \quad (21) $$

The notation is chosen to respect the condition $\rho_1 \ge \rho_2$. The dispersion relation (20) is then described using 7 parameters: $(U, \rho_1, \rho_2, B_0, P, k_x, k)$. We can then reduce the physical space to only 4 dimensionless parameters: $(A, M_c, M_{Ac}, \phi)$, with $\cos \phi = k_x/k$, in (22). This computation is detailed in Appendix A. No assumption on the magnetic field vector has yet been imposed, leading to this general formulation. We find $\mathcal{D}$, the dispersion relation function

$$ \mathcal{D}(A, M_c, M_{Ac}, \phi) = 0, \quad (22) $$

with

$$ \begin{aligned} \mathcal{D}(A, M_c, M_{Ac}, \phi) = & \frac{1+A}{1-A} \cdot \frac{\mathcal{D}_1 \sqrt{\mathcal{D}_2}}{\sqrt{\mathcal{D}_3}} (A, M_c, M_{Ac}, \phi, W_-) \\ & + \frac{\mathcal{D}_1 \sqrt{\mathcal{D}_2}}{\sqrt{\mathcal{D}_3}} (-A, M_c, M_{Ac}, \phi, W_+), \end{aligned} \quad (23) $$

with $W_\pm = \omega/\gamma_U \pm 1$ and $\mathcal{D}_{i=1..3}$ defined in Appendix A. In the case where $k_z = 0$, the dispersion relation can be expressed from only two variables, $\mathcal{D}(A, M_c)$, for the purely hydrodynamic case, and from three variables, $\mathcal{D}(A, M_c, M_{Ac})$, for a magnetized fluid. The novelty of the latter formulation lies in the introduction of a convective Alfvénic Mach number.

Furthermore, one may use the following scaling law to convert easily from $B$ and $M_{Ac}$ to $\beta$-plasma

$$ \beta = \frac{P_0}{P_{mag}} = \frac{2}{\Gamma} \left( \frac{M_{Ac}}{M_c} \right)^2, \quad (24) $$

with $P_{mag} = B_0^2/2$, the magnetic pressure. $(A, M_c, \beta)$ is another set of independent parameters fully characterizing KHI development.

## E. Magnetized homogeneous plasma ($A=0$)

### 1. Analytical growth rate

From the dispersion relation (20), a generalized solution for the specific case of a magnetized homogeneous plasma ($A=0$, $k_z=0$ and $B \neq 0$) can be derived. The details of the computations as well as the identification of stability regions are given in Appendix B. The normalized growth rate can be written as

$$ \frac{\gamma_U^2}{\gamma_U^2} = - \left(1 + \frac{1}{M_c^2}\right) + \frac{2}{M_c} \sqrt{1 + \frac{1}{4M_c^2} - \frac{1}{2M_{c,tot}^2}}, \quad (25) $$

where $M_{c,tot}^2 = M_c^2 + M_{Ac}^2$ and $\gamma/\gamma_U$ is given by the principal square root of (25). To the authors' knowledge, the simple expression of Eq. (25) has not previously been published. Again, the novelty and consequent simplicity comes from the global convective Mach number, $M_{c,tot}$. Note that $M_{c,tot}$ can be written as
$$M_{c,tot}^2 = U^2 \left( \frac{1}{c_{A,m}^2} + \frac{1}{c_{s,m}^2} \right), \quad (26)$$

with the mean Alfvénic speed, $c_{A,m} = (v_{A1} + v_{A2})/2$, the mean sound speed, $c_{s,m} = (a_{01} + a_{02})/2$, the magnetoacoustic speed, $c_{mgs} = \sqrt{c_{A,m}^2 + c_{s,m}^2}$, so that the total Mach number can be written as $M_{c,tot} = U/(c_{A,m}c_{s,m}/c_{mgs})$. We see the velocity $U$ in Eq. (26) is normalized by the harmonic mean of the magnetoacoustic speed rather than directly with the magnetoacoustic speed, $c_{mgs}$.

The square root term in Eq. (25) is always positive for $0 \le M_{Ac}^{-1} \le \sqrt{1/(1-\sqrt{3}/2)} \approx 2.73$. As shown in the subsequent section, stabilization is reached for $0 \le M_{Ac}^{-1} \lesssim 1.1$, thus, equation (25) is fully defined.

In the limit of a zero $M_{Ac}$, the solution converges to the expression of Miles$^{24}$ and Choudhury$^{25}$ (presented previously in Eq. (13)). The growth rate appears to be simply correlated to the convective Alfvénic Mach number. This result also supports our choice to use a convective Alfvénic Mach number $M_{Ac}$ in the general dispersion relation for any Atwood number.

## 2. KHI stability thresholds

We can find the limiting $M_{Ac}^{-1}$ above which the flow is stable. To do so, we set $\gamma/\gamma_U = 0$ in Eq. (25). $M_{Ac}^{-1}$ then becomes

$$M_{Ac,lim}^{-1}(A=0, M_c) = \sqrt{\frac{1 - M_c^2/2}{1 - M_c^2 + M_c^4/2}} \quad (27)$$

The domain of validity of this equation is given by $M_c \in [0, \sqrt{2}]$ signifying that any flow traveling at Mach number greater than $\sqrt{2}$ will not generate KHI. At lower Mach number, to completely stabilize the flow, the applied magnetic field varies depending on $M_c$. Studying the (27), we find the threshold starts increasing for the low-$M_c$ regimes and then decreases after a peak. The peak occurs for $M_c = \sqrt{2-\sqrt{2}} \approx 0.765$, which, when inserted into (27), leads to $M_{Ac,lim}^{-1} \approx 1.1$. Therefore, KHI magnetic stabilization occurs for $0 \le M_{Ac}^{-1} \lesssim 1.1$.

Using the dispersion relation (11), for $M_c \in [0, \sqrt{2}]$, we find the $\beta \le \beta_{lim}$ threshold for which KHI is stabilized to be

$$\beta_{lim}(A=0, M_c) = \frac{2(1-M_c^2+M_c^4/2)}{\Gamma M_c^2 - M_c^4/2} \quad (28)$$

Finally, the stabilization can also be formulated in terms of a magnetic field limit $B_{lim}$ for which the flow is stabilized for $B \ge B_{lim}$, with $B_{lim}$ defined as

$$B_{lim}(A=0, M_c) = M_c \sqrt{\Gamma P_0} \sqrt{\frac{1 - M_c^2/2}{1 - M_c^2 + M_c^4/2}} \quad (29)$$

$B_{lim}$ varies with a bell shape with respect to $M_c$ with two notable values: $B_{lim}(A=0, M_c) = 0$ and $B_{lim}(A=0, \sqrt{2}) = 0$. For $P_0 = 1.0$, the bell-shaped curve starts from $B_{lim} = 0.0$ at $M_c = 0.0$, reaches a maximum $B_{lim} = \sqrt{\Gamma P_0} \approx 1.29$ at $M_c = 1.0$ and ends at $M_c = \sqrt{2}$.

The complete set of dimensionless parameters $(A, M_c, M_{Ac}, \phi)$ allows us to find a series of scaling laws that ease the understanding of KHI dynamics for future laser experiments or numerical simulations. In the next section, we present a global analysis using the new set of parameters to clarify the various KHI regimes for nonhomogeneous, magnetized, compressible flows.

## III. KHI GROWTH RATES FOR NONHOMOGENEOUS MAGNETIZED COMPRESSIBLE FLOWS

### A. Method

The dispersion relation $\mathcal{D} = 0$ in Eq. (22) is solved numerically which results in a single value of $\gamma$ for each quadruplet $(A, M_c, M_{Ac}, \phi)$. For $A=0$, we use the analytical solution given in (25). For nonzero Atwood numbers, all the solutions are overstable so that we need to solve numerically $\omega_r \neq 0^{16,17}$. There are two calculation steps in order to find $\gamma$.

The first step is to estimate a couple $(\omega_r, \gamma)$, and then specify the value of $\gamma$ with a more precise calculation result. The couple $(\omega_{r,1st}, \gamma_{1st})$ is found by minimizing $\mathcal{D}$ by iterating in a set of $512^2 (\omega_r, \gamma)$. The second step is made by setting up $\omega_{r,2nd} = \omega_{r,1st}$ and iterating over $\gamma$ until finding the minimum for $|\mathcal{D}|$ at $\gamma_{2nd}$. The method is illustrated in Fig. 2 where we plot $|\mathcal{D}|$ versus the normalized growth rate. The solution of $|\mathcal{D}| = 0$ corresponding to a given $\gamma/\gamma_U$ is the positive (i.e. unstable) solution. This is shown for five quadruplets $(A, M_c, M_{Ac}, \phi)$ with $A = 0.8, M_c = 0.8, \phi = 0$ and $M_{Ac}^{-1} \in \{0.0, 0.25, 0.5, 0.75, 1.0\}$.

### B. Stability diagrams of KHI for $\phi = 0$

We quantify the behavior of KHI for configurations where $\phi = 0$. This corresponds to the most interesting case currently for laser experiments (cf. section V) for which the KH waves propagate only in the direction of the flow.

### 1. Stability diagrams for fixed ($M_{Ac}^{-1}$, A)

The normalized growth rate, $\gamma/\gamma_U$, of the dispersion relation (20), is solved numerically using the technique described in section III A to obtain the diagrams presented in Figs. 3 and 4. The normalized growth rates are plotted as a function of the convective Mach number for a constant Atwood number. Each line plot corresponds either to a different $M_{Ac}^{-1}$ in Fig. 3 or a different magnetic field B in Fig. 4.
<figure><img src="image_2.png" /><figcaption>FIG. 2. Solutions $\gamma/\gamma_U$ of the dispersion relation (Eq. (22)) for $A = 0.8$, $M_c = 0.8$, $M_{Ac}^{-1} \in \{0.0, 0.25, 0.5, 0.75, 1.0\}$ and $\phi = 0$.</figcaption></figure>

Dimensionless growth rates are determined numerically starting from $M_c = 0.05$, due to the imprecision of the calculation for lower convective Mach numbers, up to $M_c^+ \le \sqrt{2}$, the maximal convective Mach number above which stability is reached. In the incompressible regime, the solution of Chandrasekhar¹ (Eq. (11)) is used and plotted from $M_c = 0$ up to $M_c = 0.01$ in Fig. 3 as red tick marks for all Atwood numbers. Figs. 3 and 4 show that the incompressible limit of Chandrasekhar¹ is verified by the numerical solutions of the compressible dispersion relation (20). Moreover, this solution appears to be accurate up to $M_c \approx 0.05$.

Diagrams for $A = 0$ are generated using the analytical relation of (25). Note that for an unmagnetized flow, the $M_{Ac}^{-1} = 0$ level line corresponds to Miles²⁴, equation (5.9)'s solution. The dimensionless growth rate decreases with increasing compressibility, i.e. increasing $M_c$, up to complete stabilization at $M_c = \sqrt{2}$. $\gamma/\gamma_U$ reaches its maximum 1.0 for $M_c = 0.0$ and $A = 0$.

For nonzero $A$, $\gamma/\gamma_U$ decreases with increasing Atwood numbers, as shown in Shimony et al.²⁶ for unmagnetized flows. This behavior also applies for magnetized flows, see Fig. 3. For $A < 0.6$, the maximal value of $\gamma/\gamma_U$ is reached at $M_c = 0$ for unmagnetized flows. Nonzero $A$ also reshapes the stability curves, causing the growth rate to increase with increasing convective Mach number for some regions. Adding a constant magnetic field shrinks $\gamma/\gamma_U$ with respect to $M_c$ for all Atwood numbers. Decreasing $\gamma/\gamma_U$ with increasing $M_c$ occurs only up to $A \sim 0.6$ and $M_{Ac}^{-1} \le 0.7$. Otherwise the normalized growth rate can increase with increasing $M_c$, but only for $M_c < M_{c,max} < 1$, with $M_{c,max}$ defined as the convective Mach number for which the growth rate is maximum for a given $A$.

For each iso-$M_{Ac}^{-1}$ and $A$ curves in Fig. 3, the growth

rates are defined for a single continuous range of $M_c$ and are zero outside this interval. This means secondary KH (SKH) waves (González and Gratton ¹⁶,¹⁷), which would be defined in a second disjoint domain in ($M_c$, $A$), do not exist in the configuration with a parallel magnetic field. Only the main KH (MKH) waves are a solution of the compressible dispersion relation.

## 2. Stability diagrams for fixed ($A, A$)

In Fig. 4, the iso-values of $B$ are calculated for a constant pressure $P_0 = 1.0$. Using (24) and assuming two constant pressures $P_{M1}$ and $P_{M2}$ while conserving the couple ($M_c, M_{Ac}$) such that $\beta_{M1} = \beta_{M2}$, we can therefore replace any $B$-field value by another. Therefore, Fig. 4 remains universal and useful in experimental considerations for specific initially constant $\beta$-values.

Contrary the profiles plotted versus $M_{Ac}^{-1}$, the normalized growth rates have a bell shape as soon as $B \ne 0$. When increasing the magnetic field, the minimum convective Mach number above which the flow destabilizes, noted $M_c^-$, and the maximum convective Mach number before stabilization, noted $M_c^+$, vary less and less as the Atwood number increases.

In the low-$M_c$ regime, defined as $M_c < M_{c,max}$, complete KHI stabilization is only possible thanks to the magnetic field. In the high-$M_c$ regime ($M_c > M_{c,max}$), three effects can each lead to complete stabilization: the magnetic field effect, the Atwood number and the compressibility effect. Compressibility also causes the iso-$B$ lines to be closer to each other, thus diminishing the magnetic stabilization.

For $B = 0$ and for $M_c > M_{c,max}$, the flow is completely stabilized when $M_c^+ = \sqrt{2}$ for $A = 0$. $M_c^+$ then decreases down to $M_c^+ \approx 1.25$ for increasing $A$.

## C. Stability diagrams for fixed ($A, M_c, M_{Ac}^{-1}$)

Fig. 5 presents two-dimensional plots of $\gamma/\gamma_U$ in the 2D ($A, M_c$)-parameter space for various values of $M_{Ac}^{-1}$. The color scale indicates the magnitude of the dimensionless growth rate, which goes from dark blue (0) to yellow (1). The white curves are iso-$\gamma/\gamma_U$ lines.

The upper left panel of Fig. 5 corresponds to the unmagnetized case ($M_{Ac}^{-1} = 0$). Shimony et al.²⁶ have performed a series of 15 simulations for the unmagnetized case, which they then fit. They obtained similar 2D profiles (their figure 2b) as those presented here. Here, by solving the theoretical dispersion formulation for all parameters, we can reach highly-resolved profiles without fits. The $\gamma/\gamma_U$ iso-lines are almost horizontal at low $A$, showing that the growth rate, while highly dependent on $M_c$, only weakly depends on $A$ in the range $A \le 0.2$ for $M_c \le 0.4$ up to $A \le 0.8$ for $M_c \ge 1.0$. Conversely, for high $A$ ($A \ge 0.8$), the trend is reversed with contours of constant $\gamma/\gamma_U$ which are now vertical, staying quasi-constant for varying $M_c$, but changing greatly with respect to $A$. These trends are particularly
<figure><img src="image_3.png" /><figcaption>FIG. 3. Dimensionless growth rate, $\gamma/\gamma_U$, versus $M_c$ for iso-values of the convective Alfvénic Mach number $M_{Ac}$ and (from top left to bottom right) $A = 0.0, 0.2, 0.4, 0.6, 0.8, 0.95$. (red ticks) The incompressible solution of Chandrasekhar$^1$ for $M_c = 0$.</figcaption></figure>

<figure><img src="image_4.png" /><figcaption>FIG. 4. Dimensionless growth rate, $\gamma/\gamma_U$, versus $M_c$ for iso-values of the magnetic field $B$ for $P_0 = 1.0$ and (from top left to bottom right) $A = 0.0, 0.2, 0.4, 0.6, 0.8, 0.95$.</figcaption></figure>
marked for $\gamma/\gamma_U < 0.7$. The contour lines bend more and more for $\gamma/\gamma_U > 0.7$. The strongest $\gamma/\gamma_U > 0.9$ are confined to low A and low $M_c$.

In the magnetized case, the diagrams for $M_{Ac}^{-1} > 0$ show that decreasing the Alfvénic Mach number stabilizes the system. The corresponding contour plots of $\gamma/\gamma_U$ remain qualitatively the same as their pure hydrodynamics counterpart for all values of $M_{Ac}^{-1}$ up to 0.8. For $M_{Ac}^{-1} = 0.8$ and above, the dimensionless growth rate level lines are no longer vertical for $M_c \lesssim 1$ and large A $\gtrsim 0.6$.

Furthermore, the limit $M_c = M_c^+$ for which the flow is stabilized, decreases with respect to increasing $M_{Ac}^{-1}$. The stability zone $M_c > M_c^+$ extends more toward $M_c = 1$ as $M_{Ac}^{-1}$ increases so that the $M_c \in [1, \sqrt{2}]$ region becomes more stable. The lower stability zone $M_c < M_c^-$ appears only for $M_{Ac}^{-1} = 1$, which forms an island of instability. Thus, for all Atwood numbers, in the $M_{Ac}^{-1} = 1$ limit, the convective Mach number for which the maximum dimensionless growth rate $\max(\gamma/\gamma_U)$ is reached departs from $M_c = 0$ (where the maximum is otherwise) and increases (up to $M_c \approx 0.6$ as soon as $A > 0.15$). Otherwise, for $M_{Ac}^{-1} \le 0.8$, $\max(\gamma/\gamma_U)$ is always reached for $M_c = 0.0$ and $A = 0.0$.

For sufficiently large $M_{Ac}^{-1} \approx 1$, the island structure also creates a zone of stability for A close to zero ($A < 0.15$). This means a slight inhomogeneity in the plasma ($A \neq 0$) is required to destabilize the flow.

## 1. Minimal stabilizing magnetic field

In this section, for each $(A, M_c)$ couple we find the minimal magnetic field threshold, $B_{lim}$, above which the flow becomes stable. To do so, we assume that the flow is stable for growth rates below an arbitrary threshold set as $\gamma/\gamma_U = 0.01$.

This leads to the stability diagrams shown in Fig. 6, representing the threshold in terms of either $B_{lim}$ (left), $M_{Ac,lim}^{-1}$ (center), or $\beta_{lim}$ (right). We used a polynomial interpolation of high order to better draw the limit of $\gamma/\gamma_U = 0$. The red area corresponds to regimes for which the flow is unstable for all A below 0.95. The blue region is stable for all A due to the presence of a sufficiently strong magnetic field. The dashed level lines represent the boundary of the stability region for a given A, noted $A_{th}$. For instance, all flows with A below $A_{th} = 0.8$, which are unstable due to a too low magnetic field, lie under the green dashed line in the yellow or red regions. The yellow area corresponds to a buffer region for which the stability is conditioned to an A range whose maximum is strictly above 0 and below 0.95.

For low convective Mach ($M_c < 0.25$), we compute the incompressible solution of Chandrasekhar¹ (Eq. (12)). We verify these solutions, plotted in solid lines in Fig. 6, which connect well with the A threshold limit (dashed line). This also shows that the incompressible limit is well approximated up to $M_c = 0.25$ for $B_{lim}$ and $\beta_{lim}$ while it applies only up to $M_c = 0.15$ for the $M_{Ac,lim}^{-1}$ stability diagram.

In Fig. 6(center), the $A=0$ level line corresponds to Eq. (27) and to the $\mathcal{C}_5$ curve in Figure 13 of Appendix B

(in which the ordinate axis is $M_{Ac}$ instead of $M_{Ac}^{-1}$). For $A=0$ in Fig. 6(center) and $M_c > 0.1$, the $M_{Ac,lim}^{-1}$ threshold increases up to a peak, the position of which can be theoretically found at $M_c = M_{c,crit}$ ($A=0 = \sqrt{2-\sqrt{2}} \approx 0.765$). Evaluating Eq. (27) at 0.765 gives the position of the maximum $M_{Ac,lim}^{-1}$ ($A=0 \approx 1.1$). The $M_{Ac,lim}^{-1}$ threshold decreases from its maximum reached at $M_c = 0.765$ and vanish when approaching $M_c = \sqrt{2}$. This decay shows again that an increase of the compressibility effect results in the reduction of $\gamma/\gamma_U$ as $M_c$ increases. That is, smaller $M_{Ac}^{-1}$ becomes sufficient to stabilize the flow for all Atwood numbers.

Furthermore, all iso-A threshold curves intersect at the same point. The coordinate of the intersection can be determined for $A=0$ by setting $M_c = 1.0$ in Eq. (27). This gives an intersection located at $(M_c = 1.0, M_{Ac,lim} = 1.0)$. Note that this intersection is never where the maximum value of the $M_{Ac,lim}^{-1}$ is. Similarly, the $B_{lim}$ diagram (Fig. 6(left)) is generated from Fig. 6(left or center) setting $P_0 = 1$. All iso-A threshold curves intersect at $M_c = 1.0$ for $B \approx 1.29$. In addition, the maximum value of the $B_{lim}$ is now located at $M_c = 1.0$, which makes the $M_{c,crit} = 1.0$ value to become a critical Mach for which the strongest magnetic field is required to damp KHI and above which compressible effects come into play. The fact that for all A, the threshold curves are the same at $M_c = 1.0$ shows that the B-field stabilization effect dominates among the three effects near $M_c = 1.0$. This behavior is also shown in the $\beta$-diagram in Fig. 6(right).

## IV. NUMERICAL COMPARISONS WITH THE 3D MHD GAMERA CODE

This section aims to compare the analytical results of the previous sections with 2D MHD simulations using the 3D MHD GAMERA code (Grid Agnostic MHD for Extended Research). GAMERA is a global compressible MHD code, which is used for magnetospheric physics applications and space weather⁵¹⁻⁵⁵. GAMERA is based on the LFM (Lyon-Fedder-Mobarry) code Lyon, Fedder, and Mobarry⁶⁰ which has been used by the space physics community for over four decades⁶¹. GAMERA uses an arbitrary non-orthogonal curvilinear grid, new memory access patterns, hybrid parallelism, and high order Gaussian quadrature techniques and numerical schemes⁵⁰. The GAMERA code has also been used for the simulation of magnetized KHI⁵⁶,⁵⁷. In Briard et al.⁵⁸, we characterized GAMERA's accuracy for KHI by the reproduction of the seminal KHI stability diagrams of Miura and Pritchett²² for homogeneous (zero Atwood number), compressible, shear layers with finite width. Stability diagrams were further numerically extended to low Atwood numbers ($A \le 0.3$). The simulations of the next section will complement these results for $A = 0.2$ and $A = 0.8$.

### A. The vortex sheet limit with GAMERA

We reach the vortex sheet limit by imposing a very small value of the shear width interface, $\sigma$, in the code. We con-
<figure><img src="image_5.png" /><figcaption>FIG. 5. 2D plots of the dispersion relation solution. $M_c$ is plotted versus A for iso-values of the convective Alfvénic Mach number $M_{Ac}$ and with $\gamma/\gamma_U$ as the colorbar. The plots are ranged from $M_{Ac}^{-1} = 0.0$ to $M_{Ac}^{-1} = 1.0$. This figure was highly inspired by the work of Shimony et al.²⁶.</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>FIG. 6. KHI stability diagrams of (left) the magnetic field threshold $B_{lim}$, (center) the convective Alfvénic Mach number threshold, and (right) the plasma-$\beta$ threshold for which the flow is either (blue) stable or (red) unstable, plotted as function of the convective Mach number $M_c$. $B_{lim}$ and $\beta_{lim}$ are computed here for $P_0 = 1.0$. The red, green and blue dashed lines correspond to the A threshold, $A_{th} \in \{0.0, 0.8, 0.95\}$ so that all $A \in [0, A_{th}]$ are unstable below that line. (plain lines) The incompressible solution of Chandrasekhar¹ is plotted from $M_c = 0$ up to $M_c = 0.25$.</figcaption></figure>

sider a 2D monoatomic (adiabatic coefficient $\Gamma = 5/3$) flow in a box size of length $L = 1.0$, with a pressure $P_0 = 1.0$. The eighth order centered numerical scheme of GAMERA is used with a 1024² grid. The KH configuration for simulating the vortex sheet limit is inspired by the work of McNally, Lyra, and Passy⁶², for which the initial profiles of both the density and the velocity are given by

$$ \rho(t=0) = \begin{cases} \rho_2 - \rho_m e^{(y-1/4)/\sigma}, & y \in [0; 1/4] \\ \rho_1 + \rho_m e^{(-y+1/4)/\sigma}, & y \in [1/4; 1/2] \\ \rho_1 + \rho_m e^{(y-3/4)/\sigma}, & y \in [1/2; 3/4] \\ \rho_2 - \rho_m e^{(-y+3/4)/\sigma}, & y \in [3/4; 1] \end{cases} \quad (30a) $$
<figure><img src="image_7.png" /><figcaption>FIG. 7. ($\gamma/\gamma_U, M_c$) stability diagrams for (top) A = 0.2 and (bottom) A = 0.8. The parameters of the 18 simulations performed using GAMERA are represented by blue and orange numbered circles.</figcaption></figure>

$$U_x(t=0) = \begin{cases} U_2 - U_m e^{(y-1/4)/\sigma}, & y \in [0; 1/4] \\ U_1 + U_m e^{(-y+1/4)/\sigma}, & y \in [1/4; 1/2] \\ U_1 + U_m e^{(y-3/4)/\sigma}, & y \in [1/2; 3/4] \\ U_2 - U_m e^{(-y+3/4)/\sigma}, & y \in [3/4; 1] \end{cases}, \quad (30b)$$

$$\text{where } \rho_m = (\rho_2 - \rho_1)/2, U_m = (U_2 - U_1)/2.$$

The density and velocity profiles are then $\mathcal{C}^0$ with smooth transitions at $y = 1/4$ and $y = 3/4$ where KHI occurs and with a symmetric plane at $y = 1/2$. Furthermore, we use periodic boundary conditions for the x and y directions. These periodic boundary conditions were similarly used in Briard et al.⁵⁸ in which we showed a good accuracy with former theoretical growth rates for unbounded plasmas²². We set a single-mode perturbation with the initial vertical velocity profile defined by $v_y = v_{0,y} \sin(k_x x)$, where the $k_x = 4\pi/L$ mode is defined such that $k_x \sigma \ll 1$. The perturbation amplitude, $v_{0,y}$, and $\sigma$ were adjusted in the simulations to ensure the linear regime is satisfied.

## B. Growth rate computation with GAMERA

In the vortex sheet limit, we numerically compute the kinetic energy density in the y direction of the $k^{th}$ mode (with $k = k_x$), $\mathcal{E}_{y,k}$, and calculate the numerical growth rate $\gamma_{num}$ with the method developed in Briard et al.⁵⁸.

This is illustrated in Fig. 8 in which we show the density profile at $t = 1, 2, 3s$ (Fig. 8 left) and the associated kinetic energy density in the y direction versus time (Fig. 8 right) for one of the simulations discussed below (Case 10). We identify the linear phase of the growth of $\mathcal{E}_{y,k}$ in Fig. 8 (right) and apply a linear regression to find the slope of the fitting line (red line). The slope corresponds to $2\gamma_{num}$⁵⁸, leading to the numerical growth rate $\gamma_{num}$. The uncertainty of the fit, $\Delta\gamma_{num}$, is also given.

## C. Numerical results

A total of 18 simulations were carried out with GAMERA for different plasma flow conditions defined in terms of $(A, B, M_c)$ triplets (setting $\phi = 0$) numbered and shown in Fig. 7. The A number is set to either 0.2 in Fig. 7 (top) or 0.8 in Fig. 7 (bottom). Calculations are done for $k_x \sigma \in [0.025, 0.113]$ with a constant $k_x = 4\pi \text{rad.m}^{-1}$. All parameters, as well as the interface width $\sigma$ and the perturbation velocity amplitude $v_{0,y}$, are reported in Table I. Table I also contains the computed numerical and theoretical dimensionless growth rates, $\gamma_{num}$ and $\gamma_{th}$ as well as the absolute and relative errors. The absolute error between the theoretical and the numerical solution is calculated by $\varepsilon_{abs}^{\pm} = (\gamma_{th} - \gamma_{num} \pm \Delta\gamma_{num})/\gamma_U$. We select the largest absolute error by taking $\varepsilon_{abs} = \max(|\varepsilon_{abs}^{+}|, |\varepsilon_{abs}^{-}|)$. The relative error is defined as $\varepsilon_{rel}^{\pm} = \varepsilon_{abs}^{\pm}/(\gamma_{th}/\gamma_U) = (\gamma_{th} - \gamma_{num} \pm \Delta\gamma_{num})/\gamma_{th}$. Similarly, we select the largest relative error by taking $\varepsilon_{rel} = \max(|\varepsilon_{rel}^{+}|, |\varepsilon_{rel}^{-}|)$. Thus, errors reported in Table I contains the error bounds and are the largest possible errors.

The largest relative error in Table I is between 2.50% and 34.20% while the largest absolute error is between 0.5% and 15.6%. The mean value of $\varepsilon_{abs}$ computed over all cases is 5.14%. The mean value of $\varepsilon_{rel}$ computed over all cases is 16.5% showing the numerical growth rates are on average close to the theoretical value. The relative error is not correlated with either the plasma parameters A or $M_c$ in Table I. $M_{Ac}^{-1}$ appears to be the critical parameter driving the precision of the numerical growth rate, $\gamma_{num}$. The magnetic field smooths the small scales preventing mixing in the shear layer⁵⁸, which explains why $\varepsilon_{rel}$ diminishes from $M_{Ac}^{-1} = 0.1$ to $M_{Ac}^{-1} \sim 0.5$ where the minimal differences are found. Below $M_{Ac}^{-1} = 0.4$, $\varepsilon_{rel}$ remains below 20%. Increasing further $M_{Ac}^{-1} > 0.5$ does now increase the relative error as the flow gets closer to the stability limit.

Although a given set of $(A, M_c, M_{Ac}^{-1})$ always corresponds to a single normalized growth rate, this may not be the case of the associated error since there exist a multitude of different initial parameters $(U, \rho_1, \rho_2, B, P, k_x)$ leading to this $(A, M_c, M_{Ac}^{-1})$ set. We therefore assess the variation of the
<table><caption>TABLE I. (lines) Parameters and results of the 18 cases studied by simulations and theory. Columns correspond to (from left to right) the A number, B-field value, convective Mach number, shear layer width, perturbation amplitude, theoretical dimensionless growth rate, numerical dimensionless growth rate and the error $\Delta\gamma_{num}/\gamma_U$ of the linear regression, the largest difference between the theoretical and the numerical dimensionless growth rate ($\epsilon_{abs}$), the largest relative difference between the theoretical and the numerical dimensionless growth rate ($\epsilon_{rel}$). Computations made for $\rho_1 = 1.0$ and $P_0 = 1.0$ using a $1024^2$ grid.</caption><thead><tr><th>#</th><th>A</th><th>B</th><th>M<sub>c</sub></th><th>M<sub>Ac</sub><sup>-1</sup></th><th>σ</th><th>v<sub>0,y</sub></th><th>γ<sub>th</sub>/γ<sub>U</sub></th><th>(γ<sub>num</sub> ± Δγ<sub>num</sub>)/γ<sub>U</sub></th><th>ε<sub>abs</sub> [%]</th><th>ε<sub>rel</sub> [%]</th></tr></thead><tbody><tr><td>1</td><td>0.2</td><td>0.2</td><td>0.170</td><td>0.911</td><td>0.002</td><td>10<sup>-3</sup></td><td>0.400</td><td>(4.5 ± 0.1) × 10<sup>-1</sup></td><td>6.0</td><td>15</td></tr><tr><td>2</td><td>0.2</td><td>0.2</td><td>0.392</td><td>0.395</td><td>0.004</td><td>10<sup>-4</sup></td><td>0.792</td><td>(7.7 ± 0.3) × 10<sup>-1</sup></td><td>5.2</td><td>6.6</td></tr><tr><td>3</td><td>0.2</td><td>0.2</td><td>1.10</td><td>0.141</td><td>0.004</td><td>10<sup>-4</sup></td><td>0.400</td><td>(3.20 ± 0.01) × 10<sup>-1</sup></td><td>8.10</td><td>20.2</td></tr><tr><td>4</td><td>0.2</td><td>0.8</td><td>0.601</td><td>1.03</td><td>0.004</td><td>10<sup>-4</sup></td><td>0.200</td><td>(2.0 ± 0.2) × 10<sup>-1</sup></td><td>2.0</td><td>10</td></tr><tr><td>5</td><td>0.2</td><td>0.8</td><td>0.810</td><td>0.765</td><td>0.002</td><td>4.54 × 10<sup>-6</sup></td><td>0.380</td><td>(3.20 ± 0.02) × 10<sup>-1</sup></td><td>6.20</td><td>16.3</td></tr><tr><td>6</td><td>0.2</td><td>0.8</td><td>1.19</td><td>0.521</td><td>0.002</td><td>4.54 × 10<sup>-6</sup></td><td>0.200</td><td>(2.03 ± 0.02) × 10<sup>-1</sup></td><td>0.500</td><td>2.50</td></tr><tr><td>7</td><td>0.2</td><td>1.1</td><td>0.790</td><td>1.08</td><td>0.002</td><td>4.54 × 10<sup>-6</sup></td><td>0.100</td><td>(7.80 ± 0.07) × 10<sup>-2</sup></td><td>2.27</td><td>22.7</td></tr><tr><td>8</td><td>0.2</td><td>1.1</td><td>0.935</td><td>0.911</td><td>0.004</td><td>4.54 × 10<sup>-6</sup></td><td>0.207</td><td>1.55 × 10<sup>-1</sup></td><td>5.20</td><td>25.1</td></tr><tr><td>9</td><td>0.2</td><td>1.1</td><td>1.14</td><td>0.747</td><td>0.004</td><td>4.54 × 10<sup>-6</sup></td><td>0.100</td><td>(1.149 ± 0.007) × 10<sup>-1</sup></td><td>15.6</td><td>15.60</td></tr><tr><td>10</td><td>0.8</td><td>0.2</td><td>0.229</td><td>0.677</td><td>0.008</td><td>10<sup>-3</sup></td><td>0.400</td><td>(3.73 ± 0.01) × 10<sup>-1</sup></td><td>2.80</td><td>7.00</td></tr><tr><td>11</td><td>0.8</td><td>0.2</td><td>0.515</td><td>0.301</td><td>0.003</td><td>10<sup>-3</sup></td><td>0.564</td><td>(4.6 ± 0.2) × 10<sup>-1</sup></td><td>12</td><td>22</td></tr><tr><td>12</td><td>0.8</td><td>0.2</td><td>0.997</td><td>0.155</td><td>0.002</td><td>4.54 × 10<sup>-6</sup></td><td>0.400</td><td>(3.19 ± 0.06) × 10<sup>-1</sup></td><td>8.70</td><td>21.8</td></tr><tr><td>13</td><td>0.8</td><td>0.8</td><td>0.658</td><td>0.942</td><td>0.008</td><td>10<sup>-3</sup></td><td>0.200</td><td>(2.32 ± 0.04) × 10<sup>-1</sup></td><td>3.60</td><td>18.0</td></tr><tr><td>14</td><td>0.8</td><td>0.8</td><td>0.844</td><td>0.734</td><td>0.005</td><td>10<sup>-3</sup></td><td>0.317</td><td>(2.650 ± 0.002) × 10<sup>-1</sup></td><td>5.22</td><td>16.50</td></tr><tr><td>15</td><td>0.8</td><td>0.8</td><td>1.13</td><td>0.548</td><td>0.0068</td><td>4.54 × 10<sup>-6</sup></td><td>0.200</td><td>(1.938 ± 0.003) × 10<sup>-1</sup></td><td>0.650</td><td>3.250</td></tr><tr><td>16</td><td>0.8</td><td>1.1</td><td>0.823</td><td>1.04</td><td>0.004</td><td>10<sup>-3</sup></td><td>0.100</td><td>(1.20 ± 0.08) × 10<sup>-1</sup></td><td>2.80</td><td>28.0</td></tr><tr><td>17</td><td>0.8</td><td>1.1</td><td>0.948</td><td>0.899</td><td>0.004</td><td>4.54 × 10<sup>-6</sup></td><td>0.178</td><td>(1.553 ± 0.001) × 10<sup>-1</sup></td><td>2.28</td><td>12.80</td></tr><tr><td>18</td><td>0.8</td><td>1.1</td><td>1.11</td><td>0.768</td><td>0.004</td><td>4.54 × 10<sup>-6</sup></td><td>0.100</td><td>(1.341 ± 0.001) × 10<sup>-1</sup></td><td>3.42</td><td>34.20</td></tr></tbody></table>

<figure><img src="image_8.png" /><figcaption>FIG. 8. (Left) 2D plots of the plasma density using GAMERA for time (from left to right) 1.0, 2.0 and 3.0 s for Case 10. (Right) Kinetic energy density of the velocity perturbation in the y direction, $\epsilon_{y,k}$, associated to the second mode $k = 4\pi/L$ for Case 10.</figcaption></figure>

error with respect to varying $(U, \rho_1, \rho_2, B, P, k_x)$ for the fixed $(A, M_c, M_{Ac}^{-1})$ of case numbers 3 and 18 in Table I. This is done the following way. A given set of $(A, M_c, M_{Ac}^{-1})$ leads to a constant $\gamma_{num}/\gamma_U = \alpha$. We solve for $\gamma_{num} = \alpha\gamma_U$, where $\alpha$ is 0.4 for case 3 and 0.1 for case 18. We find $\gamma_U$ by recasting $U$ in terms of the convective Mach number. This gives

$$ \gamma_U = k_x U = \frac{k_x}{2} M_c a_{01} \left( 1 + \sqrt{\frac{1+A}{1-A}} \right), \quad (31) $$

such that $\gamma_U$ is parameterized by $(A, M_c)$ and the sound speed $a_{01} = \sqrt{\Gamma P_0 / \rho_1}$ (i.e. $P_0$ and $\rho_1$). Fixing $P_0$, we

vary only $\rho_1$ for cases 3 and 18 and find $\gamma_{num}$ from $\gamma_U$ for fixed $(A, M_c, M_{Ac}^{-1})$. Seven simulations for cases 3 and 18 are made with $\rho_1$ varying between 0.438 and 63.0 for case 3 and between 0.0907 and 13.1 for case 18. For each, we plot the relative error $\epsilon_{rel}$ with respect to $\gamma_{num}$ in figure 9. For case 3, $\epsilon_{rel}$ varies from 10.5% to 20.2%, without any special trend. For case 18, the relative error is rather large and constant (≈ 34%).

We believe there are two main reasons for the differences between theory and simulation. First, diminishing the shear layer width generates a numerical instability with
higher modes that increase the error. As mentioned by McNally, Lyra, and Passy$^{62}$, small sound waves can propagate through the 2D periodic boundary conditions creating spurious noise and errors. This noise further amplifies when reducing the shear layer width. Therefore, parasite modes can develop numerically at the thin shear interface due to the coupling between the spatio-temporal resolution and the periodic boundary conditions. These modes can perturb the flow and directly deteriorate the growth of the principal imposed mode. We believe this is the case for case 3 for which the low $M_{Ac}^{-1}$ does not damp enough the growth of small scales, which further interact with the large scales at the interface. This explains why the relative error varies with $\gamma$ for case 3 in Fig. 9 and depends on the set parameter $(U, \rho_1, \rho_2, B, P, k_x)$ although we fixed $(A, M_c, M_{Ac}^{-1})$. One way to limit the scale interaction is to increase the initial perturbation amplitude $v_{0,y}$ (as done in Table I). By giving more energy to the mode $k_x = 4\pi/L$ at the beginning, the high frequencies become negligible. However, nonlinear effects may appear if the condition $v_{0,y}k_x/(2\pi) \ll U/L$ is not fulfilled. The so-called errors simply reflect the departure from a purely linear regime in a vortex sheet. Thus, the difficulty lies in numerically solving the vortex sheet approximation. These scale interactions are the main reasons for the presence of errors for $M_{Ac}^{-1} \lesssim 0.5$, which increases as $M_{Ac}^{-1}$ tends toward zero. It is worth noting that the appearance of secondary KHI can be avoided by adding a dissipation term in the MHD equations. However, this is beyond the scope of the present work.

Second, the increase of $M_{Ac}^{-1}$ also results in potential errors since increasing the magnetic field ($M_{Ac}^{-1} \gtrsim 0.5$) causes the flow to approach the stability limit. This is the case for case 18 for which $M_{Ac}^{-1} = 0.768$. The error does not depend anymore on the parameter set $(U, \rho_1, \rho_2, B, P, k_x)$ for a fixed triplet $(A, M_c, M_{Ac}^{-1})$ in Fig. 9. The error remains high and constant because $\gamma$ is harder to fit close to the stability limit. The high $\gamma$ on the x-axis is compensated with $\gamma_U$ as we fix $\gamma/\gamma_U = 0.100$ such that the $(M_c, M_{Ac}^{-1})$ parameters are close to $(M_c, M_{Ac,lim}^{-1})$ even if $\gamma$ increases.

Besides small numerical differences (close to stabilization or due to the vortex sheet limit), this work shows the universality of a dimensionless growth rate variable $\gamma/\gamma_U$ for a given set of $(A, M_c, M_{Ac})$.

## V. APPLICATION TO LASER EXPERIMENTS

### A. The Hurricane et al. 43's experiment at low-$M_c$

Hurricane et al. 43's experiment is a purely hydrodynamic experiment (unmagnetized configuration) in the high density regime (HEDP). The experiment of laser shots (number 51086, 51090 and 51097) was carried out at the OMEGA laser facility at the Laboratory for Laser Energetics$^{63}$. The target was designed$^{49}$ to develop ablative KHI between two solid materials of different densities ($C_{50}H_{47}I_3$, iodinated polystyrene, at $\rho_1 = 1.43 \text{ g/cm}^3$ and CRF, carbon resor-

cinol foam, at $\rho_2 = 0.1 \text{ g/cm}^3$). The interface between the two materials is initially given by a sinusoidal shape $h(x) = a\sin(x/\lambda)$, where the parameters $a = 30 \mu\text{m}$ and $\lambda = 400 \mu\text{m}$ are chosen to obtain several eddy turn-overs in the duration of the experiment and be large enough to be diagnosable using x-ray radiography. In a 1 nanosecond pulse, the laser deposits $\sim 4$ of kilojoules on one end of the target at the intensity of $10^{15} \text{ W/cm}^2$. This generates a laser blast-wave (i.e. a shock) that then traverses the length of the target. The shock wave then transforms the materials in a plasma flow reaching megabar pressures. Vorticity is created at the plastic interface by baroclinic effect$^{49}$ that generates the well known KH roll-up structures. Despite incompressible KH growth being maximized with $A = 0$, a jump in density, with a gradient perpendicular to the pressure gradient of the shock, is required in order to maximize the baroclinic effect in this experiment design$^{49}$. During the experiment, the plasma flow conditions change in time and space within the post-shock region, between the ablator and the propagating laser blast-wave. Therefore, only the early stage of the experiment can be compared with the KH linear analysis as it becomes nonlinear later in the experiment. The linear regime is observed for instants $10 < t < 40 \text{ ns}$ after the passage of the shock and before the interaction of shock-reflection according to Hurricane et al. $^{43,64}$. Moreover, an estimate of the kinematic viscosity of the plasma in the vicinity of the KH vortex$^{65}$ is $0.0014 \text{ cm}^2/\text{s}$ so that the Reynolds number is of the order of $Re \sim 10^6$, which justifies the use of the Euler equations in the modeling$^{43,64}$. Transition to turbulence and turbulence may develop at later times and break-up of the vortex structures$^{43}$, which is not seek to be captured by the simulations, which are confined to the early times. From now on, we refer to Hurricane et al. 43's experiment as the H09 experiment.

<figure><img src="image_9.png" /><figcaption>FIG. 9. Relative errors for Case 3 and Case 18 versus $\gamma_{num}$ for $P_0 = 1.0$ with varying density $\rho_1$.</figcaption></figure>

777

778

779

780

781

782

783

784

785

786

787

788

789

790

791

792

793

794

795

796

797

798

799

800

801

802

803

804

805

806

807

808

809

810

811

812

813

814

815

816

817

818

820

821

822

823

824

825

826

827

828

829

830

831

832

833

834

835

836

837

838

839

840

841

842

843

844

845

846

847

848

849

850

851

852

853

854

855

856

857

858

859

860

861

862
The flow conditions of the H09 experiment are reported in table II (first row). Since the linear regime begins for $\Delta t = 10$ ns after the passage of the shock, we choose the conditions of the flow at $\Delta t = 12$ ns. This experiment was conducted at the low-$M_c$ values, which vary from $M_c = 0.0$ to $M_c = 1.0$ and at high-A numbers, which range from A = 0.5 to A = 0.95 for the early stage. The rolls that develop at $\Delta t = 12$ ns are estimated at $M_c = 0.312$ and A = 0.9 according to the numerical simulation of the H09 experiment presented in figure 8 of Hurricane et al. $^{43}$.

In Table III, we present the calculated growth rates for H09 with the conditions reported in table II.

## B. Experimental growth rate at $\Delta t = 12$ ns

The H09 experimental growth rate, $\gamma_{exp}$, is defined in the linear regime by $h(t) = h_0 e^{\gamma_{exp} t}$, where $h(t)$ measures the size of the vortex roll at the instant $t$ initiated by a defect of size $h_0$. Since the instability begins only after $\Delta t = 10$ ns, we estimate the growth rate between the mean instants $t_1 = 11.8$ ns and $t_2 = 24.3$ ns with measurement uncertainties of $\delta t = 5.0$ ns. The measurements show that $h_1 = 42.0 \mu$m and $h_2 = 72.2 \mu$m$^{43}$ (figure 9) with negligible uncertainties, yielding a $\gamma_{exp,12} = 4.4 \times 10^{-2}$ ns$^{-1}$. We also calculate the value of the experimental growth rate, $\gamma_{exp,13}$, between $t_1$ and $t_3 = 32.1$ ns and $\gamma_{exp,14}$ then between $t_1$ and $t_4 = 36.8$ ns. Computing the average of the 3 experimental growth rates ($\gamma_{exp,12}$, $\gamma_{exp,13}$, $\gamma_{exp,14}$), we obtain a final average estimate $\langle \gamma_{exp} \rangle = \gamma_{exp} = (5 \pm 1.4) \times 10^{-2}$ ns$^{-1}$. The uncertainty of $\gamma_{exp}$ is defined as $\delta \gamma_{exp,1j} = \gamma_{exp,1j} \delta t / (t_j - t_1)$ for $j \in \{2, 3, 4\}$.

$\gamma_{exp}$ is reported in Table III(first row) and compared with both the theoretical growth rate, $\gamma_{th}$ (second row), and the numerical growth rate, $\gamma_{num}$ (third row). Both computations are described in the following subsections.

## C. Numerical growth rate

The GAMERA simulation is performed by setting the shear layer's thickness, $\sigma$, to be the smallest allowable value on a $1024^2$ grid by using an initial perturbation amplitude $v_{0,y} = 10^{-2}$ for which the linear approximation is satisfied. The numerical growth rate is computed using GAMERA as done in section IV by extracting the slope of the temporal evolution of the $k^{th}$ mode of the kinematic energy $E_{y,k}^r$ (see also Briard et al. $^{58}$).

Although the numerical simulation flow parameters of table II are normalized, they are set to similar convective Mach and Atwood numbers as those used in H09's experiment. In order to compare the numerical growth rate with the experimental one, we use the H09's experimental growth rate of $\gamma_U = k_x U \approx 0.1045$ ns$^{-1}$ at $\Delta t = 12$ ns ($k_x = 2\pi/400 \mu$m$^{-1}$) to convert the numerical dimensionless growth rate, $(\gamma/\gamma_U)_{num}$, to physical units via $\gamma_{num} = \gamma_U \cdot (\gamma/\gamma_U)_{num}$.

## D. Theoretical growth rate at $\Delta t = 12$ ns.

The first calculation of the theoretical growth rate is obtained by exploiting the results of Fig. 5's diagrams for the flow conditions of Table II at $\Delta t = 12$ ns. The result is reported in Table III(second row).

In Table III, we compute the average absolute errors $\langle \epsilon_{abs} \rangle = (\epsilon_{abs}^+ + \epsilon_{abs}^-)/2$, and the average relative errors $\langle \epsilon_{rel} \rangle = (\epsilon_{rel}^+ + \epsilon_{rel}^-)/2$, between either the experimental or numerical growth rates and the theoretical value. The experimental growth rate has a small average relative error of 10% while the mean relative error of the growth rate from the numerical simulation is also quite small, with an error of 5.3%. We note that the 10% mean relative error on the growth rate we estimate for the experiment is essentially due to the measurement uncertainties which are $\pm 1.4 \times 10^{-2}$ ns$^{-1}$.

One could argue that the calculation of the theoretical growth rate is a rough approximation since it does not take into account the time variations of the flow conditions from $\Delta t = 0$ to $\Delta t = 12$ ns. Hurricane et al. $^{43}$ indeed deplored that linear theory supposes a constant velocity and density shear values which may change the growth rate at earlier times. We study the possible changes of conditions and instability growth in the next section and will interestingly confirm our first growth rate estimate is already accurate.

## E. Comparison with simulated growth rates at $t = 25$ ns

Hurricane et al. $^{43}$, figure 8 perform simulations of the H09 flow with the Lawrence Livermore National Laboratory LLNL code C-based arbitrary Lagrangian Eulerian CALE$^{66}$. Their simulations provide the profile of density, velocity, and pressure versus the position $x$ at 25 ns, which is still within in the linear regime (valid till $\approx 40$ ns). We convert these profiles into $M_c$ and A values and plot them in Fig. 10(left). The corresponding dimensionless growth rate $\gamma/\gamma_U$ is calculated theoretically using (20) and plotted versus time. This provides all growth rates behind the shock from $x = 0$ to $x = 1.82$ mm at $t = 25$ ns. We see in Fig. 10(left) that $\gamma/\gamma_U$ starts from $\approx 0.5$ and increases up to $\approx 0.6$ decreasing after to $\approx 0.4$ so that the variations are $\pm 0.1$ around 0.5. This confirms the value of Table III taken at 12 ns after the shock and discussed in the previous section is accurate. The reason why such a large variation of both the Atwood number (from 0.5 to 1.0) and the convective Mach number (from 1.0 to 0.0) does not lead to large variation of the dimensionless growth rates are explained in Fig. 10(center) in which we plot the 2D ($M_c$, A, $\beta = \infty$) stability diagram of $\gamma/\gamma_U$ and report the (A, $M_c$)-path of the H09 simulation. One can see the path follows the 0.5 white level line of the dimensionless theoretical growth rate. The latter remains close to 0.5 although there is a strong decay of the Mach number. At the end of the path, the growth rate is controlled by the Atwood number. A small increase in A from 0.90 to 0.94 at $x = 1.2$ mm and $x = 0.65$ mm causes a significant decrease of $\gamma/\gamma_U$ by 37% (almost independently of the Mach number varying from 0 to 1.0). Overall, this analysis shows the
<table><caption>TABLE II. Plasma flow conditions of the H09' experiment for (Exp.) the experimental point at $\Delta t$ = 12 ns after the passage of the shock, (Th.) KHI linear theory, and (Num.) the GAMERA simulation. Simulation are done for $\sigma = 0.002$ and $v_{0,y} = 10^{-2}$).</caption><thead><tr><th>Case (Δt = 12 ns)</th><th>B</th><th>ΔU</th><th>ρ<sub>1</sub></th><th>M<sub>c</sub></th><th>A</th><th>λ</th><th>P</th></tr></thead><tbody><tr><td>Exp.</td><td>0</td><td>13.32 [μm/ns]</td><td>1.43 [g/cm<sup>3</sup>]</td><td>0.312</td><td>0.9</td><td>400 [μm]</td><td>0.55 [Mbar]</td></tr><tr><td>Th.</td><td>0</td><td>13.32 [μm/ns]</td><td>1.43 [g/cm<sup>3</sup>]</td><td>0.312</td><td>0.9</td><td>400 [μm]</td><td>0.55 [Mbar]</td></tr><tr><td>Num.</td><td>0</td><td>2.16</td><td>1.0</td><td>0.312</td><td>0.9</td><td>0.5</td><td>1.0</td></tr></tbody></table>

growth rate at 25 ns before the shock does not change significantly, also remaining close to the experimental growth rate at $\Delta t$ = 12 ns.

### F. Growth rates estimated at fixed position from $t_{AS} = 0$ to $t_{AS} = 17$ ns

If we now assume that the shock travels at constant speed $U_S$, as in a conventional pressure-driven shock tube, and the flow density, velocity, and temperature profiles are roughly unchanged during the shock propagation, an observer located at a fixed position will see those flow profiles pass by. The shock velocity is $U_S = X_S/t_S \approx 72.8$ µm/ns in the H09 simulations. We choose the observer position to be at the shock front $X_S = 1.82$ mm at $t_S = 25$ ns and again use the profiles of Fig. 10(left). For instance, the duration to see passing the flow conditions at $x_0 = 1.0$ mm from $x_0 = 1.0$ mm to $x = X_S = 1.82$ mm is $t_S + (X_S - x_0)/U_S = 36.3$ ns. For simplicity, we define the time after shock, a shifted time $t_{AS}$, by

$$t_{AS}(x) = \frac{X_S - x}{U_S} \quad \text{with} \quad t_{AS}(x = X_S) = 0. \quad (32)$$

Profiles at $x = 0$ in Fig. 10(left) will be seen by an observer located at $X_S = 1.82$ mm after $t_{AS} = 25$ ns. Thus, we are able to convert the previous growth rates in $x$ in profiles over time from $t = 0$, ($t_{AS} = 25$ ns) to $t = 25$ ns ($t_{AS} = 0$). $t_{AS}$ is reported as a top side axis in Fig. 10(left) from $t_{AS} = 0$ to $t_{AS} = 17$ ns. In Fig. 10(center), the H09 (A, $M_c$)-path of the simulation can be labeled by the time $t_{AS}$ in units of [ns] (red crosses).

In order to compute a global growth rate prior to the nonlinear regime, we define a weighted temporal average by

$$\overline{\mathcal{F}}(t_{AS}) = \frac{1}{t_{AS}} \int_{0}^{t_{AS}} \mathcal{F}(t')dt', \quad (33)$$

where $\mathcal{F}$ is $M_c$, A or $\gamma/\gamma_U$. For $t_{AS} = 12$ ns, we find $\overline{M_c}(12 \text{ ns}) = 0.49$, $\overline{A}(12 \text{ ns}) = 0.88$ leading to $\gamma/\gamma_U = 0.505$. Furthermore, if we calculate directly the average $\overline{\gamma/\gamma_U}$, we find a similar value of $\overline{\gamma/\gamma_U}(12 \text{ ns}) = 0.48$. Averaged values are consistent with the first evaluation at $t = 12$ ns reported in Table III. Large variations in (A, $M_c$)-path have little influence on the dimensionless growth rate, which interestingly remains in a narrow zone varying by $\pm 0.1$ around 0.5. This result is in good agreement with the experimental estimation $\gamma_{exp}$ and the first computation method of $\gamma_{th}$

of Table III. If we now assume the shock follows the self-similar solution of a 1D blast wave and travels at speed $U_S \sim d((E/\rho_2)^{1/3}t^{2/3})/dt$, with E being the energy deposited by the laser per unit area and $\rho_2$ is the lower density$^{67}$, we find $\overline{\gamma/\gamma_U} = 0.55$. Since the shock is slower, it changes the timing in Fig. 10(left and center) from 0 to 23.5 ns (instead of 0 to 17). The ticks of the top axis in Fig. 10(left) would be {0, 16.2, 20.1, 22.9} instead of {0, 5, 10, 15}. The $\gamma/\gamma_U$ stays quasi-constant around 0.55 even longer, during 15 ns instead of 4 ns with constant speed. Thus, varying the shock speed does not change the conclusions significantly.

A fundamental result is that the variations of the growth rate with the flow conditions are driven instead by $\gamma_U$ evaluated for a homogeneous and unmagnetized flow while $\gamma/\gamma_U \simeq 0.505$. This is shown by plotting $\gamma_U/\gamma_{U,X_S}$ in Fig. 10(left), with $\gamma_{U,X_S}$ the maximum value of $\gamma_U$ which is equal to the growth rate evaluated at the shock front, $\gamma_{U,X_S} = 0.268 \text{ ns}^{-1}$ for $k_x = 2\pi/400 \text{ µm}^{-1}$. $\gamma_U/\gamma_{U,X_S}$ varies greatly from 0 to 1, following closely the variation of the convective Mach number. This can also be seen by plugging $M_c = 2U/(a_{01} + a_{02})$ into $\gamma_U = k_x U$ recast in terms of $M_c$ and given by Eq. 31. Therefore, the large variations of KHI are carried by $\gamma_U$ so that $\gamma \simeq 0.5\gamma_U \simeq 0.5M_c f(a_{01}, A)$. An important consequence is the need to dimension the stability diagrams with respect to $\gamma_U$ for experiment design and to directly work with $\gamma$.

### G. Consequences for an experiment with a magnetic field

Magnetic stabilization for H09 is shown in Fig. 10(right) for $\beta$-plasma ranging from $1/\beta = 0.02$ to $1/\beta = 0.405$. The growth rate calculated for $\overline{A}(12 \text{ ns}) = 0.88$ and $\overline{M_c}(12 \text{ ns}) = 0.49$ diminishes for increasing $1/\beta$ and becomes zero for $1/\beta > 0.405$. At a given $\beta$, damping is stronger for low $M_c$, in agreement with the diagrams of Fig. 4. $\gamma/\gamma_U$ rapidly vanishes for small $M_c$.

Since the variations are driven by the incompressible growth rate (for a homogeneous and unmagnetized flow $\gamma_U$), the growth rate is now dimensionned for the condition of Table II and plotted versus $M_c$ at iso-values of $M_{Ac}^{-1}$ in Fig. 11(left). Using the pressure $P_0 = 0.55$ Mbar, it is possible to represent the stability diagram for fixed magnetic field using both the diagram of Fig. 11(left) and the scaling law in Eq. (24). This is plotted in Fig. 11(center). This shows that the dimensionned growth rate has a quite different shape
<table><caption>TABLE III. Calculated growth rates using the plasma flow conditions of Table II. Comparison of growth rates computed by (Th.) KH linear theory, (Exp.) from H09's$^{43}$ experimental point taken at $\Delta t = 12$ ns, and (Num.) GAMERA numerical simulations performed for $\sigma = 0.002$ and $v_{0,y} = 10^{-2}$. We compute the mean absolute error, $\langle \epsilon_{abs} \rangle = |\gamma - \gamma_{th}|/\gamma_U$ and the mean relative error, $\langle \epsilon_{rel} \rangle = \langle \epsilon_{abs} \rangle/(\gamma_{th}/\gamma_U)$.</caption><thead><tr><th>Case(Δt = 12 ns)</th><th>γ [1/ns]</th><th>γ<sub>U</sub> [1/ns]</th><th>(γ ± Δγ)/γ<sub>U</sub></th><th>⟨ε<sub>abs</sub>⟩ [%]</th><th>⟨ε<sub>rel</sub>⟩ [%]</th></tr></thead><tbody><tr><td>Exp.</td><td>(5 ± 1.4) × 10<sup>-2</sup></td><td>0.1045</td><td>0.5 ± 0.13</td><td>4.6</td><td>10</td></tr><tr><td>Th.</td><td>4.75 × 10<sup>-2</sup></td><td>0.1045</td><td>0.454</td><td>-</td><td>-</td></tr><tr><td>Num.</td><td>(4.5 ± 0.1) × 10<sup>-2</sup></td><td>0.1045</td><td>0.43 ± 0.01</td><td>2.4</td><td>5.3</td></tr></tbody></table>

<figure><img src="image_10.png" /><figcaption>FIG. 10. (Left) Conversion of 1D curves of Hurricane et al. $^{43}$, figure 8's (H09) simulations in terms of (red line) $M_c$ and (blue line) A. (black line) The theoretical $\gamma/\gamma_U$ (Eq. (20)) computed for the H09 varying ($M_c$, A) plotted versus (bottom axis) the position x in [mm] and (top axis) the time after shock ($t_{AS}$) in [ns]. (green line) $\gamma_U$ normalized by its maximum found at the shock front, $\gamma_{U,X_s} = 0.268 \text{ ns}^{-1}$. The shock front is set at $t_{AS} = 0$. (Center) 2D ($M_c$, A, $\beta^{-1}$) stability diagrams of KHI showing the color coded $\gamma/\gamma_U$ with white level lines ranging from 0.1 to 1.0. (red line) The path of $\gamma/\gamma_U$ for the first 17 ns of the H09 experiments. The growth rate is estimated using ($M_c$, A) computed from the simulated dynamic flow conditions reported in Hurricane et al. $^{43}$, figure 8 and plotted in 10(left). Annotations of the time after the shock in [ns] are reported in red. (Right) The ($M_c$, A, $\beta^{-1}$) stability diagrams of $\gamma/\gamma_U$ plotted for nonzero $1/\beta$ values ranging from $1/\beta = 0.02$ to $1/\beta = 0.405$.</figcaption></figure>

from the dimensionless one in Fig. 11(left) and is required for experiment design. The experimental growth rate of H09 with its error bars and the growth rate of the GAMERA simulation are reported in red and green respectively. For analysis purpose, we represent the difference between the unmagnetized growth rate and the magnetized one at constant $\tilde{B}$ for the same pressure $P_0$, ($\gamma(0) - \gamma(\tilde{B})$) in Fig. 11(right). Maximal growth rate is reached for $M_c$ close to 1. stabilization from compressible effects occurs for convective Mach above $M_c = \sqrt{2}$ without magnetic field (see theoretical results in section II), a limiting Mach which slightly decreases towards $M_c \sim 1$ as B increases. The experimental, numerical, and theoretical growth rates are all in good agreement on average. If we take the Fig. 11(center), at H09 flow conditions for $\Delta t = 12$ ns, by moving down vertically, the flow is completely stabilized by applying a constant $\tilde{B}$-field greater than 90 T. This is further discussed in the next paragraphs.

Fig. 12 shows the associated dimensionned 2D stability

diagrams for different B-field for all A numbers. At fixed A, the growth rate maxima are reached in a restricted area near $M_c = 1.0$. For $\tilde{B} = 90.0$ T, if we take the minima and the maxima of the convective Mach number for which the flow is unstable (i.e. the frontiers of the dark blue regions), we see a quasi-symmetric distribution of the growth rate on both side of $M_c = 1.0$ for all Atwood numbers and all magnetic field values. We observe that the $M_c = 0.312$ regime of the H09 experiment is completely stabilized for $\tilde{B} > 90.0$ T, as discussed previously.

Bouffetier$^{68}$ presented results of experiments performed at the OMEGA laser facility as well as resistive MHD simulations of the experiments performed with the FLASH code$^{69-71}$ (see Bouffetier$^{68}$ [Figure 5.22]). In these experiments and simulations, the same initial conditions as H09 (Table II) are taken, except for a slight difference in the Atwood number set at A = 0.93 (instead of 0.90 here). FLASH simulations showed that applying an external $\tilde{B} = 40$ T only
1084

1085

1086

1087

1088

1089

1090

1091

1092

1093

partially stabilizes KHI. In these simulations, the external magnetic field reduces the growth rate by 20% compared with an experiment without external magnetic field. When these values are compared with the theoretical growth rates (for times up to 25 ns after the passage of the shock), there are a 6% relative error for the unmagnetized case and a 5% relative error for the magnetized case. Furthermore, we can complement the analysis and establish that the flow would be completely stabilized using about twice this value, $\vec{B} = 90$ T according to the diagram of Figure 11(Center).

1139

1140

1141

1142

1143

1144

1145

1146

1147

1148

1149

1150

1151

1152

1153

1154

1155

1156

1157

1158

1159

1160

1161

1162

1163

1164

1165

1166

1094

### H. $M_c \ge 1$ regime for future NIF experiments

1095

1096

1097

1098

1099

1100

1101

1102

1103

1104

1105

1106

1107

For all Atwood numbers and magnetic field values, stability curves peak at $M_c \sim 1$, the regime of highest instability. The $M_c \approx 1$ regime is a pivoting regime between the low-$M_c$ ($M_c \le 1$) regime and high-$M_c$ ($M_c \ge 1$) regime. Increasing $M_c$ up to $\sim 1$ destabilizes the flow in the $M_c \le 1$-regime. However, an increase of the compressibility in the high-$M_c$ regime from, from $M_c = 1$ to $M_c = \sqrt{2}$, has a stabilizing effect. The growth rate level lines decay sharply with increasing $M_c$. The stabilization by the magnetic field is weaker. Magnetic field level lines become closer to each other as if they were compressed by the shock. In comparison with the low-$M_c$ regimes ($M_c < 1$), higher magnetic fields are needed to reach similar growth rates.

1108

1109

1110

1111

1112

1113

1114

1115

1116

1117

1118

1119

1120

1121

1122

1123

The peaked value reached by $\gamma$ close to $M_c = 1.0$ intensifies and sharpens for increasing A. Thus, in a laser experiment devoted to study KH growth, one would think to increase the A number to increase the growth rate. However, as A increases, the growth rate levels become sharper leading to a narrower range in convective Mach number of the probability to reach the maximal growth rate. Conversely, by decreasing the Atwood number, the range in convective Mach number of the probability to reach the maximal growth rate increases, though for a lower maximal growth rate. There is thus a choice to make between either targeting a maximal growth rate spread over a narrow range of convective Mach numbers around 1 or targeting a lower growth rate but covering a wider convective Mach number range around 1 for which the instability has a higher probability to develop at high rate (but not at the highest rate).

1124

1125

1126

1127

1128

1129

1130

1131

1132

1133

1134

1135

1136

1137

1138

Our study shows that, even if linear stability theory is an approximation, it remains applicable during the first instants of KHI development and can guide us to determine whether the flow is stable or not and to provide an estimation of the growth rate. Thus, these analyses can be used to help design high-A, fully compressible experiments of KHI produced with the NIF laser. Though, this does not remove the complexity of the physics involved in the experiment. As the blast wave moves parallel to the interface between a low-density foam and a high-density plastic, an instability-driven turbulent mixing layer is created in the postshock flow due to surface roughness$^{44}$. In addition, Hurricane 49 specified that there may be self-generation of a magnetic field in the experiment by the Biermann battery effect. Moreover, a magnetic Reynold's number of about $\sim 6$ is estimated for

1167

the plasma generated in the experiment, so that any magnetic field may resistively diffuse out in about 6 eddy-turnover times. The stabilizing effect of the magnetic field-line bending is reduced in the presence of magnetic diffusion. Shiva-moggi$^{72}$ and Chhajlani and Purohit$^{73}$ show an effect coming from resistivity for linear KHI growth in the incompressible regime. In their configuration, resistivity induces an overstable flow for velocities exceeding a threshold, beyond which the ideal MHD flow is stabilized. Resistive effects can alter the spectral behavior of the linear KHI and the growth rate becomes a function of $\eta_m^{1/3}$, where $\eta_m$ represents resistivity. However, as mentioned in section V, the estimation of the linear KHI using ideal MHD seems correct when compared with resistive FLASH simulations$^{68}$. At early times in the linear regime, resistive effects seem to be negligible for the early stage, whereas they become more significant in the nonlinear regime. For instance, in Praturi and Girimaji$^{74}$, numerical simulations of KHI with resistive MHD were performed in which $R_{em} = 250$. They find resistivity activates tearing modes in the nonlinear regime, leading to disruption of KH vortices. Praturi and Girimaji$^{74}$ showed that the resistive instability is most pronounced for weak magnetic fields ($M_A^{-1} < 0.1$), is delayed with increasing $M_A$, and does not develop for strong magnetic fields ($M_A^{-1} > 0.5$). Future works should take into account resistivity effects to clarify its importance for KHI. Both turbulent and resistive effects, not included in our study, constitute other interesting and potentially important aspects to investigate further.

## VI. CONCLUSIONS

1168

1169

1170

1171

1172

1173

1174

1175

1176

1177

1178

1179

1180

1181

1182

1183

1184

1185

1186

1187

1188

1189

1190

1191

1192

In this article, we studied the development of the KH instability for a compressible and nonhomogeneous magnetized flow in the vortex sheet approximation. The magnetic field was set to be constant and parallel to the flow. We derived a new analytical expression of KHI's growth rates for a homogeneous ($A = 0$) magnetized compressible plasmas. A novel aspect is the introduction of a convective Alfvénic Mach number, $M_{Ac} = 2U/(v_{A1} + v_{A2})$, analogous to the convective Mach number introduced by Choudhury$^{25}$ for unmagnetized compressible flows, $M_c = 2U/(a_{01} + a_{02})$, as well as the sum of both these convective Mach numbers, defining a total Mach number as $M_{c,tot}^2 = M_c^2 + M_{Ac}^2$. This expression tends both toward the analytical expression of Chandrasekhar$^{1}$'s growth rate for nonhomogeneous ($A \ne 0$) magnetized flows in the incompressible limit and toward Miles$^{24}$'s expression for homogeneous ($A = 0$) unmagnetized compressible flows (though rewritten with a convective Mach number following Choudhury$^{25}$). For nonhomogeneous, magnetized, compressible flows, the dispersion relation is written in terms of only four dimensionless parameters that completely describe the development of the KHI growth: the Atwood number A, the convective Mach number $M_c$, the convective Alfvénic Mach number $M_{Ac}$, and the angle $\phi$ formed by the wavevector $k$ with the direction of the flow. Since $\phi = 0$ for a parallel magnetic field, the dispersion relation is then solved numerically for all $(A, M_c, M_{Ac})$ in or-
<figure><img src="image_11.png" /></figure>

FIG. 11. (Left) Growth rate $\gamma$[1/ns] versus $M_c$ for varying $M_{Ac}^{-1}$ and for $A = 0.9$. (Left insert) For $M_c \le 0.05$, the Chandrasekhar $^1$ solution is plotted in red dotted curves as the incompressible limit, converging toward the full solution presented up to $M_c = 0.1$. (Center) Growth rate $\gamma$[1/ns] versus $M_c$ for varying $\tilde{B}$[T] and for $A = 0.9$. H09 growth rate is represented as a red cross with a vertical error bar. The numerical point of GAMERA is plotted as a green dot. The case $M_{Ac}^{-1} = 0$ and $B = 0$ are plotted in a red solid curve to better see the differences with Hurricane et al. $^{43}$'s point. (Right) The gap between the unmagnetized and magnetized growth rates ($\gamma(0) - \gamma(\tilde{B})$) versus $M_c$ plotted for different iso-values of $\tilde{B}$. The three diagrams are done for $\rho_1 = 1.43 \text{ g/cm}^3$ and $P_0 = 0.55 \text{ Mbar}$ given by the extracted flow conditions on H09 numerical simulation at $\Delta t = 12 \text{ ns}$ after the shocked has passed$^{43}$.

<figure><img src="image_12.png" /></figure>

FIG. 12. 2D ($M_c$, $A$, $\tilde{B}$) stability diagrams for the dimensionned growth rate $\gamma$[1/ns] with contours of constant $\gamma$. $\tilde{B}$ is computed from $\beta$ using $P_0 = 0.55 \text{ Mbar}$ and $\rho_1 = 1.43 \text{ g/cm}^3$, the conditions of the H09 experimental point at $\Delta t = 12 \text{ ns}$.

<footnote n="1194">1194 der to generate a series of stability diagrams for the all range of parameters (e.g. Figs. 3 and 4). The numerical solution also tends toward Chandrasekhar $^1$'s growth rate for nonhomogeneous ($A \neq 0$) magnetized flows in the incompressible limit.</footnote>

<footnote n="1199">1199 Instability diagrams were presented in the</footnote>

<footnote n="1200">1200 ($\gamma/\gamma_U$, $A$, $M_c$, $M_{Ac}^{-1}$) space parameter. With the vortex sheet approximation, the normalization by $\gamma_U = k_x U$ is a powerful concept, which carries a linear dependence in $M_c$ and simplifies the shape of the diagrams. The linear normalization with respect to the wavenumber $k_x$ and the convective velocity $U$ is specific to the vortex sheet</footnote>
1206

1207

1208

1209

1210

1211

1212

1213

1214

1215

1216

1217

1218

1219

1220

1221

1222

1223

1224

1225

1226

1227

1228

1229

1230

1231

1232

1233

1234

1235

1236

1237

1238

1239

1240

1241

1242

1243

1244

1245

1246

1247

1248

1249

1250

1251

1252

1254

1255

1256

1257

1258

1259

1260

approximation. For shear layers with finite widths, $\gamma_U$ varies with a bell-shape structure with respect to $k_x^{22}$.

We showed that KHI develops differently depending on the low-$M_c$ regime ($M_c < M_{c,max}$) or on the high-$M_c$ ($M_c > M_{c,max}$) regime (with $M_{c,max}$ the convective Mach number for which the dimensionless growth rate is maximum). For $M_c < M_{c,max}$, KHI is mostly stabilized by the magnetic field up to sufficiently high Atwood numbers for which the dimensionless growth rate is greatly attenuated. For high $M_c$, however, three effects can contribute to stabilize the flow: the magnetic field as in the low-$M_c$ regime, the compressibility effect and the A effect. Compressibility diminishes the magnetic field stabilizing contribution of KHI. The magnetic field effect is stronger for $M_c < 1$ than for $M_c > 1$ regions.

We also presented a series of 3 diagrams (varying either $B_{lim}$, $M_{Ac,lim}$ and $\beta_{lim}$) showing which magnetic field strength is required for stabilizing a given flow (Fig. 6).

A total of 18 simulations were carried out with the global GAMERA MHD code for different flow conditions defined from the stability diagrams in terms of $(A, B, M_c)$ triplets selected to cover a wide range of parameters. Numerical simulations confirmed the validity of the solution in the linear regime for sufficiently small initial perturbations. The mean value of $\epsilon_{abs}$ and $\epsilon_{rel}$ computed over the 18 simulations were 5.14% and 16.5% respectively, signifying the GAMERA growth rates were on average close to the theoretical value. Small numerical errors were mainly due to two aspects. The vortex sheet approximation is a hard limit to reach numerically. Small scales developed at the shear layer interface and easily interacted with the large scales. The numerical growth rate thus departed slightly from the purely linear growth rate. Increasing $M_c^{-1}$ smoothed the flow and allowed only the development of the large scales, improving accuracy. However, large $M_c^{-1}$ values led to stabilization, which in turn was source of error since the growth rate was small and harder to numerically extract from the simulations.

In section V, we proposed a method to dimension and analyze KH laser experiments involving nonhomogeneous compressible flows, which we tested on Hurricane et al. 43's unmagnetized experiment and associated simulations. This required us to work with dimensionned parameters and directly with the growth rate, $\gamma$. We showed good agreement between the newly derived theoretical growth rate and the experimental growth rate, with an error between 5% and 10% at early times ($t < 25ns$). For the same conditions, the GAMERA simulated growth rate had an error between 2.5% and 5%. We further explained why we believe there were small temporal and spatial variations of the KHI normalized growth rate during this experiment. We also theoretically estimated that complete magnetic stabilization of this experiment would have been possible for large values of $\tilde{B}$ above 90T.

The KHI growth rate does always increase with $M_c$ increasing for values $M_c < 1$, then peak close to $M_c = 1$ and

1261

1262

1263

1264

1265

1266

1267

1268

1269

1270

1271

1272

1273

1274

1275

1276

1277

1278

1279

1280

1281

vanish for $M_c > 1$ at $M_c = \sqrt{2}$. Diagrams showed how sharp $\gamma$ evolves with respect to $(A, M_c)$ in order to improve the flow parameters and reduce uncertainties. We can thus predict the height of the KH roll-up structures for the first instances of the growth and the precise time of the experiment radiography (itself requiring a second unique intense laser shot).

Even though experimental conditions are not ideal and present some alterations compared with the academical approach of this article, we believe it is possible to estimate the global effect of the magnetic field from the linear stability approach at early times. Our study enables an overall view of KHI for all possible parameters and eases the choice of the magnetic field in experiments for KHI stabilization. Nevertheless, a complete study of the resistive MHD focusing on laser experiments will be necessary to fully understand their effect. These results have guided us for designing compressible high Atwood number magnetized/unmagnetized plasma experiments at NIF (experiment number N240825-001 and N240825-002) and conducting the analysis of NIF experiment number N231113-001, N231113-002, N240825-001 and N240825-002. They will be the subject of a future article.

## ACKNOWLEDGMENTS

1283

The author(s) declare financial support was received for the research, authorship, and/or publication of this article. The NIF facility time for project P-000794 - Calming Kelvin-Helmholtz Instability with pre-imposed B field: bringing Solar Wind and Magnetosphere Physics into the laboratory - was allocated through the NIF Discovery Science Program. AM and VM were supported by the NASA DRIVE Science Center for Geospace Storms (CGS) under award 80NSSC22M0163. This manuscript has been co-authored by Lawrence Livermore National Security, LLC under Contract No. DE-AC52-07NA2 7344 with the US. Department of Energy. The United States Government retains, and the publisher, by accepting the article for publication, acknowledges that the United States Government retains a non-exclusive, paid-up, irrevocable, world-wide license to publish or reproduce the published form of this manuscript, or allow others to do so, for United States Government purposes.

1284

1285

1286

1287

1288

1289

1290

1291

1292

1293

1294

1295

1296

1297

1298

1299

1300

## Appendix A

1302

We write here the dispersion relation $\mathcal{D}$ (22) using the four dimensionless parameters discussed in section II: the Atwood number $A = (\rho_1 - \rho_2)/(\rho_1 + \rho_2)$, the convective Mach number $M_c = 2U/(a_{01} + a_{02})$, the convective Alfvénic Mach number $M_{Ac} = 2U/(v_{A1} + v_{A2})$, and the angle, $\phi$, made between the flow and the wavevector directions. We note $W_{\pm} = \omega/\gamma_U \pm 1$, with $\omega$ the complex angular frequency and $\gamma_U = k_x U$, we have

1303

1304

1305

1306

1307

1308

1309

1310
$$
\mathcal{D}_1(A, M_c, M_{Ac}, W_-) = \left( 1 - \frac{1}{W_-^2} \left( \frac{2}{M_{Ac}} \right)^2 \left( \frac{1}{1 + \sqrt{\frac{1+A}{1-A}}} \right)^2 \right) W_-^2, \quad (A1)
$$

1311

$$
\mathcal{D}_2(A, M_c, M_{Ac}, W_-) = 1 + \left(\frac{M_c}{M_{Ac}}\right)^2 - \frac{1}{W_-^2} \left(\frac{2}{M_{Ac}} \frac{1}{1 + \sqrt{\frac{1+A}{1-A}}}\right)^2, \quad (A2)
$$

1312

$$
\begin{equation}
\begin{split}
\mathcal{D}_3(A, M_c, M_{Ac}, \phi, W_-) = {}& 1 - W_-^2 \left( \frac{M_c}{2} \left( 1 + \sqrt{\frac{1+A}{1-A}} \right) \right)^2 (\cos \phi)^2 \\
& + \left( \frac{M_c}{M_{Ac}} \right)^2 - \frac{1}{W_-^2} \left( \frac{2}{M_{Ac}} \frac{1}{1 + \sqrt{\frac{1+A}{1-A}}} \right)^2 .
\end{split}
\tag{A3}
\end{equation}
$$

1313

The dispersion relation $\mathcal{D}$ can then be rewritten as

$$
\frac{1+A}{1-A} \cdot \frac{\mathcal{D}_1 \sqrt{\mathcal{D}_2}}{\sqrt{\mathcal{D}_3}} (A, M_c, M_{Ac}, \phi, W_-) = - \frac{\mathcal{D}_1 \sqrt{\mathcal{D}_2}}{\sqrt{\mathcal{D}_3}} (-A, M_c, M_{Ac}, \phi, W_+) . \quad (A4)
$$

1314

1315

1316

1317

1318

1319

It can be seen in Eq. (A4) that the solution of the gen-
eral dispersion relation, $\omega = \gamma_U(W_- + 1)$, depends only on
four dimensionless parameters: $(A, M_c, M_{Ac}, \phi)$. The nov-
elty of this expression is the use of a convective Alfvénic
Mach number, similarly to $M_c$, which allows a convenient
simple writing.

1338

1339

1340

1341

1342

1343

1320

Appendix B

1344

1321

1322

1323

1324

1325

1326

The dispersion relation (20) can be solved analytically for
the specific case of a magnetized homogeneous plasma (A =
0, B ≠ 0) with k_z = 0. To do that, it is convenient to deal
with nondimensional variables; accordingly, we define φᵢ ≡
Ωᵢ/kₓU and use the dimensionless parameters {M_c, M_{Ac}} to
write the dispersion relation in the following compact form

1344

1345

1346

1347

1348

1349

1350

1351

1352

1353

$$
\mathcal{D}(\bar{\omega}) \equiv \sum_{m=1}^{2} \frac{(M_{Ac}^2 \phi_m^2 - 1)^{1/2} (1 - M_{c,tot}^2 \phi_m^2)^{1/2}}{(M_c^2 \phi_m^2 - 1)^{1/2}} = 0, \quad (B1)
$$

1327

1328

1329

1330

1331

1332

1333

1334

1335

where $M_{c,tot}^2 = M_c^2 + M_{Ac}^2$. To ensure the amplitude of the perturbation decays away from the interface between the fluids, the signs of the square roots in (B1) must be chosen so that their real parts are positive. This is equivalent to define the arguments of complex numbers between $-\pi$ and $\pi$. Equation (B1) is then rationalized by squaring both sides and rearranging to yield a sextic polynomial equation. After factoring out the trivial solution $\bar{\omega} = 0$, it reduces to a biquadratic equation

1354

$$
M_c^2 M_{c,tot}^2 \bar{\omega}^4 - 2 M_{c,tot}^2 (M_c^2 + 1) \bar{\omega}^2 + M_{c,tot}^2 M_c^2 - 2 M_{c,tot}^2 + 2 = 0, \quad (B2)
$$

where $\bar{\omega} (= \omega/k_x U)$. It is important to note that in writing the dispersion relation (B1) in the form given by (B2), the parent equation has been squared so that not all admissible roots of (B2) are genuine. One should therefore be careful to discard any spurious solution by checking whether it satisfies the original dispersion relation.

We know that the general biquadratic equation

$$
\frac{1}{4} \bar{\omega}^4 + p \bar{\omega}^2 + q = 0 \quad (B3)
$$

possesses four zeros whose characteristic can be determined
from the polynomial coefficients $p$ and $q$. Figure 13(left)
shows the nature and analytical expressions of these roots
for all possible values of parameters $p$ and $q$. As can be seen,
the $(p,q)$-plane is divided into four regions that are bounded
by the $p$-axis and the parabola $q = p^2$. Equation (B3) can
have four genuinely complex zeros ($q > p^2$), four real zeros
($p \le 0 \le q \le p^2$), four purely imaginary zeros ($p > 0$ and
$0 \le q \le p^2$) or two real and two purely imaginary zeros ($q <
0$). With reference to (B2), we have

$$
p = -\frac{M_c^2 + 1}{2M_c^2} \quad \text{and} \quad q = \frac{M_{c,\text{tot}}^2 M_c^2 - 2M_{c,\text{tot}}^2 + 2}{4M_c^2 M_{c,\text{tot}}^2}, \quad (B4)
$$

so that

$$
p^2 - q = \frac{4M_c^2 M_{c,tot}^2 + M_{c,tot}^2 - 2M_c^2}{4M_{c,tot}^2 M_c^4}. \quad (B5)
$$
<figure><img src="image_13.png" /><figcaption>FIG. 13. (Left) The general characteristic of the four roots to the quartic Eq. (B3). (Right) Corresponding solutions for the couple ($M_c, M_{Ac}$).</figcaption></figure>

Equations (B4) define the mapping $\mathcal{M}$ between the $(p,q)$-space (Figs. 13(left)) and the $(M_c, M_{Ac})$-space (Figs. 13(right)). The color shaded areas in the $(p,q)$-plane are mapped onto the same-color shaded areas in the $(M_c, M_{Ac})$-plane. The limit defined by the parabola $q = p^2$ is mapped into the $(M_c, M_{Ac})$-plane onto the curve $\mathcal{C}_1$ in Figure 13(right). The line $q = 0$ is mapped onto the line $\mathcal{C}_5$ in Figure 13(right), which will be shown to be the stability limit. Furthermore, since $M_c$ and $M_{Ac}$ are positive, the domain of $\mathcal{M}$ is restricted to the set of $(p,q)$ that verify $p < -1/2$ and $(4q - 12p - 8p^2 - 5)(4q - 4p - 3) \ge 0$. The latter inequality gives the boundary of the hatched area, plotted as a dashed-line in Figure 13(left). The set of points lying outside the domain corresponding to the hatched area in Fig. 13(left) does not have a counterpart in the $(M_c, M_{Ac})$-plane.

We then check for the possible existence of spurious roots in the unhatched area (Fig. 13(left)). First and foremost, it is readily verified that $\mathcal{D}(\bar{\omega}) = \mathcal{D}(-\bar{\omega})$ and $[\mathcal{D}(\bar{\omega})]^* = \mathcal{D}(\bar{\omega}^*)$ where the star denotes the complex conjugate. It follows that the solutions for $\bar{\omega}$ of $\mathcal{D}(\bar{\omega}) = 0$ are either real or purely imaginary, which implies that only spurious roots are found in the blue-shaded area.

Suppose now that $\bar{\omega}$ is real. If the square roots appearing in (B1) were then real, they would have to be positive to satisfy the requirement that the disturbances remain bounded at large distances from the interface. But then (B1) could not be satisfied. Therefore, for the dispersion relation to have real solutions, some of the square roots in $\mathcal{D}(\bar{\omega})$ must be purely imaginary. We first consider the roots found in the green-shaded area. In this case, the relative signs of the radicands of the square roots in (B1) define five regions in the $(M_c, M_{Ac})$-plane, which are denoted by Roman numbers in Fig. 13(right). The regions are bounded by the lines $\mathcal{C}_2, \mathcal{C}_3$, and $\mathcal{C}_4$ (shown as dot-dashed lines), defined respectively by

the equations

$$M_{Ac} + (2M_{Ac} - 1)M_{c,tot} = 0, \qquad (B6)$$

$$M_{Ac} - (2M_{Ac} - 1)M_{c,tot} = 0, \qquad (B7)$$

$$M_{Ac} - M_c = 0. \qquad (B8)$$

Since $\mathcal{D}(\bar{\omega})$ is even in $\bar{\omega}$, only the roots $\bar{\omega}_{\pm} \equiv \sqrt{-p + \sqrt{q} \pm \sqrt{-p - \sqrt{q}}}$ need be considered. The signs of the radicands in $\mathcal{D}(\bar{\omega}_{\pm})$ for the five regions are reported in Table IV. Consider first $\bar{\omega} = \bar{\omega}_{+}$. We see that, in regions I and II, all the radicands are positive, and hence $\bar{\omega}_{+}$ does not satisfy Eq. (B1). In regions III and IV, only $M_c^2\phi_2^2 - 1$ and $M_{Ac}^2\phi_2^2 - 1$ are negative. The corresponding square roots are thus purely imaginary, but since the imaginary unit appears both in the numerator and denominator of the second term of (B1) they cancel out and $\bar{\omega}_{+}$ is also spurious. Finally, in region $\bigvee$, only $M_c^2\phi_1^2 - 1$ and $1 - M_{c,tot}^2\phi_2^2$ are positive. This implies that the first term of (B1) is now negative, while the second is positive. Region $\bigvee$ therefore allows real roots of the dispersion relation (B1). It is a matter of simple algebra to show that indeed $\bar{\omega}_{+}$ is a genuine solution in this region. By the same token, we can show that $\bar{\omega}_{-}$ is a genuine zero only for $(M_c, M_{Ac})$ in regions $\bigvee$ and $\bigwedge$. It is also straightforward to verify that for the other real zeros in the red-shaded area, namely $\bar{\omega} = \pm(2\sqrt{p^2-q-2p})^{1/2}$, all the radicands in $\mathcal{D}(\bar{\omega})$ are positive and hence are not solutions to (B1). Next, we verify that $\bar{\omega}_3 = i(2\sqrt{p^2-q+2p})^{1/2}$ is solution to (B1). To show that this root satisfies the dispersion relation (B1), we need to prove that one of the terms (say the first) in the sum (B1) is purely imaginary, because the other one is the complex conjugate. It is suitable here to use the polar form of complex numbers and write $M_c^2\phi_1^2 - 1, M_{Ac}^2\phi_1^2 - 1$, and $1 - M_{c,tot}^2\phi_1^2$ as $Z_1 = r_1 \exp i\phi_1, Z_2 = r_2 \exp i\phi_2$, and $Z_3 = r_3 \exp i\phi_3$, respectively. We then have
<table><caption>TABLE IV. Signs of the radicands in the dispersion relation (B1) in the green-shaded area of Fig. 13(right), for $\bar{\omega} = \sqrt{-p+\sqrt{q}} \pm \sqrt{-p-\sqrt{q}}$.</caption><thead><tr><th rowspan="2">Region</th><th colspan="2"><i>M</i><sub>c</sub><sup>2</sup><i>&phi;</i><sub>1</sub><sup>2</sup> &minus; 1</th><th colspan="2"><i>M</i><sub>Ac</sub><sup>2</sup><i>&phi;</i><sub>1</sub><sup>2</sup> &minus; 1</th><th colspan="2"><i>M</i><sub>c,tot</sub><sup>2</sup><i>&phi;</i><sub>1</sub><sup>2</sup> &minus; 1</th><th colspan="2"><i>M</i><sub>c</sub><sup>2</sup><i>&phi;</i><sub>2</sub><sup>2</sup> &minus; 1</th><th colspan="2"><i>M</i><sub>Ac</sub><sup>2</sup><i>&phi;</i><sub>2</sub><sup>2</sup> &minus; 1</th><th colspan="2"><i>M</i><sub>c,tot</sub><sup>2</sup><i>&phi;</i><sub>2</sub><sup>2</sup> &minus; 1</th></tr><tr><th><i>&omega;</i><sub>+</sub></th><th><i>&omega;</i><sub>&minus;</sub></th><th><i>&omega;</i><sub>+</sub></th><th><i>&omega;</i><sub>&minus;</sub></th><th><i>&omega;</i><sub>+</sub></th><th><i>&omega;</i><sub>&minus;</sub></th><th><i>&omega;</i><sub>+</sub></th><th><i>&omega;</i><sub>&minus;</sub></th><th><i>&omega;</i><sub>+</sub></th><th><i>&omega;</i><sub>&minus;</sub></th><th><i>&omega;</i><sub>+</sub></th><th><i>&omega;</i><sub>&minus;</sub></th></tr></thead><tbody><tr><td>I</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td></tr><tr><td>II</td><td>&gt; 0</td><td>&lt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&lt; 0</td><td>&gt; 0</td><td>&lt; 0</td><td>&gt; 0</td><td>&lt; 0</td></tr><tr><td>III</td><td>&gt; 0</td><td>&lt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&gt; 0</td><td>&lt; 0</td></tr><tr><td>IV</td><td>&gt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&lt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&gt; 0</td><td>&lt; 0</td></tr><tr><td>V</td><td>&gt; 0</td><td>&gt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&gt; 0</td><td>&gt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&lt; 0</td><td>&gt; 0</td></tr></tbody></table>

$$ \frac{(M_{Ac}^2 \phi_1^2 - 1)^{1/2} (1 - M_{c,tot}^2 \phi_1^2)^{1/2}}{(M_c^2 \phi_1^2 - 1)^{1/2}} = \left(\frac{r_2 r_3}{r_1}\right)^{1/2} e^{\frac{i}{2}(\phi_2 + \phi_3 - \phi_1)}. \quad (B9) $$

1419

1420

1421

1422

To compute $\varphi_2 + \varphi_3 - \varphi_1$, an analysis of the signs of the real and imaginary parts of $Z_1, Z_2$, and $Z_3$ must be performed in the three regions ($\bigvee\mathrm{I}, \bigvee\mathrm{II}$ and $\bigvee\mathrm{III}$) shown in Fig. 13(right). These regions are bounded by $\mathscr{C}_5$ and the curves

$$ \mathrm{Re}\{M_c^2 \phi_1^2 - 1\} = 0, \quad (\mathrm{B}10) $$

$$ \mathrm{Re}\{M_{Ac}^2 \phi_1^2 - 1\} = 0, \quad (\mathrm{B}11) $$

which are denoted by $\mathscr{C}_6$ and $\mathscr{C}_7$, respectively. $\mathrm{Re}\{\cdot\}$ refers to the real part of complex numbers. For $(M_c, M_{Ac})$ belonging to the region $\bigvee\mathrm{I}$, we can show that $Z_1$ and $Z_2$ lie in the first quadrant of the complex plane, while they are in the second quadrant when $(M_c, M_{Ac})$ is in the region $\bigvee\mathrm{III}$. If $(M_c, M_{Ac})$ is in the region $\bigvee\mathrm{III}$, $Z_1$ and $Z_2$ lie in the second and the first quadrant, respectively. The complex number $Z_3$ is found in the third quadrant for $(M_c, M_{Ac})$ in any of the

1423

1424

1425

1426

three regions. These imply that

$$ \varphi_3 = -\pi + \arctan \left( i \frac{Z_3^* - Z_3}{Z_3^* + Z_3} \right), $$

whereas

$$ \varphi_1 = \begin{cases} \arctan \left( i \frac{Z_1^* - Z_1}{Z_1^* + Z_1} \right), & \text{in } \bigvee\mathrm{I}, \\ \pi + \arctan \left( i \frac{Z_1^* - Z_1}{Z_1^* + Z_1} \right), & \text{in } \bigvee\mathrm{II} \text{ and } \bigvee\mathrm{III}, \end{cases} $$

and

$$ \varphi_2 = \begin{cases} \arctan \left( i \frac{Z_2^* - Z_2}{Z_2^* + Z_2} \right), & \text{in } \bigvee\mathrm{I} \text{ and } \bigvee\mathrm{II}, \\ \pi + \arctan \left( i \frac{Z_2^* - Z_2}{Z_2^* + Z_2} \right), & \text{in } \bigvee\mathrm{III}. \end{cases} $$

By $\arctan(\cdot)$ we denote the principal branch of the inverse tangent function; its values are restricted to the interval $(-\frac{\pi}{2}, \frac{\pi}{2})$.

Therefore,

$$ \varphi_2 + \varphi_3 - \varphi_1 = \begin{cases} -\pi + \sum_{m=0}^{3} (-1)^{m+1} \arctan \left( i \frac{Z_m - Z_m^*}{Z_m + Z_m^*} \right), & \text{in } \bigvee\mathrm{I} \text{ and } \bigvee\mathrm{III}, \\ -2\pi + \sum_{m=0}^{3} (-1)^{m+1} \arctan \left( i \frac{Z_m - Z_m^*}{Z_m + Z_m^*} \right), & \text{in } \bigvee\mathrm{II}. \end{cases} \quad (\mathrm{B}12) $$

1427

1428

To compute the sum of arctangent functions in (B12), we introduce the complex number

$$ \mathscr{P} \equiv \frac{8Z_1^* Z_2 Z_3^*}{(Z_1 + Z_1^*)(Z_2 + Z_2^*)(Z_3 + Z_3^*)}. \quad (\mathrm{B}13) $$

1429

1430

1431

This number has the nice property that its argument is equal to the arctangent sum in (B12). Expanding $\mathscr{P}$ and resolving the result into real and imaginary parts, we obtain

$$ \mathrm{Re}(\mathscr{P}) = \frac{(M_c^2 - M_{Ac}^2)(M_{c,tot}^2 - 2M_{Ac}^2 M_{c,tot}(1 + \bar{\omega}^2) + M_{Ac}^2)}{(M_c^2(1 + \bar{\omega}^2) - 1)(M_{Ac}^2(1 + \bar{\omega}^2) - 1)(1 - M_{c,tot}^2(1 + \bar{\omega}^2))}, \quad (\mathrm{B}14) $$

$$ \mathrm{Im}(\mathscr{P}) = 2i\bar{\omega}M_{Ac}^2 \frac{M_c^2 M_{c,tot}^2 \bar{\omega}^4 - 2M_{c,tot}^2(M_c^2 + 1)\bar{\omega}^2 + M_{c,tot}^2 M_c^2 - 2M_{c,tot}^2 + 2}{(M_c^2(1 + \bar{\omega}^2) - 1)(M_{Ac}^2(1 + \bar{\omega}^2) - 1)(1 - M_{c,tot}^2(1 + \bar{\omega}^2))}. \quad (\mathrm{B}15) $$
We observe that the biquadratic polynomial in the numer- ator of Im($\mathcal{P}$) is (B2). Consequently, Im($\mathcal{P}$) = 0 for $\bar{\omega} = i(2\sqrt{p^2-q+2p})^{1/2}$. Thus, depending on the sign of the real part of $\mathcal{P}$, the sum is equal to either 0 or $\pi$. In (B14), $(M_c^2 - M_{Ac}^2)$ is negative because $M_{Ac} > M_c$ in the red-shaded area. The second factor in the numerator is positive, as can be shown by substituting $\bar{\omega}$ for $i(2\sqrt{p^2-q+2p})^{1/2}$. Hence, the sign of Re($\mathcal{P}$) is determined by the opposite sign of its denominator. As we already know from the above discus- sion, $(1 - M_{c,tot}^2(1 + \bar{\omega}^2))$ is negative in regions $\mathcal{V}\bar{\mathcal{V}}$, and $\mathcal{V}\bar{\mathcal{V}}$ while the product $(M_c^2(1 + \bar{\omega}^2) - 1)(M_{Ac}^2(1 + \bar{\omega}^2) - 1)$ is positive in regions $\mathcal{V}\bar{\mathcal{V}}$ and $\mathcal{V}\bar{\mathcal{V}}$ and negative in $\mathcal{V}\bar{\mathcal{V}}$. Re($\mathcal{P}$) is thus positive in regions $\mathcal{V}\bar{\mathcal{V}}$ and $\mathcal{V}\bar{\mathcal{V}}$ and negative in $\mathcal{V}\bar{\mathcal{V}}$. As a result, by (B12), we find that $\varphi_2 + \varphi_3 - \varphi_1 = -\pi$ for $(M_c, M_{Ac})$ belonging to the three regions $\mathcal{V}\bar{\mathcal{V}}$, $\mathcal{V}\bar{\mathcal{V}}$, and $\mathcal{V}\bar{\mathcal{V}}$. This proves that $Z_1^{-1/2}Z_2^{1/2}Z_3^{1/2}$ is purely imaginary and con- sequently $\bar{\omega} = i(2\sqrt{p^2-q+2p})^{1/2}$ is a genuine solution to the dispersion relation (B1). Finally, the parity of $\mathcal{P}$ in $\bar{\omega}$ implies that $-i(2\sqrt{p^2-q+2p})^{1/2}$ is also a solution to (B1).

<footnote n="1">¹S. Chandrasekhar, *Hydrodynamic and Hydromagnetic Stability* (Oxford University Press, 1961).</footnote>

<footnote n="2">²X. Garbet, C. Fenzi, H. Capes, P. Devynck, and G. Antar, “Kelvin–Helmholtz instabilities in tokamak edge plasmas,” Physics of Plasmas **6**, 3955–3965 (1999).</footnote>

<footnote n="3">³J. R. Myra, D. A. D'Ippolito, D. A. Russell, M. V. Umansky, and D. A. Baver, “Analytical and numerical study of the transverse kelvin–helmholtz instability in tokamak edge plasmas,” Journal of Plasma Physics **82**, 905820210 (2016).</footnote>

<footnote n="4">⁴H. Hasegawa, M. Fujimoto, T.-D. Phan, H. Reme, A. Balogh, M. Dunlop, C. Hashimoto, and R. TanDokoro, “Transport of solar wind into Earth’s magnetosphere through rolled-up Kelvin-Helmholtz vortices,” Nature **430**, 755–758 (2004).</footnote>

<footnote n="5">⁵G. Gnavi, F. Gratton, C. Farrugia, and L. Bilbao, “Supersonic mixing layers: Stability of magnetospheric flanks models,” in *Journal of Physics: Conference Series*, Vol. 166 (IOP Publishing, 2009) p. 012022.</footnote>

<footnote n="6">⁶S. Haaland, A. Runov, A. Artemyev, and V. Angelopoulos, “Characteristics of the flank magnetopause: THEMIS observations,” Journal of Geophysical Research: Space Physics **124**, 3421–3435 (2019).</footnote>

<footnote n="7">⁷B. Walsh, E. Thomas, K.-J. Hwang, J. Baker, J. Ruohoniemi, and J. Bonnell, “Dense plasma and Kelvin-Helmholtz waves at Earth’s dayside magnetopause,” Journal of Geophysical Research: Space Physics **120**, 5560–5573 (2015).</footnote>

<footnote n="8">⁸C. Foullon, E. Verwichte, V. Nakariakov, K. Nykyri, and C. Farrugia, “Magnetic Kelvin-Helmholtz instability at the sun,” The Astrophysical Journal Letters **729**, L8 (2011).</footnote>

<footnote n="9">⁹A. Ershkovich, “Kelvin-Helmholtz instability in type-1 comet tails and associated phenomena,” Space Science Reviews **25**, 3–34 (1980).</footnote>

<footnote n="10">¹⁰X. Li, J. Zhang, S. Yang, Y. Hou, and R. Erdélyi, “Observing Kelvin-Helmholtz instability in solar blowout jet,” Scientific Reports **8**, 8136 (2018).</footnote>

<footnote n="11">¹¹H. Lee and J. Kim, “Two-dimensional Kelvin-Helmholtz instabilities of multi-component fluids,” European Journal of Mechanics-B/Fluids **49**, 77–88 (2015).</footnote>

<footnote n="12">¹²A. Lobanov and J. Zensus, “A cosmic double helix in the archetypical quasar 3c273,” Science **294**, 128–131 (2001).</footnote>

<footnote n="13">¹³M. Faganello and F. Califano, “Magnetized Kelvin-Helmholtz instability: theory and simulations in the Earth’s magnetosphere context,” Journal of Plasma Physics **83**, 535830601 (2017).</footnote>

<footnote n="14">¹⁴D. Southwood, “The hydromagnetic stability of the magnetospheric boundary,” Planetary and Space Science **16**, 587–605 (1968).</footnote>

<footnote n="15">¹⁵Z.-Y. Pu and M. Kivelson, “Kelvin-Helmholtz instability at the magnetopause: Solution for compressible plasmas,” Journal of Geophysical Research: Space Physics **88**, 841–852 (1983).</footnote>

<footnote n="16">¹⁶A. González and J. Gratton, “The Kelvin-Helmholtz instability in a com-</footnote>

<footnote>pressible plasma: the role of the orientation of the magnetic field with respect to the flow,” Journal of plasma physics **51**, 43–60 (1994).</footnote>

<footnote n="17">¹⁷A. González and J. Gratton, “The role of a density jump in the Kelvin–Helmholtz instability of compressible plasma,” Journal of plasma physics **52**, 223–244 (1994).</footnote>

<footnote n="18">¹⁸I. Mann, A. Wright, K. Mills, and V. Nakariakov, “Excitation of magnetospheric waveguide modes by magnetosheath flows,” Journal of Geophysical Research: Space Physics **104**, 333–353 (1999).</footnote>

<footnote n="19">¹⁹K. Mills, A. Wright, and I. Mann, “Kelvin-Helmholtz driven modes of the magnetosphere,” Physics of Plasmas **6**, 4070–4087 (1999).</footnote>

<footnote n="20">²⁰H. Turkakin, R. Rankin, and I. Mann, “Primary and secondary compressible Kelvin-Helmholtz surface wave instabilities on the earth’s magnetopause,” Journal of Geophysical Research: Space Physics **118**, 4161–4175 (2013).</footnote>

<footnote n="21">²¹A. Walker, “The Kelvin-Helmholtz instability in the low-latitude boundary layer,” Planetary and Space Science **29**, 1119–1133 (1981).</footnote>

<footnote n="22">²²A. Miura and P. Pritchett, “Nonlocal stability analysis of the MHD Kelvin-Helmholtz instability in a compressible plasma,” Journal of Geophysical Research: Space Physics **87**, 7431–7444 (1982).</footnote>

<footnote n="23">²³H. Turkakin, R. Rankin, and I. Mann, “Emission of magnetosound from MHD-unstable shear flow boundaries,” Journal of Geophysical Research: Space Physics **121**, 8740–8754 (2016).</footnote>

<footnote n="24">²⁴J. Miles, “On the disturbed motion of a plane vortex sheet,” Journal of Fluid Mechanics **4**, 538–552 (1958).</footnote>

<footnote n="25">²⁵S. Choudhury, “Nonlinear evolution of the Kelvin–Helmholtz instability of supersonic tangential velocity discontinuities,” Journal of Mathematical Analysis and Applications **214**, 561–586 (1997).</footnote>

<footnote n="26">²⁶A. Shimony, W. Wan, S. Klein, C. Kuranz, R. Drake, D. Shvarts, and G. Malamud, “Construction and validation of a statistical model for the nonlinear Kelvin-Helmholtz instability under compressible, multimode conditions,” Physics of Plasmas **25** (2018).</footnote>

<footnote n="27">²⁷R. Ong and N. Roderick, “On the Kelvin-Helmholtz instability of the Earth’s magnetopause,” Planetary and Space Science **20**, 1–10 (1972).</footnote>

<footnote n="28">²⁸T. Jones, J. Gaalaas, D. Ryu, and A. Frank, “The MHD Kelvin-Helmholtz instability. II. The roles of weak and oblique fields in planar flows,” The Astrophysical Journal **482**, 230 (1997).</footnote>

<footnote n="29">²⁹A. Otto and D. Fairfield, “Kelvin-Helmholtz instability at the magnetotail boundary: MHD simulation and comparison with geotail observations,” Journal of Geophysical Research: Space Physics **105**, 21175–21190 (2000).</footnote>

<footnote n="30">³⁰M. Faganello, F. Califano, and F. Pegoraro, “Being on time in magnetic reconnection,” New Journal of Physics **11**, 063008 (2009).</footnote>

<footnote n="31">³¹S. Baek and S.-I. Sohn, “Motion of a current-vortex sheet in the magnetic Kelvin-Helmholtz instability,” Physical Review E **108**, 035107 (2023).</footnote>

<footnote n="32">³²Q.-C. Shi, H.-H. Zhang, Z.-J. Zhao, Z.-H. Chen, and C. Zheng, “Effect of the transverse magnetic field on the Kelvin-Helmholtz instability of the supersonic mixing layer,” Physics of Plasmas **30** (2023).</footnote>

<footnote n="33">³³K. G. Brown and S. R. Choudhury, “Kelvin-Helmholtz instabilities of high-velocity magnetized shear layers with generalized polytrope laws,” Quarterly of Applied Mathematics **58**, 401–423 (2000).</footnote>

<footnote n="34">³⁴V. D. Kuznetsov and N. S. Dzhalilov, “Sixteen-moment approximation for a collisionless space plasma: Waves and instabilities,” Plasma Physics Reports **35**, 962–975 (2009).</footnote>

<footnote n="35">³⁵O. A. Hurricane, D. A. Callahan, D. T. Casey, A. R. Christopherson, A. L. Kritcher, O. L. Landen, S. A. Maclaren, R. Nora, P. K. Patel, J. Ralph, D. Schlossberg, P. T. Springer, C. V. Young, and A. B. Zylstra, “Energy Principles of Scientific Breakeven in an Inertial Fusion Experiment,” Phys. Rev. Lett. **132**, 065103 (2024).</footnote>

<footnote n="36">³⁶D. J. Strozzi, H. Sio, G. B. Zimmerman, J. D. Moody, C. R. Weber, B. Z. Djordjević, C. A. Walsh, B. A. Hammel, B. B. Pollock, A. Povilus, et al., “Design and modeling of indirectly driven magnetized implosions on the NIF,” Physics of Plasmas **31** (2024).</footnote>

<footnote n="37">³⁷M. J.-E. Manuel, B. Khiar, G. Rigon, B. Albertazzi, S. Klein, F. Kroll, F.-E. Brack, T. Michel, P. Mabey, S. Pikuz, et al., “On the study of hydrodynam- ical instabilities in the presence of background magnetic fields in high- energy-density plasmas,” Matter and Radiation at Extremes **6** (2021).</footnote>

<footnote n="38">³⁸J. D. Sadler, S. Green, S. Li, Y. Zhou, K. A. Flippo, and H. Li, “Faster ab- lative Kelvin–Helmholtz instability growth in a magnetic field,” Physics of Plasmas **29** (2022).</footnote>

<footnote n="39">³⁹A. Casner, L. Masse, S. Liberatore, P. Loiseau, P. E. Masson-Laborde,</footnote>

<footnote>1496</footnote>

<footnote>1497</footnote>

<footnote>1498</footnote>

<footnote>1499</footnote>

<footnote>1500</footnote>

<footnote>1501</footnote>

<footnote>1502</footnote>

<footnote>1503</footnote>

<footnote>1504</footnote>

<footnote>1505</footnote>

<footnote>1506</footnote>

<footnote>1507</footnote>

<footnote>1508</footnote>

<footnote>1509</footnote>

<footnote>1510</footnote>

<footnote>1511</footnote>

<footnote>1512</footnote>

<footnote>1513</footnote>

<footnote>1514</footnote>

<footnote>1515</footnote>

<footnote>1516</footnote>

<footnote>1517</footnote>

<footnote>1518</footnote>

<footnote>1519</footnote>

<footnote>1520</footnote>

<footnote>1521</footnote>

<footnote>1522</footnote>

<footnote>1523</footnote>

<footnote>1524</footnote>

<footnote>1525</footnote>

<footnote>1526</footnote>

<footnote>1527</footnote>

<footnote>1528</footnote>

<footnote>1529</footnote>

<footnote>1530</footnote>

<footnote>1531</footnote>

<footnote>1532</footnote>

<footnote>1533</footnote>

<footnote>1534</footnote>

<footnote>1535</footnote>

<footnote>1536</footnote>

<footnote>1537</footnote>

<footnote>1538</footnote>

<footnote>1539</footnote>

<footnote>1540</footnote>

<footnote>1541</footnote>

<footnote>1542</footnote>

<footnote>1543</footnote>

<footnote>1544</footnote>

<footnote>1545</footnote>

<footnote>1546</footnote>

<footnote>1547</footnote>

<footnote>1548</footnote>

<footnote>1549</footnote>

<footnote>1550</footnote>

<footnote>1551</footnote>

<footnote>1552</footnote>

<footnote>1553</footnote>

<footnote>1554</footnote>

<footnote>1555</footnote>

<footnote>1556</footnote>

<footnote>1557</footnote>

<footnote>1558</footnote>

<footnote>1559</footnote>

<footnote>1560</footnote>

<footnote>1561</footnote>

<footnote>1562</footnote>

<footnote>1563</footnote>

<footnote>1564</footnote>

<footnote>1565</footnote>

<footnote>1566</footnote>

<footnote>1567</footnote>

<footnote>1568</footnote>

<footnote>1569</footnote>

<footnote>1570</footnote>

<footnote>1571</footnote>

<footnote>1572</footnote>