Relativistic Two-fluid Simulations of Guide Field Magnetic Reconnection

Snapshots of runs 3b with a guide field $B_{G}/B_{0}=0.5$ are presented in Figure 1. Shown in Figure 1a is an outflow structure of the plasma mean number flow,

We see that a fast outflow jet travels to the magnetic island (plasmoid) in the downstream region. Since there are ambient current sheet plasmas, the plasmoid exhibits a crab claw structure in the downstream, similar to large-scale MHD simulations of nonrelativistic magnetic reconnection (e.g., Abe & Hoshino (2001)). Inside the outflow channel, the flow speed $\langle v_{x}\rangle$ is $0.8\sim 0.85c$. The maximum 4-velocity of each species is $u_{xs}\sim 2.1c$. These speeds are slightly slower than the antiparallel counterpart. The backward flows around $x\sim 70L$ (blue regions; $\langle v_{x}\rangle<0$) are reverse flow from the plasmoid. Since the downstream plasmoid expands, the surrounding magnetic fields are compressed and then its magnetic pressure expels the plasma along the field lines. Sharp boundaries at $x\sim 60L$ are the fronts of such backward flows.

Shown in Figure 1b is the out-of-plane magnetic field $B_{y}/B_{0}$. We see that the guide field is compressed at the edge of the plasmoid and in the narrow reconnection outflow region, because the incoming reconnection flows transport and pile up the $B_{y}$ magnetic flux from the upstream region. Such $B_{y}$-compression changes the composition of outgoing energy flux in relativistic magnetic reconnection.

Figure 1c shows the charge separation, $(\gamma_{p}n_{p}-\gamma_{e}n_{e})/(\gamma_{p}n_{p}+\gamma_{e}n_{e})$. The charge density $\rho_{c}$ looks similar. We see that positrons dominate in the upper of the outflow region, and electrons in the lower. The charge separation generates the electric fields in the $x$–$z$ plane: $E_{x}<0$ in the upper half of the inflow region, $E_{x}>0$ in the lower half of the inflow region, and $E_{z}<0$ in the rightward outflow region. The $\bm{E}\times\bm{B}$ condition with the guide field $B_{y}$ is consistent with the reconnection flow pattern. It is impressive to see that the non-neutral layers extend over the large spatial domain and that the charge separation is strong $\sim$0.5. Such a non-neutral structure is qualitatively consistent with positron-rich or electron-rich regions around the reconnecting diffusion region in previous kinetic work (Zenitani & Hoshino, 2008).¹¹1Charge signs are opposite from Figure 2(c) in Zenitani & Hoshino (2008), because the guide field is set oppositely. Stripes around $x\lesssim 60L$ are transient oscillations invoked by the reverse flow. As the system evolves, the plasmoid and the reverse flows move further away. We also see strong perturbations around $x\sim 120L$, where the plasmoid hits the ambient plasmas.

Later, the system approaches to a quasi-steady stage. Shown in Figure 1d is a late-time structure of $B_{y}$ at $t=400\tau_{c}$ in run 3b. We see that the outflow structure in Figure 1b straightforwardly extends toward the $+x$ direction. The backward flow fronts are outside the presented domain at this time. The in-plane magnetic fields show a narrow Petschek-type structure. Interestingly, the $B_{y}$-structure is weakly bifurcated, and then the two peaks separate the outflow channel into three layers. We examine the outflow channel structure in Section 3.4. Around $x\sim 115L$, the bifurcated peaks exhibit a very weak wave structure, probably due to a velocity–shear driven instability. This introduces a minor oscillation inside the central channel in the further downstream; however, it does not change the global outflow structure. We also confirmed that the global outflow structure does not depend on the domain size by comparing runs 3a and 3b. Therefore, the outflow boundary condition does not change the physics.

We believe that the large distance to the upstream boundary ($z=60L$) is virtually eliminates boundary effects on the simulation. The plasmoid and reconnected field lines are well confined in $|z|<10$–$20L$ in panels in Figure 1. In the steady stage of $t=400\tau_{c}$, the slope angle of the field lines is only $\sim 0.1$ around the outflow region, and it becomes even smaller further upstream. For example, in the vicinity of the top boundary, the field line at $(0,+56.2L)$ is connected to $(120L,+60L)$.

The reconnection electric field $E_{y}$ at the $X$-point presents the transfer rate of the antiparallel magnetic fields from the upstream region. Since the surrounding conditions change in time, the normalized rate or “reconnection rate” gives a good measure of the system evolution. Here we define the rate as

where the subscript $in^{\prime}$ denotes the upstream properties at $x=0,z=20L$. The time evolution of the rates is shown in Figure 2a, and we see that the (normalized) reconnection becomes constant after the initial ramp up for $t\gtrsim 50\tau_{c}$. For comparison, the rate for run 3a is presented, too. Runs 3a and 3b are in excellent agreement except for around $t\sim 400\tau_{c}$ (indicated by an arrow) when the boundary effect temporally comes back to the reconnection point in the smaller run 3a. In Table 1, we present the peak values of the quasi-steady reconnection rates for all other runs.

So, what sustains such quasi-steady evolution of magnetic reconnection? In other words, what is responsible for the reconnection electric field $E_{y}$ around the $X$-point in our simulations? We study the composition of $E_{y}$ near the $X$-point using the Ohm’s law. Combining the momentum equations (Equation 2) of the two species, we obtain the following relation:

where $h_{s}=w_{s}/(n_{s}mc^{2})$ is the specific enthalpy. Note that we drop the time derivatives considering the quasi-steady condition ($\partial_{t}\sim 0$). In addition, along the inflow line ($x=0$), we find $n_{p}=n_{e}$, $\gamma_{p}=\gamma_{e}$, $v_{z,p}=v_{z,e}=\langle{v_{z}}\rangle$, $h_{p}u_{y,p}=-h_{e}u_{y,e}$, and $B_{z}=0$, and so the equation can further be simplified. We can also approximate that the artificial viscosity works for $h_{s}\bm{u}_{s}$ in a form of $-\nu_{z}{\partial_{zz}}$ with an effective viscous coefficient $\nu_{z}$. The $x$-derivative terms (${\partial}_{x}$ and ${\partial_{xx}}$) are negligible in this case. As a result, the $y$-component of Equation 10 along the inflow line yields

where $\eta_{eff}=(\tau_{fr}/q^{2}_{p})$ is the effective resistivity by the inter-species friction. Terms in the right-hand side represent the field convection, the advective inertia of the fluid, the frictional resistivity, and the viscous inertia, respectively.

Based on Equation 3.2, we study the composition of the reconnection electric field $E_{y}$ in Figure 2b. Here, $\nu_{z}$ is represented by the typical value around the neutral point. Outside the central diffusion region, the convection term mostly explains the electric field. The frictional resistivity term works inside the diffusion region where the current is strong. However, it explains only a quarter of $E_{y}$. We find that the local Lorentz factor $\gamma_{p}{\sim}1.6$ partially suppresses this term in Equation 3.2. Instead, the other two inertial effects are important. The advective term is the biggest contributor. In addition to the acceleration in the ${\pm}y$ directions ($u_{y,s}$), the specific enthalpy $h_{s}\sim 1+4p_{s}/(n_{s}mc^{2})$ increases because dissipation process heats incoming plasmas. In a very vicinity of the neutral plane ($z=0$), the viscous inertial effect replaces the advection. Although we introduced it for the numerical stability, we believe it reasonably represent the physics because the kinetic effect plays a quasi-viscous role near the neutral plane (Hesse et al., 2004).

We carried out a similar analysis in the antiparallel case (run 1), too. The frictional resistivity plays a bigger role, but it explains only one-third of the electric field. The local Lorentz factor $\gamma_{p}{\sim}1.6$ similarly suppresses the frictional resistivity. Regardless of the presence of the guide field, the inertial effects play a role to sustain magnetic reconnection in the relativistic two-fluid system.

Next, we investigate how the guide field changes the energy conversion process near the reconnection region in the quasi-steady stage. We consider the composition of the energy flux $\bm{F}$ in the following way:

In Equation 13, we decompose the Poynting flux into the contribution by the guide field ($B_{y}$) and the one by the reconnecting magnetic fields ($B_{x}$ and $B_{z}$). The plasma energy flux is decomposed into three parts. The pressure term transports the gas internal energy and it also contains the work by the gas. Since it recovers the classical enthalpy flux in the nonrelativistic limit ($\rightarrow\sum_{s}\frac{5}{2}p_{s}\bm{v}_{s}$ with $\Gamma=5/3$), we call this term “enthalpy flux” in this work. Note that the rest-mass contributions are considered separately from this enthalpy. The last two stand for the bulk kinetic energy ($\rightarrow\sum_{s}\frac{1}{2}mn_{s}{v}_{s}^{2}\bm{v}_{s}$) and the matter flow ($\rightarrow mc^{2}\sum_{s}n_{s}\bm{v}_{s}$), respectively.

We consider the energy budget around the reconnecting square region of $|z|<10L,|x|<40L$. The top panel in Figure 3 presents the composition of the incoming energy flux ($-\bm{F}\cdot\bm{\hat{z}}$) per cross section at the inflow boundary ($z=10L$) as a function of the guide field. They are normalized by the typical Poynting flux $c(B^{2}_{x}/4\pi)$ at the upstream border $x=0,z=10L$. The bottom panel in Figure 3 presents the outgoing energy flux ($-\bm{F}\cdot\bm{\hat{x}}$) at the outflow boundary ($x=40L$), evaluated by the same unit as the upstream values. They are evaluated at the specific time steps in different runs, but they are carefully selected so that the energy budget becomes nearly steady in the rectangle region. We see that both two fluxes are well balanced. Although it may be difficult to recognize in the figure, the matter fluxes are also well balanced.

For comparison, rescaled values of reconnection rates (Table 1) are overplotted on the top panel in Figure 3. It is reasonable to see that they are proportional to the incoming energy flux except for $B_{y}$-Poynting flux. The reconnection rate is relevant to the inflow speed $v_{in}$, which transports the upstream antiparallel magnetic flux $B_{x}$ into the reconnected field $B_{z}$. In these cases, under the similar upstream conditions, the $B_{x}$-Poyting flux decreases as ${\sim}\mathcal{R}$ and the accompanying plasma flux also decreases accordingly. We also see that the reconnection rate decreases but it weakly depends on the guide field amplitude $B_{G}/B_{0}$. This trend is consistent with a lot of previous nonrelativistic surveys (e.g., Huba (2005)).

We see that reconnection well dissipates the upstream antiparallel magnetic energy. When the guide field is weak, the Poynting energy is mostly converted into the plasma energy. Importantly, the energy is converted to the enthalpy flux rather than the bulk kinetic energy. The plasma pressure becomes strong in the outflow region in order to balance the strong upstream magnetic pressure. In magnetically dominated cases, plasma temperature becomes relativistic, and then such relativistic pressure substantially enhances the enthalpy flow. Therefore, the enthalpy flux dominates the bulk kinetic energy.

As the guide field increases, the relevant Poynting flux becomes a major component of the incoming energy flux. Simultaneously, the outgoing energy flux is dominated by the Poynting flux of the guide field. In the case of $B_{G}/B_{0}=0.5$ (run 3), as presented in Section 3.1, the guide field is compressed by a factor of ${\sim}3$. Consequently, the outgoing guide field Poynting flux increases and then it exceeds the enthalpy flux. When the guide field is stronger, the outflow is dominated by the Poynting flux by $B_{y}$. In the outflow region, the ratio of the $B_{y}$-Poynting flux to the plasma energy flux is $\sigma_{y,out}=0,0.32,1.21,4.29$, and $8.99$, respectively. It is interesting to see that these ratios are similar to the initial upstream $\sigma_{y,in}$ parameter (Table 1).

Under the single-fluid ideal RMHD approximation, we have the following relations from the polytropic law, the continuity, and the out-of-plane flux conservation:

where the prime sign (^′) denotes the single-fluid properties. They immediately suggest that $\sigma_{y}$ increases through the weak compression by magnetic reconnection,

In the case of run 3b, the typical compressional ratio of $\sim 2$ suggests that $\sigma_{y}$ increases by a factor of $\sim 1.6$. However, since our two-fluid system contains non-ideal effects, $\frac{d}{dt}(p^{\prime}/n^{\prime\Gamma})>0$ tends to reduce $\sigma_{y}$. In the stronger guide field cases, it is known that plasmas behave incompressibly and so the system is more likely to preserve $\sigma_{y,in}\sim\sigma_{y,out}$.

We look at the structure of the outflow region. Shown in Figure 4 are one-dimensional profiles across the outflow channel, at $x=90L$ at $t=400\tau_{c}$. The cut line is indicated by a white line in Figure 1d. From these profiles, we see that the outflow region consists of three characteristic layers, separated by two peaks of $B_{y}$: (1) a positron-rich boundary layer on the upper side ($2L\lesssim z\lesssim 4L$), (2) a central channel with high plasma pressure, and (3) an electron-rich boundary layer on the lower side.

The central channel contains a minor oscillating structure in plasma pressure (Figure 4a). Inside the central channel, the plasma temperature is hot and the outgoing fluid velocities are roughly similar $v_{s,x}\sim 0.8c$ across the channel. Since the Lorentz force $q_{s}(v_{s,x}B_{y})$ works differently in the presence of the guide field $B_{y}$, the positron and electron outflow channels are oppositely displaced to the ${\pm}z$-directions. The structures of the density, the pressure, the outflow speed $v_{x}$, and the inflow speed $v_{z}$ (magnified 10 times larger in Figure 4c) are all consistent with the displacement of the flow channels. The central channel becomes wider in the further downstream and we confirm that the system keeps the three-layer structure at least $x\sim 180L$.

We see that currents in the boundary layers are responsible for the global magnetic field topology. For example, $j_{x}$ explains the compression of the guide field $B_{y}$, and a “bifurcated” $j_{y}$-structure is consistent with the Petschek-type structure of in-plane magnetic fields (Figure 1d). We also see minor reverse currents near the center, but they are less important. Note that relative motion between the two species is significant in the boundary layers (Figure 4c), because plasmas sustain the currents in low-dense regions. Interestingly, in the upper positron-rich boundary layer, we see a fast electron flow in the $y$-direction and indeed electrons are major contributors to $j_{y}$ there. This is because the Lorentz force $j_{x}B_{z}$ modulates the local bulk outflow to the $-y$-direction, and so the electron flow looks overemphasized in our simulation frame.

In addition, in order to see the difference between our two-fluid model and the single-fluid RMHD model, we evaluate an equivalent single-fluid properties in the following way. We translate the conserved properties with the single-fluid properties (with the prime sign),

and then we calculate the single-fluid primitive variables by solving a quartic equation (Zenitani et al., 2009). The single-fluid pressure is assumed to be isotropic. Combining two fluid pressure into a single-fluid pressure does not lead to an isotropic pressure unless the relative velocity vanishes. The obtained pressure $p^{\prime}$ is overplotted in Figure 4a. We see that $p^{\prime}$ is in excellent agreement with the total two-fluid pressure $(p_{p}+p_{e})$ in the central channel. On the other hand, the discrepancy is significant $p^{\prime}/(p_{p}+p_{e}){\sim}1.4$ in the two boundary layers and so it may widen the outflow channel in single-fluid simulations. It is reasonable that the discrepancy is observed where the relative motion between the species is significant (Figure 4c).

In the boundary layers the local frozen-in condition also breaks down. We similarly analyze the composition of the reconnection electric field $E_{y}$ based on Equation 10. Figure 4d presents the profile of $E_{y}$, along with the field convection and the $z$-convection of fluid inertia terms,

In Figure 4d, the inertia terms are magnified 10 times larger to emphasize them. We found that the non-ideal contribution (the difference between $E_{y}$ and the $\langle\bm{v}\rangle\times\bm{B}$ term) is mostly explained by these terms. The positron inertial term appears in the lower electron-rich boundary layer, and electron term in the lower electron-rich layer. For example, in the positron-rich boundary region, both the $z$-component of the positron velocity and the $y$-momentum of the positrons are much smaller than that of the electrons. As a result, the electron inertial term dominates even though the positron density is larger. Around $z\sim\pm 3L$, the inertial contributions from both two fluids sustain 12% of $E_{y}$. Because of the large Reynolds number $S\sim 3000$, the frictional resistivity effect is negligible and we confirmed that the viscosity is negligible, too. Compared with the reconnection site where the nonideal terms are 100% responsible for the electric field, the nonideal effects are smaller, but they still play a role here.

Shown in Figure 5 is the $x$-dependence of the ratio of the single-fluid pressure ($p^{\prime}$) to the total fluid pressure ($p_{p}+p_{e}$), which is a good measure of the two-fluid effects. As the central outflow channel expands downstream, the two-fluid effects are always localized in the relevant boundary regions. Their amplitude and spatial width are similar until the backward flows around the plasmoid hit the boundary layers.