Can Multi-threaded Flux Tubes in Coronal Arcades Support a Magnetohydrodynamic Avalanche?

Abstract

Magnetohydrodynamic (MHD) instabilities allow energy to be released from stressed magnetic fields, commonly modelled in cylindrical flux tubes linking parallel planes, but, more recently, also in curved arcades containing flux tubes with both footpoints in the same photospheric plane. Uncurved cylindrical flux tubes containing multiple individual threads have been shown to be capable of sustaining an MHD avalanche, whereby a single unstable thread can destabilise many. We examine the properties of multi-threaded coronal loops, wherein each thread is created by photospheric driving in a realistic, curved coronal arcade structure (with both footpoints of each thread in the same plane). We use three-dimensional MHD simulations to study the evolution of single- and multi-threaded coronal loops, which become unstable and reconnect, while varying the driving velocity of individual threads. Experiments containing a single thread destabilise in a manner indicative of an ideal MHD instability and consistent with previous examples in the literature. The introduction of additional threads modifies this picture, with aspects of the model geometry and relative driving speeds of individual threads affecting the ability of any thread to destabilise others. In both single- and multi-threaded cases, continuous driving of the remnants of disrupted threads produces secondary, aperiodic bursts of energetic release.

1 Introduction

Despite incredible advances, in both observations and modelling, our understanding of the exact means through which magnetic energy is released in the solar corona remains limited. Solutions to the coronal heating problem likely rely upon a combination of processes (Parnell and De Moortel, 2012). One of the strongest candidates for dissipating magnetic energy in the solar corona is reconnection (Priest and Forbes, 2000), which releases energy stored in coronal magnetic structures. This idea lies at the heart of the “nanoflare” scenario (Parker, 1988), in which the corona is heated by many small, yet frequent, impulsive releases of magnetic energy in localised reconnection events. Such nanoflares are extremely challenging to probe, or, indeed, even to observe on such small scales, leaving the community reliant on models underpinned by limited observational evidence.

Magnetic fields that permeate the active solar corona often emerge as, or subsequently form, distinct, cylindrical magnetic structures: flux tubes. Twisted flux tubes are commonly associated with bursts of energy released during flares, which attest to their destabilisation and to energy transfer to local heating, waves, and particle acceleration. Such effects produce distinct, observable signatures across the electromagnetic spectrum.

Insight into aspects of the behaviour of flux tubes has been gained through a medley of observational evidence, analytical deduction, and numerical modelling. A common approach in modelling flux tubes is to construct straight magnetic cylinders linking two parallel photospheric planes (representing the solar surface, following Parker, 1972). Ideal MHD kink-mode instabilities, which presuppose a large degree of twist, result in a “kinking” of the axis and destabilisation. Advances in understanding coronal heating using this hypothesis include those of Browning and Van der Linden (2003), Browning et al. (2008), Hood, Browning, and van der Linden (2009), Bareford, Browning, and van der Linden (2010), Bareford, Hood, and Browning (2013). Such models have enabled forward modelling and comparison with observation (Haynes and Arber, 2007; Botha, Arber, and Srivastava, 2012; Pinto, Vilmer, and Brun, 2015; Snow et al., 2017), as well as assessing the role played by certain physical effects, or combinations thereof, known to be present in the solar corona (such as by Botha, Arber, and Hood, 2011; Reale et al., 2016). Straight cylindrical models neglect curvature of flux tubes, which must arc between two photospheric footpoints and whose cross-sectional area may expand from footpoint to apex. Whether, how, and when instabilities may occur can vary: straight cylinders and the semi-tori representing curved coronal loops are merely topologically, and not geometrically, equivalent.

Of particular interest here is the work of Bareford et al. (2016) examining the influence of geometry upon certain properties of a single, highly twisted flux tube. Comparisons of straight cylindrical flux tubes and curved models suggest that curvature tends to reduce the twist necessary for a kink instability, causing current density to concentrate near the apex (in contrast with a more uniform distribution in straight cylindrical tubes).

While models of single flux tubes have continued to evolve, so too have observational techniques. One recent focus has been to probe the substructure within the tubular architecture readily apparent in the active corona. Torsional flows along remarkably twisted atmospheric formations have been reported (Cirtain et al., 2013; De Pontieu et al., 2014). In turn, this has affected modelling, with a recent focus upon loops containing multiple, interacting threads. Such models have, to date, treated flux tubes between two separate, parallel planes (per Parker, 1972). Examples include threads initialised close to a stability threshold and capable of destabilising other threads (Tam et al., 2015; Hussain, Browning, and Hood, 2017), which can thereafter release a great deal of energy through the acceleration of particles (as demonstrated by Threlfall, Hood, and Browning, 2018), and which, more importantly, can trigger a runaway “avalanche” process of successive destabilizations in, and the associated releases of energy from, many threads (Hood et al., 2016). On the other hand, photospheric driving, in addition to causing existing flux tubes to interact and reconnect by relative footpoint motion (O’Hara and De Moortel, 2016), can also create threads that subsequently interact and mutually destabilise (Reid et al., 2018, 2020).

The present work aims to explore how curvature affects the process and properties of energetic release, caused by driven photospheric motions in a multi-threaded model, and, in particular, whether (and, if so, how) an MHD avalanche process can take place in such a three-dimensional (3D) MHD model of an arcade, wherein all the threads start and end in the same plane. Similarly, how far the results of previous models of MHD instabilities and avalanche processes (with straight cylindrical flux tubes) hold for this toroidal geometry remains to be evaluated. Key to the avalanche process is that the magnetic field is able to store excess energy and then quickly release it once a threshold is reached. The rapid release of magnetic energy occurs when an ideal MHD instability is triggered, and the subsequent restructuring of the field creates more current sheets and an avalanche of releases of energy.

Issues to be addressed include how the field itself can affect the onset, evolution, and energetic output of instabilities among many threads. In order to address these and other questions, Section 2 describes our model, Section 3 reports the instability of a single thread in a magnetic arcade, and Section 4 and Section 5 extend the model, displaying results from the simulation of seven threads subject to different initial driving profiles. Section 6 presents an analysis of these results, before Section 7 outlines our conclusions.

2 Model

The MHD equations are solved in a Lagrangian remap scheme, described by Arber et al. (2001). In non-dimensionalised form, Lare3d solves:

$$\begin{aligned} \frac{\partial {\rho }}{\partial {t}} + \nabla \cdot (\rho \,{\boldsymbol {v}}) =& 0 \end{aligned}$$

(1a)

$$\begin{aligned} \rho \frac{\partial {{\boldsymbol {v}}}}{\partial {t}} + \rho {\boldsymbol {v}} \cdot \nabla {\boldsymbol {v}} =& - \nabla p + {\boldsymbol {j}}\times {\boldsymbol {B}} + {\boldsymbol {F}}_{ \mathrm{visc.}} \end{aligned}$$

(1b)

$$\begin{aligned} \frac{\partial {{\boldsymbol {B}}}}{\partial {t}} =& \nabla \times ({\boldsymbol {v}} \times {\boldsymbol {B}}) + \eta \nabla ^{2} {\boldsymbol {B}} \end{aligned}$$

(1c)

$$\begin{aligned} \frac{\partial {P}}{\partial {t}} + {\boldsymbol {v}}\cdot \nabla P =& - \gamma P \,\nabla \cdot {\boldsymbol {v}} + \eta j^{2} + q_{\mathrm{visc.}} \end{aligned}$$

(1d)

for plasma density $\rho$, velocity ${\boldsymbol {v}}$, thermal pressure $P$, magnetic field ${\boldsymbol {B}}$, and associated current density ${{\boldsymbol {j}}}=\nabla \times {\boldsymbol {B}}$. (Vacuum permeability is $\mu _{0}=4\,\pi \times 10^{-7}~\mbox{H}\,\mbox{m}^{-1}$, which, in dimensionless terms, becomes unity, and $\sigma$ is the usual electrical conductivity.) Here, the ratio of specific-heat capacities is $\gamma =\frac{5}{3}$, and the code ensures divergence-free solutions for the magnetic field [$\nabla \cdot {\boldsymbol {B}}=0$].

The dimensionless variables, as evolved in time in Lare3d, are calculated assuming a magnetic field strength $B_{0}$, length $L_{0}$, and mass density $\rho _{0}$. Normalising quantities are here chosen as $B_{0}=10^{-3}~\mbox{T}=10~\mbox{G}$, $L_{0}=10^{7}~\mbox{m}=10~\mbox{Mm}$, and $\rho _{0}=1.67\times 10^{-12}~\mbox{kg}\,\mbox{m}^{-3}$, leading to the typical and normalising values in Table 1. Henceforth, lengths are quoted normalised with reference to $L_{0}$, and times with reference to Alfvén times: $\tau _{\mathrm{A}}$.

Table 1 Normalising quantities. A specific physical quantity listed in the first column is denoted by the symbol shown in the second and normalised with reference to the value in the third column.

Our computational domain is $x \in \left [-x_{\mathrm{max.}},x_{\mathrm{max.}}\right ], y \in \left [-y_{\mathrm{max.}},y_{\mathrm{max.}}\right ], z \in \left [0,z_{ \mathrm{max.}}\right ]$. Boundary conditions prescribed in $x$ are periodic (representing neighbouring arcades, which are common in active regions). Those in $y$ and $z$ are static, perfectly conducting, and with zero normal derivatives. Identical conditions exist on the lower $z$-boundary (our photosphere) for all variables except velocity, in which the driver is imposed according to Section 2.3. Our results and conclusions do not depend on the actual choice of boundary conditions (other than the form of the driver).

2.1 Dissipation

Shock viscosities are used to capture discontinuities and other steep gradients, which otherwise would not be resolved on the numerical grid. Such viscosities contribute to the viscous dissipation in the system; further details have been given by Reid (2020). Viscosity enters as a force term [${\boldsymbol {F}}_{\mathrm{visc.}}$] in Equation 1b and a heating term [$q_{\mathrm{visc.}}$] in Equation 1d.

We apply an anomalous magnetic diffusivity [$\eta _{a}$]; this is active where some measure of current exceeds a critical threshold. Such an approach reflects the enhanced resistivity (above classical values) believed to occur in extreme conditions in the corona, such as in very strong current layers.

Slow photospheric motions create each thread. Over time, these motions also generate internal currents within each thread. Accordingly, the magnetic energy of the loop system will increase above the initial potential value, and it is this “free” energy that can be tapped when an instability is eventually triggered. (Theoretically, all energy above that in the potential field is available to be dissipated, but only the energy above that in the linear, force-free, constant-$\alpha$ field of equal helicity is readily accessible.) At the onset of instability, a thin, laminar current sheet rapidly forms along part of the thread; such current is generally faster growing, and ultimately far stronger, than the internal current of the thread. Typically, in previous investigations (as in that of Tam et al., 2015), a threshold is chosen based upon the magnitude of the current [$\left |{{\boldsymbol {j}}}\right |$]. This threshold would activate anomalous resistivity, having been breached by the current sheet of the instability, but not by internal currents. Two aspects of the present model warrant a different approach.

Firstly, the curved initial magnetic field, outlined in Section 2.2, whose strength diminishes with height, naturally gives rise to significant currents near the footpoints of any thread. Twisting creates a toroidal current; in this geometry that component is predominantly $j_{z}$ near the boundaries. Secondly, an additional problem for cases with more than one thread is variation of field strength with height, causing threads closer to the lower boundary to exhibit stronger magnetic field. Such threads possess, in general, stronger internal currents, compared with outer (higher) threads. Therefore, the criterion for enhanced resistivity must be triggered by specific current sheets and account for both internal thread currents and variation of $\left |{\boldsymbol {B}}\right |$ with height in the model. To address these issues, we define a variable

$$ \zeta = \frac{\sqrt{j_{x}^{2}+j_{y}^{2}}}{\left |{\boldsymbol {B}}\right |} $$

(2)

and a threshold $\zeta _{\mathrm{crit.}}$. In this way, magnetic diffusivity $\eta$ is $\eta _{\mathrm{a}}=10^{-4}$ where $\zeta \ge \zeta _{\mathrm{crit.}}$, otherwise $\eta =0$. This allows strong current sheets related to an instability to be targeted with anomalous resistivity (while ignoring strong currents associated with twisting in low-lying regions of greater field strength). Although resistivity depends upon $\zeta$ (rather than the more usual $\left |{{\boldsymbol {j}}}\right |$), the effects of resistivity upon the evolution of the magnetic field, in Equation 1c, and in the energy equation, Equation 1d, remain unchanged. Ohmic heating, in particular, remains the familiar $\eta j^{2}$ term in the energy equation.

2.2 Background Magnetic Field

An arcade-like background magnetic field could be constructed in several ways. In the light of our intention to drive the field at the base, two features in particular are desirable: firstly, a vertical photospheric field that reverses sign at the centre of the domain and, secondly, sufficiently large regions of (near-)uniform field strength on either side of a polarity inversion line (PIL). To that end, we modify a hyperbolic trigonometric field structure (as used by Howson, De Moortel, and Fyfe, 2020), adding further Fourier modes to create wide regions of near-uniform photospheric $B_{z}$ at the edges of the arcade, while smoothly reversing sign at the domain centre. The general form of the field is:

$$\begin{aligned} B_{y}\left (y,z\right ) = \sum _{j=1}^{N}{ a_{j} \frac{\cosh {\left (jk\left (z-z_{\mathrm{max.}}\right )\right )}}{\sinh {\left (-jkz_{\mathrm{max.}}\right )}} \cos {\left (jky\right )} } \end{aligned}$$

(3a)

$$\begin{aligned} B_{z}\left (y,z\right ) = \sum _{j=1}^{N}{ a_{j} \frac{\sinh {\left (jk\left (z-z_{\mathrm{max.}}\right )\right )}}{\sinh {\left (-jkz_{\mathrm{max.}}\right )}} \sin {\left (jky\right )} } \end{aligned}$$

(3b)

in which $k=\frac{\pi }{2y_{\mathrm{max.}}}$. Coefficients $a_{j}$ are constructed for the several Fourier modes as:

$$ a_{j} = \left \{ \textstyle\begin{array}{c@{\quad }c} \frac{3\,\pi +8}{18\,\pi } & j=3 \\ \frac{1}{9\,\pi } & j=6 \\ 18 \frac{\left (64c^{6}+32c^{5}-96c^{4}-40c^{3}+36c^{2}-2\right )\left (j^{2}-9\right )+9c\left (j^{2}-6\right )}{j\pi \left (j^{2}-36\right )\left (j^{2}-9\right )} & j\ne 3,6 \end{array}\displaystyle \right . $$

with $c=\cos {\left (\frac{j\pi }{6}\right )}$. Preliminary tests suggested that fixing $N=20$ yields a sufficiently wide region of approximately uniform magnetic field strength, which could contain several threads. Figure 1 illustrates the magnetic field profile seen in Equations 3a – 3b, comparing cases where $N=1$ and $N=20$. The size of the regions of near-uniform magnetic field strength also increases with the number of Fourier modes.

Figure 1

Components of the modelled magnetic field. Vertical ($B_{z}$, seen in blue) and horizontal ($B_{y}$, in red) field components, at $z=0, t=0$, for a model incorporating one Fourier mode or twenty modes in Equations 3a – 3b (shown with dashed or solid lines, respectively); $x$ is the invariant direction in the model.

In tests, simulations with $N=1$ resulted in excessively strong currents at the footpoints of each thread in the domain, which grew rapidly over time in response to an imposed rotation. Variation of magnetic field strength over the radius was responsible; these difficulties are avoided if $N\geq 20$.

2.3 Photospheric Driver

Within the magnetic arcade, a vortical driving motion is applied at the footpoints of the thread(s) to be created, twisting the magnetic field. In order to aid comparison with previous cases of photospherically driven flux tubes, the form of the driver matches that imposed by Reid et al. (2018):

$$\begin{aligned} v_{\phi }\left (r,t\right ) =& v_{0}f\left (r\right )D\left (t\right ) \end{aligned}$$

(4a)

$$\begin{aligned} f\left (r\right ) =&\left \{ \textstyle\begin{array}{l@{\quad }l} \frac{r}{a}\left (1-\frac{r^{2}}{a^{2}}\right )^{3} & \qquad \qquad \qquad r\le a \\ 0 & \qquad \qquad \qquad r>a \end{array}\displaystyle \right . \end{aligned}$$

(4b)

$$\begin{aligned} D\left (t\right ) =&\left \{ \textstyle\begin{array}{l@{\quad }l} 0 & \qquad t \leq t_{\mathrm{s}} \\ \frac{1}{2}\left (1-\cos {\left ( \frac{(\pi \left (t-t_{\mathrm{s}}\right )}{t_{\mathrm{e}}-t_{\mathrm{s}}} \right )}\right ) & \qquad t_{\mathrm{s}}< t< t_{\mathrm{e}} \\ 1 & \qquad t \geq t_{\mathrm{e}} \end{array}\displaystyle \right . \end{aligned}$$

(4c)

Driving is imposed over a radius $a$, the minor toroidal radius of the flux tube. Velocity, with amplitude governed by $v_{0}$, is gradually introduced through a function $D\left (t\right )$. In the straight cylindrical case, the spatial form of the driver produced zero net axial current. However, the fact that the axial field here is not completely uniform at the photospheric boundary means that the net current is small but non-zero. While Reid et al. (2018) use a driver beginning at $t=0$, we postpone the start of rotation in order that the magnetic arcade relaxes towards potential. The ramp-up phase begins at time $t_{\mathrm{s}}$ and ends at $t_{\mathrm{e}}$, after which $v_{\phi }\left (r,t\right )$ remains constant. Initial tests have determined that a delay of $50\,\tau _{\mathrm{A}}$ allows the model to reach sufficiently close to minimal energy (hence $t_{\mathrm{s}}=50\,\tau _{\mathrm{A}}$), and that a duration extent of $10\,\tau _{\mathrm{A}}$ is sufficient to avoid unnecessary shocks (and hence $t_{\mathrm{e}}=60\,\tau _{\mathrm{A}}$).

In Section 2.2, we noted that our model aimed to fit many threads in the near-uniform field region; the size of this region places constraints on the radius $a$. The radius must be large enough to resolve the formation of associated current sheets, while being small enough to fit many threads in the near-uniform field region. Similarly, our choice of velocity amplitude $v_{0}$ is also constrained by the need to be faster than slow numerical diffusion at the base, yet slower than the coronal Alfvén speed. This issue of timescale separation is addressed by Bowness, Hood, and Parnell (2013).

Previous investigations have shown that the relative amplitude of the photospheric driver in different threads plays a key role in determining the sequence in which straight cylindrical multi-threaded flux tubes destabilise (Reid et al., 2018). We aim to extend this to consider a set of a toroidal loops, and to study the impact of threads of different length, curvature, and driving speed.

2.4 Instability

The onset of an MHD avalanche, in a flux tube containing many threads, requires a trigger. An exponential rise in the kinetic energy accompanies the first disruption (seen later in, for example, Figure 3). This rise is due to an ideal MHD instability, which begins the avalanche process. It is possible that the ideal kink instability, widely used in straight cylindrical cases (e.g. Tam et al., 2015; Hussain, Browning, and Hood, 2017), is responsible. A sufficiently twisted thread will “kink” (i.e. displace radially outward), with this mode characterised by the rapid development of a helical current sheet surrounding the unstable thread. The twist [$\Phi$] along field lines necessary for this instability is a critical parameter. However, critical thresholds (and, indeed, how one quantifies twist, discussed, for example, by Threlfall, Hood, and Priest, 2018; Threlfall, Wright, and Hood, 2020) are subject to considerable uncertainty, and may deviate significantly from theoretical values, which are usually specific to particular field configurations; for example, $\Phi \gtrsim 3.3\,\pi$ for uniformly twisted tubes (Hood and Priest, 1979). For a loop with radially varying speed, Gerrard, Hood, and Brown (2004) determine an average twist after a given time [$t$]. When applied to our driving velocity, the same approach gives the average twist at time $t$ here:

$$ \langle \Phi \rangle =\frac{v_{0}}{2a}\left [t-\frac{1}{2}\left (t_{ \mathrm{s}}+t_{\mathrm{e}}\right )\right ] $$

(5)

(for $t>t_{\mathrm{e}}$; for our profile, twist peaks on the axis at four times this value).

Our introduction of toroidal geometry into our multi-threaded arcade makes the torus instability (Bateman, 1978) a plausible trigger mechanism for the avalanche, although our flux tubes do not establish significant non-zero net current in each thread. The torus instability is analogous to the kink instability, but it occurs in a different set of conditions, including a large and negative rate of change of magnetic field along the major toroidal axis [$R$]. Such a rate of change is commonly quantified by the decay index $n=-{\mathrm{d}\log {\left |{\boldsymbol {B}}_{\mathrm{ex.}}\right |}}/{ \mathrm{d}\log {R}}$, where ${\boldsymbol {B}}_{\mathrm{ex.}}$ denotes the external magnetic field. A field is unstable to the torus instability where $n>n_{\mathrm{crit.}}$. The critical threshold $n_{\mathrm{crit.}}$ varies depending on the configuration. For a purely poloidal external field, Bateman (1978) derives $n_{\mathrm{crit.}} = 1.5$, but several effects, such as expansion, can lower such values, or raise them as high as $n_{\mathrm{crit.}}=2$ (Titov and Démoulin, 1999; Kliem and Török, 2006; Aulanier et al., 2010; Zuccarello, Aulanier, and Gilchrist, 2015). However, the original derivation of the torus instability featured a non-zero net toroidal current, surrounded by a purely poloidal magnetic field. This is not the case here and the condition for instability will be modified. Thus, the value and importance of $n$ is subject to considerable uncertainty, and so we consider a plausible range of critical $n$ (after Syntelis, Archontis, and Tsinganos, 2017, who also consider a flux tube with zero total axial current in the initial state).

We will provide diagnostic information (i.e. the twist and decay index at each thread destabilisation) for additional context. While it is interesting to determine if the critical conditions at the onset of the ideal MHD instability correspond to either a kink or torus mode, the important point, for the avalanche, is that there is an ideal MHD instability. This instability must grow on an Alfvén timescale.

3 Single-Threaded Case

Before considering the behaviour of a loop containing many threads in the arcade field, described in Section 2.2, we put this in context by first considering the behaviour of a single, continuously driven thread. Within a domain with $x_{\mathrm{max.}}=y_{\mathrm{max.}}=1.0,z_{\mathrm{max.}}=2.0$, the prescribed form of velocity, Equations 4a – 4c, is implemented on both sides of the PIL (found at $y=0$). The velocity profile is centred on $\left (x,y\right )=\left (0,\pm 0.65\right )$, with a radius $a=0.2$ and with $v_{0}=0.008$. The general 3D structure of the system can be seen in Figure 2. Blue field lines in Figure 2a highlight the part of the domain that will be continuously driven, subject to the driving profile seen in Figure 2b. The driven magnetic field region becomes increasingly twisted as time progresses, until such time as an ideal instability forms a current sheet in the flux tube. The volume-integrated magnetic, kinetic, and internal energies and the instantaneous components of heating are seen in Figure 3. As the bottom boundary is constantly driven, Poynting flux is continuously being injected into the system. After the avalanche process has started, the accumulated heating (consisting of Ohmic and viscous heating) eventually settles down to around $35\%$ of the time-integrated injected Poynting flux.

Figure 2

Single thread: field and driver. (a) Field structure at $t=0$, showing field lines traced throughout the domain (red) and specifically chosen field lines traced from the driven region (blue), above a magnetogram showing vertical field strength at the base ($B_{z}\left (z=0\right )$). (b) Contours of the local vertical field strength at the base of the simulation, with planar velocity vectors ($v_{x}$, $v_{y}$) overlaid, illustrating photospheric driving at the base at full speed ($t=t_{\mathrm{e}}$).

Figure 3

Single-threaded case, energies: Magnetic (solid, black), internal (dash–dotted, green), and kinetic (dotted, blue) energy components with their initial values subtracted, shown together with Ohmic (thick, red) and viscous (dashed, red) instantaneous heating in the system. The two dashed, vertical lines indicate times at which strong current sheets form, illustrated in contours of current in Figure 5, and hence there are large releases of energy.

Components of energy and heating neatly describe distinct phases of system behaviour. Prior to photospheric driving (at $t=50\,\tau _{\mathrm{A}}$), there is a slight decrease in magnetic energy and a sharper reduction in kinetic energy as the system relaxes towards minimum energy. At the onset of driving, all energies increase sharply; while kinetic energy appears to reach a quasi-steady state faster, magnetic and internal energies continue to increase still further, and then gradually begin levelling off. Viscous heating also follows the start of driving, but it remains very small. The magnetic field is gradually twisted, while internal currents associated with the thread form (as seen in Figure 5a).

The second phase of the experiment begins when a current sheet forms along the twisted thread at approximately $t=330\,\tau _{\mathrm{A}}$. Current in the sheet eventually satisfies the condition $\zeta >\zeta _{\mathrm{crit.}}$, triggering anomalous resistive effects. An example of the 3D structure of such a current sheet can be seen in Figure 4a; locations which satisfy $\zeta >\zeta _{\mathrm{crit.}}$ are difficult to see in contour images (Figure 5). Resistivity, where activated, causes Ohmic heating.

Figure 4

Single-threaded case, 3D configuration: (a) illustrates the 3D configuration of the magnetic field during the initial energy release event (at $t=330\,\tau _{\mathrm{A}}$), while (b) illustrates a secondary (less energetic) eruption ($t=430\,\tau _{\mathrm{A}}$). Ambient field is shown in red, while field lines traced from the photospheric driving region at the base are shown in blue, with green isosurfaces of current where $\zeta \geq 0.75\,\zeta _{\mathrm{crit.}}$ (lowered below $\zeta _{\mathrm{crit.}}$ for illustrative purposes).

Figure 5

Contours of toroidal current [$j_{y}$] above the polarity inversion line [$y=0$], in a single-threaded loop.

Following this event, the final phase of behaviour begins, as the tube remnants undergo restructuring and gradual expansion, while the magnetic field continues to be driven at the base. This phase is characterised by additional, aperiodic formation of small, secondary current sheets, scattered throughout the tube remnants. Fast Fourier Transforms (FFT) of the energy and heating components (not shown) detect no distinct frequency or period, hence our description of the formation of these current sheets as “aperiodic”. Examples of such current sheets are present at $t=430\,\tau _{\mathrm{A}}$, and are visible in Figure 4b as thin 3D isosurfaces, or as highly concentrated regions of current in Figure 5c throughout the loop. The heating from such current sheets is weaker than, but comparable to, that from the first. We end the numerical experiment at $t=500\,\tau _{\mathrm{A}}$, however, there is nothing to suggest that this final phase (of aperiodic formation and dissipation of current sheets) would not continue indefinitely.

Regarding the onset of the initial instability, the thread has an average twist of $1.75\,\pi$ (and maximum $7\,\pi$) at $t=330\,\tau _{\mathrm{A}}$. Similarly, we recover a decay index $n\approx 1.57$ on the axis at the apex; some variation occurs with radius, but $n$ is around, or just below, 1.5 across most of the cross-section.

During the first major instability, the flux tube is more twisted than in later disruptions (Figure 4b). Much of the initial twist dissipates in the first reconnection event: the field is more relaxed, and heating arises from the continual conversion of injected magnetic energy to thermal/internal energy. The Poynting flux injected through the driven boundary dissipates at a largely constant rate. No steady state in energy is achieved: energy releases are aperiodic and vary in size and location within the loop remnants. Later energy releases are far smaller than, but similar to, the first. These aspects are readily apparent in the energies displayed in Figure 3.

4 Seven Threads: Identical Drivers

Our investigation of a single thread puts into context the behaviour observed when additional threads are included in our arcade model. We now consider seven threads, all created by photospheric driving as before. The footpoints are hexagonally packed, with a row of three threads between two rows of two threads; the general arrangement can be seen accompanying our later discussion (in Figure 11a). The central thread is centred on $\left (x,y\right )=\left (0,\pm 1.3\right )$. To accommodate these additional threads, the length of the simulation domain is doubled in $y$ (i.e. $y_{\mathrm{max.}}=2, y\in \left [-2,2\right ]$). All thread footpoints are rotated with $v_{0}=0.02$ over a radius $a=0.1$ (i.e. half that used in the single-threaded case).

The 3D magnetic configuration and its evolution are illustrated in Figure 6; all seven threads are initially untwisted (Figure 6a), but they also travel higher into the domain than in the single-threaded case. All threads are twisted at the same rate, but a large reconnection event (at the time shown in Figure 6b) causes the lowest row of threads to destabilise. The times of disruption of each thread, and values of twist and decay index at those times, are recorded in Table 2, for comparison with later experiments. The destabilisation and fragmentation of the bottom row are followed very soon afterwards by the same in the middle row, with a short pause before the uppermost row of threads is disrupted (Figure 6c). Another window on the evolution is provided by the energy and heating components of the system shown in Figure 7, coupled with contours of the current structures present at the apex of each of the threads (Figure 8).

Figure 6

Seven threads, identical drivers: 3D magnetic configuration at different times during the experiment, where specific field lines have been coloured corresponding to their initial position in one of seven regions of photospheric driving, and green isosurfaces (if present) identify regions where $\zeta >0.75\,\zeta _{\mathrm{crit.}}$ (lowered from $\zeta _{\mathrm{crit.}}$ for illustrative purposes). For key to colours and destabilisation times, see Table 2.

Figure 7

Energies of seven threads, with identical driving speeds: Magnetic (solid, black), internal (dash–dotted, green), and kinetic (dotted, blue) energy components with their initial values subtracted, shown together with Ohmic (thick, red) and viscous (dashed, red) instantaneous heating of the system. The six dashed, vertical lines indicate times associated with mergers of threads, illustrated as contours of current in Figure 8.

Figure 8

Seven threads, identical drivers: Contours of toroidal current [$j_{y}$] above the polarity inversion line [$y=0$], at various times throughout the simulation (the key to contour levels is seen in colour bar).

Table 2 Seven threads: times of disruption, average twist $\langle \Phi \rangle$, and decay index $n$ for each thread in two experiments: one with all threads driven at same rate (Section 4) and one where the central thread is driven three times faster than surrounding threads (Section 5). Each thread is identified by a number and the assigned colour in 3D images (Figure 6) or field-line end-point locations (Figure 11). Times of disruption are discerned by eye from the contours of current in the $x$–$z$-plane at $y=0$. (As a consequence of the limited output cadence of simulation data, times of disruption are limited to specific discrete values. Other diagnostic information (for example the current in Figures 8 and 10) is only available at this lower cadence. Disruption times could alternatively be identified using (higher cadence) energy data, shown in Figures 7 and 9, distinguished by an exponential rise in the kinetic energy, associated with each instability. However, such rises are difficult to identify for subsequent disruptions in the overall energies of the system.) Twist is determined from field lines twisting helically around the axis of each thread, at the times of disruption. Decay indices quoted are evaluated at the intersection of the axis of each thread with the plane $y=0$ (following the practice of Zuccarello, Aulanier, and Gilchrist, 2015).

As with the single-threaded case described in Section 3, the system initially tends towards a potential state before the rotation commences. All three components of energy rapidly increase following the onset of the driver, quickly forming seven threads. The twist in field lines reaches the apex of the arcade before $t=140\,\tau _{\mathrm{A}}$: seven distinct current structures are visible in Figure 8a. The lower pair of threads become destabilised at/before $t=150\,\tau _{\mathrm{A}}$ (seen in Figure 8b); this destabilisation occurs much earlier than in the single-threaded case. Very shortly after this, threads in the middle row are all disrupted together (Figure 8c), followed a short time later by the uppermost row of threads ($t=167.5 \,\tau _{\mathrm{A}}$; Figure 8d). The time taken for all three rows of threads to destabilise is less than $20\,\tau _{\mathrm{A}}$ (compared with the complete experiment duration of $500\,\tau _{\mathrm{A}}$). This concentrated and related series of disruptions is responsible for the largest bursts of Ohmic heating seen in Figure 7 at $t=150\,\mbox{--}\,170\,\tau _{\mathrm{A}}$. However, these initial bursts significantly raise the internal energy of the system, to the extent that the change in internal energy becomes much larger than the change in magnetic energy. One should note, however, that these trends are found upon subtracting the initial energy values; the plasma-$\beta$ is not materially affected. This contrasts with Section 3, where changes in magnetic and internal energy remained closely matched, but with a slightly larger change in the magnetic component throughout the single-threaded case.

Remnants of all seven threads now sporadically create thin current sheets, scattered throughout the volume and prone to reconnect; in particular, many of these later current sheets are associated with the lower and middle rows of threads. At much later times (such as $t=380\,\tau _{\mathrm{A}}$ in Figure 8f), the residue of current from all threads has expanded significantly, but is still pervaded by small current sheets, which occasionally, and irregularly, trigger anomalous resistivity. Notably, the ongoing motions cause some threads apparently (and partially) to reform, as can be seen towards the bottom left of Figure 8f.

Regarding the identification of a specific instability responsible for each disruption, Figure 8 appears to suggest thread mergers (and not specific instabilities in individual threads) are responsible for many of the larger reconnection events and associated energy releases in this model. We will further analyse the times, twist, and decay indices recorded in Table 2 in Section 6, following a second experiment.

As in the single-threaded case, photospheric driving introduces a continuous Poynting flux from below and, in this case, marginally more of the injected energy ($\approx 58\%$) is ultimately converted to heat. Differences in number of threads, driving speeds, and footpoint areas between the cases further complicate a comparison of Poynting-flux conversion between these cases.

5 Seven Threads: Faster Central Thread

Our final investigation considers the effect of driving speed in individual threads; by enhancing the driving speed of the central thread, we focus on the role (if any) of variations in speed in the cascade process. In the same configuration as in Section 4, we repeat the simulation with the central thread twisted three times faster ($v_{0}=0.03$) than the others ($v_{0}=0.01$). Other parameters, including the relative positions of threads and resistivity conditions, remain unchanged. The 3D evolution of the system is similar to that seen in Section 4; we will describe the evolution of this experiment using energetics and heating (in Figure 9) and through contours of current through the apex of the multi-threaded loop (in Figure 10).

Figure 9

Energies of seven threads, with faster central thread: Single-threaded case, energies: Magnetic (solid, black), internal (dash–dotted, green), and kinetic (dotted, blue) energy components with their initial values subtracted, shown together with Ohmic (thick, red) and viscous (dashed, red) instantaneous heating of the system. The six dashed, vertical lines indicate times of interest examined using contours of current in Figure 10.

Figure 10

Seven threads, central thread fastest: Contours of toroidal current [$j_{y}$,] above the polarity inversion line [$y=0$] at various times throughout the simulation (the key to contour levels is seen in colour bar).

The early evolution of the system, before $t\approx 100\,\tau _{\mathrm{A}}$, largely matches that described in Section 4 for threads driven at the same speed; a sufficiently potential state is realised before photospheric drivers create each thread. Differences between multi-threaded experiments emerge soon after; current associated with the central thread builds up much faster than in surrounding threads (as shown, e.g., in Figure 10a). The central thread is the first to destabilise (Figure 10b). This destabilisation is accompanied by a burst of Ohmic heating, enhanced internal energy, and decreased magnetic energy ($t\approx 130\,\tau _{\mathrm{A}}$, Figure 9). Two more energy releases soon follow, as the remaining two threads in the middle row merge with the remnants of the central thread (Figure 10c). Once the middle row is entirely disrupted, there is a short pause ($\approx 20\,\tau _{\mathrm{A}}$) in which little energy is dissipated or heating is seen. At $t\approx 160\,\tau _{\mathrm{A}}$, the remnants of the central row of threads begin to disrupt the lower row (Figure 10d). Two corresponding bursts of Ohmic heating are seen in Figure 9 at approximately $t=165\,\tau _{\mathrm{A}}$ and $175\,\tau _{\mathrm{A}}$, respectively. The two uppermost threads remain largely unaffected until $\approx 220\,\tau _{\mathrm{A}}$, when they also begin to merge with the central row of threads (Figure 10e), leading to yet another burst of Ohmic heating. Between the disruption of this final pair of threads and the end of the simulation, further sporadic bursts of Ohmic heating continue to occur. As before, these are driven by the continual photospheric driving, forming thin, fragmented current sheets throughout the expanding remnants of the threads. Unlike in Section 4, these bursts appear more concentrated in time, with long intervals between larger bursts of Ohmic heating. Although often insufficiently large (and hence unable to trigger anomalous resistive effects), several strong, thin, fragmented current sheets are clearly visible in the final state of the experiment, seen in Figure 10f.

Compared with the previous multi-threaded case (Section 4), the variation in driving speeds generates significant differences in the magnitude of Poynting flux injected. Despite further differences in the avalanche process, a very similar overall proportion ($\approx 54\%$) of Poynting flux is ultimately converted into heat. Once more, we record the approximate times of disruption of each thread, the associated mean twist, and the relevant values of the decay index, in Table 2.

6 Analysis

Photospheric boundary motions, acting upon a curved magnetic arcade, have been shown to be capable of creating and supporting both single- and multi-threaded flux tubes. The results of Sections 4 and 5 clearly illustrate that disruption in one thread may trigger an avalanche-like process, wherein all the threads in a multi-threaded tube are destabilised. Our discussion will focus on a comparison of our multi-threaded cases and the nature of the initial instability.

One issue that warrants clarification is the relative timing of instability. In the single-threaded experiment, the initial instability occurs much later ($t=330\,\tau _{\mathrm{A}}$) than in either multi-threaded case ($t=150\,\tau _{\mathrm{A}}$ in Section 4, $t=127.5\,\tau _{\mathrm{A}}$ in Section 5). This is the result of different driving speeds: $v_{0}$ is much smaller in the single-threaded case ($v_{0}=0.008$) than in the multi-threaded cases (in Section 4, $v_{0}=0.02$; in Section 5, the central thread has $v_{0}=0.03$, the outer threads $v_{0}=0.01$). As mentioned in the introduction, and by Bowness, Hood, and Parnell (2013), successful MHD simulations require compressed timescales; while the slower speed here suffices for the first case, greater driving speeds were necessary to model the more complex, multi-threaded behaviour. Faster speeds relative to the single-threaded case, by factors of 2.5 and 3.75, respectively, bring forward the time of the initial instability proportionately. Hence, it occurs approximately three or four times sooner (after the start of driving) in each multi-threaded case than in the single-threaded case (respectively, after $100\,\tau _{\mathrm{A}}$ or $75\,\tau _{\mathrm{A}}$, compared with $280\,\tau _{\mathrm{A}}$).

Comparing multi-threaded cases, both configurations ultimately undergo a very similar process and achieve end states that contain common traits. One visual representation of this is the distribution of end positions of field lines starting within specific threads seen in Figure 11. The field lines associated with each thread become mixed with those from neighbouring threads by the end of the experiment in Figure 11. The specific pattern in each experiment differs; this is to be expected when the threads destabilise and interact in a different order. However, the remnants of all threads have expanded significantly, coalesced in the lanes created between photospheric driving regions, and ultimately fill much more of the volume of the simulation than when the driver was initiated. Similar mixing is observed in straight cylindrical cases (e.g. Figure 7 of Reid et al., 2018); with the arcade geometry, there is comparable, but perhaps slightly less widespread, entangling of flux. However, seven threads in three rows may be insufficient to draw out any such influence: experiments containing more threads, or more rows of threads, may bring more evidence of geometrical effects.

Figure 11

End-points of specific field lines within seven-threaded loop experiments. Field lines are traced from footpoints at positive polarity to those at negative, coloured according to their initial position in a specific driving region or thread number (described in Table 2) and overlaid onto the vortical driving profile at the base of the domain at the (a) beginning and (b) end of the uniform driven amplitude experiment, with (c) showing the final state of the faster central thread experiment.

Our experiments recover a tendency of entire rows of threads (at the same height) to destabilise simultaneously. Furthermore, the expansion of the field above the PIL also causes a significant gap between the uppermost and middle rows. Separation attenuates the upward spread of instability and destabilisation. The arcade geometry of our model and variations in field strength with height lead to a preferential upward expansion near the apex of the loop. These are clear effects imparted on the system evolution by our arcade geometry.

The preferential destabilisation of specific rows (with height above the PIL) may be influenced by other factors. Varying driving speed allows for more direct comparisons with previous multi-threaded models (e.g. Hood et al., 2016), wherein a single thread typically triggers further disruptions. For a central thread driven three times faster than surrounding threads, the resulting avalanche proceeds differently from the case where all threads are driven at the same rate. In Section 5, we recover a very similar picture to other multi-threaded models (e.g. those of Reid et al., 2020), reporting that individual threads destabilise their neighbours in turn. Despite some interaction between the middle and top rows, the top row destabilises largely as in the previous case. The driving speed of threads on the top row is half that seen in Section 4, consistent with the merger of the top row occurring twice as long after the twisting starts. Thus, it is possible to obtain at least a partial cascade of thread eruptions, but the avalanche process can be attenuated (or, indeed, suppressed entirely) by the field geometry.

The evidence for the nature of the initial instability (found in Table 2) is mixed. In the single-threaded case, the average twist is $\langle \Phi \rangle =1.68\,\pi$ when the thread destabilises. This value is lower than in comparable investigations; our value is smaller than in the equivalent straight cylindrical case (Reid et al., 2018), and is nearly half as small as threshold values in models containing curvature and expansion ($\Phi \approx 3.26\,\pi$ in Gerrard, Hood, and Brown, 2004) or line-tied loops ($\Phi \approx 3.3\,\pi$, determined by Hood and Priest, 1979). Small values of twist at instability are less common, but not unheard of (e.g. Rappazzo et al., 2019, finding a mean twist $\approx \frac {\pi }{2}$). Bareford et al. (2016) assert that expansion and curvature of flux tubes significantly accelerate the onset of instability. Our model contains more curvature and expansion than previous models, in order that thread footpoints are anchored in near-vertical field regions of the photosphere. We also note that average twist is an ensemble measure. Considering the peak twist, some field lines have twist around $7\,\pi$ at the onset of instability: this value far exceeds equivalent threshold values in the literature. Turning to our multi-threaded configurations, in Section 4 the lowest row of threads destabilise at $\langle \Phi \rangle =2.67\,\mbox{--}\,2.74\,\pi$, while the central thread in Section 5 achieves $\langle \Phi \rangle =4.12\,\pi$ at the point of instability. These values are comparable with threshold values of Hood and Priest (1979) and Gerrard, Hood, and Brown (2004), but more than more than a factor of two smaller than those of Reid et al. (2018). The twist of subsequent threads at destabilisation depends on the experiment. If all threads are driven at the same speed, threads that become unstable later are twisted more than the first. If the central thread is driven faster, it becomes unstable first and then affects others; hence subsequent thread disruptions exhibit less twist. The recorded values of twist imply that the kink instability is unlikely to be responsible for any later disruptions.

Turning to the torus instability, the value of $n$ at the apex in the single-threaded case ($n=1.572$) and the range over the cross-section are consistent with values commonly associated with the torus instability. In the multi-threaded cases, values of $n$ are generally lower: the bottom-most threads in Section 4 show values $n=1.04\,\mbox{--}\,1.06$, lower than such typical values. The faster central thread in Section 5 achieves $n=1.272$ at disruption, again less indicative of instability. Values of $n$ increase with height, making the uppermost threads more likely than the lower to be susceptible to the torus instability. As a consequence of the superposition of Fourier modes in order to achieve near-vertical field at the footpoints, the loop here has a flatter, off-circular shape, evident in the red field lines of Figure 2a. This is less favourable to such an instability, and the rise of the threads before instability is neither far nor sustained. As such, and given the small net current, it is unlikely that the instability here is a conventional torus instability (Démoulin and Aulanier, 2010; Rees-Crockford et al., 2020).

Here, as elsewhere, it is seen that separating conjectured instabilities can be challenging. For the spread of the avalanche, the most important facet is how an ideal instability in a single thread triggers the destabilisation of additional, stable threads.

One final issue with models of curved flux tubes relies upon our chosen velocity scales. In Section 3, our model assumes mean velocities $\approx 8.12\times 10^{-4}\,v_{\mathrm{A}}$. Given our normalisation (Table 1), this corresponds physically to $561~\mbox{m}\,\mathrm{s}^{-1}$. Typical velocity scales in simulations of photospheric velocity patterns often approach several tens or hundreds of metres per second (Gizon and Birch, 2005; Rieutord and Rincon, 2010). We approach realistic driving timescales in this model, in addition to our realistic magnetic structure. Our MHD approach affords some flexibility in normalisation, allowing us to consider speeds even closer to those observed.

7 Conclusions and Future Work

Our investigation means that we are now able to answer the question posed in the title of this work: multi-threaded curved flux tubes are indeed capable of supporting a magnetohydrodynamic avalanche. Our magnetic-arcade field has been shown to be capable of supporting flux tubes, formed through rotational photospheric driving. As with straight cylindrical examples (e.g. Reid et al., 2018), continuous driving of individual threads leads to instability and the destabilisation of other stable threads, whereupon anomalous resistive effects release and redistribute the field lines of the entire tube. The initial instabilities show characteristics common to both kink and torus instabilities. Continual driving of (post-disruption) flux tube remnants leads to the formation of secondary current sheets, which cause aperiodic bursts of Ohmic heating, releasing more energy from the field. Thus an arcade field may continue to yield bursty reconnection events when continuously driven from below.

Our findings diverge from similar models based on straight cylindrical geometry; our cases containing more than one thread demonstrate that the geometry of the magnetic field itself (here, a curved arcade structure) and the imposed driving speed of individual threads both play a significant role in the formation and evolution of an MHD avalanche. In the case where the speed is the same across all footpoints, threads preferentially destabilise through mergers with other threads of the same height and field strength, separate from the others; twisting all threads at the same rate in a single flux tube leads to a gradual destabilisation of threads, from strong to weak field strength. However, driving a single thread unstable more quickly in a region of strong magnetic field can lead to a cascade of destabilisation of other threads, which more closely resembles previous examples between parallel planes; individual threads can destabilise neighbours, regardless of field strength. The geometric configuration (and resulting variation in field strength) may also attenuate the cascade in some parts of the loop.

Several opportunities for further investigation readily present themselves. Our curved model geometry provides an important first step towards a more realistic model of a magnetic arcade, and yet relies on a finite number of Fourier modes. Each additional mode affects curvature and will hence potentially affect the mode of instability recovered: this merits further consideration. Similarly, our geometry requires a novel critical parameter [$\zeta$] in order to trigger anomalous resistive effects upon specific current sheets: would the use of such a critical threshold alter/affect the findings seen in straight cylindrical cases (e.g. Reid et al., 2020)?

Looking ahead still further, this model still omits several relevant (yet potentially complex) physical factors that warrant consideration: these include the effect of gravity; a stratified atmosphere (in particular, stratifying density and temperature); and thermodynamic effects, such as thermal conduction and radiation. Our model has generated many complex, highly time-dependent current sheets, distributed throughout the volume of the flux tube; further analysis of the distribution, lifetimes, and strengths of currents is warranted, together with examination of the impact of model parameters (including field strength and resistivity). From these current sheets, heating is continuously being produced, the magnitude, composition, and distribution of which will be discussed in greater detail. How do the electric fields associated with these current sheets compare with those observed or inferred solar observations, and might these electric fields be capable of accelerating particles? Finally, our initial photospheric driving pattern is relatively simplistic: could a more complex photospheric driver create and/or drive flux tubes in this model?