# Annihilation of Magnetic Islands at the Top of Solar Flare Loops

Yulei Wang School of Astronomy and Space Science, Nanjing University, Nanjing
210023, People’s Republic of China  Key Laboratory for Modern Astronomy and
Astrophysics (Nanjing University), Ministry of Education, Nanjing 210023,
People’s Republic of China  Xin Cheng School of Astronomy and Space Science,
Nanjing University, Nanjing 210023, People’s Republic of China  Key Laboratory
for Modern Astronomy and Astrophysics (Nanjing University), Ministry of
Education, Nanjing 210023, People’s Republic of China  Mingde Ding School of
Astronomy and Space Science, Nanjing University, Nanjing 210023, People’s
Republic of China  Key Laboratory for Modern Astronomy and Astrophysics
(Nanjing University), Ministry of Education, Nanjing 210023, People’s Republic
of China  Quanming Lu CAS Key Laboratory of Geospace Environment, Department
of Geophysics and Planetary Science, University of Science and Technology of
China, Hefei 230026, People’s Republic of China  CAS Center for Excellence in
Comparative Planetology, People’s Republic of China

###### Abstract

The dynamics of magnetic reconnection in the solar current sheet (CS) is
studied by high-resolution 2.5-dimensional MHD simulation. With the commence
of magnetic reconnection, a number of magnetic islands are formed
intermittently and move quickly upward and downward along the CS. When
colliding with the semi-closed flux of flare loops, the downflow islands cause
a second reconnection with a rate even comparable with that in the main CS.
Though the time-integrated magnetic energy release is still dominated by the
reconnection in main CS, the second reconnection can release substantial
magnetic energy, annihilating the main islands and generating secondary
islands with various scales at the flare loop top. The distribution function
of the flux of the second islands is found to follow a power-law varying from
$f\left(\psi\right)\sim\psi^{-1}$ (small scale) to $\psi^{-2}$ (large scale),
which seems to be independent with background plasma $\beta$ and if including
thermal conduction. However, the spatial scale and the strength of the
termination shocks driven by main reconnection outflows or islands decrease if
$\beta$ increases or thermal conduction is included. We suggest that the
annihilation of magnetic islands at the flare loop top, which is not included
in the standard flare model, plays a non-negligible role in releasing magnetic
energy to heat flare plasma and accelerate particles.

##  1 Introduction

Magnetic reconnection, a fundamental energy-releasing process in magnetized
plasma, is believed to be the core mechanism driving solar eruptions including
solar flares, coronal mass ejections (CMEs), and other solar activities.
According to the standard flare model, the sheared magnetic field is stretched
and forms a large-scale thin current sheet (CS) in the wake of the eruption
(Kopp & Pneuman, 1976; Lin & Forbes, 2000; Lin et al., 2015). The reconnection
in the CS in turn drives the CME eruption and the formation of the flare
loops. Theoretically, the squeezed thin CS is unstable as the tearing mode
instability grows (Priest & Forbes, 2000; Pucci & Velli, 2013; Loureiro &
Uzdensky, 2015). Especially, in the corona of high Lundquist-number, the
plasmoid instability will dominate the dynamics of magnetic reconnection
(Samtaney et al., 2009; Bhattacharjee et al., 2009; Huang & Bhattacharjee,
2010; Huang et al., 2011; Huang & Bhattacharjee, 2013, 2012; Huang et al.,
2017; Loureiro & Uzdensky, 2015; Shen et al., 2011; Zhao & Keppens, 2020). The
classical CS predicted by Sweet-Parker is thus split and the reconnection
enters into the fast scheme and generates a number of magnetic islands.

Magnetic islands are believed to be closely related to fine structures in the
CS (Shen et al., 2011; Lin et al., 2015). The coalescence of magnetic islands
can further produce secondary CSs and islands (Bárta et al., 2011). Various
researchers have analyzed the cascading law of magnetic island flux
theoretically or numerically, which predicts that the distribution function of
magnetic island flux $f\left(\psi\right)$ in the CS follows a power-law
between $\psi^{-1}$ to $\psi^{-2}$ (e.g., Uzdensky et al., 2010; Huang &
Bhattacharjee, 2012; Shen et al., 2013; Lynch et al., 2016; Ye et al., 2019).
The fragmented and turbulent characteristics of the CS have also been
suggested by recent observations (e.g., Cheng et al., 2018).

The flare loop top, interacting with the reconnection downflows, also shows
complex dynamics. The turbulent characteristics of the loop top have been
studied by high-resolution simulations (e.g., Cai et al., 2019; Ye et al.,
2020; Cai et al., 2021). The generation of flare quasi-periodic pulsations
(QPPs) and the supra-arcade downflows (SADs) above the flare loop top is
interpreted to be closely related to the fast reconnection downflows (Takasao
& Shibata, 2016; Guo et al., 2014). Abundant MHD shocks have been predicted by
simulations (Takasao et al., 2015). Especially, as an important prediction of
the standard solar flare model, the fast mode termination shocks (TSs), formed
when the downflows hit the relatively steady high-density structure at the
loop top, are believed to be an effective mechanism to accelerate particles
(Tsuneta & Naito, 1998; Shen et al., 2018; Kong et al., 2019, 2020; Ruan et
al., 2020), which are supported by radio imaging observations (Aurass & Mann,
2004; Chen et al., 2015). The downflow magnetic islands also affect the
dynamics of the flare loops (Ye et al., 2020; Cai et al., 2019), the formation
of TSs (Shen et al., 2018), and even the acceleration of energetic electrons
(Kong et al., 2020).

Nevertheless, mainly owing to the limitation of observations and simulations,
magnetic reconnection in the CS and its relation to the dynamics at the flare
loop top are still far from being fully understood. In this paper, we perform
high-resolution 2.5-dimensional MHD simulations of the CS reconnection in a
high-Lundquist-number and low-$\beta$ coronal environment and focus on the
dynamics of magnetic islands evolving towards flare loops. We find that the
downflow magnetic islands quickly annihilate after colliding with the tip of
cusp-shaped flare loops. Such a second reconnection process is characterized
by horizontal CSs, smaller-scale islands, and a comparable reconnection rate
with that in the main CS. In Section 2, we describe our method. Section 3
presents the main results, which are followed by a summary and discussion.

##  2 Method

In this work, the MHD equation we solve is as follows:

$\displaystyle\frac{\partial\rho}{\partial t}+\nabla\cdot\left(\rho{\bf
u}\right)$ $\displaystyle=$ $\displaystyle 0\,,$
$\displaystyle\frac{\partial\left(\rho{\bf u}\right)}{\partial
t}+\nabla\cdot\left(\rho{\bf u}{\bf u}-{\bf BB}+P^{*}\right)$ $\displaystyle=$
$\displaystyle 0\,,$ $\displaystyle\frac{\partial e}{\partial
t}+\nabla\cdot\left[\left(e+P^{*}\right){\bf u}-{\bf B}\left(\bf{B}\cdot{\bf
u}\right)\right]$ $\displaystyle=$
$\displaystyle\nabla\cdot\left(\kappa_{\parallel}\hat{\bf{b}}\hat{\bf{b}}\cdot\nabla
T\right)\,,$ (1) $\displaystyle\frac{\partial{\bf B}}{\partial
t}-\nabla\times\left(\bf{u}\times\bf{B}\right)$ $\displaystyle=$
$\displaystyle-\nabla\times\left(\eta{\bf J}\right)\,,$ $\displaystyle{\bf J}$
$\displaystyle=$ $\displaystyle\nabla\times{\bf B}\,,$

where, $P^{*}=p+B^{2}/2$, $e=p/\left(\gamma-1\right)+\rho u^{2}/2+B^{2}/2$,
$\kappa_{\parallel}=\kappa_{0}T^{2.5}$, and standard notations of variables
are used. In our model, only the thermal conduction (TC) along magnetic field
is considered which is much more important than perpendicular components and
$\kappa_{0}$ is set as $10^{6}\,\mathrm{erg\cdot s^{-1}\,cm^{-1}\,K^{-3.5}}$
(Yokoyama & Shibata, 2001). Here, all variables have been normalized according
to constant units. The unit of space is chosen as the scale of simulation
region $L_{0}=5\times 10^{9}\,\mathrm{cm}$ that is comparable with the value
for a typical observed flare event (see, Cheng et al., 2018). We also assume
that the coronal plasma is composed of fully ionized hydrogen with the
electron density $n_{e}=10^{10}\,\mathrm{cm^{-3}}$ and the averaged particle
mass $\bar{m}=0.5m_{p}=8.36\times 10^{-25}\,\mathrm{g}$, which gives the unit
of mass density $\rho_{0}=2n_{e}\bar{m}=1.67\times
10^{-14}\,\mathrm{g/cm^{3}}$ (Priest & Forbes, 2000). Here, $m_{p}$ is the
mass of proton. The unit of magnetic strength is set as a typical coronal
value, namely, $B_{0}=20\,\mathrm{Gauss}$ (also see Chen et al., 1999; Ye et
al., 2020). Based on $L_{0}$, $\bar{m}$, $\rho_{0}$, and $B_{0}$, normalized
units of other variables are deduced as
$u_{0}=B_{0}/\sqrt{\mu_{0}\rho_{0}}=4.36\times 10^{7}\,\mathrm{cm/s}$,
$t_{0}=L_{0}/u_{0}=114.61\,\mathrm{s}$,
$p_{0}=\rho_{0}u_{0}^{2}=3.18\,\mathrm{Pa}$,
$T_{0}=\bar{m}u_{0}^{2}/k_{B}=11.52\,\mathrm{MK}$, $\kappa_{\parallel
0}=k_{B}\rho_{0}u_{0}L_{0}/\bar{m}=6.02\times 10^{11}\,\mathrm{erg/s\cdot
cm\cdot K}$, and $\eta_{0}=L_{0}u_{0}=2.18\times 10^{17}\,\mathrm{cm^{2}/s}$,
where $\mu_{0}$ is the magnetic permeability in vacuum and $k_{B}$ denotes the
Boltzmann constant. The ratio of specific heat is $\gamma=5/3$. Under this
configuration, the unit time, velocity, and temperature are also on the same
order as the observational ones given by Cheng et al. (2018).

The initial magnetic field is set according to the CSHKP model (Chen et al.,
1999; Ye et al., 2020)

$\displaystyle B_{x}$ $\displaystyle=$ $\displaystyle 0\,,$ $\displaystyle
B_{y}$ $\displaystyle=$ $\displaystyle\begin{cases}\sin\left(\frac{\pi
x}{2\lambda}\right)\,,&\left|x\right|\leqslant\lambda\,,\\\
1\,,&x>\lambda\,,\\\ -1\,,&x<-\lambda\,,\end{cases}$ (2) $\displaystyle B_{z}$
$\displaystyle=$ $\displaystyle\sqrt{1-B_{x}^{2}-B_{y}^{2}}\,,$

where $\lambda$ denotes the half-width of the CS, which is set as 0.1. To
approximate the chromosphere, transition region, and corona, the initial
density distribution follows Takasao et al. (2015), namely,

$\rho\left(y\right)=\rho_{chr}+\frac{\rho_{cor}-\rho_{chr}}{2}\left[\mathrm{tanh}\left(\frac{y-h_{chr}}{l_{tr}}+1\right)\right]\,,$
(3)

where, $\rho_{cor}=1$, $\rho_{chr}=10^{5}$, $h_{chr}=0.1$, and $l_{tr}=0.02$.
Fast reconnection is initially triggered by a localized anomalous resistivity
and thus the resistivity distribution is

$\eta=\begin{cases}\eta_{b}+\eta_{a}\mathrm{exp}\left[-\frac{x^{2}+\left(y-h_{\eta}\right)^{2}}{l_{\eta}^{2}}\right]\,,&t\leq
t_{\eta}\,,\\\ \eta_{b}\,,&t>t_{\eta}\,,\end{cases}$ (4)

where, $\eta_{a}=5\times 10^{-4}$, $h_{\eta}=0.5$, $l_{\eta}=0.03$,
$t_{\eta}=5$, and the background resistivity is uniformly set as
$\eta_{b}=5\times 10^{-6}$ to simulate a high-Lundquist-number environment.
Physically, the localized anomalous resistivity can be caused by microscopic
instabilities in the current sheet such as the lower hybrid drift and/or ion
acoustic instabilities which can boost the process of reconnection (Yokoyama &
Shibata, 2001). Because $\eta_{a}\gg\eta_{b}$, the anomalous resistivity
dominates the evolution during $t<t_{\eta}$, which provides similar initial
states at $t=t_{\eta}$ for the cases listed in Tab. 1. The background pressure
$p_{b}$ is uniform to obtain a static initial state and is assigned as shown
in Tab. 1. The boundary conditions are arranged as follows. The left ($x=-1$)
and right ($x=1$) are free boundaries, the top ($y=4$) is no-inflow boundary,
and the bottom ($y=0$) is symmetric boundary. To reduce the influences of
numerical boundaries, we only use data in the region defined by
$x\in\left[-0.5,0.5\right]$ and $y\in\left[0,2\right]$ for analysis.

Cases | $p_{b}$ | $\beta$ | TC
---|---|---|---
1 | 0.02 | 0.04 | No
2 | 0.05 | 0.10 | No
3 | 0.10 | 0.20 | No
4 | 0.02 | 0.04 | Yes
Table 1: Configurations of all simulation cases.

The above system is numerically solved with Athena 4.2 (Stone et al., 2008).
We use the HLLD Riemann solver (Miyoshi & Kusano, 2005), 3-order piecewise
parabolic flux reconstruction algorithm (Stone et al., 2008), and the Corner
Transport Upwind (CTU) method (Gardiner & Stone, 2008) to solve the
conservation part of Eq. 1. The resistivity and TC are calculated by the
explicit operator splitting method. To suppress numerical dissipation and get
a uniformly high resolution (Shen et al., 2011; Ye et al., 2020), we set high-
precision uniform Cartesian grids, namely, 3840 and 7680 grids on $x$ and $y$
directions, respectively. The corresponding pixel scales in two directions are
$\Delta x=\Delta y=26\,\mathrm{km}$. The maximum simulation time in our
simulation is $t_{max}=15$, which corresponds to $28.65$ minutes in physical
time.

##  3 Results

###  3.1 Overview

Figure 1: Snapshots of mass density $\rho$, temperature $T$, and current
density $J_{z}$ at four typical moments for Case 4. The 1st column at $t=5$
shows the initial status when the localized anomalous resistivity is shut
down. The erupting plasmoid, the stretched main CS, and the flare loop “seed”
are clearly displayed by the 2nd column at $t=7.55$. The 3rd column records
the moment when a relatively large magnetic island starts to collide with the
loop-top flux at $t=11$. The 4th column exhibits the final state of the flare
loop systems at $t=15$. The x-marker and the dashed box in (b2) denote
respectively the principal X-point in the main CS and the inflow region for
calculating the averaged magnetic strength and Alfvén speed in Eq. 5. The
dashed box in (b3) marks the loop-top region analyzed in Fig. 3. The three
dashed horizontal lines in (b4), namely, $y=0.25$, $y=y_{lt}\left(t\right)$,
and $y=1$, define the loop-top region, $R_{1}$, and the principal reconnection
site in the main CS, $R_{2}$. Here, $y_{lt}\left(t\right)$ denotes the height
of the top of the cusp-shaped loops, which dynamically changes with the
evolution of the system.

Figure 1 shows four typical snapshots of density $\rho$, temperature $T$, and
out-of-plane current $J_{z}$ for Case 4. Before $t=5$, the reconnection is
dominated by the anomalous resistivity at $y=0.5$ where the main CS is
squeezed. After $t=5$, the anomalous resistivity is set to zero and the
reconnection is determined by the uniform background resistivity $\eta_{b}$.
At $t=5$, a plasmoid grows above the X-point and starts to move upward, while
the lower part of the main CS forms a ”seed” of flare loops with two foot-
points line-tied in the high-density chromosphere (1st column of Fig. 1). As
the ejection of the upper plasmoid, the main CS is rapidly stretched (2nd
column of Fig. 1). When its aspect ratio is large enough, the main CS is
quickly split by the plasmoid instability (also see Shen et al., 2011). The
downflow magnetic islands carrying substantial kinetic and thermal energies
are intermittently formed and then collide with the cusp-shaped flare loops.
At $t=11$, a relatively large magnetic island starts to interact with the
loop-top flux (3rd column of Fig. 1) and disappears at $t=11.65$. Fed by the
downflow islands, the flare loops also manifest an obvious expansion. At
$t=15$, the simulation stops and the flare loops evolve into their final
states (4th column of Fig. 1). The distributions of density and temperature
are consistent with the results of Ye et al. (2020) and highly resemble
observed flare loops (e.g., Sun et al., 2014).

###  3.2 Reconnection in the Main CS

Figure 2: Evolutions of the magnetic reconnection rate in the main CS. The
blue, orange, and yellow solid curves show $\mathcal{R}_{cs}\left(t\right)$
for Cases 1, 2, and 3, respectively, while the purple dashed curve presents
Case 4 in which TC is considered. Here we only exhibit the reconnection rate
after $t=5$, because it is dominated by the anomalous resistivity before
$t=5$.

We first estimate the reconnection rate in the main CS as shown by the area
$R_{2}$ in Fig. 1(b4)). It is represented by the reconnection electric field
at the principal X-point (see also, Yokoyama & Shibata, 2001), namely,

$\mathcal{R}\left(t\right)=\frac{\max\left(\eta_{b}\left|J_{zxp}\left(t\right)\right|\right)}{B_{in}\left(t\right)u_{in}\left(t\right)}\,,$
(5)

where, $J_{zxp}$ denotes the out-of-plane current density at X-points,
$\max\left(\cdot\right)$ means taking the maximum value in the target CS
region, and $\mathcal{R}$ is normalized by the product of inflow magnetic
strength $B_{in}$ and Alfvén speed $u_{in}$ which are averaged in the adjacent
region of the principal X-point defined by $x\in\left[-0.05,0.05\right]$ and
$y\in\left[y_{xp}-0.05,y_{xp}+0.05\right]$ (see, Fig. 1(b2)). Here, $y_{xp}$
denotes the y-coordinate of the principal X-point. The null points are
determined by the method developed by Parnell et al. (1996).

For each simulation case, the reconnection rate in the main CS,
$\mathcal{R}_{cs}$, shows a similar evolution trend (Fig. 2). To be specific,
the average value of $\mathcal{R}_{cs}$ first rises quickly and then reaches
its peak value which is on the order of $0.01$. We can see significant
oscillations in $\mathcal{R}_{cs}$ curves which correspond to the formation of
magnetic islands (Yokoyama & Shibata, 2001). The overall evolution of the
reconnection rate in the main CS is consistent with previous numerical results
of 2D reconnection (e.g., Ni et al., 2012; Huang & Bhattacharjee, 2016;
Zenitani & Miyoshi, 2020). Differing $\beta$ mainly affects the initial rising
stage before $t=7$. With the increase of $\beta$, the reconnection rate at
this stage becomes smaller. However, when the reconnection is fully developed
($t>10$), the evolution of $\mathcal{R}_{cs}$ is less affected by $\beta$.
After including TC, the reconnection rate before $t=6$ gets smaller and the
average value also decreases slightly after $t=10$ (Fig. 2).

###  3.3 Second Reconnection of Magnetic Islands

Figure 3: Two typical magnetic island annihilation events in Case 4 and Case
1. Panels (a), (b), and (c) draw the current, the horizontal velocity, and the
vertical velocity around the peak time of $\mathcal{R}_{lt}$. In panel (d),
the blue and orange solid curves plot the reconnection rate $\mathcal{R}_{lt}$
and the height of CSs, respectively. And the dashed blue curve depicts the
reconnection rate $\mathcal{R}_{cs}$ in the main CS. The box and ”x” symbol in
panels (a) mark the horizontal CS and the principal X-point, respectively. The
principal TS is marked by the box in panel (c) and the spatial scale of TS is
approximately represented by its scale along the $x$-direction. Panels (e)-(f)
display a similar event in Case 1 where TC is not included.

To study the loop-top annihilation of magnetic islands in detail, we analyze a
typical event that appears at $t\in\left[11,11.65\right]$ in Case 4 (3rd
column of Fig. 1). As the magnetic island collides with the loop-top flux, a
horizontal current sheet forms (see, Fig. 1(c3) and Fig. 3(a)).
Correspondingly, a pair of horizontal jets appear and quickly move along two
opposite directions (Fig. 3(b)). The peak speed of these jets reaches $0.8$,
which is of a similar order as the fastest downflow in the main CS. When these
horizontal jets encounter the edge of the flare loop, they start to move
downward, leading to the increase of downward speed $u_{y}$ (Fig. 3(c)), as
observed by Tian et al. (2014). It is thus clearly revealed that the second
reconnection enables the magnetic island annihilation at the flare loop top.
During the second reconnection of islands, new small-scale magnetic islands
are also generated in the loop-top CSs (see, Fig. 3(a) and (b)), similar to
the “fragmenting coalescence” picture proposed by Bárta et al. (2011).

To estimate the reconnection rate of the horizontal CS, we trace it using a
box with the height of $0.008$ and width of $0.05$ (see, Fig. 3(a)). After
locating all X-points in this region, we use Eq. 5 to estimate
$\mathcal{R}_{lt}$, the reconnection rate in the horizontal CS. It is shown
that, as the magnetic island annihilates, the height of the horizontal CS
keeps decreasing (Fig. 3(d)). While for the reconnection rate within, it first
increases, then reaches its peak value, and finally decreases (Fig. 3(d)).
Obviously, the peak reconnection rate is comparable with that in the main
vertical CS. We also analyze a similar event in Case 1 where TC is not
considered (see, Fig. 3(e)-(h)). It is found that the evolution of horizontal
CS, outjets and the reconnection rate are similar to Case 4.

Figure 4:  Illustrations of the energy conversion and magnetic island
distribution of the second reconnection at the flare loop top. Panels (a),
(b), and (c) show respectively the evolutions of the effective changing rate
of magnetic, kinetic, and internal energy in the main CS and at the flare loop
top. The orange (blue) curves show the results with (without) TC. For each
case, the solid and dashed curves depict the evolution of $\tilde{P}$ in
regions $R_{1}$ and $R_{2}$, respectively. The physical unit of $\tilde{P}$ is
$6.94\times 10^{18}\,\mathrm{erg/\left(cm\cdot s\right)}$. Panel (d) depicts
the CDF of magnetic island flux $N\left(\psi,t\right)$ at loop-top region
$R_{1}$ in Case 4. To obtain $N\left(\psi,t\right)$, we set $30$ bins spaced
uniformly in $\log\left(\psi\right)$ over the range
$\psi\in\left[10^{-8},0.0215\right]$ in units of
$\psi_{0}=10^{11}\,\mathrm{Mx/cm}$.

We further study the energy conversion during the second reconnection at the
loop-top region and compare it with that in the main CS. In an area $S$, the
effective changing rates of magnetic, kinetic, and internal energies are
respectively defined by (Forbes, 1988)

$\displaystyle\tilde{P}_{m}$ $\displaystyle=$
$\displaystyle\frac{\mathrm{d}}{\mathrm{d}t}\int_{S}\frac{B^{2}}{2}\mathrm{d}x\mathrm{d}y+\oint_{\partial
S}\left({\bf E}\times{\bf B}\right)\cdot{\bf n}\mathrm{d}l\,,$ (6)
$\displaystyle\tilde{P}_{k}$ $\displaystyle=$
$\displaystyle\frac{\mathrm{d}}{\mathrm{d}t}\int_{S}\frac{\rho
u^{2}}{2}\mathrm{d}x\mathrm{d}y+\oint_{\partial S}\frac{\rho u^{2}}{2}{\bf
u}\cdot{\bf n}\mathrm{d}l\,,$ (7) $\displaystyle\tilde{P}_{i}$
$\displaystyle=$
$\displaystyle\frac{\mathrm{d}}{\mathrm{d}t}\int_{S}\frac{p}{\gamma-1}\mathrm{d}x\mathrm{d}y+\oint_{\partial
S}\frac{\gamma p}{\gamma-1}{\bf u}\cdot{\bf n}\mathrm{d}l\,,$ (8)

where, $\partial S$ denotes the boundary, ${\bf E}=\eta{\bf J}-{\bf
u}\times{\bf B}$ is the electric field, ${\bf n}$ is the unit normal vector on
the boundary, and $dl$ denotes the line element. The negative value of
$\tilde{P}_{m}$ implies the release of magnetic energy. We calculate the
evolutions of $\tilde{P}_{m}$, $\tilde{P}_{k}$, and $\tilde{P}_{i}$ at the
loop-top region, $R_{1}$, and that in the main CS, $R_{2}$ (see, Fig. 1(b4)).
It’s found that, during most of the time, the magnetic energy release is
dominated by the main CS reconnection, especially during the early stage (Fig.
4(a)). However, the second reconnection at the loop-top region is able to
further release magnetic energy (see, negative spikes of solid curves in Fig.
4(a)). For example, in Case 4, during $t\in\left[11,11.65\right]$, the
effective magnetic energy release rate in $R_{1}$ reaches about half of that
in the main CS (Fig. 4(a)). In $R_{2}$, the evolutions and magnitudes of
$\tilde{P}_{k}$ and $\tilde{P}_{i}$ are similar (see, dashed curves in Figs.
4(b) and (c)), which shows that magnetic energy keeps converting to kinetic
and internal energies. However, in $R_{1}$, the energy conversion process is
relatively complicated. When a magnetic island enters $R_{1}$, its velocity is
damped by the magnetic tension force of the flare loop. As a result, during
the second reconnection, $\tilde{P}_{k}$ is mostly negative (see, solid curves
in Fig. 4(b)), the magnetic energy and kinetic energy significantly increase
the internal energy (see positive peaks of $\tilde{P}_{i}$ in Fig. 4 (c)).
Furthermore, the time-integrated effective changes of magnetic, kinetic, and
internal energies during $t\in\left[5,15\right]$ are calculated by
$\tilde{W}_{m}=\int_{5}^{15}\tilde{P}_{m}\mathrm{d}t$,
$\tilde{W}_{k}=\int_{5}^{15}\tilde{P}_{k}\mathrm{d}t$, and
$\tilde{W}_{i}=\int_{5}^{15}\tilde{P}_{i}\mathrm{d}t$, respectively. For Case
4, $\tilde{W}_{m}$ released in $R_{1}$ is $13.7\%$ of that in $R_{2}$ (see,
Tab. 2), which shows that the total release of magnetic energy during solar
flares is dominated by the reconnection in the main CS region and further
replenished by that at the loop-top region. In $R_{2}$, about half of released
magnetic energy is converted to internal energy and half becomes kinetic
energy. However, in $R_{1}$, the main increment of the internal energy (about
$72.5\%$) is due to magnetic energy release and the rest is from the damping
of kinetic energy.

Region | $R_{1}$ | $R_{2}$
---|---|---
$\tilde{W}_{m}$ | -0.0153 | -0.1115
$\tilde{W}_{k}$ | -0.0091 | 0.0527
$\tilde{W}_{i}$ | 0.0211 | 0.0582
Table 2: Time-integrated increments of magnetic, kinetic, and internal
energies in regions $R_{1}$ and $R_{2}$ for Case 4. The normalizing unit of
$\tilde{W}_{m}$, $\tilde{W}_{k}$, and $\tilde{W}_{i}$ is $W_{0}=7.96\times
10^{20}\,\mathrm{erg/cm}$. Summations of $\tilde{W}_{m}$, $\tilde{W}_{k}$, and
$\tilde{W}_{i}$ in both areas should be zero theoretically but small
deviations may be caused by numerical errors. Figure 5: PDFs of magnetic
island area (a) and flux (b) in the loop-top region $R_{1}$. To obtain the CDF
of area, we set $30$ bins spaced uniformly in $\log\left(S\right)$ over the
range $S\in\left[1,2.1\times 10^{4}\right]$ in units of $\Delta L=\Delta
x\Delta y$. The configuration for calculating CDF of flux is the same as Fig.
4(b).

The magnetic energy cascading at the loop-top region is closely related to the
distribution of corresponding magnetic islands. The loop-top magnetic islands
generated during the second reconnection can be characterized by the evolution
of cumulative distribution function (CDF) of island flux
$N\left(\psi,t\right)$, which measures the number of islands with flux
$\geq\psi$ at moment $t$ (Fig. 4(d)). When the downflow magnetic islands occur
the second reconnection as marked by negative peaks of $\tilde{P}_{R1}$, the
CDF of island flux is significantly enhanced (see duration
$t\in\left[11,11.65\right]$ in Fig. 4(d) as an example). In Fig. 5, we analyze
the probability distribution functions (PDFs) of magnetic island area and flux
at the loop-top region $R_{1}$ based on the same method widely used for
calculating plasmoid PDFs in CS (e.g., Shen et al., 2013; Lynch et al., 2016;
Ye et al., 2019). The PDF of magnetic island flux, $f\left(\psi\right)$, is
determined by
$f\left(\psi\right)=\mathrm{d}N_{\tau}\left(\psi\right)/\mathrm{d}\psi$, where
$N_{\tau}\left(\psi\right)=\int N\left(\psi,t\right)\mathrm{d}t$ and the
integration is taken over duration $t\in\left[5,15\right]$. We also normalize
the PDF to satisfy $\int f\left(\psi\right)\mathrm{d}\psi=1$. The procedure is
similar for calculating $f\left(S\right)$, the PDF of the island area. For
most of the sample domain, namely, $S<10^{3}$ and $\psi<10^{-3}$, both
$f\left(S\right)$ and $f\left(\psi\right)$ follow a power-law independent with
$\beta$ and TC. The slopes of PDFs vary from $-1$ to $-2$ as the island scale
increases. For the large-scale end, the PDFs slightly vary for different
cases, but may be unreliable because the number of large-scale islands we
collected is very limited.

###  3.4 Formation of termination shocks

Figure 6: The histograms of spatial scale and maximum Mach number of TSs at
the loop-top for all simulation cases. As shown in Fig. 3(c) and (g), at each
moment, we recognize the TS at the loop-top by a box. Then we measure the
scale along the $x$-direction, $L_{TS}$, as an approximation of the spatial
scale of TS. The strength of the TS is approximated by the maximum fast mode
Mach number, $M_{TS}$, near the front of TS. Furthermore, we uniformly sample
200 snapshots in the time period of $t\in\left[5,15\right]$ and count TS in
each snapshot. The bin width for panel (A) is $\delta L=0.00175$ and that for
in Panel (b) is $\delta M=0.055$.

As predicted by the standard model, the reconnection downflows may drive a TS
(Forbes & Priest, 1983; Forbes, 1986; Forbes & Malherbe, 1986, 1991; Yokoyama
& Shibata, 1997, 1998, 2001; Shen et al., 2018; Kong et al., 2019, 2020).
Here, we further study the influences of $\beta$ and TC on the formation of
TSs. Because the TS is highly dynamic, in particular, it could be fragmented
as the magnetic islands pass (see, Shen et al., 2018), we mainly focus on the
principal TS, namely the largest and strongest segment of the TS front (see,
Fig. 3(c) and (g)). Without TC, the scale and strength of the principal TSs
both decrease with the increase of $\beta$ (Fig. 6). On the other hand, after
TC is considered, although the histogram of $M_{TS}$ is similar to that
without TC (see the blue and purple curves in Fig. 6(b)), the distribution of
$L_{TS}$ concentrates towards the lower end. It shows that the increase of
$\beta$ and TC can suppress the scale and strength of TSs.

##  4 Summary and Discussion

In this paper, we study the loop-top annihilation of magnetic islands formed
in main vertical CS in detail. Besides the reconnection in main eruption CS,
we find that the reconnection also occurs at the flare loop top, as
illuminated by Bárta et al. (2011). The corresponding reconnection rate is
found to be even of the same magnitude as that in the main CS. Interestingly,
this process is highly similar to the magnetic reconnection taking place in
the earth magnetotail (Lu et al., 2015; Wang et al., 2010, 2016). We thus
suggest that the magnetic island annihilation at the flare loop top, as
revealed here but neglected in the standard flare model, play a non-negligible
role in further releasing magnetic energy so as to heat flare plasma and
accelerate particles.

The fast magnetic reconnection in the main CS is dominated by the generation
of plasmoids, which subsequently move quickly upward and downward along the
CS. This is largely in agreement with previous results (e.g., Bhattacharjee et
al., 2009; Shen et al., 2011, 2018; Kong et al., 2019, 2020) and implies the
existence of fine structures in the CS during solar eruptions. The second
reconnection characterized by horizontal CSs is initiated when magnetic
islands collide with the flux above the flare loops. As the magnetic islands
moving towards the flare loop, the reconnection rate first grows and then
decreases with the peak being comparable with that in the main CS. The second
reconnection also generates smaller-scale secondary magnetic islands which
enhance the turbulence and thus further facilitate the release of magnetic
energy at the flare loop top. Moreover, the whole reconnection process seems
to be intermittent because of the presence of magnetic islands, it is thus
expected that the acceleration of non-thermal particles is intermittent. This
might be used to explain the QPPs observed during flares.

The cascading law of magnetic islands at the loop top seems to be less
affected by background $\beta$ or TC. The PDF of the magnetic island flux
generated during the second reconnection follows a power-law varying from
$f\left(\psi\right)\sim\psi^{-1}$ to $\psi^{-2}$, similar to the PDF of the
magnetic island flux in the main CS obtained by Loureiro et al. (2012), Shen
et al. (2013) and Lynch et al. (2016). However, for the spatial scale and
strength of the TSs formed by the downward outflows or plasmoids, it is found
that both of them decrease with the increase of plasma beta or including TC.

It is worthy of mentioning that simulations of magnetic reconnection depend on
the spatial resolution. The results presented here are obtained under a high-
resolution spatial grid ($3840\times 7680$), which can realize the simulations
with background $S=10^{6}$ as discussed by Ye et al. (2020). The thickness
scales of CSs in both $R_{1}$ and $R_{2}$ are on the order of $0.01$ which is
$0.5\%$ of the $x$-domain (see, Fig. 1 and Fig. 3). Under the current
resolution, the CS thickness can be resolved by about 20 grids. Furthermore,
the scale of internal singular layer in classic tearing mode theory is
$\delta\sim aS_{a}^{-1/4}$ (Pucci et al., 2017), where $a$ denotes the width
of a current sheet, $S_{a}=av_{A}/\eta$ is the locale Lundquist, and $v_{A}$
is the Alfvén speed. The number of grids resolving the inner singular layer
can be estimated by

$\frac{\delta}{\Delta L}\sim\frac{aS_{a}^{-1/4}}{\Delta
L}=\frac{\eta^{1/4}}{\Delta L}\frac{a^{3/4}}{v_{A}^{1/4}}\,,$ (9)

where $\Delta L$ denotes the dimensionless size of mesh grid. Setting $\Delta
L=2/3840$, $\eta=5\times 10^{-6}$, $v_{A}=1$, and $a=0.01$, we have
$\delta/\Delta L\sim 3$. It’s shown that the small-scale inner singular layer
can still be roughly resolved. We have also studied the influence of
resolutions on our conclusions by running cases with different grid
configurations. It is shown that the valid conclusions are only for relatively
high-resolution simulations.

In our numerical model, the fine dynamics of the chromosphere might be unreal,
because the effects of partially ionized plasma (PIP) is not included. It is
true that partial ionization has influences on the properties of shocks, the
CS structure, and the fractal reconnection (Imada et al., 2011; Hillier et
al., 2016). However, in this paper, we mainly focus on the reconnection
dynamics in the high-temperature corona. Thus, the effects of partially
ionized chromosphere on our results are very limited. The only influence is on
the structure and dynamics in the post-flare loops where the chromospheric
plasma is evaporated into.

We would like to thank the anonymous referee for valuable suggestions. This
research is supported by Natural Science Foundation of China grants 11722325,
11733003, 11790303, 11790300, and 11805203.

## References

  * Aurass & Mann (2004) Aurass, H., & Mann, G. 2004, ApJ, 615, 526, doi: 10.1086/424374
  * Bárta et al. (2011) Bárta, M., Büchner, J., Karlický, M., & Skála, J. 2011, ApJ, 737, 24, doi: 10.1088/0004-637x/737/1/24
  * Bhattacharjee et al. (2009) Bhattacharjee, A., Huang, Y.-M., Yang, H., & Rogers, B. 2009, Phys. Plasmas, 16, 112102, doi: 10.1063/1.3264103
  * Cai et al. (2021) Cai, Q., Feng, H., Ye, J., & Shen, C. 2021, ApJ, 912, 79, doi: 10.3847/1538-4357/abee27
  * Cai et al. (2019) Cai, Q., Shen, C., Raymond, J. C., et al. 2019, Monthly Notices of the Royal Astronomical Society, 489, 3183, doi: 10.1093/mnras/stz2167
  * Chen et al. (2015) Chen, B., Bastian, T. S., Shen, C., et al. 2015, Science, 350, 1238, doi: 10.1126/science.aac8467
  * Chen et al. (1999) Chen, P. F., Fang, C., Tang, Y. H., & Ding, M. D. 1999, ApJ, 513, 516, doi: 10.1086/306823
  * Cheng et al. (2018) Cheng, X., Li, Y., Wan, L. F., et al. 2018, ApJ, 866, 64, doi: 10.3847/1538-4357/aadd16
  * Forbes (1986) Forbes, T. G. 1986, ApJ, 305, 553, doi: 10.1086/164268
  * Forbes (1988) —. 1988, Sol. Phys., 117, 97, doi: 10.1007/BF00148576
  * Forbes & Malherbe (1986) Forbes, T. G., & Malherbe, J. M. 1986, ApJ, 302, L67, doi: 10.1086/184639
  * Forbes & Malherbe (1991) —. 1991, Sol. Phys., 135, 361, doi: 10.1007/BF00147508
  * Forbes & Priest (1983) Forbes, T. G., & Priest, E. R. 1983, Sol. Phys., 84, 169, doi: 10.1007/BF00157455
  * Gardiner & Stone (2008) Gardiner, T. A., & Stone, J. M. 2008, Journal of Computational Physics, 227, 4123, doi: https://doi.org/10.1016/j.jcp.2007.12.017
  * Guo et al. (2014) Guo, L. J., Huang, Y. M., Bhattacharjee, A., & Innes, D. E. 2014, ApJ, 796, L29, doi: 10.1088/2041-8205/796/2/l29
  * Hillier et al. (2016) Hillier, A., Takasao, S., & Nakamura, N. 2016, A&A, 591.  https://doi.org/10.1051/0004-6361/201628215
  * Huang & Bhattacharjee (2010) Huang, Y.-M., & Bhattacharjee, A. 2010, Phys. Plasmas, 17, 062104, doi: 10.1063/1.3420208
  * Huang & Bhattacharjee (2012) —. 2012, Phys. Rev. Lett., 109, 265002, doi: 10.1103/PhysRevLett.109.265002
  * Huang & Bhattacharjee (2013) —. 2013, Phys. Plasmas, 20, 055702, doi: 10.1063/1.4802941
  * Huang & Bhattacharjee (2016) —. 2016, ApJ, 818, 20, doi: 10.3847/0004-637x/818/1/20
  * Huang et al. (2011) Huang, Y.-M., Bhattacharjee, A., & Sullivan, B. P. 2011, Phys. Plasmas, 18, 072109, doi: 10.1063/1.3606363
  * Huang et al. (2017) Huang, Y.-M., Comisso, L., & Bhattacharjee, A. 2017, ApJ, 849, 75, doi: 10.3847/1538-4357/aa906d
  * Imada et al. (2011) Imada, S., Murakami, I., Watanabe, T., Hara, H., & Shimizu, T. 2011, ApJ, 742, 70, doi: 10.1088/0004-637x/742/2/70
  * Kong et al. (2020) Kong, X., Guo, F., Shen, C., et al. 2020, ApJ, 905, L16, doi: 10.3847/2041-8213/abcbf5
  * Kong et al. (2019) —. 2019, ApJ, 887, L37, doi: 10.3847/2041-8213/ab5f67
  * Kopp & Pneuman (1976) Kopp, R. A., & Pneuman, G. W. 1976, Sol. Phys., 50, 85, doi: 10.1007/BF00206193
  * Lin & Forbes (2000) Lin, J., & Forbes, T. G. 2000, Journal of Geophysical Research: Space Physics, 105, 2375, doi: https://doi.org/10.1029/1999JA900477
  * Lin et al. (2015) Lin, J., Murphy, N. A., Shen, C., et al. 2015, Space Sci. Rev., 194, 237, doi: 10.1007/s11214-015-0209-0
  * Loureiro et al. (2012) Loureiro, N. F., Samtaney, R., Schekochihin, A. A., & Uzdensky, D. A. 2012, Phys. Plasmas, 19, 042303, doi: 10.1063/1.3703318
  * Loureiro & Uzdensky (2015) Loureiro, N. F., & Uzdensky, D. A. 2015, Plasma Phys. Controlled Fusion, 58, 014021, doi: 10.1088/0741-3335/58/1/014021
  * Lu et al. (2015) Lu, S., Lu, Q., Lin, Y., et al. 2015, Journal of Geophysical Research: Space Physics, 120, 6286, doi: https://doi.org/10.1002/2015JA021213
  * Lynch et al. (2016) Lynch, B. J., Edmondson, J. K., Kazachenko, M. D., & Guidoni, S. E. 2016, ApJ, 826, 43, doi: 10.3847/0004-637x/826/1/43
  * Miyoshi & Kusano (2005) Miyoshi, T., & Kusano, K. 2005, Journal of Computational Physics, 208, 315, doi: https://doi.org/10.1016/j.jcp.2005.02.017
  * Ni et al. (2012) Ni, L., Ziegler, U., Huang, Y.-M., Lin, J., & Mei, Z. 2012, Phys. Plasmas, 19, 072902, doi: 10.1063/1.4736993
  * Parnell et al. (1996) Parnell, C. E., Smith, J. M., Neukirch, T., & Priest, E. R. 1996, Phys. Plasmas, 3, 759, doi: 10.1063/1.871810
  * Priest & Forbes (2000) Priest, E., & Forbes, T. 2000, Magnetic Reconnection: MHD Theory and Applications (Cambridge: Cambridge University Press), doi: DOI: 10.1017/CBO9780511525087
  * Pucci & Velli (2013) Pucci, F., & Velli, M. 2013, ApJ, 780, L19, doi: 10.1088/2041-8205/780/2/l19
  * Pucci et al. (2017) Pucci, F., Velli, M., & Tenerani, A. 2017, ApJ, 845, 25, doi: 10.3847/1538-4357/aa7b82
  * Ruan et al. (2020) Ruan, W., Xia, C., & Keppens, R. 2020, ApJ, 896, 97, doi: 10.3847/1538-4357/ab93db
  * Samtaney et al. (2009) Samtaney, R., Loureiro, N. F., Uzdensky, D. A., Schekochihin, A. A., & Cowley, S. C. 2009, Phys. Rev. Lett., 103, 105004, doi: 10.1103/PhysRevLett.103.105004
  * Shen et al. (2018) Shen, C., Kong, X., Guo, F., Raymond, J. C., & Chen, B. 2018, ApJ, 869, 116, doi: 10.3847/1538-4357/aaeed3
  * Shen et al. (2011) Shen, C., Lin, J., & Murphy, N. A. 2011, ApJ, 737, 14, doi: 10.1088/0004-637X/737/1/14
  * Shen et al. (2013) Shen, C., Lin, J., Murphy, N. A., & Raymond, J. C. 2013, Phys. Plasmas, 20, 072114, doi: 10.1063/1.4816711
  * Stone et al. (2008) Stone, J. M., Gardiner, T. A., Teuben, P., Hawley, J. F., & Simon, J. B. 2008, The Astrophysical Journal Supplement Series, 178, 137, doi: 10.1086/588755
  * Sun et al. (2014) Sun, J. Q., Cheng, X., & Ding, M. D. 2014, ApJ, 786, 73, doi: 10.1088/0004-637x/786/1/73
  * Takasao et al. (2015) Takasao, S., Matsumoto, T., Nakamura, N., & Shibata, K. 2015, ApJ, 805, 135, doi: 10.1088/0004-637x/805/2/135
  * Takasao & Shibata (2016) Takasao, S., & Shibata, K. 2016, ApJ, 823, 150, doi: 10.3847/0004-637x/823/2/150
  * Tian et al. (2014) Tian, H., Li, G., Reeves, K. K., et al. 2014, ApJ, 797, L14, doi: 10.1088/2041-8205/797/2/l14
  * Tsuneta & Naito (1998) Tsuneta, S., & Naito, T. 1998, ApJ, 495, L67, doi: 10.1086/311207
  * Uzdensky et al. (2010) Uzdensky, D. A., Loureiro, N. F., & Schekochihin, A. A. 2010, Phys. Rev. Lett., 105, 235002, doi: 10.1103/PhysRevLett.105.235002
  * Wang et al. (2010) Wang, R., Lu, Q., Du, A., & Wang, S. 2010, Phys. Rev. Lett., 104, 175003, doi: 10.1103/PhysRevLett.104.175003
  * Wang et al. (2016) Wang, R., Lu, Q., Nakamura, R., et al. 2016, Nat. Phys., 12, 263, doi: 10.1038/nphys3578
  * Ye et al. (2020) Ye, J., Cai, Q., Shen, C., et al. 2020, ApJ, 897, 64, doi: 10.3847/1538-4357/ab93b5
  * Ye et al. (2019) Ye, J., Shen, C., Raymond, J. C., Lin, J., & Ziegler, U. 2019, Monthly Notices of the Royal Astronomical Society, 482, 588, doi: 10.1093/mnras/sty2716
  * Yokoyama & Shibata (1997) Yokoyama, T., & Shibata, K. 1997, ApJ, 474, L61, doi: 10.1086/310429
  * Yokoyama & Shibata (1998) —. 1998, ApJ, 494, L113, doi: 10.1086/311174
  * Yokoyama & Shibata (2001) —. 2001, ApJ, 549, 1160, doi: 10.1086/319440
  * Zenitani & Miyoshi (2020) Zenitani, S., & Miyoshi, T. 2020, ApJ, 894, L7, doi: 10.3847/2041-8213/ab8b5d
  * Zhao & Keppens (2020) Zhao, X., & Keppens, R. 2020, ApJ, 898, 90, doi: 10.3847/1538-4357/ab9a31