Simulation Models for Exploring Magnetic Reconnection

Abstract

Simulations have played a critical role in the advancement of our knowledge of magnetic reconnection. However, due to the inherently multiscale nature of reconnection, it is impossible to simulate all physics at all scales. For this reason, a wide range of simulation methods have been crafted to study particular aspects and consequences of magnetic reconnection. This article reviews many of these methods, laying out critical assumptions, numerical techniques, and giving examples of scientific results. Plasma models described include magnetohydrodynamics (MHD), Hall MHD, Hybrid, kinetic particle-in-cell (PIC), kinetic Vlasov, Fluid models with embedded PIC, Fluid models with direct feedback from energetic populations, and the Rice Convection Model (RCM).

Keywords: Magnetic reconnection; Magnetosphere; Numerical methods; Plasma physics; Plasma simulation; Solar corona; Turbulence.

Fig. 1

Physics of the Hall effect in reversing magnetic fields. A reversing anti-parallel magnetic field in the $\pm \hat{\mathbf{x}}$ direction is drawn using black arrows, with the neutral line as the dashed line. A uniform electric field in the $\hat{\mathbf{z}}$ direction is drawn in green. The smaller blue trajectory is that of an electron $E\times B$ drifting towards the neutral line with bulk velocity ${\mathbf{u}}_e$. An ion $E\times B$ drifting far from the neutral line has bulk ion velocity ${\mathbf{u}}_i$ identical to ${\mathbf{u}}_e$, so the electric field in this region is given by $\mathbf{E}=-\mathbf{u}\times \mathbf{B}/c$, where $\mathbf{u}$ is the single fluid (MHD) bulk flow velocity. Within an ion gyroradius $r_{Li}$ of the neutral line, the electron continues with bulk flow ${\mathbf{u}}_e$, while the ion demagnetizes (in the red trajectory) and ${\mathbf{u}}_i$ becomes small. In this region, there is a non-zero current density $\mathbf{J}$ (orange), and the electric field in this region is predominantly given by the Hall electric field $\mathbf{E}=\mathbf{J}\times \mathbf{B}/n_{e}ec$

Fig. 2

Sketch of how the Hall effect impacts linear waves in MHD. The wave structure for (a) shear Alfvén waves, (b) parallel propagating whistler waves, and (c) kinetic Alfvén waves. The equilibrium magnetic field $B_{y0}$ is the dashed black arrow. The perturbed magnetic field ${\mathbf{B}}_1$ are the black arrows. The perturbed ion bulk flow ${\mathbf{u}}_i$ are the red arrows, which occur in (a) because the wavelength is much larger than the ion gyroradius. In (b) and (c), the wavelength is at or below ion gyroscales, so the bulk flow is due to electron motion in the blue arrows. In (c), there is a large out-of-plane equilibrium magnetic field $B_{z0}$, so the magnetic perturbation $B_{z1}$ changes the magnetic pressure which to first order and requires a change to the gas pressure $p$

Fig. 3

“GEM Challenge” result showing Hall-MHD simulations faithfully obtain the rate of change of reconnected flux obtained in kinetic models. The legend describes the simulation approach for each curve; resistive-MHD is far slower than the other models. Adapted from Birn et al. (2001)

Fig. 4

Pressure and out-of-plane magnetic field in Hall MHD reconnection simulations. The legend describes the simulation approach for each curve; resistive-MHD is far slower than the other models. Adapted from Rogers et al. (2003)

Fig. 5

Evolution of (a) out-of-plane electric field $E_y$ at the X-point ($x=0$, $z=0$) and (b) root-mean-square value of the normal component of magnetic field $B_z$ evaluated over the length of the reconnecting current sheet (between two outflow regions) Vertical dashed lines mark the times at which $E_y$ attains its peak value $E_y^{peak}$ for different values of the current sheet half-thickness $ϵ$. (c) Scaling of $E_y^{peak}$ with $ϵ$ (in log -scale) from simulations (blue circles) and a fit $E_y^{peak}=0.05/{ϵ}^{1.15}$ (red line). Adapted from Jain and Sharma (2015b)

Fig. 6

Yee lattice and electric and magnetic fields in a cell in a 3D case, where the length of each side is $\mathrm{\Delta }x$. Electric fields $\mathit{E}$ are defined at the midpoint of each side of the cube, while magnetic fields $\mathit{B}$ are defined at the center of each face of the cube. In a 2D case, all the quantities are defined in the $x$-$y$ plane, projecting each position onto the cell in the $x$-$y$ plane

Fig. 7

2D PIC simulation result for magnetotail reconnection. (a) Electron fluid velocity $V_{ex}$, (b) electron temperature $T_e$, (c) zoom-in view of $T_e$, and (d) the heating term in Eq. (65). Adapted from Hesse et al. (2018a, 2019)

Fig. 8

Simulation results overplotted with magnetic field lines. (a) The $\mathbf{J}\mathbf{\cdot }\mathbf{E}$ term from Poynting’s theorem; (b) In-plane electron flow field; (c) $E_N$, the normal component of the electric field; (d) $E_{\parallel }$, the component of the electric field parallel to the magnetic field; (e) $B_M-B_{M,0}$, the change in the out-of-plane component of the magnetic field from its (spatially constant) initial value; (f) $S_L$, the horizontal component of the Poynting flux. Adapted from Swisdak et al. (2018)

Fig. 9

3D PIC simulation results. (a) 3D view of reconnection with a limited X-line extent, where the thin current sheet region extends $30d_i$ in $y$ (Huang et al. 2020). The mass ratio in the simulation is 25. (b) The current density on the $x=0$ plane (left) and magnetic field $B_z$ on the $z=0$ plane (middle and right) (Liu et al. 2019). The mass ratio is 75. The gray shaded area represents the “suppressed reconnecting region”. Adapted from Huang et al. (2020) and Liu et al. (2019)

Fig. 10

An example of particle trajectory analysis. Adapted from Oka et al. (2010)

Fig. 11

(a) An example of guiding-center drift analysis. Particle energization due to different drift currents for electrons (top) and ions (bottom). $j_c$ is due to particle curvature drift. $j_g$ is due to particle grad-$B$ drift. $j_m$ is due to magnetization. $j\prime \prime =j_c+j_g+j_m$. $\dot{K_e}$ and $\dot{K_i}$ are the energy change rates for electrons and ions, respectively. They are all normalized by $m_e{c}^2{\omega }_{pe}$. (b) Similar to (a) but shows energy dependent values (Li et al. 2019b) at a given time. Adapted from Li et al. (2017)

Fig. 12

PIC simulations of reconnection in shocks. (a) 2D simulation domain and the current density $J_z$. (b) Electron fluid velocity $V_{ex}$. (c) Ion fluid velocity $V_{ix}$. (d) 3D simulation domain and the current density $J_z$. Adapted from Bessho et al. (2019) and Ng et al. (2022)

Fig. 13

A schematic shows the improvement of the MHD-AEPIC (right) model from the MHD-EPIC (left) model. Adapted from Chen et al. (2023)

Fig. 14

Evolution of FTEs. Viewed from the Sun, a series of snapshots are shown with magnetic field lines colored by ion velocity $u_{iz}[\mathrm{k}\mathrm{m}/\mathrm{s}]$. Adapted from Chen et al. (2017)

Fig. 15

Energetic electron spectra from kglobal. A log-log plot of the electron differential density F(W) versus energy ($W$) at multiple times from a reconnection simulation with a guide field $B_g/B_0=0.25$. A power-law develops after $t/{\tau }_A\sim 3-5$. Inset: The late-time F(W) for several guide fields, illustrating the dependence on the ratio of the guide-to-ambient magnetic field. Adapted from Arnold et al. (2021)

Fig. 16

Sketch of the typical sampling of the plasma distribution function $f(\mathbf{x},\mathbf{v})$ adopted in PIC (left) and Eulerian (methods). Adapted from Finelli (2022)

Fig. 17

Magnetic reconnection modeled using three different codes (HVM, HVLF and PIC). Left: normalized reconnection rate $R/B_0{v}_{\mathrm{A}}$. The inset shows the time interval corresponding to the ticker curves in the main plot. To ease the comparison, the curves in the inset are shifted in time. Right: (first row) out-of-reconnection-plane magnetic field $\mathrm{\Delta }{B}_z/B_0=(B_z-B_z(t=0))/B_0$ showing the expected Hall quadrupolar pattern; (second row) electric field parallel to the ambient magnetic field $E_{\parallel }/\overline{E_{\parallel }}$, where $\overline{E_{\parallel }}$ is the root mean square of $E_{\parallel }$ in the shown region; (third row) current density in the out-of-plane direction $J_z$; (fourth row) electron current density in the plane $J_{\mathrm{e}}^{(\mathrm{i}\mathrm{n}-\mathrm{p}\mathrm{l}\mathrm{a}\mathrm{n}\mathrm{e})}$. The superposed black or white curves are the magnetic field lines. The three columns show results from the three different models. The left column show results from HVM at the simulation time $t=237.5{\mathrm{\Omega }}_{\mathrm{cp}}^{-1}$; the center column show results from the HVLF code at the simulation time $t=232.5{\mathrm{\Omega }}_{\mathrm{cp}}^{-1}$; the right column show results from the PIC code at the simulation time $t=235.0{\mathrm{\Omega }}_{\mathrm{cp}}^{-1}$. ${\mathrm{\Omega }}_{\mathrm{cp}}$ is the proton cyclotron frequency. Adapted from Finelli et al. (2021)

Fig. 18

(a) Plasma $\beta$; (b) proton $V_z$. The black lines show the magnetic field lines. (c) Magnetic field strength (black solid) and plasma density (red dashed) fluctuations from the virtual spacecraft location indicated with the black dot in panel (a). The anticorrelation between magnetic field and density fluctuations is compatible with mirror mode waves. Adapted from Hoilijoki et al. (2017)

Fig. 19

Evolution of the magnetotail current sheet in a 3D-3V Vlasiator simulation. The panels show the current sheet surface (defined as $B_r=0$) at different times, $t=1300\text{s}$ (a), $t=1400\text{s}$ (c), $t=1470\text{s}$ (g). The color of the surface corresponds to the current density $J$. The yellow line indicates the flow reversal between the Earthward and tailward reconnection outflow. The magenta and green lines are locations where $B_r=0$ and $B_z=0$ and correspond to X-lines and O-lines (differentiated using the sign of $\partial B_z/\partial r$, which is positive at the X-lines and negative at the O-lines). The primary reconnection line is where the X-line (magenta) and flow reversal (yellow) contours are approximately co-located. The background grid shows the coordinates but also the magnetic-field topology: the black grid shows areas where the magnetic field is directed northward, and the white grid shows the areas where it is southward-directed. Adapted from Palmroth et al. (2023)