Reconnection along a Separator in Shock Turbulence

Numerous structures conducive to magnetic reconnection are frequently observed in the turbulent regions at quasi-parallel shocks. In this work, we use a particle-in-cell simulation to study 3D magnetic reconnection in shock turbulence. We identify and characterize magnetic null points, and focus on reconnection along the separator between them. We identify a reconnection region with strong parallel current, a finite parallel potential, and counterrotating electron flows. Electrons are shown to be accelerated by the parallel electric field before being scattered at the null.


Introduction
Numerous reconnecting current sheets have been observed at and downstream of the transition region at Earth's bow shock (Phan et al. 2018;Wang et al. 2019;Gingell et al. 2020;Stawarz et al. 2022).In the quasi-parallel regions, these are the result of reflected particles exciting ion-ion instabilities (Wang et al. 2019;Bessho et al. 2020;Gingell et al. 2023) leading to the formation of intense current sheets, and further secondary instabilities causing electron-scale structures to form.
Kinetic and hybrid simulations have been a useful tool in studying the formation of reconnection regions and their importance in these transition regions.Previous investigations (Bessho et al. 2022) have demonstrated how electrons are accelerated, with the trapping of electrons in islands being a key mechanism, while Gingell et al. (2023) has studied a wide range of parameters to determine which are favorable for reconnection at quasi-parallel shocks.3D simulations have shown that a wider range of reconnection regimes can be accessed (Ng et al. 2022) as the additional degree of freedom allows weak guide field reconnection and differently oriented reconnection planes.
The simulation studies listed above have focused on 2D or quasi-2D reconnection regions.In these quasi-2D systems, the "out-of-plane" direction parallel to the current is treated as slowly varying, with instabilities allowing the development of 3D structure.Generalizations of magnetic reconnection to 3D systems (Schindler et al. 1988) have linked a finite value of U ∥ = −∫E ∥ ds integrated along a field line to global topology change.3D structures are expected in the solar corona, at the magnetopause, and different regions of the magnetosphere (Li et al. 2021a and references therein).Based on theoretical considerations, there are different forms of reconnection that can take place at magnetic nulls (Priest & Pontin 2009), or along magnetic separators (Parnell et al. 2010a;Stevenson & Parnell 2015)."Slippage" reconnection can also occur in sheared flux tubes (Kuniyoshi et al. 2021).
Though many of the studies of 3D reconnection are in the context of solar physics (e.g., Priest & Pontin 2009;Parnell et al. 2010a;Pontin & Wyper 2015;Cheng et al. 2023) and take place at MHD scales, there have been recent simulations and observations that discuss kinetic scales.The Magnetospheric Multiscale (MMS) mission has observed the evolution of a null in the flank of the magnetopause (Ekawati & Cai 2023), while reconstruction using the magnetic field data of Cluster has allowed the identification of magnetic nulls and separators including events in turbulent regions (Guo et al. 2022).Kinetic simulations have been used to study nulls in different environments, where it is shown that energy dissipation is stronger in regions close to spiral nulls (Olshevsky et al. 2016).
In this paper we study 3D reconnection using a kinetic simulation of a quasi-parallel shock.We focus on a single event where we analyze reconnection along a separator between two magnetic nulls.We compare the properties of the plasma in this region and relate it to prior studies of separator reconnection, and analyze the acceleration of electrons in the reconnection region.We show that parallel acceleration can increase the electron energy by an amount comparable to the electron temperature, but is likely to be overshadowed by the acceleration of electrons trapped in flux ropes.

Simulation Setup
We perform 3D simulations of a quasi-parallel shock using the fully kinetic particle-in-cell (PIC) code VPIC (Bowers et al. 2008a(Bowers et al. , 2008b)).The initial condition consists of a uniform plasma and electromagnetic fields, with q = B B cos , where θ Bn is the angle between the magnetic field and shock normal ( x).The initial plasma moves in the negative x-direction with velocity −V flow .The lower x boundary uses conducting walls for fields and reflecting walls for particles, while plasma is injected and the z-component of the electric field is imposed at the upper x boundary with the initial field and flow values.The y and z boundaries are periodic.electron plasma frequency, Ω ce the electron cyclotron frequency, β the ratio between thermal pressure and magnetic pressure for either species, and M A = V flow /v A the Alfvén Mach number of the injected plasma.As the simulation develops, the shock front propagates from the lower x boundary in the positive x-direction.Unless otherwise mentioned in the text, length scales in the paper are normalized to d e , and velocities to c, and number densities to the initial upstream density.Aside from the smaller θ Bn , the physical conditions are similar to Ng et al. (2022).

Results
An overview of the simulation is shown in Figure 1.In the quasi-parallel shock geometry, the interaction between incident and reflected ions leads to the generation of electromagnetic waves in the foreshock.Consequently, numerous current sheets form in the transition region and downstream of the shock, as shown in Figure 1(a).
In Figure 1(b), a 1D cut of the electron density along x with y = 86, z = 94 is shown.At this time, the shock is undergoing reformation and a new shock front is forming around x = 500.For the analysis of magnetic topology, we focus on a specific ( ) d 120 e 3 volume in the simulation domain, centered at (349,86,94).The x extent of this region is shaded in Figure 1.
The region of interest in this work is the magnetic topology around the null points shown in Figure 2.Here the nulls have been located using the trilinear interpolation method of Haynes & Parnell (2007), and their character has been determined by using the eigenvalues and eigenvectors of the matrix at the null.These are not the only nulls in the volume, but the analysis is confined to this specific region due to how it illustrates reconnection and electron acceleration.
In this figure the two null points are marked by spheres.The upper null (red) is a radial null with the matrix M ij having three real eigenvalues, while the lower null (green) is a spiral null with M ij having one real and two complex eigenvalues.The eigenvalues of the radial and spiral nulls are (0.17, −0.11, −0.066) and (0.042, −0.022+0.019i,−0.022−0.019i),respectively.The spiral structure is clearly visible in the fan plane of the lower null, while the straighter field lines associated with the radial null are clearest just below the upper null, though the changing magnetic field in the vicinity makes them challenging to see.At each null, the spine direction is shown by a white line.Field lines seeded in the vicinity of the upper null are shown in the darker pink, with the spine and fan planes marked by arrows, while the field topology around the lower null is shown by the lighter pink lines.The separator, which connects the two nulls, is shown by the thicker cyan line.
An analysis of the plasma parameters in this region is shown in Figure 3.In this figure, field lines are traced starting from the vicinity of the separator and shaded with different quantities.Particularly relevant to reconnection are the parallel electric field E ∥ , parallel current J ∥ , and the electron vorticity ω.Below the field line panels are J ∥ and U ∥ = −∫E ∥ ds evaluated along the separator.Figure 3(a) shows that there is a region of intense parallel current just below the upper (red) null point.E ∥ in this region is highly nonuniform, which will be better quantified below.Toward the lower null point, the field lines are twisted and are beginning to form a flux rope where both J ∥ and E ∥ reverse sign.
The dynamics between the two nulls may be interpreted in terms of reconnection along a separator (Parnell et al. 2010a;Stevenson & Parnell 2015).As mentioned earlier, there are two regions with intense current density-the vertical region just below the red null toward the right, and the twisted region toward the left.Unlike the systems studied in Parnell et al. (2010a) and Stevenson & Parnell (2015) where there are also multiple regions with an E ∥ signature, these regions have oppositely signed E ∥ .A more detailed view of the electron flows is shown in Figure 4.Here two planes perpendicular to B corresponding to the local maximum and minimum of vorticity along the separator are shown.For each plane, the electron velocity is shown after subtracting the local bulk velocity, and in-plane streamlines are plotted.Using the local coordinate systems as defined in the figure, the streamlines are shaded according to the x-components of the velocity, while the planes are shaded according to the y-components of the velocity.Viewing from above, the upper plane shows a clockwise rotation, while the lower plane shows a counterclockwise rotation.
The parallel electron temperature and parallel electron velocity are also shown in Figures 3(e) and (f).Here, the parallel electron velocity is measured in the electron frame at the red null point.The region where the parallel velocity is large is consistent with the parallel current.On approaching the red null point, there is a reduction in the parallel velocity, with a corresponding increase in the parallel pressure.This is consistent with parallel momentum balance.The parallel temperature increases toward the end of the region with negative J ∥ close to the red null, but is reduced moving along the field lines past the null, suggesting that scattering is taking place.
In generalized magnetic reconnection, a finite value of the quantity U ∥ = −∫E ∥ ds can be used to determine if reconnection has global consequences (Schindler et al. 1988).However, there is an assumption of an external ideal region, which is challenging to identify and may not exist in the turbulent environment.For completeness, we have evaluated this quantity, as shown in Figure 3(g).Large parallel potential variations can be seen in the s < 20 and s > 50 regions, where s is the distance along the field line, corresponding to the strong E ∥ in Figure 3(c).
In prior MHD studies of reconnection between two nulls (Parnell et al. 2010a(Parnell et al. , 2010b)), there are cases with multiple regions of enhanced reconnection along different parts of the separator, though the sign of E ∥ in these studies remains the same.As mentioned before, in Figures 3(c) and (g) there are two main regions where the potential variation is large (s < 20 and s > 50), with the parallel electric field being oppositely signed on average in these regions.We interpret this as two reconnection regions where the field lines are twisting in opposite directions.
In the rest of the paper, we focus on the vertical region with the strong negative parallel current just below the red null point.The average density and magnetic field in this region are 2.9 and 0.87 respectively, giving Bv A = 0.044, while the average E ∥ is −0.003, giving 〈E ∥ 〉/(Bv A ) ≈ −0.06.Interestingly, this is similar to the typical values in 2D and 3D  reconnection sites with slow variation in the direction of the current (Birn et al. 2001;Cassak et al. 2017), though whether this is generally true requires further study.
Figure 3(d) shows the parallel current density and the inplane structure of the magnetic field perpendicular to the separator.Here the field lines show an O-point-like structure, which is a possibility during separator reconnection (Parnell et al. 2010a;Stevenson & Parnell 2015).This has implications for the detection of reconnection events, which we will discuss later.
We now consider electron acceleration along the separator and close to the null, with a focus on the current layer.Figure 5 shows examples of trajectories of self-consistent tracer particles in this region.The shading of the trajectories indicates the simulation time, and the overplotted field lines use data from tΩ ci = 17, where the red and black parts of the trajectories meet.As such, the particles are in the region of strong J ∥ and mostly negative E ∥ during red parts of the trajectory, which takes place during the interval from tΩ ci = 17 to 17.125.The electrons travel toward the null, where they are scattered in different directions.
The work done by the parallel and total electric fields is shown in the lower panels, where W = q∫E • v e dt and W ∥ = q∫E ∥ v ∥ dt, shown by the blue and orange lines respectively.Initially, the acceleration is primarily due to the parallel electric field, after which the trajectories and acceleration mechanisms differ.
During the period between tΩ ci = 17 to tΩ ci = 17.125, the parallel energy gain of the particles is approximately 0.1m e c 2 , which is consistent with the maximum potential difference in Figure 3. Particles do show acceleration before and after this interval, with #1 and #2 showing larger energy increases due to work done by the parallel electric field of approximately 0.25m e c 2 .
Because of the spatial and temporal scales of the system, the system is evolving throughout the period we analyze, while the parallel potential is calculated at a fixed time.Also, the particles do not travel along the separator for the entire duration.We investigate the potential difference for nearby field lines at tΩ ci = 17 in Figure 6 to understand its spatial variation.The field lines are seeded just after the twist in the separator close to the lower null.There are common features, such as the initial potential increase around s = 20, though the magnitude varies with position.At larger s, the field lines diverge, with the yellow line showing the largest potential difference of slightly less than 0.2m e c 2 .This is more consistent with the work done by the parallel electric field for particle #1 in particular, which travels toward the right (of the figure), close to this field line as it passes the null.Particles #2 and #3 initially travel toward the left of Figure 5 when close to the null (around tΩ ci = 17.125 when the trajectories transition from red to blue) where the field lines show a smaller potential difference, with both showing deceleration after being accelerated along the  separator.Finally, it should also be reemphasized that the field structure is drifting in the negative z-and positive x-direction (toward the left of the figure).
To evaluate the importance of the parallel acceleration, we focus on the trajectory of particle #1. Figure 7 shows the fate of this particle as it continues on its journey in the shock.There is further parallel acceleration as it moves past the null, followed by another jump in the work done due to both parallel and perpendicular acceleration at another reconnection site (shown in panel (a) at tΩ ci = 17.5).After this, the particle enters a fluxrope-like structure (panel (b) at tΩ ci = 18), where the work is primarily done by perpendicular electric fields, as can be seen in Figure 5, where the parallel work remains approximately constant after tΩ ci = 17.5.The energy increases at the second reconnection site and in the flux rope are larger than the increase close to the null, but are still of the same order of magnitude.

Discussion and Summary
In this work we analyze 3D PIC simulation results of shock turbulence in the light of prior studies on separator reconnection.There are multiple differences between the configuration studied in this work and previous studies of separator reconnection.Many of these studies are at the MHD scale with a focus on solar flares (e.g., Priest & Pontin 2009;Parnell et al. 2010aParnell et al. , 2010b;;Stevenson & Parnell 2015;Threlfall et al. 2015), and take place in idealized systems, while this work focuses on an event at kinetic scales within the shock transition region.Although the topological structures are the same, they are formed in a turbulent environment rather than from initial conditions constructed for the formation of nulls and separators, which are then driven to reconnect.Nevertheless, it is still instructive to use these studies to help understand our results.
From the electron vorticity shown in Figure 3, the field topology changes are driven by counterrotating plasma flows, similar to what is seen in Parnell et al. (2010aParnell et al. ( , 2010b) ) and Stevenson & Parnell (2015), though the field lines here follow the electron, rather than the ion flow.The presence of multiple actively reconnecting regions is similar to Parnell et al. (2010b), where there are multiple regions with strong E ∥ .The small current and electric field at the null points suggests that reconnection is not taking place at these null points (Priest & Pontin 2009).
With respect to acceleration, Threlfall et al. (2015) used test particles to study electron acceleration along MHD-scale separator current layers, showing that electrons close to the reconnection regions are accelerated with most of the work done by the parallel electric field, similar to what is seen from the first part of the sample trajectories in Figure 5.The gain in energy depends on the parallel potential, which is controlled by the length scale of the reconnection regions and the parallel electric field.In the solar corona, this provides a means to accelerate electrons to high energies due to the large length scales involved (∼10 6 m scale).In the fluctuating magnetosheath, it is difficult to provide an estimate of how large reconnection regions could be, but one expects correlation length scales of ∼10d i close to reconnection regions from statistical studies (Stawarz et al. 2022).For the event in this work, the reconnection region is approximately 5d i long.Because of the numerical parameters used in this simulation such as the reduced m i /m e and ω pe /Ω ce , the energy gain shown in the earlier analysis is larger than would be expected at the bow shock.If we assume E ∥ /(Bv A ) remains the same, and use typical upstream values of n 0 = 5 cm −3 , B 0 = 5 nT for the solar wind at Earth (Russell 2001) and scaling B and n in the reconnection region appropriately, the energy gain would be eE ∥ L ≈ 31 eV, compared to the electron temperature of ≈18 eV using these parameters and b = 2 e .We may also consider more extreme upstream parameters such as those during solar flares (Tsurutani et al. 2006).For example, using B 0 = 50 nT and n 0 = 15 cm −3 , we find an energy gain of approximately 1 keV.Understanding the parameters for which this type of acceleration is significant will require further study.
Another difference between studies at MHD scales and kinetic scales is that the timescale of the electron transit along the separator is comparable to the timescale of the magnetic field evolution.This may contribute to the discrepancy between the potential difference and the energy gain, in addition to the spatial variation discussed earlier.
Other mechanisms of electron acceleration by reconnection in shocks have also been studied.In Bessho et al. (2023), it has been shown that electrons are accelerated at both electron-scale and ion-scale reconnection sites.During interactions with multiple electron-scale reconnection sites, the electrons are Fermi-reflected, but the most energetic electrons are those that are trapped in ion-scale magnetic islands.Both these mechanisms involve the perpendicular electric field, and the only significant parallel acceleration seen in that work is due to the pseudopotential across the shock.In this work the parallel electric field provides an additional mechanism for electron acceleration, with the energy gain comparable to the electron temperature.Trapping in a flux rope still provides a larger energy increase than the parallel acceleration, similar to the previous shock work (Bessho et al. 2022), or current sheet studies (Dahlin et al. 2017;Li et al. 2021b).Further work will be required to determine if the additional acceleration due to trapping in a moving island can still be achieved in 3D.
With respect to observations, this work is relevant to the quasi-parallel regions of Earth's bow shock.Current studies of reconnecting current sheets use various diagnostics to identify candidate events, including |J|, B reversals and electron flow perturbations (Gingell et al. 2020;Stawarz et al. 2022).Further analysis is then used to identify if these signatures are consistent with the reconnection geometry.However, in 3D reconnection, the magnetic topology may not have an X-point like structure in the 2D plane perpendicular to the separator (Parnell et al. 2010a(Parnell et al. , 2010b;;Stevenson & Parnell 2015).Figure 3(d) shows the magnetic field lines in a plane perpendicular to the separator within the reconnection region.Here the structure is more similar to an O-point, which may not be recognized by existing searches for reconnection which generally assume X-point-like structures.Although magnetic nulls have been studied by in situ measurements, their detection depends on accurate reconstruction of the magnetic field (Guo et al. 2022).Even when identifying strong J ∥ , E ∥ and counterrotating flows, it may be challenging to identify this form of reconnection if the null points are sufficiently far from the reconnection regions.
To summarize, we have identified 3D reconnection along a separator in the quasi-parallel shock transition region.We have identified signatures such as J ∥ , E ∥ and counterrotating flows, which show similarities to separator reconnection studied at MHD scales.In the actively reconnecting region, electrons are accelerated by the parallel electric field, with the energy increase comparable to the work done by the reconnection electric field, before being scattered in different directions at the magnetic null.The work done by the parallel electric field is smaller than the later energy increase after trapping in a flux rope, and much smaller than the energy gain for electrons trapped in islands as seen in 2D simulations (Bessho et al. 2023).Finally, we have discussed how these results relate to MHD-scale separator reconnection, and the potential of seeing this form of reconnection in observations such as that from NASA's MMS mission.

Figure 1 .
Figure 1.Overview of the shock simulation at tΩ ci = 17.The upstream region is to the right and plasma is flowing in the negative x-direction.(a) Magnitude of the current density showing multiple current structures.(b) A cut at y = 86, z = 94 of the electron density showing the reforming shock.The shaded area illustrates the x extent of the volume used in the later figures.

Figure 2 .
Figure 2. Topology around the null points in the quasi-parallel shock simulation.The nulls are marked by red and green spheres.Field lines showing the spine and fan for the upper and lower nulls are darker and lighter respectively.The separator is the thicker cyan line.The white lines show the direction of the spines at each null, and are 20d e long.The total length along the separator is 98d e .

Figure 3 .
Figure 3. Field lines in the vicinity of the null points showing (a) J ∥ , (b) ω e,∥ , and (c) E ∥ .(d) Parallel current density and magnetic field lines in a plane perpendicular to the separator (looking down from the right (red) null point) showing an O-point-like structure.(e) Parallel electron temperature.(f) Parallel electron velocity in the electron frame at the right (red) null point.(g) Parallel current and parallel potential along the separator, starting from the left (green) null.(h) |B| and ¢ |E where ¢ = + É E u B s s s

Figure 4 .
Figure 4. Rotating electron flows in planes perpendicular to the separator where the vorticity shows its local maximum and minimum.The color scale for vorticity along the separator is the same as Figure 3(c).Electron velocities are calculated after subtracting the mean flow from ( ) d 6 e 3 regions along the local maxima and minima.Red and yellow lines show local coordinate systems used to plot the in-plane velocity components.Streamlines are colored by the xcomponents; planes are colored by the y-components.

Figure 5 .
Figure 5. (Top) Particles accelerated along the separator with trajectories shaded by simulation time.The distances (s) corresponding to Figure 3 are marked by the boxes.(Bottom) Work done by total and parallel electric fields.Shading corresponds to the colors above representing time.Figure 6. (Top) Field lines in the vicinity of the separator along with the sample electron trajectories.(Bottom) Integral of the parallel electric field along the field lines with colors corresponding to the field lines above.

Figure 6 .
Figure 5. (Top) Particles accelerated along the separator with trajectories shaded by simulation time.The distances (s) corresponding to Figure 3 are marked by the boxes.(Bottom) Work done by total and parallel electric fields.Shading corresponds to the colors above representing time.Figure 6. (Top) Field lines in the vicinity of the separator along with the sample electron trajectories.(Bottom) Integral of the parallel electric field along the field lines with colors corresponding to the field lines above.

Figure 7 .
Figure 7. Extended trajectory of particle #1 showing field topology at (a) tΩ ci = 17.5 with a reconnection region marked by the arrow and (b) tΩ ci = 18 with a flux rope.
The simulation domain is 1500 × 500 × 200 ( ) d e 3 covered by 3000 × 1000 × 400 cells, and is initialized with 150 particles per species per cell.Physical parameters used in the simulation are ω pe /Ω ce = 4, m i /m e = e i, θ Bn = 20°, and M A = 10.Here ω pe is the