© The Authors 2026

1. Introduction

The Kelvin–Helmholtz instability (KHI) can develop when the velocity shear exceeds the Alfvén speed associated with the magnetic-field component parallel to the shear (Chandrasekhar 1961). At Earth’s low-latitude magnetopause, this condition is frequently met during intervals of a strongly northward or southward interplanetary magnetic field, and signatures of surface waves and rolled-up vortices have been widely observed (Hasegawa et al. 2004; Eriksson et al. 2016; Stawarz et al. 2016; Settino et al. 2024). Kelvin–Helmholtz (KH) activity has also been reported in the magnetospheres of other planets, including Mercury (Sundberg et al. 2012; Aizawa et al. 2020), Mars (Ruhunusiri et al. 2016; Koh et al. 2025), Saturn (Ma et al. 2015; Dialynas 2018), and Jupiter (Ranquist et al. 2019; Montgomery et al. 2023), raising the question of how magnetospheric size influences the multi-scale development of the instability.

The KHI is also ubiquitous in heliospheric and astrophysical shear flows. In the heliosheath and heliotail, global magnetohydrodynamic (MHD) simulations have shown that velocity shear can drive KH and related instabilities that contribute to large-scale turbulence and plasma transport (e.g., Kieokaew et al. 2021; Ma et al. 2025). Similar shear-driven processes are expected to operate in magnetized astrophysical jets and other collisionless boundary layers, where the nonlinear evolution of KH vortices can influence magnetic-field restructuring and energy dissipation (Chow et al. 2023a,b). Large-scale studies have highlighted the global impact of shear-driven instabilities but do not resolve the kinetic-scale processes that regulate particle interpenetration within vortices.

After saturation, KH vortices generate strong gradients in the velocity, density, and magnetic field down to ion and electron scales, triggering secondary instabilities that broaden the mixing layer and enhance plasma transport across the solar wind–magnetosphere boundary (Nykyri & Otto 2001; Hasegawa et al. 2004; Nakamura et al. 2017a). Kinetic studies and in situ observations have shown that this transport is mediated by processes such as reconnection, wave–particle interactions, and vortex-induced reconnection (e.g., Nakamura et al. 2011), enabling solar-wind plasma entry even under northward interplanetary magnetic field conditions (Nykyri et al. 2006; Stawarz et al. 2016; Nakamura et al. 2017b, 2020). While MHD simulations capture the global KH dynamics and reconnection-driven transport (Faganello et al. 2012; Ma et al. 2017; Ferro et al. 2024), they do not resolve the kinetic mechanisms governing particle diffusion and heating. Despite these advances, quantitatively characterizing plasma mixing in fully kinetic KHI remains challenging, particularly in distinguishing genuine plasma interpenetration from large-scale advection (Nakamura et al. 2013; Karimabadi et al. 2013).

Several studies have addressed plasma transport and mixing in KH dynamics using both simulations and observations. In numerical models, distinguishing plasma populations originating on opposite sides of the shear layer has been achieved through population tagging and passive tracers, and this has allowed boundary-layer broadening and reconnection-driven transport to be quantified (e.g., Nakamura & Daughton 2014; Ma et al. 2017). In spacecraft observations, where particle origin cannot be directly tracked, mixing parameters have been introduced based on energy-separated magnetospheric and magnetosheath populations (Settino et al. 2022, 2024; Radhakrishnan et al. 2024). These observational diagnostics provide valuable classification tools for KH-driven particle crossing but do not allow for direct measurements of the true particle interpenetration or disentangling of the mixing from heating and spectral evolution.

Observational studies have reported mixed magnetosheath and magnetospheric particle populations within KH vortices and boundary layers together with signatures of reconnection and vortex-induced reconnection (e.g., Hasegawa et al. 2004; Stawarz et al. 2016; Nakamura et al. 2017b; Settino et al. 2022; Radhakrishnan et al. 2024). These measurements provide compelling evidence that KH activity facilitates plasma entry. However, in situ observations are inherently local and do not directly quantify the degree of irreversible cross-field interpenetration over global magnetopause scales. Fully kinetic simulations can therefore help clarify the mechanisms and spatial localization of mixing within KH structures. These simulations provide the possibility of labeling particles according to their initial spatial origin, thus enabling a direct quantification of cross-boundary transport. Label-based mixing diagnostics have been introduced in previous particle-in-cell (PIC) studies to investigate plasma entry during KH evolution. In particular, Nakamura & Daughton (2014) performed 3D PIC simulations and quantified plasma transport using a mixing fraction based on particle origin. In that configuration, secondary interchange-type instabilities developed and contributed to turbulent plasma transport across the boundary layer.

In this study, we examined a different kinetic regime characterized by an in-plane magnetic field and a symmetric-shear configuration. We focused on how mixing develops during the nonlinear evolution of KH vortices under these conditions, and we quantified the mixing spatial localization, species dependence, and temporal evolution, relating it to kinetic reconnection diagnostics. Rather than estimating global transport rates, our goal is to determine where and how genuine particle interpenetration occurs and whether mixing remains localized at kinetic scales.

To this end, we combined particle labeling with a complementary density-based tracer and reconnection diagnostics to quantify the spatial structure of plasma interpenetration during the nonlinear stage of the KHI using fully kinetic PIC simulations of an electron–ion plasma. This framework enables a direct comparison between particle-origin mixing metrics and reconnection activity, allowing us to assess whether and to what extent magnetic reconnection is associated with localized cross-field transport. In contrast to previous studies that primarily inferred plasma transport from large-scale topology changes or global mixing fractions, we focused on the spatial localization and species dependence of genuine particle interpenetration at kinetic scales. By combining particle-origin labeling, a complementary density-based tracer, and reconnection diagnostics, we provide a quantitative assessment of how and where magnetic reconnection activity coincides with enhanced cross-field mixing during the nonlinear evolution of the KHI.

2. Methodology and simulation setup

We performed 2D (PIC) simulations using semi-implicit methods (the implicit moment method and the energy conserving semi-implicit method) implemented in the iPIC3D and ECsim codes (Markidis & Lapenta 2010; Lapenta 2017; Bacchini 2023; Croonen et al. 2024). Details of the numerical methods are given in Appendix A.

We modeled a collisionless magnetized plasma with a uniform per-species density (n₀) and an initial magnetic field composed of an in-plane component along x and a guide field B_z, 0 = 10B_x, 0. All quantities were normalized using ion (mass m_i and charge q_i) reference scales: the ion gyro-frequency Ω_c, i = q_iB_tot, 0/(m_ic), where $B_{\mathrm{tot},0}=\sqrt{B_{x,0}^2+B_{z,0}^2}$ is the total magnetic field magnitude at initialization; the ion plasma frequency ${\omega }_{p,i}=\sqrt{4\pi n_0{q}_{\mathrm{i}}^2/m_{\mathrm{i}}}$; the ion inertial length d_i = c/ω_p, i; and the Alfvén speed $v_{\mathrm{A}}=B_{\mathrm{tot},0}/\sqrt{4\pi n_0{m}_{\mathrm{i}}}$. At initialization, the ion cyclotron-to-plasma frequency ratio is Ω_c, i/ω_p, i = 0.05.

The initial equilibrium configuration was implemented following the approach of Cerri et al. (2013), ensuring a self-consistent shear-velocity profile within a finite-Larmor-radius (FLR) framework. Although this is not an exact kinetic equilibrium, it retains first-order FLR corrections, which significantly reduce the spurious fluctuations that arise from pressure-tensor responses at ion kinetic scales when using fluid-like MHD initial conditions, thereby yielding a more accurate representation of shear layers only a few ion Larmor radii thick (Henri et al. 2013). Ions and electrons were initialized with drifting Maxwellian distributions with an initial drift velocity V_drift = 0.28 c. Our boundary conditions are periodic in all directions. We imposed an initial velocity profile with a double shear:

$$\begin{array}{c}\hfill V_x(y)=V_{\mathrm{drift}}[tanh\left(\frac{y-y_{\mathrm{sh},1}}{\delta }\right)-tanh\left(\frac{y-y_{\mathrm{sh},2}}{\delta }\right)-1],\end{array}$$(1)

where δ = 3d_i represents the initial shear-layer thickness, and the shears are located at y_sh, 1 = L_y/4 and y_sh, 2 = 3L_y/4. The shear layers were chosen to study configurations where vorticity and the guide magnetic field are parallel (at y_sh, 1) or antiparallel (at y_sh, 2). The computational domain was a 2D box of size 150d_i × 400d_i with 2304 × 6144 grid points, giving a spatial resolution of Δx = 0.0651 d_i = 0.5208 d_e = 5.21 λ_D, with λ_D = v_th, e/ω_p, e as the Debye length. The time step was $\mathrm{\Delta }t=0.1{\omega }_{\mathrm{p},\mathrm{i}}^{-1}$. This choice resolved the electron gyration with approximately 17.6 time steps per cyclotron rotation. The ion-to-electron mass ratio was m_i/m_e = 64, with 640 particles per cell. The temperature ratio was T_e/T_i = 0.2, and the plasma-β = n₀kT/(B₀²/(8π)) for ions and electrons was β_i = 0.5 and β_e = 0.1, respectively. The most unstable KH mode had a wavelength of L_x/3, leading to the formation of three vortices in each shear layer as the system entered the nonlinear phase. We identified four distinct stages in the evolution: a linear phase (t = 0–$250{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$) with exponential growth, an early nonlinear phase (t = 250–$450{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$), a merging phase (t = 450–$650{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$) dominated by vortex coalescence, and a late nonlinear phase ($t>650{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$) characterized by increasingly turbulent behavior. In this study, we focused on the nonlinear evolution, including the early, merging, and late nonlinear phases.

We tracked the plasma mixing by assigning each particle a label (ℓ) based on its initial y-position, with ℓ = 0 for particles initialized in y < y_sh, 1 or y > y_sh, 2 and ℓ = 1 for those initialized in y_sh, 1 < y < y_sh, 2. This labeling allowed us to distinguish particles originating from different regions and to track their movement across shear layers, providing a new measure of plasma transport and mixing for both ions and electrons that isolates genuine interpenetration from coherent vortex advection and systematically quantifies its spatial and temporal evolution at kinetic scales.

3. Results: Plasma mixing and magnetic reconnection

In Fig. 1 we show snapshots representative of the evolution of the vortices at the lower shear during the nonlinear stage of the KHI. The quantity shown is the electron mixing fraction (Nakamura & Daughton 2014),

$$\begin{array}{c}\hfill F_{\mathrm{e}}=\frac{(n_{\mathrm{e},\ell =1}-n_{\mathrm{e},\ell =0})}{n_{\mathrm{e}}},\end{array}$$(2)

Fig. 1. Electron mixing fraction (Fe) in the lower shear layer during the nonlinear phase. From left to right, the rolled-up vortices at $t=389{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$, their merging at $t=501{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$, and the late nonlinear stage at $t=946{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$ are shown. White lines trace the in-plane magnetic field. In the right panel, two ejected plasma parcels (P1: red; P2: blue) are highlighted.

where n_{e, ℓ = 1} and n_{e, ℓ = 0} are the electron number densities initially on either side of the shear. Regions with F_e = −1, F_e = 1, and |F_e|< 1 correspond to unmixed plasma on either side and mixed plasma, respectively. The ion mixing fraction and the upper shear layer for both species behave similarly and are not shown.

The left panel of Fig. 1 shows three rolled-up vortices at the early nonlinear stage. As they merge (middle panel), a dynamic mixing layer forms, where plasma parcels can temporarily cross regions before returning to their origin (Fig. 1, middle and left panels). Magnetic tension limits permanent cross-field transport: although reconnection allows plasma parcels to detach from their original field lines, the strong in-plane field pulls them back. Consequently, magnetic tension counteracts plasma transport and forces macroscopic parcels to rebound toward their region of origin, limiting their ability to permanently migrate into the opposite region. In Appendix B we further elaborate on this conjecture, showing how the presence of an in-plane magnetic field influences the development of the KHI and prevents plasma from mixing efficiently.

While the mixing fraction, F_e, provides the qualitative picture described above, it does not allow us to adequately quantify actual mixing because it does not capture the movement of small groups of particles potentially crossing between regions. To quantify interpenetration more accurately, we used a new tracer, $\stackrel{\sim }{n}$, based on the labeling scheme described in Sect. 2. For each particle species, it is defined as the absolute change in particle density on either side of the shear,

$$\begin{array}{c}\hfill \stackrel{\sim }{n}(\mathit{x},t)={\displaystyle \sum _{\ell }}|n_{\ell }(\mathit{x},t)-n_{\ell }(\mathit{x},t=0)|,\end{array}$$(3)

summing over all labeled populations. This new tracer highlights regions where plasma from different initial regions mixes, providing a clearer visualization than the mixing fraction. In essence, this tracer identifies any cell of the domain in which particles with a certain label (ℓ; belonging to a specific side of the shear at t = 0) are present at time t. In regions where particles with a certain ℓ were not initially present, the tracer assumes a value > 0.

In Fig. 2 we compare the spatial distribution of $\stackrel{\sim }{n}/n_0$ around both shear layers at different stages of the KHI for ions and electrons. In the top two rows, we observe that ions mix more efficiently both close to and farther away from the frontier. In the bottom two rows, we observe that electrons mix less efficiently, especially at early times. At later stages, there are localized regions around plasma parcels and the frontier where the electron and ion mixing level increases, reaching values of up to 30%. These peaks are very localized, and in general plasma mixing remains only mildly efficient, with most particles penetrating the opposite plasma region only marginally. For ions, higher mixing levels are observed, particularly within plasma parcels that have fully crossed from one region to the other. In the right panel of Fig. 1, two such parcels (P1 and P2) are highlighted and also shown in Figs. 2h and p, where they exhibit higher internal mixing, particularly for the ions in Fig. 2h. This suggests that while plasma elements are advected coherently by the KH dynamics, intra-parcel diffusion promotes additional mixing. Other regions farther away from the frontier also show elevated mixing (see all panels in Fig. 2 at times $t=723,946{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$).

Fig. 2. Morphological evolution of the lower (y = ysh, 1) and upper (y = ysh, 2) shear layers at t = 445, 501, 723, and $946{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$ during the nonlinear stage of the KHI. Each panel shows a region around each shear that is 100 × 150 di2. The color bar shows $\stackrel{\sim }{n}/n_0$ in percentages, which describes the mixing of the two plasma regions normalized to the initial value (see text). The top two rows show the ion mixing (a–d upper shear, e–h lower shear), and the bottom two rows show the electron mixing (i–l upper shear, m–p lower shear).

To quantify plasma mixing, we computed the volume-averaged $\stackrel{\sim }{n}/n_0$ within the 100 × 150 d_i² regions highlighted in Fig. 2. The top panel of Fig. 3 shows its temporal evolution in both shear layers during the nonlinear phase. Although mixing can locally reach 30%, its strong localization keeps the volume-averaged fraction low. Ions cross the interface more efficiently than electrons, reaching averages of ∼7% compared to ∼3% for electrons. Around $t\simeq 300{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$, when fully rolled-up vortices form, mixing increases for both species, especially ions, then plateaus after $t\simeq 850{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$. As Fig. 2 shows, low average mixing is expected since mixing mostly occurs near the interface and within plasma parcels crossing it coherently. We also show here that mixing is similar for different layers, implying that the anti-alignment between the vorticity and magnetic field plays no evident role in this dynamic.

Fig. 3. Temporal evolution of plasma mixing, current density, and X-points. The bottom and top shear layers (ysh, 1, solid; ysh, 2, dashed) are indicated in all panels. Top panel: percentage of mixed plasma, $\stackrel{\sim }{n}/n_0$, for electrons (blue) and ions (red). Middle panel: maximum absolute out-of-plane current density, |Jz/(n0qic)|,  for electrons (blue) and ions (red). Bottom panel: number of X-points in each shear layer.

To understand the link between magnetic-field dynamics and mixing, we analyzed magnetic-reconnection activity by measuring the out-of-plane current density and X-point formation in the shear layers. Magnetic reconnection has previously been identified as a key process in the nonlinear evolution of KH vortices and suggested as a potential driver of plasma transport (e.g., Hasegawa et al. 2004; Nakamura et al. 2011). Here, we directly quantified its association with particle interpenetration using our mixing tracer. The middle panel of Fig. 3 shows the temporal evolution of the maximum out-of-plane current density, J_z, for each species and shear layer. In our 2D geometry, the reversal of the in-plane magnetic field is sustained by a J_z sharply localized at X-points, making it an effective proxy for magnetic reconnection. We observe that J_z evolves similarly in both shears, with electron currents generally larger in magnitude, as expected from the mass difference between species. For $t\gtrsim 200{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$, the current density begins to increase in both species, reaching saturation around $t\simeq 400{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$. This coincides with the onset of the early nonlinear KHI and the beginning of the vortex-merging stage (left and middle panels of Fig. 1) as well as the previously noted rise in the mixing percentage. The bottom panel of Fig. 3 shows the number of X-points over time in each shear layer. The X- and O-points are identified as nulls of the magnetic-flux gradient, ∇ψ = 0, and distinguished by the sign of the Hessian determinant, following standard 2D reconnection diagnostics (Servidio et al. 2009; Wan et al. 2013; Olshevsky et al. 2016). The X-point count increases from $t\gtrsim 300{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$, which corresponds to the early nonlinear stage of the KHI. These diagnostics indicate that enhanced reconnection activity – higher max|J_z| and more X-points – coincides with increased plasma mixing ($\stackrel{\sim }{n}/n_0$). As reconnection breaks and reconnects field lines, particles from opposite sides of the shear are transported along newly connected lines, promoting cross-field mixing in frontier regions and detached plasma parcels. The spatial correspondence between X-point clusters and regions of elevated mixing, together with their synchronized temporal evolution, provides strong evidence that magnetic reconnection facilitates the localized mixing observed during the nonlinear development of the KHI.

4. Conclusions

We have investigated the efficiency of plasma mixing driven by the nonlinear development of the magnetized KHI using high-resolution 2D fully kinetic simulations with FLR effects (Cerri et al. 2013) and diagnostics, including particle labeling, a new mixing tracer, and reconnection proxies. Unlike previously introduced observational mixing parameters (Settino et al. 2022, 2024), which rely on the energy-space separation of particle populations, our diagnostics are based on direct origin tracking in a fully kinetic framework. This enabled a quantitative assessment of the true particle interpenetration and its localization at electron scales, providing a complementary perspective to spacecraft-based mixing measures. Our results show that in this setup, mixing between magnetosheath and magnetospheric plasmas is present but relatively inefficient, even after KH vortices form and the nonlinear stage develops. Mixing is highly localized, reaching 30% only in sparse locations within a narrow interface layer and isolated advected parcels, and it is generally absent farther from the shear layers. On average, ions achieve somewhat higher mixing levels than electrons, and electron mixing remains weak and largely confined to localized structures. We can explain this behavior by considering that ions are significantly less magnetized than electrons, a characteristic that allows them to partially cross the shear interface, whereas electrons remain tightly bound to field lines. With a reduced mass ratio (m_i/m_e = 64), the actual electron magnetization in physical plasmas would be even higher, and their mixing correspondingly lower, meaning that electron mixing here likely represents an upper bound. Enhanced mixing occurs mainly within plasma parcels that detach and reattach across the interface, suggesting a significant contribution from magnetic reconnection and individual particle dynamics. Regions located between the interface and these parcels exhibited an intermittent increase in mixing coinciding with enhanced magnetic reconnection signatures, such as an elevated out-of-plane current and higher X-point counts. We find that enhanced mixing preferentially occurs in regions where reconnection signatures are present, suggesting that magnetic reconnection facilitates localized cross-field particle interpenetration during the nonlinear development of the KHI.

While these results provide new insights into localized plasma transport during the nonlinear KHI, several limitations should be noted. The simulations are 2D and use a reduced ion-to-electron mass ratio, which may overestimate electron mobility and restrict the range of fully 3D reconnection dynamics. In addition, the present simulations adopt a symmetric configuration across the shear layer, whereas the Earth’s magnetopause is typically asymmetric in density, temperature, and magnetic-field strength. Asymmetric configurations can trigger additional interchange-type instabilities and enhance turbulent transport (Nakamura & Daughton 2014; Dargent et al. 2019). The symmetric setup considered here isolates the roles of velocity shear, magnetic geometry, and magnetic reconnection in regulating cross-field transport. Our results should therefore be interpreted as characterizing plasma mixing in a magnetized, symmetry-controlled regime and providing a reference case against which more realistic asymmetric configurations can be assessed. While strongly asymmetric conditions can enhance turbulent transport and increase the net amount of cross-boundary plasma exchange, the mechanisms identified here – localized mixing near reconnection sites and within advected plasma parcels – are expected to remain key ingredients of kinetic KHI evolution in asymmetric configurations. The finite horizontal extent of the simulation domain limits the ultimate growth of individual vortices due to periodic boundary conditions. In the reference configuration, vortex merging leads to structures approaching the box size in x, which prevents further large-scale horizontal expansion. However, we ran additional simulations with larger domains, confirming that increasing the number of vortices does not qualitatively alter the localized nature of mixing at kinetic scales. Although longer evolution times could increase the cumulative transport, the strong in-plane magnetic field constrains the vertical broadening of the mixing layer, and interpenetration remains confined to narrow regions near reconnection sites and detached plasma parcels. These results should therefore be interpreted as characterizing the kinetic regulation of mixing within KH vortices rather than the long-term global evolution of the magnetopause boundary layer.

Despite these limitations, this study provides clear, fully kinetic evidence that plasma mixing driven by the nonlinear KHI can remain intrinsically localized and strongly influenced by in-plane magnetic-field components. By introducing a complementary mixing tracer and combining it with particle labeling and reconnection diagnostics, we identified the specific structures – detached plasma parcels and reconnection sites – where cross-field transport is enhanced. Our results demonstrate that enhanced mixing preferentially coincides with magnetic-reconnection signatures and thus establish a quantitative connection between kinetic reconnection activity and localized plasma interpenetration in nonlinear KH dynamics.

In future work, we will extend these studies to fully 3D setups and explore realistic mass ratios. We will also conduct a more detailed investigation of reconnection physics at the particle level, which could further clarify its role in facilitating cross-field transport.

Acknowledgments

S.F. is supported by the FWO PhD fellowship “Investigating Magnetospheric Plasma Dynamics with Large–scale Fully Kinetic Simulations” (grant no. 1126325N). F.B. acknowledges support from the FED-tWIN programme (profile Prf-2020-004, project “ENERGY”) issued by BELSPO, and from the FWO Junior Research Project G020224N granted by the Research Foundation – Flanders (FWO). The resources and services used in this work were provided by the VSC (Flemish Supercomputer Center), funded by the Research Foundation – Flanders (FWO) and the Flemish Government.

References

Appendix A: Semi-implicit PIC methods

The PIC codes typically use explicit time integration, which requires resolving the Debye length and the electron plasma frequency and is constrained by the Courant–Friedrichs–Lewy condition. Semi-implicit methods relax these constraints by linearizing or exactly treating the coupling between particles and fields, allowing simulations to retain small temporal and spatial scales, including electron-scale dynamics, while employing large domains. In this work, the first part of the simulation uses the implicit moment method (Brackbill & Forslund 1982), which approximates the particle–field coupling and provides enhanced stability compared to explicit methods. The subsequent stage employs the energy conserving semi-implicit method (ECSIM; Lapenta 2017). The ECSIM replaces the nonlinear particle–field coupling with an exact relation, conserves energy exactly, retains the stability of fully implicit methods, and allows larger time steps and grid sizes while accurately capturing electron-scale dynamics.

The results shown in this work correspond to the nonlinear stage of the KHI. We started the simulation with the implicit moment method until the KHI develops and then switched to ECSIM at the early stages of the nonlinear phase so that all results presented here are obtained using the ECSIM method.

Appendix B: Effect of the in-plane magnetic field on plasma mixing

To further support the interpretation proposed in the main text, we compare in Fig. B.1 the nonlinear evolution of the lower shear layer in simulations performed with and without an in-plane magnetic field B_x. These simulations employ the same physical and numerical setup as the reference run discussed in Sect. 2, but with a spatial resolution reduced by a factor of four. Despite the reduced resolution, the simulations capture the large-scale nonlinear dynamics of the KHI and enable a clear assessment of the role of the in-plane magnetic field in regulating plasma mixing.

Fig. B.1. Panels (a) and (c): Percentage of mixed plasma, $\stackrel{\sim }{n}/n_0$, for ions and electrons in the lower shear layer with a finite Bx at $t=334{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$. Panels (b) and (d): Corresponding ion and electron distributions in the lower shear layer without Bx at the same time. Each panel shows a region around the shear layer of size 100 × 150 di2. A power-law normalization is used to enhance the visibility of low-amplitude structures.

Figures B.1a and B.1c show the spatial distribution of the percentage of mixed plasma, $\stackrel{\sim }{n}/n_0$, for ions and electrons, respectively, in the presence of a finite B_x. In this case, as KH vortices develop, plasma mixing remains weak and is largely confined to the narrow interface between the two plasma regions. Although vortical structures form and advect plasma across the shear layer, the in-plane magnetic field limits cross-field transport, resulting in well-defined vortices and localized mixing (up to 30%) primarily concentrated near the frontier, in line with the discussion presented in the main text. In contrast, Figs. B.1b and B.1d show the corresponding ion and electron distributions, always for the lower shear, this time in the absence of an in-plane magnetic field. In this case, the nonlinear evolution of the KHI is accompanied by significantly enhanced plasma mixing (up to 50%), with particles transported far from the initial interface and mixing regions extending well into both plasma domains. The absence of magnetic tension allows vortical motions to more efficiently redistribute plasma across the shear layer, leading to a broader and more spatially extended mixing region.

This comparison indicates that the limited and localized mixing observed in the main simulations is a robust physical effect associated with the presence of the in-plane magnetic field, rather than a consequence of numerical resolution. The in-plane field thus plays a key role in constraining cross-field transport and limiting the efficiency of plasma mixing in the nonlinear phase of the instability.

All Figures

Fig. 1. Electron mixing fraction (Fe) in the lower shear layer during the nonlinear phase. From left to right, the rolled-up vortices at $t=389{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$, their merging at $t=501{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$, and the late nonlinear stage at $t=946{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$ are shown. White lines trace the in-plane magnetic field. In the right panel, two ejected plasma parcels (P1: red; P2: blue) are highlighted.

In the text

Fig. 2. Morphological evolution of the lower (y = ysh, 1) and upper (y = ysh, 2) shear layers at t = 445, 501, 723, and $946{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$ during the nonlinear stage of the KHI. Each panel shows a region around each shear that is 100 × 150 di2. The color bar shows $\stackrel{\sim }{n}/n_0$ in percentages, which describes the mixing of the two plasma regions normalized to the initial value (see text). The top two rows show the ion mixing (a–d upper shear, e–h lower shear), and the bottom two rows show the electron mixing (i–l upper shear, m–p lower shear).

In the text

Fig. 3. Temporal evolution of plasma mixing, current density, and X-points. The bottom and top shear layers (ysh, 1, solid; ysh, 2, dashed) are indicated in all panels. Top panel: percentage of mixed plasma, $\stackrel{\sim }{n}/n_0$, for electrons (blue) and ions (red). Middle panel: maximum absolute out-of-plane current density, |Jz/(n0qic)|,  for electrons (blue) and ions (red). Bottom panel: number of X-points in each shear layer.

In the text

Fig. B.1. Panels (a) and (c): Percentage of mixed plasma, $\stackrel{\sim }{n}/n_0$, for ions and electrons in the lower shear layer with a finite Bx at $t=334{\mathrm{\Omega }}_{\mathrm{c},\mathrm{i}}^{-1}$. Panels (b) and (d): Corresponding ion and electron distributions in the lower shear layer without Bx at the same time. Each panel shows a region around the shear layer of size 100 × 150 di2. A power-law normalization is used to enhance the visibility of low-amplitude structures.

In the text