{'fraction_non_alphanumeric': 0.04095175464817361, 'fraction_numerical': 0.012733903270162014, 'mean_word_length': 4.307789941297795, 'pattern_counts': {'":': 0, '<': 1, '<?xml version=': 0, '>': 5, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 3, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
The Weibel Instability inside the Electron-Positron Harris Sheet

Recent full-particle simulations of electron-positron reconnection have revealed that the Weibel instability plays an active role in controlling the dynamics of the current layer and maintaining fast reconnection. A four-beam model is developed to explore the development of the instability within a narrow current layer characteristic of reconnection. The problem is reduced to two coupled second-order differential equations, whose growing eigenmodes are obtained via both asymptotic approximations and finite difference methods. Full particle simulations confirm the linear theory and help probe the nonlinear development of the instability. The current layer broadening in the reconnection outflow jet is linked to the scattering of high-velocity streaming particles in the Weibel-generated, out-of-plane magnetic field.


I. INTRODUCTION
The temperature anisotropy driven Weibel instability [1] is thought to play an important role in several astrophysical systems. For instance, Weibel-mediated collisionless shocks in relativistic jets, pulsar winds, and gamma-ray bursts have been suggested ( [2] [3] [4]) as a possible particle acceleration mechanism. The Weibelgenerated magnetic field scatters particles, enabling them to bounce back and forth across the shock front, leading to acceleration via the first-order Fermi mechanism. [5] Recently the role of the Weibel instability in electronpositron (pair) reconnection has begun to receive notice. Magnetic reconnection is a fundamental problem in plasma physics, and is ubiquitous in astrophysical phenomena involving magnetic fields, where it is the principal mechanism for transforming magnetic energy into kinetic energy and heat. Historically, the greatest difficulty in modeling reconnection has been to demonstrate that it is fast enough to match observations of energy release in, for instance, solar flares. By comparing multiple simulation models (e.g., two-fluid, hybrid, and full particle-incell (PIC)) the GEM challenge [6] demonstrated that inclusion of the Hall term in the generalized Ohm's law was sufficient to produce fast reconnection. However, recent studies of electron-positron reconnection (in contrast to the usual electron-proton case) raised serious questions about the necessity of the Hall term for producing fast reconnection. In contrast with electron-proton plasmas, the mass symmetry in pair plasmas eliminates the Hall term. Yet, simulations suggest pair reconnection is still fast. Bessho & Bhattacharjee [7] attribute this fact to contributions from the off-diagonal components of the pressure tensor. Daughton & Karimabadi [8] later discussed the role of island formation along the reconnection layer. However Swisdak et al. [9] recently proposed that the Weibel instability, driven by an temperature anisotropy arising as inflowing plasma mixes with outflow from the x-line, localizes the reconnection layer, and leads to fast pair reconnection.
As an example of the importance of the Weibel instability in pair reconnection, we show, in Fig. 1, that suppressing the Weibel instability strongly influences the morphology of the current layer. In Fig. 1(a) we show a standard pair reconnection simulation where the Weibel instability causes the current layer to become turbulent and broaden, which opens the outflow exhaust as in Petschek's model [10]. In Fig. 1(b), we show the results of a simulation in which we suppress the instability by forcing the out-of-plane component of the magnetic field to zero. It is evident that without the turbulence provided by the Weibel instability the narrow current layer extends to the system size. The longer current layer reduces the reconnection rate in the simulation of Fig. 1(b) by one third. Since the Weibel instability plays such an important role in maintaining fast reconnection in pair plasma, a thorough understanding of its development in current layers is crucial.
Although Swisdak et al. [9] proposed that the Weibel instability strongly influences reconnection, they only briefly considered the effects of the current layer environment on the instability's development. In pair reconnection the outflow layer is typically narrow (on the order of a few electron inertial lengths) and confined on both sides by regions of strong magnetic field. In this work we take a closer look at the role of the current layer and, in particular, how it slows (or prevents) the onset of the Weibel instability. We also try to understand how the unstable Weibel mode is able to open (broaden) the current layer.
In section II of this paper we introduce our analytic model and its assumptions. Section III includes the derivation of the homogeneous dispersion relation and a comparison with kinetic theory. In Section IV we introduce the inhomogeneity arising from the reconnection geometry, and then numerically compute the eigenmodes and benchmark them with asymptotic theory in the limits of large and small current layer widths. In section V, we report on particle simulations that produce Weibel modes inside a Harris sheet. In section VI the implications for pair plasma reconnection are discussed and the downstream turbulent structure is compared with that of unstable Weibel modes.

II. THE GOVERNING EQUATIONS
A cartoon of the basic Weibel instability is shown in Fig. 2. Consider a uniform unmagnetized plasma with beams counter-propagating in the x-direction (more generally, possessing an anisotropic temperature with T x > T y , T z ). If a sinusoidal magnetic field component B z (y) arises from noise the positively charged particles with velocities v x > 0 will converge towards the x-axis because of the V × B Lorentz force while those with v x < 0 will diverge. The combination leads to a current density J x (y) of the correct sign to amplify B z , thus driving the mode unstable. (Negatively charged particles move in the opposite directions but have the same net effect.) To explore the structure of the Weibel instability in reconnection generated current layers, we use a fluid description. Since we are considering a pair plasma we include four species, α ∈ {p+, p−, e+, e−}, in our model: species of positrons and electrons with bulk velocities V and −V. In standard notation the governing equations are: We assume the pressure tensor can be written in the diagonal form and that the temperature components do not vary in space or time.
To particularize our coordinates we take the counter streaming velocities to be parallel tox. The perturbed magnetic field of the Weibel mode can then, without loss of generality, be taken to be parallel toẑ and the wavevector parallel toŷ. All physical quantities are assumed uniform in both the x and z directions, ∂/∂ z = ∂/∂ x = 0.
Our initial state is characterized by T yy,α = T zz,α = T.
The effective x-direction temperature is T x = mV 2 x + T xx . Note that by assuming ∂/∂x = 0 we eliminate T xx from the equations. The effective temperature anisotropy is then determined only by the streaming velocity.
Due to the mass and charge symmetry between electrons and positrons, we can collapse these eighteen equations, Eq. (7)-Eq. (12), into two second-order differential equations (see Appendix A for details): with the following definitions:χ ≡ n(Ṽ y,p+ −Ṽ y,p− − V y,e+ +Ṽ y,e− ); Ω is the gyrofrequency based on B x ; C s ≡ T /m is the sound speed; and d ≡ mc 2 /8πne 2 ≡ c/ω p is the skin depth with ω p the plasma frequency.

III. THE WEIBEL INSTABILITY IN A HOMOGENEOUS PLASMA
In a uniform magnetized plasma a dispersion relation can be found by combining Eqs. (13) and (14) and letting ∂/∂ y → ik y : Clearly while the streaming temperature (the first term of the right hand side) serves as the instability driver, thermal effects in the y-direction (the second term) and the background magnetic field (the third term) stabilize the mode.
In the even simpler case of a strongly anisotropic unmagnetized plasma, the kinetic growth rate of the Weibel instability is [9], [11]: where v th,x = T x /m is the thermal velocity in the x direction. It is evident that this matches the first term of Eq. (15) with v th,x → V x . The validity of using counter streaming cold plasma beams to analyze a single warm plasma was demonstrated by Davidson et al. [12], who showed that this instability is not affected by the detailed shape of the plasma distribution function, but only by the effective temperature.

IV. THE WEIBEL INSTABILITY INSIDE A HARRIS SHEET
The results of Section III were derived in the context of a homogeneous plasma. However for reconnection simulations it is important to study how the development of the Weibel instability proceeds inside a narrow current layer. In this section we seek to understand whether, and to what degree, the inhomogeneities in the plasma density and background magnetic field affect the instability.
Our equilibrium is taken to be the usual Harris layer with a slight modification that incorporates an anisotropic plasma temperature B x = B x0 tanh(y/δ); n p,e = 2n = n 0 sech 2 (y/δ); where B x0 , n 0 , and V x0 are constants, δ measures the width of current sheet and the subscripts p/e stand for positron/electron. We will refer to this setup as Profile A. The profiles are shown in Fig. 3(a). For the sake of comparing with full particle pair reconnection simulations, all physics quantities are presented in the same normalized units as those in Swisdak et al. [9]: the magnetic field to the asymptotic value of the reversed field B x0 , the density to n 0 which is the value at the center of the current sheet minus a possible uniform background density, velocities to the electron Alfvén speed V A,e , lengths to the electron inertial length d e ≡ mc 2 /4πn 0 e 2 , times to the inverse electron cyclotron frequency Ω −1 ce , and temperatures to m e V 2 A,e . We use T = 0.25m e V A,e in the analysis of this section, which is the initial set up of the pair reconnection simulation of [9], except for the latter's inclusion of a uniform background density of 0.2n 0 .
To find the modes of Eqs. (13) and (14) we discretize the governing equations in the y direction (imposing zero derivative boundary conditions) and find the eigenvalues of the resulting matrix. We use a grid size of δ/100 and domain size of 20δ, both of which are sufficient to ensure covergence. Before describing the numerical results, however, we investigate the behavior of the equation analytically.
A. kyde 1, δ/de 1 In this limit the current layer thickness and the wavelength of the instability are much smaller than the inertial length d e . From general considerations we expect the mode to be harder to excite in such circumstances, unless V x is large. It is straightforward to combine Eqs. (13) and (14) into a single second-order ordinary differential equation by eliminating the second term in Eq. (14), which is small since k y d e 1. The result is For the parameter regime of interest, d/d e ∼ Ω/Ω ce ∼ O(1). Therefore for every term in Eq. (19) to be of the same order (except, perhaps, the Ω 2 term), the following scaling rules must apply: After substituting for the functional form of the inhomogeneity and changing variables to u = cosh(y/δ)χ, we have where and Ω 0 ≡ B x0 e/mc. For V 2 x0 > C 2 s d 2 e /2δ 2 , Q is negative near y = 0 and positive at large y and Eq. (21) therefore has bounded solutions.
We further simplify Q by Taylor expanding for small y/δ, and neglecting the C 2 s tanh 2 and Ω 2 0 tanh 2 terms in Eq. (22), which are small in the ordering given in Eq. (20), This equation can be solved in terms of Hermite polynomials with eigenvalues that give a maximum growth rate of Clearly the density shear term, C 2 s /δ 2 , which arises from the second term of Eq. (19), has the strongest stabilizing effect. By letting γ max = 0, we obtain the critical temperature anisotropy required for the unstable mode in the small δ/d e limit, where T x0 ≡ (T x ) y=0 = mV 2 x0 . We plot γ versus anisotropy in Fig. 4(a) for δ = 0.1d e and in Fig. 4(b) for δ = 2d e . In the δ = 0.1d e case, the full equations (Eqs. (13), (14)) and the reduced equation (Eq. (21)) result in the same curve because the wavenumber of the growing mode is large enough to validate the approximation. The analytical solution from Eq. (24) follows the correct trend and matches the numerical results, particularly in the large anisotropy limit. We observe that γ is proportional to V x0 (i.e., V x0 ∼ (T x0 /T ) 1/2 ) as the scaling rule (Eq. (20)) suggests. By comparing the growth rate of the unmagnetized homogeneous plasma from Eq. (15), it is also clear that the instability is severely suppressed by the inhomogeneity. In contrast, the δ = 2d e case results in a somewhat closer match because of the increase in the inhomogeneity scale length.
Although in this limit the search for bounded modes is rather complicated, it is possible to gain some insight by expanding the homogeneous dispersion relation, Eq. (15) in the small k y limit: If we only keep terms of O(k 2 y ) we find that which is always positive and has no bounded modes. It is only by keeping the next term (i.e., O(k 4 y )) in the expansion that bounded modes exist, a fact that will guide our treatment of the full system.
We begin by neglecting theẼ x term in Eq. (14), solving for E x , and then substituting the result back into the equation, ultimately giving We then use this approximation in Eq. (13) in the small y/δ limit and assume that we are close to marginal stability. The resulting equation is: The parameter regime of interest to us is V x0 /V A,e ∼ C s /V A,e ∼ Ω 0 /Ω ce ∼ O(1). By requiring each term in this equation to be the same order we arrive at the ordering 1.
(30) To proceed, we transform the equation to Fourier space The quantity inside the square brackets has a maximum at Since we are looking for the maximally growing mode, we Taylor expand this quantity in s ≡ k y − k y0 around s = 0: χ.
(32) The solution of this equation can again be written in the form of Hermite polynomials with a maximal eigenvalue of Without the second term, which arises from the background magnetic field, the growth rate scales as , which is essentially the same result as the unmagnetized homogeneous relation, Eq. (15), in the k y d e 1 limit. From Eq. (33), we derive the marginal criterion in the large δ/d e limit, Finally we plot the threshold of marginal instability for different values of δ in Fig. 5. The numerical solution of our model with Profile A is carefully benchmarked by these asymptotic theories in both the small and large δ (or k y ) limits. The (blue) dashed and (red) dot-dashed lines are discussed later.

V. THE SMALL BOX PIC SIMULATION
In order to confirm our linear theory and study the nonlinear development of the Weibel instability, we conduct several simulations with the particle-in-cell (PIC) code p3d [13]. The electromagnetic fields are defined on gridpoints and advanced in time with an explicit trapezoidal-leapfrog method using second-order spatial derivatives. The Lorentz equation of motion for each particle is evolved by a Boris algorithm where the velocity v is accelerated by E for half a timestep, rotated by B, and accelerated by E for the final half timestep. To ensure that ∇ · E = 4πρ a correction electric field is calculated by inverting Poisson's equation with a multigrid algorithm. Although the code permits other choices, we work with fully periodic boundary conditions. We consider a system periodic in the x − y plane. The simulations presented here are two-dimensional, i.e., ∂/∂z = 0. The initial equilibrium consists of a Harris current sheet superimposed on an ambient population of uniform density n b , where B x,h , n h , n b are constants, subscript h stands for Harris sheet and δ is the half width of the current sheet. We input the Harris plasma with an initial temperature T xx,h > T yy = T zz ≡ T . The background plasma has an isotropic temperature T . Therefore, where T and T xx,h are constants. This equilibrium is very similar to Profile A defined in Eqs. (17)-(18) except for the inclusion of a background density. In the following we take the parameters n h = 0.8n 0 , n b = 0.2n 0 , B 2 x,h = 0.8B 2 x0 , δ = 2d e in order to compare later to our full simulations of pair reconnection, and refer to the initial condition as Profile B (see Fig. 3(b)).
In order to prevent the simultaneous growth of the tearing mode [14] we consider a domain with L x = 4d e L y = 50d e . We let the Weibel mode grow from noise for different values of the temperature anisotropy, measure the growth rate, and compare the results with the theoretical values from the full four-beam model in Fig. 6(b). Kinetic effects only slightly reduce the growth rates of our four-beam results and the eigenfunction predicted by the model is in good agreement with the simulations (see Fig. 6(a)). The anisotropy threshold for the four-beam model with Profile B is shown as a (blue) dashed line in Fig. 5.
The evolution of the temperature anisotropy and Weibel-generated B z from four small-box runs with different initial conditions is shown in Fig. 6(c)-(d). The anisotropies (at y = 0) decrease in time to slightly above the marginal value ∼ 2.3, while the amplitude of B z simultaneously rises. The increase of B z scatters the hot streaming plasma, reducing the central anisotropy. Only a small part of the energy released transfers to B z . The scattering increases T yy , and therefore the central pressure, causing the layer to expand. As a result, the ambient B x increases due to compression (not shown).
In homogeneous plasmas we can predict the saturation level by analyzing the particle motion in the ydirection, dV y /dt = −(e/mc)V x B z . Roughly speaking, the Weibel instability saturates when the magnetic field grows to a value such that the particles become magnetically trapped and can no longer amplify the field. Trapping occurs when the particle excursion alongŷ, ∆y, is comparable to the wavelength during the mode growth time, Equivalently, The saturation occurs when the magnetic bounce frequency, ω B , is comparable to the fastest linear growth rate. This is essentially the empirical result from Davidson et al. [12]. In the strong anisotropy limit, γ max ∼ V x k y,max (Eq. (15)), we obtain the simple saturation criterion k y,max ρ z ∼ 1, where ρ z is the gyro-radius in the B z field. In our inhomogeneous plasmas the same principle guides the saturation of the Weibel mode, although the predicted saturation value of B z is smaller due to the inhomogeneity.

VI. IMPLICATIONS FOR PAIR RECONNECTION
A. ky structure The basic features of a pair plasma reconnection simulation are shown in Fig. 7. The Weibel instability manifests itself as a chess-board-like structure in the downstream out-of-plane magnetic field (see Fig. 7 (b)). The structure travels downstream with the outflow speed from the x-line, which implies that the Weibel is a purely growing mode in the frame of the outflowing plasma. Swisdak et al. [9] proposed that the anisotropy (with typical magnitude T x /T y ∼ 2−4) driving the instability arises from cold inflowing plasma mixing with outflow from the x-line. Here we apply our analytical results to these observations. After fitting the inhomogeneity (Fig. 7(d)) seen in the reconnection simulations by Eqs. (35)-(36) with parameters n h = 0.25n 0 , n b = 0.2n 0 , B 2 x,h = 0.8B 2

x0
(denoted Profile C) we can find the necessary temperature anisotropy for Weibel to be unstable within the current layer. The (red) dot-dashed curve in Fig. 5 is the marginal criterion for the Weibel instability, and it gives us the minimum temperature anisotropy T x /T ∼ 3.3 for a typical reconnection layer of width 4d e (i.e., δ = 2d e ). This predicted minimum is within the observed anisotropy values (2 ∼ 4) and thus demonstrates that the reconnection current layer can support the Weibel mode. The even parity of B z produced by the Weibel mode ( Fig. 6(a)) is also seen in the reconnection simulation. Farther downstream, the nonlinear saturation of the Weibel instability stops the system from moving to even higher temperature anisotropies and keeps the system near marginal stability. The magnitude of B z saturates at ∼ 0.1 − 0.4B x0 in Fig. 6(d), which is comparable to the value shown in Fig. 7(b). Therefore by measuring the time, t s , for B z to saturate, we can estimate the half-length of the reconnection nozzle as L ≈ V A,e t s since B z develops from the instability while the plasma is convected out from the x-point at roughly the electron Alfvén velocity. As a result, anisotropy values of 2.5 ∼ 4 with saturation times from Fig. 6(d) give a predicted half-nozzle length in the range 25 ∼ 80d e , which compares favorably with the observed reconnection layer half-length of 60d e . Although these small runs use higher density plasmas than the reconnection runs, we expect the saturation behavior to be similar. Furthermore, the (red) dot-dashed curve in Fig. 6 (b), which represents the four-beam model solution with Profile C, indicates that an anisotropy of 4.0 produces a growth rate of ∼ 0.2Ω ce , which again leads to a nozzle half-length of ∼ 50d e (i.e., we expect its evolution to be similar to the solid curve in Fig. 6 (d)).
In contrast to our small-box runs where reconnection was suppressed, we expect the background to be noisier in a simulation that allows both reconnection and the Weibel instability to develop. However that will not strongly affect our results. Since the out-of-plane field B z grows exponentially from the initial noise we expect the estimated nozzle length to scale logarithmically with the initial noise level. Hence the predicted length via this estimate is insensitive to the noise in B z and we can expect our small-box runs to still give a reliable comparison to reconnection simulations.

B. kx structure
We now focus on the origin and magnitude of the finite k x in the region downstream of the x-line. Loosening our assumptions by letting ∂/∂ x = 0 in the linearized fourbeam model from Eq. (1)-(4), we can numerically solve for the two-dimensional dispersion relation in a homogeneous plasma. One new feature is the introduction of a two-stream instability (i.e., k y = 0 mode) that has previously been noted by Zenitani & Hesse [15] (see Appendix B for the dispersion relation). As shown in Fig. 8, the two-stream instability has a higher growth rate than the Weibel instability. Hence it always grows in front of the Weibel structures with wavelength ∼ 2π/3 ≈ 2d e , where the factor of 3 comes from the wavenumber for the maximum growth mode measured in Fig. 8(c). We do observe a double peaked distribution of the x-direction velocity and the E x signature of the two-stream instability (Fig. 7 (e), (c)) in the appropriate region with a wavelength comparable to the predicted value. The nonlinear development of the two-stream instability tends to merge the counter-streaming distributions, resulting in a singlehumped distribution farther downstream with T x > T y . The result is a transition from coexisting Weibel and twostream instabilities to a pure Weibel instability, as can be seen in Fig. 7(b).
The two-stream maximum growth rate predicted by the four-beam model is three times larger than that of the Weibel instability (divide the maximum growth rate of the solid curve in Fig. 8(c) by that of Fig. 8(b)), which seems inconsistent with their relatively close development in the reconnection simulations. However this disagreement can be reduced by the introduction of kinetic effects and a finite T xx , both of which are not included in our four-beam model. In particular, a finite T xx suppresses the growth of the two-stream instability. Kinetic theories predict lower growth rates for both the two-stream and Weibel instabilities, with the ratio of their fastest growth rates decreasing to about 1.6 (divide the maximum growth rate of the dashed curve in Fig. 8(c) by that of Fig. 8(b)). This lower ratio helps to explain the nearly coincident signatures of both instabilities in pair reconnection. In general, we also expect that a full kinetic treatment of Harris sheet inhomogeneous plasmas will produce slightly lower Weibel growth rates than those of our four-beam model, a feature that has already been observed in Fig. 6(b).
The two-stream instability can not explain the longer x-direction variance of the chess-board-like structure farther downstream. In that region, B x decreases from the asymptotic value of the reversed field (≈ B x0 ) at the nozzle edge to zero at the nozzle symmetry line while the temperature anisotropy remains large. Although Swisdak et al. [9] use a parity argument to argue against the possibility that the firehose instability plays a role, it is perhaps plausible that it could couple to the stronger Weibel mode and provide the finite k x . However this possibility is again ruled out by directly solving for the homogeneous dispersion relation of the firehose instability via a one fluid double-adiabatic model with finite Larmor radius corrections [16]. The firehose instability in a homogenous plasma with n p,e = 0.45n 0 , T x = 4.0T has the strongest growth rate ∼ 0.035Ω ce at B x ∼ 0.44B x0 and k x d e ∼ 0.07, which is an order of magnitude smaller than the growth rates of the Weibel and two-stream instabilities. Moreover, the predicted wavelength is too large to explain the observed k x .
A possible mechanism for the observed k x is proposed here. First imagine that the instability is confined in the y direction between reflecting walls (see Fig. 9). The reflecting walls mimic the strong B x field at the boundaries of the current layer. In this geometry there is an intrinsic scale length ∆ x associated with the trajectory of a particle away from the current layer, its reflection from the wall, and its motion back towards current layer. For a uniform B z , ∆ x = 2 ρ 2 z − (ρ z − ∆) 2 , where ρ z = V x /Ω c,z is the gyro-radius based on B z and ∆ is the half distance between the walls. The particle trajectory is sketched in the dashed box of Fig. 9(a). We suggest that this intrinsic scale length controls the x-dependence of B z seen in the simulations. We define the scale ∆ x shown as a dashed box in Fig. 9(a) as a Weibel-unit and consider the interaction between two Weibel-units. When all of the converging streaming particles (dashed curves with rightward arrows) leave the center of the righthand Weibel-unit in Fig. 9(a), this unit will need a replenishment of converging streaming particles (solid curves with rightward arrows) from its leftward neighbor to maintain its central current. Hence this self-consistent arrangement can arise without x-variation.
However, suppose that due to the initial random noise the lefthand Weibel-unit acquires the opposite polarity magnetic field, as seen in Fig. 9(b). Then the intrinsic ∆ x will arise such that the source for replenishing the converging streaming particles of the righthand Weibelunit is the diverging streaming particles (solid curves with rightward arrows) of a leftward neighbor after they bounce (perhaps multiple times) against the walls. The lefthand Weibel-unit can be reinforced from its rightward neighbor in the same manner. In this configuration we can estimate the distance between two neighboring Weibel-units to be N ∆ x , where N is the number of times a particle reflects from the walls.
This qualitative explanation is supported by PIC simulations, as can be seen in Fig. 10. We confine an plasma with an anisotropic temperature within a magnetic trough with thickness 5d e . Specifically, B x = 0, T x = 4.0T , and n = 1.2n 0 between y = −2.5 ∼ 2.5d e and B x = 1.0B x0 , T x = T , and n = 0.2n 0 outside this region. The simulation is performed in a domain of size 200d e × 25d e with periodic boundary conditions in both the x and y directions. Here the B x trough serves as the reflecting walls and the trough thickness of 5d e is comparable to the reconnection nozzle thickness. Two distant Weibel-units are not able to communicate with each other, and therefore it is not surprising to see that the magnetic fields at x = −50d e and 50d e in Fig. 10(a) have opposite signs. As time evolves, the Weibel-generated field gets stronger while the long Weibel-unit (x = −30d e ∼ 30d e at Fig. 10(b)) breaks up into smaller Weibel units of opposite polarity at (c). Five test positrons with initial velocity T x /m = 1.0V A,e (i.e., T x /T = 4) are randomly placed near (10d e , 0d e ). Their trajectories are shown as white curves in Fig. 10(c) and are blown up in Fig. 10(d). They are qualitatively similar to those described in Fig. 9(b). The converging motions of those trajectories at the trough center between x = 10 ∼ 30d e and x = 45 ∼ 65d e of Fig. 10(d) contribute to the current for the out-of-plane magnetic field. These particles reflect from B x in between x = 30 ∼ 45d e , where they deposit their momenta. At late time, tΩ ce = 37.5, the entire channel relaxes to a chain of Weibel units of alternating polarity, whose length scale is approximately determined by N ∆ x with N 1.
In order to apply this idea to the downstream regions in reconnection simulations, we approximate the characteristic streaming velocity as 1.0V A,e and the averaged out-of-plane magnetic field as 0.2B x0 . This implies a gyro-radius of 5d e . The layer thickness ∆ is 2d e along the reconnection nozzle, and thus ∆ x has a scale of approximately 2 5 2 − (5 − 2) 2 = 8d e , about half of the size of the observed structure in Fig. 7(c). This mechanism roughly explains the scale size of the observed variation. Furthermore, the gyromotion of the reconnection outflow helps explain how the reconnection nozzle broadens downstream. The out-of-plane magnetic field bends the outflow momentum from the x to the y-direction. As a result this flow pushes the background Harris magnetic field away from the symmetry line, broadening the current layer.

VII. SUMMARY AND DISCUSSION
We have developed a four-beam model to study the effects on the Weibel instability of a spatial inhomogeneity arising from a current layer. We have shown that the Weibel instability is able to grow within narrow Harris sheets, and its growth rate (saturation time), saturation magnitude, and mode structure fit those values observed in pair reconnection.
This further suggests that the Weibel instability might control the current layer dynamics in pair reconnection, where the Hall term is absent. Other candidate instabilities, such as the two-stream and firehose instabilities have rather minor effects, particularly since it is not clear if the firehose instability even appears in these systems. The high-velocity outflow scatters into the transverse direction due to the Lorentz force arising from the Weibelgenerated out-of-plane magnetic field. We argue that the associated increased pressure P yy is responsible for opening the pair reconnection nozzle and shortening the current layer. As a consequence, the shorter current layer generates a higher reconnection rate.
Even though a similar temperature anisotropy also could arise in electron-proton reconnection, signatures of the Weibel instability are not seen there. The reason is that the even shorter current layer (∼ 10d e ), controlled by whistler waves [17], leaves insufficient space for the unstable Weibel mode to grow. It is an open question as to how the Weibel instability controlled reconnection transforms into whistler mediated reconnection as the electron to ion mass ratio changes.
The development of the Weibel instability in the initial state of relativistic (k B T ∼ mc 2 ) pair reconnection is described in Zenitani & Hesse [15]. There the Weibel instability is shown to grow in front of a tangential discontinuity formed by the mixing of outflowing and ambient plasma, but it is not clear if this turbulence plays a role in the steady-state development of the outflow exhaust and the associated current layer. However, by artificially suppressing the Weibel-generated out-of-plane magnetic field, they do demonstrate the ability of the Weibel mode to broaden the current layer, similar to the behavior of the present non-relativistic case. The growth rate and wave vector of the relativistic Weibel instability are smaller by a factor of γ 1/2 L , where γ L is the Lorentz factor [18]. Consequently, we expect the instability to grow more slowly during relativistic pair reconnection, not only because of the intrinsically lower growth rate but also because of the relatively larger suppressing effect of the Harris reversed field on the enlarged mode  1: (Color online) Out-of-plane electron velocity for two pair reconnection simulations. In (a) the current layer length is about ∼ 120de (i.e., |x| 60de) [9]. In (b) the current layer length of a simulation with Bz = 0 scales as the system size, ∼ 700de (i.e., |x| 350de). The structure in the large island, which may be due to a two-stream instability, does not affect the behavior near the x-line.  To derive Eqs. (13) and (14), we define the variables, χ ≡Ṽ y,p+ −Ṽ y,p− −Ṽ y,e+ +Ṽ y,e− , then combine Eqs. (7)-(12) to yield − γη = n χ + nχ , Eqs. (43)-(48) represent a set of 6 equations for the 6 variables, {η,ζ,χ,θ,B z ,Ẽ x }. Note that {n, V x , B x , T } are unperturbed quantities specifying the initial conditions. Now use Eqs.