The heating of electrons in collisionless magnetic reconnection is explored in particle-in-cell simulations with non-zero guide fields so that electrons remain magnetized. In this regime, electric fields parallel to B accelerate particles directly, while those perpendicular to B do so through gradient-*B* and curvature drifts. The curvature drift drives parallel heating through Fermi reflection, while the gradient *B* drift changes the perpendicular energy through betatron acceleration. We present simulations in which we evaluate each of these mechanisms in space and time in order to quantify their role in electron heating. For a case with a small guide field (20% of the magnitude of the reconnecting component), the curvature drift is the dominant source of electron heating. However, for a larger guide field (equal to the magnitude of the reconnecting component) electron acceleration by the curvature drift is comparable to that of the parallel electric field. In both cases, the heating by the gradient *B* drift is negligible in magnitude. It produces net cooling because the conservation of the magnetic moment and the drop of *B* during reconnection produce a decrease in the perpendicular electron energy. Heating by the curvature drift dominates in the outflow exhausts where bent field lines expand to relax their tension and is therefore distributed over a large area. In contrast, the parallel electric field is localized near X-lines. This suggests that acceleration by parallel electric fields may play a smaller role in large systems where the X-line occupies a vanishing fraction of the system. The curvature drift and the parallel electric field dominate the dynamics and drive parallel heating. A consequence is that the electron energy spectrum becomes extremely anisotropic at late time, which has important implications for quantifying the limits of electron acceleration due to synchrotron emission. An upper limit on electron energy gain that is substantially higher than earlier estimates is obtained by balancing reconnection drive with radiative loss.

## I. INTRODUCTION

Magnetic reconnection is a ubiquitous plasma process that converts magnetic energy into thermal and kinetic energy. Of particular interest is the production of non-thermal particles, in which a fraction of the plasma population is driven to energies much larger than that found in the ambient medium. Reconnection is thought to be an important driver of such particles in phenomena such as gamma ray bursts,^{13,28} stellar and solar flares,^{25} and magnetospheric storms.^{30} Recent observations of solar flares reveal the remarkable efficiency of electron acceleration: a large fraction of the electrons in the flaring region become a part of the nonthermal spectrum, with a resulting energy content comparable to that of the magnetic field.^{24,31}

Mechanisms for particle acceleration have been explored and compared in a variety of papers, e.g., Refs. 7, 14, 21, 32, 35, and 39. Previous authors^{7,15,26} examined acceleration by electric fields parallel to the local magnetic field. Drake *et al*.^{10} proposed a mechanism whereby charged particles gain energy as they reflect from the ends of contracting magnetic islands, a process analogous to first-order Fermi acceleration of cosmic rays. Similar energy gain takes place during the merging of magnetic islands.^{6,32} In reconnection in either collisionless^{4,12} or weakly collisional systems,^{1,2,22,27} current sheets break up into many magnetic islands, producing an ideal environment for this mechanism. Although such work has clearly shown that non-thermal tails can be produced in simulations of magnetic reconnection, which process or processes are the dominant drivers of these tails remains an open question.

Resistive MHD studies showed that thermal test particles injected into simulations quickly developed non-thermal features, most notably power-law tails, whose energy content was large.^{5,33} So, although kinetic simulations incur a much greater computational expense, they are necessary for studies of particle acceleration where the feedback of particles on the accelerating field is important. The current largest kinetic simulations of reconnection have domain sizes of several hundred ion inertial lengths ($di=c/\omega pi$, where $\omega pi=4\pi n0e2/mi$ and *n*_{0} is the density) in two dimensions. Three-dimensional simulations are more computationally intensive and thus smaller. By comparison, a typical solar coronal density of $109/cm3$ implies $di\u223c103$ cm, in a system with a typical scale of 10^{9 }cm. A model capable of explaining observations such as that of Krucker *et al*.^{24} will involve a significant extrapolation to physically relevant domains.

Astrophysical magnetic reconnection occurs in a wide variety of regimes. The “guide field” (magnetic field component perpendicular to the reconnecting components and parallel to the reconnection electric field) is an important parameter: anti-parallel (small guide field) reconnection is common in Earth's magnetotail, while component reconnection is common in environments such as the solar corona, where the guide field may be many times the magnitude of the reconnecting field. Anti-parallel reconnection necessarily contains regions where $|B|$ is small and particles become unmagnetized (i.e., Larmor radii exceed local scale lengths). Electrons remain magnetized during guide field reconnection even for very small values of the ratio of the guide field to the reconnecting field.^{36}

Electrons and ions, having disparate masses, can be accelerated by distinct mechanisms. Small-scale processes important for electrons may be washed out at the relatively large ion length scales. Acceleration arising from discontinuities may have a greater impact on massive particles. For example, Refs. 23 and 8 explore how heavy ions are accelerated by a pick-up process in reconnection outflows, a mechanism which does not impact electrons. We consider only electron acceleration in this paper.

Fully three-dimensional treatments are an ultimate goal in the kinetic treatment of particle acceleration. In two dimensions, magnetic islands (closed loops of magnetic flux) are efficient particle traps. However, oblique tearing modes^{4} or velocity shear instabilities^{16} in 3D may generate flux tubes. These three-dimensional equivalents of islands are porous, allowing particles to escape. 2D models dependent on the structure of the magnetic fields (e.g., the contracting islands model) may require modification in order to address acceleration in the more complex geometrical situation in 3D. Hence, a more general description of how particles are accelerated during reconnection in 2D and 3D systems is needed.

In order to establish the mechanisms for electron acceleration during reconnection, we develop a guiding-center theory that can be used to diagnose how electrons gain energy during reconnection. We identify key mechanisms, including the curvature, gradient B, and parallel electric fields, and evaluate their relative contributions in 2D kinetic simulations with a modest to strong guide field. The results are easily generalizable to 3-D configurations. This is a similar approach to that used by Ref. 18 in their treatment of particle acceleration in relativistic, antiparallel reconnection.

In Sec. II, we derive a local bulk expression for particle acceleration and discuss the physical significance of each term. In Sec. III, we delineate our simulation methods and initial conditions. We present the results from 2D simulations in Sec. IV and examine momentum spectra in Sec. V. We discuss the significance of these results in Sec. VI.

## II. ELECTRON ACCELERATION IN THE GUIDING CENTER LIMIT

In order to examine various effects contributing to electron energy evolution, we begin with a standard treatment of the guiding-center approximation given by Ref. 29. The evolution of the energy *ϵ* of a single electron in the guiding-center limit is given by

where $b=B/|B|$, $\mu =me\gamma 2v\u22a52/2B$, is the magnetic moment, *γ* is the relativistic Lorentz factor, $v\u2225=v\xb7b$, and $vc$ and $vg$ are the curvature and grad-B drifts

In Eqs. (2) and (3), the electron cyclotron frequency $\Omega ce=eB/\gamma mec$. The curvature is $\kappa =b\xb7\u2207b$. If we sum over all particles in a local region, Eq. (1) becomes

which may be rewritten as

where *U* is the total kinetic energy, * u_{E}* is the E-cross-B drift,

*c*is the electron Alfvén velocity, and $u\u2225$ is the bulk velocity parallel to the magnetic field. The parallel and perpendicular pressures are $p\u2225$ and $p\u22a5$, respectively;

_{Ae}*n*is the electron density, $\beta \u2225=8\pi p\u2225/B2$ and $\beta \u22a5=8\pi p\u22a5/B2$. It is worth noting here that the equations above are not specific to electrons, but will apply to any species for which the guiding-center approximation is valid.

The first term in Eq. (5) is the acceleration by the parallel electric field. The second term corresponds to perpendicular heating or cooling due to the conservation of *μ* (the term in parentheses is $\u223cdB/dt$). The third term drives parallel acceleration and arises from the first-order Fermi mechanism described in Refs. 6 and 10. Freshly reconnected field lines downstream from a reconnecting X-line, accelerate as a result of the tension force ($\u223cB2\kappa $) and “straighten.” Particles that reflect from this moving field line receive a Fermi “kick” and thereby gain energy. Figure 1 shows a cartoon model of this effect. A particle trapped in an elongated magnetic island whose ends are contracting due to these tension forces will repeatedly gain energy as it reflects from the ends of the island. More generally, any particle reflecting from this moving bent field line will gain energy.

Though we do not address this directly in the present paper, we note that the tension force is important for other aspects of energy conversion. For example, it drives bulk ion outflow in the reconnection exhaust and generates the counterstreaming ion distributions that have been measured in the magnetosphere and solar wind.^{17,20,34}

## III. SIMULATIONS

We explore particle heating via simulations using the particle-in-cell (PIC) code p3d.^{38} Particle trajectories are calculated using the relativistic Newton-Lorentz equations, and the electromagnetic fields are advanced using Maxwell's equations.

The initial condition consists of a uniform guide field *B _{z}* superimposed on a double-Harris equilibrium.

^{19}The magnetic field configuration is

where *B*_{0} is the asymptotic reconnecting field, *w*_{0} is the current sheet half-width, and *L _{y}* is the length of the computational domain in the

*y*-direction. The density consists of two populations, a drifting population with density

that carries the current and a uniform background with density $0.2n0$.

We use an artificial mass ratio $mp/me=25$ and speed of light $c=15cA$, where *m _{p}* and

*m*are the ion and electron masses, respectively, and

_{e}*c*is the Alfvén velocity based on

_{A}*B*

_{0}and

*n*

_{0}. These choices allow for sufficient separation of scales (between proton and electron spatial scales and electromagnetic and particle time scales, respectively) while significantly reducing the computational expense of the simulation. Lengths in our simulation are normalized to the ion skin depth $di=c/\omega pi$ (based on

*n*

_{0}) and times are normalized to the inverse ion cyclotron frequency, $\Omega ci\u22121$. The initial temperature of all species is $0.25micA2$ for both background and current sheet populations. The current sheet half-thickness is set to $w0=0.25di$ so that reconnection will onset quickly. The grid scale is $\Delta =de/4\u22480.94\lambda D$, where $de=c/\omega pe$ is the electron inertial length and $\lambda D$ is the Debye length. We use periodic boundary conditions in both directions.

The goal of the present paper is to explore the mechanisms for particle acceleration using the expression for electron energy gain in Eq. (5). Since this equation is valid only for adiabatic motion, we limit our computations to systems with a non-zero initial guide field. It was shown previously that electrons are magnetized in reconnecting systems with guide field exceeding $0.1B0$.^{36} In this paper, we therefore focus on two simulations with non-zero guide fields, both with dimensions $Lx\xd7Ly=204.8\xd7102.4$. Simulation “A” has $Bz=0.2B0$ and “B” has $Bz=1.0B0$. We note that although our simulations contain two current sheets, we will often present results from only the upper current sheet. In all cases, the other sheet exhibits similar behavior. We will then present results from a much larger simulation with dimensions $Lx\xd7Ly=819.2\xd7409.6$ with a guide field $Bz=1.0B0$.

## IV. SIMULATION RESULTS: ELECTRON HEATING

Reconnection develops rapidly from the particle noise inherent in the PIC formulation. Figure 2 shows the evolution of the electron out-of-plane current density in simulation A. The tearing instability quickly generates many magnetic islands on each current layer, which continually grow and merge due to reconnection. We halt the simulation when the islands on the two current layers begin to interact, here at $t\u2248150$. Figures 3 and 4 display the parallel and perpendicular electron temperatures in simulations A and B early and late in the simulation. In simulation A, the parallel electron temperature increases substantially within the exhausts downstream of the X-lines and within the developing magnetic islands. The perpendicular temperature increment is strongest in localized regions in the cores of magnetic islands. In contrast, $Te\u2225$ significantly exceeds $Te\u22a5$ throughout the duration of simulation B.

Figures 5 and 6 show the contributions of the various terms in Eq. (5) in the upper current sheet in simulations A and B. At a given time, each term in Eq. (5) was calculated at each grid point and then integrated over space to give the displayed curves. Some smoothing was performed to reduce the noise in the calculations but the results shown are insensitive to its details. The sum of heating terms on the right-hand side of Eq. (5) is given by the dashed black line, and should be compared to the solid black line which represents the total measured electron heating. To the extent that the two match, Eq. (5) represents a valid description of the system. The discrepancy at early time is due to the small initial size of the islands (which makes the guiding-center approach less accurate). Sharp, small-scale gradients that develop during island mergers may be a source of additional discrepancies.

In Fig. 5, which corresponds to simulation A, the curvature-drift term is the dominant source of heating and $E\u2225J\u2225$ is negligible. The grad-B and $\u2202tB$ terms are also negligible and result in net cooling. This is because magnetic reconnection releases magnetic energy and therefore reduces the magnitude of *B*. Because of *μ* conservation electrons therefore on average lose energy in the perpendicular direction. By contrast, Fig. 6 shows that in simulation B the curvature-drift and $E\u2225J\u2225$ terms are comparable, while the other terms are negligible. The increased importance of the heating from the parallel electric field in the guide field unity case is because of the long current layers that develop in this case compared with those in the case of the weak guide field. Both simulations exhibit quasi-periodic heating, which is largely due to island mergers. This can be seen by comparing times *t* = 50 and *t* = 80 in simulation A: the former exhibits only modest heating, while the latter exhibits strong heating. Figure 9 (discussed further below) reveals that at *t* = 80 two islands are merging, which causes a burst of reconnection at the rightmost X-line in the system. In contrast, reconnection is proceeding in the normal fashion at *t* = 50.

Figures 7 and 8 show the spatial distribution of the curvature and $E\u2225J\u2225$ terms for simulation A at *t* = 120 and B at *t* = 125, respectively. As expected, the curvature-driven heating is primarily located in the reconnection exhaust regions and at the ends of the islands. Heating and cooling in the island cores are due to turbulent “sloshing” of plasma inside the island. We show later that there is little net heating from this behavior. Major Contributions to the $E\u2225J\u2225$ term are localized near the X-lines in both figures. The patchy regions of alternating heating and cooling throughout the islands, which are associated with electron holes^{3,9} do not on average produce much electron heating (shown later). Note the different color scales in the two plots in Fig. 7: the maximum intensity of the heating by $E\u2225$ is much smaller than that of the curvature drift, consistent with its relatively small contribution to electron heating shown in Fig. 5.

The patchy nature of the $E\u2225J\u2225$ term makes the interpretation of these data difficult. It is not obvious, for example, whether the heating due to $E\u2225$ around the X-line or due to the electrons holes dominates. As a further diagnostic, we therefore calculate the quantity

where $U$ is a heating term, the *y*-integral is taken over the half of the box containing the current layer (varying the bounds of integration does not significantly affect the result). The slope $d\Xi (x)/dx=\u222bU(x,y)dy$ yields the heating at a given *x*.

Figure 9 shows Ξ for the curvature-drift term in simulation A at two different times, corresponding to a temporal minimum in the curvature-drift heating (*t* = 50) and a temporal maximum (*t* = 80). The merger of two islands near $X\u2248160$ drives acceleration at the X-line in the far right of the simulation. The resulting island has a larger aspect ratio (length *x* compared to width *y*) so that freshly reconnected field lines experience a greater tension force around the far right X-line. This enhances the rate of electron heating in the exhausts around this X-line. The plot of Ξ also reveals that the heating and cooling in island cores result in little net heating, as can be seen, for example, inside the island at $x\u223c165$ at *t* = 80.

Figure 10 shows Ξ for $E\u2225J\u2225$ at *t* = 100 from simulation B. The dominant heating occurs near the primary X-lines at $x\u223c30$ and 100 as well as the secondary X-lines (due to island mergers) at $x\u223c150$ and 190. Inside the islands, there is net cooling. Many of the small scale fluctuations in the $E\u2225J\u2225$ term correspond with electron holes, which are driven by electron beams generated near the X-line.^{9} Because they tend to appear as bipolar structures in the heating term, they produce little net heating.

A number of the islands exhibit dipolar heating: the curvature term makes positive and negative contributions (red and blue) at the opposite ends of an island. Figure 11 exhibits this behavior. The island on the right drives heating due to Fermi reflection at both ends, and the plot of *v _{x}* shows large inward flows indicating island contraction. By contrast, the island on the left has dipolar heating. The entire island is moving in the $\u2212x$ direction. In the simulation frame, particles see receding field lines at the left end of the island and lose energy in a reflection. Equivalently, $uE\xb7\kappa <0$. However, the magnitude of the velocity at the right end is greater than that at the left, so the cooling at the left end is more than offset by the heating at the right: Ξ shows that the total heating across the island is positive. This is ultimately an issue of frame-dependence: in the frame of the island, both ends are contracting towards the center so that $uE\xb7\kappa >0$.

## V. SIMULATION RESULTS: ELECTRON SPECTRA

During reconnection with a strong guide field, which is expected to be the generic regime in most space and astrophysical systems, the dominant mechanisms for electron acceleration are the parallel electric field and Fermi reflection associated with the curvature drift, both of which accelerate electrons parallel to the local magnetic field. An important question, therefore, is whether the energetic component of the spectrum exhibits the strong anisotropy that is reflected in the moments $T\u2225$ and $T\u22a5$ in Fig. 4. Figure 12 shows electron spectra for the momenta parallel and perpendicular to the magnetic field. These spectra are taken from a simulation with the same initial conditions as in simulation B but in a larger domain $Lx\xd7Ly=819.2\xd7409.6$ carried out to *t* = 400. The larger simulation produces much better statistics in the particle spectra compared with simulation B shown earlier. In the parallel momentum, a clear nonthermal (that is, non-Maxwellian) tail develops by *t* = 50 and continues to strengthen until the end of the simulation. The perpendicular momentum also develops a nonthermal tail, but with an intensity that is smaller by more than two orders of magnitude. We note that these energetic spectra do not form power laws, which are not expected in periodic simulations that lack a loss mechanism.^{11}

It is hence clear that the dominant nonthermal acceleration occurs in the parallel component and the anisotropy survives over long periods of time as the simulation develops. An important question is what mechanism causes the perpendicular heating of energetic electrons. If the magnetic moment were exactly preserved, such particles would not be produced because the magnetic field *B* does not reach the required values anywhere. Therefore, the increase in the perpendicular spectrum must arise from scattering either because of non-adiabatic behavior in the narrow boundary layers that develop as a result of reconnection or because of the development of an instability directly driven by the anisotropy.

The distribution of the electron magnetic moment $\mu =mv\u22a52/2B$ for both simulations is shown in Fig. 13. It is clear that *μ* is very well conserved in simulation B, especially at low energies $\mu <1$ where the electrons remain adiabatic in the presence of the strong guide field $Bz=B0$. For simulation A with $Bz=0.2B0$, there is a drop of about 20% at the lowest energies, indicating that there is scattering into higher *μ*. This further suggests that the greater perpendicular heating in simulation A is due to non-adiabatic behavior in the small guide field regime.

## VI. DISCUSSION AND CONCLUSIONS

We have presented a guiding center model to explore the heating of electrons during reconnection with modest and large guide fields. We find that for a small guide field of $0.2B0$ (with *B*_{0} the asymptotic reconnecting field) electron heating is dominated by the Fermi reflection of electrons downstream of X-lines where the tension of newly reconnected field lines drives the reconnection outflow. The electron energy gain is given by the curvature drift of electrons in the direction of the reconnection electric field. In this small guide field case, heating from the parallel electric field and that associated with betatron acceleration (which is actually an energy sink) are negligible. In the case of a stronger guide field ($1.0B0$), the heating associated with parallel electric fields and the Fermi mechanism are comparable. The greater importance of the parallel electric field is because of the elongated current layers that form during reconnection with a guide field, which is where most parallel heating by this mechanism takes place. The net electron heating from electron holes, which densely populate the separatrices and island cores, is small because positive and negative contributions cancel. For both weak and strong guide fields, island mergers lead to bursts of electron acceleration.

An important scaling question concerns the role of heating by the parallel electric field in very large systems. The acceleration by parallel electric fields is largely confined to the narrow current layers around the X-line. In contrast, the heating through Fermi reflection occurs in a broad region in the exhaust downstream of X-lines and well into the ends of magnetic islands. At early times, the sheer number of X-lines could well make parallel electric fields a significant source of heating and acceleration. However, at late time when islands may be system-size, fewer *x*-lines might remain so parallel electric fields might not produce significant acceleration. In addition, the regions in which the $E\u2225J\u2225$ term dominates have characteristic widths that scale with $de\u221dme1/2$ with *m _{e}* the electron mass. In the simulations presented here, $mp/me=25$. For a real mass ratio of $mp/me\u22481836$, the corresponding regions with $E\u2225\u22600$ are expected to be much smaller. In contrast, the curvature drift dominates electron heating on island scales, which are not expected to depend on the choice of mass ratio once islands grow to finite size.

Evidently, further simulations are required to explore how the heating mechanisms given in Eq. (5) scale with system size. One of the motivations of exploring electron acceleration in the guiding center model is to develop a generic approach for addressing acceleration mechanisms in 3D systems where simple explanations of particle acceleration in contracting islands are no longer adequate: magnetic islands will generally no longer exist because field lines in 3D systems are chaotic and therefore volume-filling. However, since the conversion of energy by the relaxation of magnetic tension is fundamental to the reconnection process, we expect that the Fermi-like acceleration mechanism will remain important in a 3D system and its role can be quantified by evaluating the heating mechanisms presented in Eq. (5).

Finally, we comment briefly about the implications of the strong anisotropy of the energetic electrons seen in the spectra in Fig. 12 for the simulation with a guide field of 1.0 *B*_{0}. Gamma-ray flares have recently been detected in the Crab Nebula with photon energies exceeding $\u2248200$ MeV. These photons exceed the upper cutoff ($\u2248160$ MeV) that is obtained by balancing energy gain from the electric field $E\u223cB$ with that from losses associated with the synchrotron radiation reaction force. One proposed solution is that electrons are accelerated to the necessary energies ($\u22481015$ eV) in a large-scale reconnecting current sheet where the magnetic field strength approaches zero and the usual synchrotron assumptions do not apply.^{37} On the other hand, constraining the electrons in a narrow layer and preventing their escape into the reconnection exhaust and downstream magnetic island is a challenge. Another possibility is that reconnection takes place in the presence of a guide field such that the acceleration of the electrons is predominantly parallel to the local magnetic field so that the anisotropic energy distribution could mitigate synchrotron losses. In such a situation, a rough upper limit on reconnection-driven energization can be obtained by balancing the Fermi drive (scaling as $\gamma /Rc$, where *R _{c}* is a typical radius of curvature of a reconnecting magnetic field) against the curvature radiation loss ($\gamma 4/Rc2$),

where $Re=e2/mec2=2.82\xd710\u221213$ cm is the classical electron radius. For the most energetic events, *R _{c}* should equal the system size. Based on the flare duration of 1 day, $Rc\u22483\xd71015$ cm and the upper limit on the electron energy is $\u03f5=\gamma mec2\u223c1015$ eV, which is in the range needed to explain the observations. Clearly, a fundamental question is whether there are scattering mechanisms that limit the degree of anisotropy of the energetic particle spectrum and therefore reduce the upper limit given in Eq. (9).

## ACKNOWLEDGMENTS

This work has been supported by NSF Grants AGS1202330 and PHY1102479, and NASA grants NNX11AQ93H, APL-975268, NNX08AV87G, NAS 5-98033, and NNX08AO83G. Simulations were carried out at the National Energy Research Scientific Computing Center.