Studying particle acceleration from driven magnetic reconnection at the termination shock of a relativistic striped wind using particle-in-cell simulations

A rotating pulsar creates a surrounding pulsar wind nebula (PWN) by steadily releasing an energetic wind into the interior of the expanding shockwave of supernova remnant or interstellar medium. At the termination shock of a PWN, the Poynting-flux-dominated relativistic striped wind is compressed. Magnetic reconnection is driven by the compression and converts magnetic energy into particle kinetic energy and accelerating particles to high energies. We carrying out particle-in-cell (PIC) simulations to study the shock structure as well as the energy conversion and particle acceleration mechanism. By analyzing particle trajectories, we find that many particles are accelerated by Fermi-type mechanism. The maximum energy for electrons and positrons can reach hundreds of TeV.


Introduction
In recent observations, high-energy emissions have been detected in both young pulsar wind nebulae (PWNe) such as the Crab nebula and mid-age PWNe such as the Geminga nebula [1]. The high-energy emissions is likely due to the high-energy electrons scattering off the cosmic microwave background (CMB) photons [2,3]. The high-energy electrons and positrons need to be produced by a particle acceleration in the nebula [4]. How pulsar winds efficiently accelerate electrons and positrons to high energies is a major puzzle and holds the key of understanding the near-earth positron anomaly [1,[4][5][6] as well as gamma rays from the Galactic Center [1,7,8]. In the synchrotron spectra of PWNe, the spectral break between the radio and the X-ray cannot be explained by synchrotron cooling and is believed to be attributed to the particle acceleration at the termination shock [9][10][11]. How the electromagnetic energy is converted into non-thermal particle acceleration at the termination shock is not fully understood.
In the case of an obliquely rotating pulsar, a radially propagating relativistic flow is continuously launched. Asymptoticly near the equatorial plane, such a flow is modeled as a striped wind containing a series of drifting Harries current sheets [12,13] with opposite polarity magnetic fields in between. Numerical simulations including magnetohydrodynamics [11,[14][15][16], particle-in-cell (PIC) [17] and test-particle simulations [18] have been used for modeling the termination shock of PWNe. Magnetic reconnection driven by the termination shock may dissipate the magnetic energy and accelerates e-mail: yclu@lanl.gov particles [17,19]. Particle acceleration in relativistic magnetic reconnection has been a recent topic of strong interests. Two candidate mechanisms of particle acceleration [20][21][22][23][24] are direct acceleration surrounding X-points and Fermi acceleration in flows generated within the reconnection layer. While several analyses [20][21][22][23] have shown that Fermi acceleration dominates particle acceleration to high energies in a spontaneous reconnection with weak guide field, this has not been studied in a shock-driven reconnection.
In this paper, we carry out two-dimensional particle-in-cell (PIC) simulations to model the relativistic striped wind interacting with the termination shock near the equatorial plane of obliquely rotating pulsars. The magnetic reconnection driven by the precursor perturbation from the shock converts the magnetic energy into particle energy and accelerates particles into a power-law energy spectrum. The maximum energy for electrons and positrons can reach hundreds of TeV if the wind has a bulk Lorentz factor γ 0 ≈ 10 6 and upstream magnetization parameter σ 0 = 10.

Simulation setup
We use the PIC code EPOCH2D [25] to study the kinetic processes in the termination shock of a relativistic striped wind. To ensure that the magnetic reconnection is driven by the physical perturbation from the shock instead of numerical instabilities such as numerical Cherenkov instability (NCI), the code is modified by the authors to use a piecewise polynomial force interpolation scheme with time-step dependency [26]. The relativistic striped wind is a steady magnetized electron-positron flow propagating along −x. The spatial profile of the electromagnetic field in the downstream rest frame S is where β 0 is the velocity of the wind normalized by the speed of light c, and λ is the wavelength of the stripes in the wind. The dimensionless parameters δ and α are such that the half thickness of the current sheet is ∆ ≈ λδ/(2π) (actually ∆ = (λδ)/(2π √ 1 − α 2 ) as shown later in this paragraph), and B y averaged over one wavelength is B y λ = B 0 [1 − 2(arccos α)/π]. We use the variable ξ = 1 δ α + cos 2π(x + β 0 ct)/λ for the phase of the electromagnetic field. The location of current sheets is determined by setting ξ = 0, i.e. x + β 0 ct = (λ/(2π)) arccos[−α]. And the location of the transitional field is determined by setting ξ = ±1, i.e. x + β 0 ct = (λ/2π) arccos[−α ± δ]. By subtracting the coordinates we get the half thickness of the current sheet ∆ = (λδ)/(2π The background cold electron/positron plasma in the wind is uniform, with constant density n cold e,p = n c0 /2 and constant temperature kT cold e,p = 0.04m e c 2 . The hot electron/positron plasma balances the magnetic pressure and keeps the steady profile of the electromagnetic field. Besides frame S which is also the downstream rest, there are several reference frames, including the center-of-mass (CM) frame S of the wind, the CM frame S e of the hot electrons, and the CM frame S p of the hot positrons. The space-time coordinate transform between S and S is [27] where γ 0 = 1/ 1 − β 2 0 is the Lorentz factor for the transformation between S and S .
The electromagnetic field in frame S can be derived by the Lorentz transform of electromagnetic field [27] E In S frame, the electric field is zero everywhere, the half thickness of the current sheet is ∆ = γ 0 ∆, and the wavelength is λ = γ 0 λ. We assume S p is moving at speed β hẑ with respect to S . The charge density-current tensor (ρ hp , 0, 0, 0) in S p frame can be transformed into the charge density-current tensor (ρ hp , j hp,x , j hp,y , j hp,z ) in S frame and (ρ hp , j hp,x , j hp,y , j hp,z ) in S frame Assuming the density of the hot electrons/positrons in the current sheet in S frame is n he,hp = n h0 /(2 cosh 2 ξ), where n h0 /n c0 = η is the overdensity relative to the cold particles outside the layer and is set to be η = 3 [13,17,20]. Then in S he,hp frame we have n he,hp = n h0 /(γ 0 γ h )/(2 cosh 2 ξ). The temperature of hot component is kT h = η T m e c 2 . In S frame, the sum magnetic pressure and hot component thermal pressure in x direction is To ensure pressure balance we need B 2 0 /(8πγ 2 0 ) = (n c0 η n η T m e c 2 )/(γ 0 γ h ), thus η T = σ 0 γ h 2η n where the magnetization parameter is σ 0 = B 2 0 4πγ 0 n c0 m e c 2 . The drift velocity β hẑ of the hot particles is determined by keeping the steady profile of electromagnetic field, i.e. in the rest frame of the wind ∇ × B = (4π/c)J is satisfied so that the electric field stays zero. Thus Thus where ω p = 4πn c0 e 2 /(γ 0 m e ) is the plasma frequency of the cold background plasma. The time in our simulation is normalized by 1/ω p , and the spatial coordinates in our simulation are normalized by c/ω p . The boundary at x = 0 is reflecting for particles and conducting for fields. The shock is selfconsistently generated by the interaction between the reflected flow and the incoming flow. The simulation is periodic in y direction. In the run we show in this paper, we have α = 0.1, ∆ = d e , λ = 640d e , γ 0 = 10 4 and σ 0 = 10. The length of the simulation box in y direction is L y = 400d e .

Results
As shown in Figure 1, the shock forms self-consistently, and propagates to the right side of the simulation box, compresses and decelerates the current sheets. The results are consistent with the previous study [17]. The shock converts the magnetic energy into the particle kinetic energy, resulting in an average particle Lorentz factor γ ≈ γ 0 (1 + σ 0 ) in the downstream, which is consistent with the jump condition of ultra-relativistic magnetized shock. The location of the shock is approximately 2), moving at speed approximately c/(2 √ 2). In the upstream, the current Figure 2. Left: particle energy distribution γdN/dγ for all the electrons and positrons in downstream region 0 < x < 354c/ω p at ω p t = 2000. Middle and right: particle trajectory for a particle in the simulation. The middle column is particle Lorentz factor (normalized by γ 0 ) vs.
x − x sh coordinate where the location of the shock jump is approximately at x sh ≈ (ct − 1000c/ω p )/(2 √ 2), and the right column is particle Lorentz factor (normalized by γ 0 ) vs. y coordinate. The red line is the trajectory where the particle is in the upstream, and the red line is the trajectory where the particle is in the downstream. The black dots are for particles at t = 0, and the red dots are for particles at ω p t = 2000.
sheets start to continuously break into a series of magnetic islands separated by X-points. The islands coalesce, grow to larger size, and further grow after passing the shock front.
The energy distribution function at downstream of the shock follows a power law f (ε) ∝ ε −p+1 with p = 1.5 as shown in Figure 2. The trajectory for a typical tracer particle that feels the selfconsistent electromagnetic field is shown in Figure 2. The trajectories of most high-energy tracer particles are similar to the one in Figure 2. Those tracer particles (i) bounce several times in the upstream and then travel into the downstream, or (ii) travel into the downstream and bounce several times in the downstream, or (iii) do (i) followed by (ii). Few particles get bounced between upstream and downstream, i.e after a particle travels into the downstream of the shock it rarely travels back into the upstream. The fact that we get a p = 1.5 power law is consistent with the picture that Fermi mechanism in the reconnection islands is much more efficient than the diffusive shock acceleration.

Conclusions and discussions
While a growing body of research using PIC methods focus on the particle acceleration in the spontaneous relativistic magnetic reconnection in the magnetically dominated regime, in this work we extended the study by setting up a simulation with shock driven magnetic reconnection at the termination shock of relativistic striped wind. The analysis shows that many particles are accelerated by Fermi-type mechanism. It is known that PIC simulations only covers a small region in the termination shock, while in reality there is enormous scale separation between the system size and the dynamical scale such as c/ω p . However, to model the particle acceleration during magnetic reconnection in a macroscopic system, it is important to determine the dominant acceleration mechanism.