Laboratory Study of the Torus Instability Threshold in Solar-relevant, Line-tied Magnetic Flux Ropes

Protrusions of the magnetic field and plasma from the solar surface often result in the formation of long, thin magnetic flux ropes (MFRs; Kuperus & Raadu 1974; Chen 1989; Rust 2003). These flux ropes are bundles of twisted magnetic field lines with foot points anchored to the solar surface through line-tying to the conductive photosphere. The ropes are often long-lived but can sometimes violently erupt, leading to coronal mass ejections (CMEs; Crooker et al. 1997; Green & Kliem 2009). Understanding the causes of these eruptions is necessary for the prediction and further understanding of space weather.

One potential mechanism that can trigger a CME from an initially stable MFR is an ideal magnetohydrodynamic (MHD) instability called the torus instability (Bateman 1978; Kliem & Török 2006). This instability occurs when the external magnetic field perpendicular to the MFR's axis, the strapping field, decays quickly enough with height. The rate of decay of the strapping field can be described by a decay index, n_s, and causes an instability when it exceeds a critical value, ${n}_{{\rm{s}}}\gt {n}_{\mathrm{cr}}^{}$. For an axisymmetric, large aspect ratio, full torus, the critical value is ${n}_{\mathrm{cr}}^{}=1.5$ (Shafranov 1966; Bateman 1978). However, extending this theory to nontoroidally symmetric line-tied systems such as MFRs has proved challenging (Isenberg & Forbes 2007; Olmedo & Zhang 2010). While some MFR simulations are consistent with ${n}_{\mathrm{cr}}^{}=1.5$ (Török & Kliem 2007; Aulanier et al. 2009; Zuccarello et al. 2015), other analytical work has found critical values in a much larger range, $0.5\lt {n}_{\mathrm{cr}}^{}\lt 2$, depending on the ratio of the apex height and foot-point half-separation (Olmedo & Zhang 2010). Some simulations have found values near the higher end of this range with $1.75\lt {n}_{\mathrm{cr}}^{}\lt 2$ (Fan & Gibson 2007; Fan 2010). This simulated range is inconsistent with some observational evidence of ${n}_{\mathrm{cr}}^{}\sim 1.3$ (Duan et al. 2019). Critical values on the lower end of the range presented in Olmedo & Zhang (2010) have been seen in recent laboratory experiments where an empirical threshold of ${n}_{\mathrm{cr}}^{}\sim 0.8$ was seen (Myers et al. 2015). This experimental value is also consistent with some solar observations (Jing et al. 2018).

One cause of the large discrepancy in the previous values of ${n}_{\mathrm{cr}}^{}$ comes from inconsistencies in definitions. For example, while some works define n_s at the rope's apex height, z_ap (Jing et al. 2018), others choose a value above the rope even up to z = 2z_ap (Chen 1996). Others still define z_ap as the location where ${n}_{{\rm{s}}}={n}_{\mathrm{cr}}^{}$ for an assumed value of ${n}_{\mathrm{cr}}^{}$ (Wang et al. 2017). Due to the difficulty of solar measurements, sometimes the decay index of the total external magnetic field is used rather than that of just the strapping field, which is used in the theory (Liu 2008). Practical measurement issues also cause some to measure n_s at a fixed height rather than at a dynamic one for each rope (Liu 2008). While the value of n_s at a fixed height is likely correlated with the true n_s(z_ap), it cannot be used to determine ${n}_{\mathrm{cr}}^{}$ or even be directly compared against it. With these considerations, it is no surprise that the literature contains many conflicting values and ranges for ${n}_{\mathrm{cr}}^{}$.

In order to reconcile the different values of ${n}_{\mathrm{cr}}^{}$, we devise a simple MFR model and use it to numerically calculate ${n}_{\mathrm{cr}}^{}$ as a function of the flux rope parameters. The relevant dimensionless parameters are the ratio of the apex height to the foot-point half-separation, the aspect ratio of the partial torus, and the normalized internal inductance. The prediction can be applied to either laboratory or solar conditions to yield different results, ${n}_{\mathrm{cr}}^{\mathrm{lab}}$ and ${n}_{\mathrm{cr}}^{\mathrm{solar}}$. These ${n}_{\mathrm{cr}}^{}$ predictions are then validated against previous (Myers et al. 2015) and more recent experimental results. Our numerical results predict critical values in the range of ${n}_{\mathrm{cr}}^{\mathrm{lab}}\sim (0.65,1.1)$ when applied to our experimental conditions. This is consistent with the previous empirical value of ${n}_{\mathrm{cr}}^{\mathrm{lab}}\sim 0.8$ but greatly increases our understanding of it and provides a detailed dependence on experimental conditions and MFR parameters. Our model is able to predict erupting ropes in our experiments with a true-positive rate (TPR) of 94%. This high detection rate indicates that our method could ultimately be used to improve space weather predictions.

We also present how these values should be translated to the conditions on the Sun, which makes the majority of our results clustered near ${n}_{\mathrm{cr}}^{\mathrm{solar}}\sim 0.9$ with a full range extending to ${n}_{\mathrm{cr}}^{\mathrm{solar}}\sim (0.7,1.2)$. The results reported here can better explain the wide range of ${n}_{\mathrm{cr}}^{\mathrm{solar}}$ values reported in the literature. This new methodology of translating laboratory MFRs to their solar counterparts paves a systematic way to quantitatively investigate the initiation of CMEs through laboratory experiments. These experiments allow for new physics insights that are required for better predictions of space weather events but difficult to obtain otherwise.

In the rest of this paper, we briefly introduce the torus instability, as well as the kink instability, another important MHD instability in arched, line-tied flux ropes. This is followed by the development of numerical models of each of the contributing forces and the resulting critical decay index applicable to our laboratory experiments. The model allows for a more accurate prediction of ${n}_{\mathrm{cr}}^{}$ than can normally be achieved though purely analytical means. The predicted values of ${n}_{\mathrm{cr}}^{}$ also lie in a range that depends on the exact parameters of the MFR. The range allows for more accurate predictions than can be achieved by a single scalar value. The results of the experiments are then translated to the conditions on the Sun to create a better prediction of MFR eruptions and subsequent CMEs.

The magnetic field, B, and electric currents, J, of an MFR can be broken up into components based on their source and direction, as shown in Figure 1(a). First, B can be separated into external and internal components based on whether they are generated by currents in the Sun or within the MFR itself. The external field is further divided into the guide field, B_g, along the axis of the MFR and the strapping field, B_s, perpendicular to this axis. The internal fields and the currents are then separated into toroidal, T, and poloidal, P, components. In a low-β plasma (where β ≡ 2μ₀P/B² is the ratio of thermal to magnetic energy in the plasma), the dominant force on an MFR is the J × B force. From the decomposition of J and B, this force can also be decomposed based on the source terms. These forces are called the hoop, strapping, and tension forces and are defined in Table 1. The hoop force caused by the toroidal current interacting with the internal poloidal field causes the MFR to expand upward, while the tension and strapping forces resist this expansion and hold it down.

Zoom In Zoom Out Reset image size

Table 1.  Breakdown of the Magnetic Forces on a Flux Rope

Force	Symbol	Source Term	Analytical Expression
Hoop force (upward)	Fh	fh = JTBPi	${F}_{{\rm{h}}}=\displaystyle \frac{{\mu }_{0}{I}_{{\rm{T}}}^{2}}{4\pi R}\left[\mathrm{ln}\left(\tfrac{8R}{a}\right)-1+\tfrac{{{\ell }}_{{\rm{i}}}}{2}\right]$ a
Strapping force (downward)	Fs	fs = −JTBs	Fs = −ITBs
Tension force (downward)	Ft	ft = −JP(Bg + BTi)	${F}_{{\rm{t}}}=\displaystyle \frac{-{\mu }_{0}{I}_{{\rm{T}}}^{2}}{8\pi R}\left[\displaystyle \frac{\left\langle {B}_{{\rm{T}}}^{2}\right\rangle -{B}_{{\rm{g}}0}^{2}}{{B}_{{\rm{P}}a}^{2}}\right]\approx -\tfrac{1}{2}\tfrac{{\mu }_{0}{I}_{{\rm{T}}}^{2}}{4\pi R}$ b

Notes. The forces are separated based on their source terms in the J × B force. The fields and currents used are shown in Figure 1(a), and the analytical expressions were derived in detail by Myers et al. (2016).

^aThis expression for F_h is different than some other sources due to the separation of F_t. ^bThis approximation is derived from the minor radius force balance.

Download table as:  ASCIITypeset image

2.1. Torus Instability Parameter

The torus instability occurs when the net force on an MFR increases with height away from equilibrium, i.e.,

$$\begin{eqnarray}&&{\left.\displaystyle \sum _{i}{F}_{i}=0\right|}_{z={z}_{\mathrm{ap}}}\quad \mathrm{and}\quad {\left.\displaystyle \sum _{i}\displaystyle \frac{\partial {F}_{i}}{\partial z}\gt 0\right|}_{z={z}_{\mathrm{ap}}},\end{eqnarray} \tag{ 1 }$$

where F_i are the constituent forces and z_ap is the equilibrium height of the rope's apex. The instability can therefore also be thought of as a loss of equilibrium, where small perturbations cause the J × B forces to push the rope away from equilibrium. Since the main downward force is caused by the external strapping field, the torus instability criterion is often cast in terms of the decay index⁷ of this field,

$$\begin{eqnarray}&&{n}_{{\rm{s}}}=-\displaystyle \frac{z}{{B}_{{\rm{s}}}}\displaystyle \frac{\partial {B}_{{\rm{s}}}}{\partial z}\gt {n}_{\mathrm{cr}}^{},\end{eqnarray} \tag{ 2 }$$

where z is the height above the photosphere and ${n}_{\mathrm{cr}}^{}=1.5$ in an axisymmetric, large aspect ratio, full torus. That is, the torus instability occurs when the strapping field decays too quickly with height above the photosphere.

2.2. Kink Instability Parameter

Another ideal MHD instability that can affect MFRs is the kink instability (Kruskal & Schwarzschild 1954; Shafranov 1956; Török et al. 2004). This instability occurs when the toroidal current, I_T, in a flux rope is too large, or, equivalently, the guide field, B_g, is too small. When the current causes the outer field lines to fully twist around the rope, the instability can set in. The instability criterion can be described in terms of the edge safety factor, q_a, or, equivalently, the twist number, N,

$$\begin{eqnarray}&&{q}_{a}\equiv \displaystyle \frac{1}{| N| }=\displaystyle \frac{2\pi a}{L}\displaystyle \frac{{B}_{{\rm{T}},a}}{{B}_{{\rm{P}},a}}\lt {q}_{\mathrm{cr}},\end{eqnarray} \tag{ 3 }$$

where a and L are the minor radius and length of the MFR, B_T,a is the edge toroidal magnetic field, B_P,a ≡ μ₀I_T/(2πa) is the edge poloidal field, and q_cr is the critical safety factor below which the instability occurs. In a toroidally symmetric system, such as a tokamak, the critical value is q_cr = 1. However, in line-tied systems where there is no toroidal symmetry, q_cr is modified (Ryutov et al. 2006). In recent experiments with arched, line-tied MFRs, the critical value was found to be near q_cr ∼ 0.8 (Myers et al. 2016). This value is also consistent with some analytical work (Hood & Priest 1981).

2.3. Torus and Kink Instabilities' Effect on MFR Eruptions

A simple model of a solar flux rope is that of a partial torus with ends anchored in the solar surface (see Figure 2). In order to approximate the path of a flux rope, we will use a "shifted-circle" model, such as in Chen (1989) and Myers et al. (2016). In this model, the axis of a flux rope is assumed to follow the path of a circle with foot points fixed on the solar surface, which is modeled as a perfectly conducting region for z < 0. By keeping the location of the foot points fixed, line-tying at the solar surface can be maintained while the flux rope height is changed.

Zoom In Zoom Out Reset image size

In this model, the radius of curvature for the MFR, R_sc, and the height of its center of curvature, z_sc, are given by

$$\begin{eqnarray}&&{R}_{\mathrm{sc}}=\displaystyle \frac{{z}_{\mathrm{ap}}^{2}+{x}_{{\rm{f}}}^{2}}{2{z}_{\mathrm{ap}}}\quad \mathrm{and}\quad {z}_{\mathrm{sc}}=\displaystyle \frac{{z}_{\mathrm{ap}}^{2}-{x}_{{\rm{f}}}^{2}}{2{z}_{\mathrm{ap}}},\end{eqnarray} \tag{ 4 }$$

where z_ap is the apex height and x_f is the (fixed) foot-point half-separation. The angle from the z-axis to the foot point, ϕ_f, which measures the fraction of the circle present above z = 0, is given by

$$\begin{eqnarray}&&{\phi }_{{\rm{f}}}=2{\tan }^{-1}\left(\displaystyle \frac{{z}_{\mathrm{ap}}}{{x}_{{\rm{f}}}}\right).\end{eqnarray} \tag{ 5 }$$

From these equations, we see that the shifted-circle model is 1D, defined by one parameter, z_ap, and one constant, x_f.

However, the shifted circle only describes the rope above the solar surface at z = 0. In order for the model to be self-consistent, the toroidal current flowing through the rope must close. Different options for this closure are shown in Figure 2. In the case of solar flux ropes, the most reasonable current closure is that of an image current in order to maintain the value of the normal magnetic field at the surface. However, in laboratory experiments, the return current flows through wires and thus should be fixed and independent of z_ap.

Another crucial part of modeling MFRs is a description of their rise and expansion. While the magnetic axis is defined by just z_ap and x_f, the hoop force depends on other parameters of the rope. These include the total toroidal current, I_T; the inverse aspect ratio, ε ≡ a/R; and the internal inductance, ${{\ell }}_{i}\equiv \langle {B}_{{\rm{P}}}^{2}\rangle /{B}_{{\rm{P}},a}^{2}$, where $\langle {B}_{{\rm{P}}}^{2}\rangle$ is the cross-section average of ${B}_{{\rm{P}}}^{2}$. The internal inductance measures the distribution of toroidal current. Specifically, ℓ_i = 0.5 for a uniform current and increases for current distributions that are more peaked around the axis. The evolution with respect to z_ap of these parameters during an instability must also be specified. The short-term evolution of I_T(z_ap) is often defined by the conservation of magnetic flux below the rope (Kliem & Török 2006; Olmedo & Zhang 2010). However, in the laboratory experiments described in Section 4, I_T is fixed on the timescale of an eruption due to the large external inductance in the capacitor banks driving the ropes, and so ${I}_{{\rm{T}}}({z}_{\mathrm{ap}})=\mathrm{Const}.$ We will also consider a self-similar expansion of the ropes. This means, in part, that $\varepsilon ({z}_{\mathrm{ap}})=\mathrm{Const}.$ Since the hoop force is only logarithmically dependent on ε, it is not very sensitive to the choice of ε(z_ap), and a self-similar expansion should be a reasonable approximation. The value of ℓ_i is determined by the radial distribution of the toroidal current. Therefore, a self-similar expansion also implies that ${{\ell }}_{i}({z}_{\mathrm{ap}})\,=\mathrm{Const}.$ However, ℓ_i does not fully describe the current distribution, and a specific distribution must be chosen in order for numerical calculations to be carried out. As a model current distribution, we have used

$$\begin{eqnarray}&&{J}_{{\bf{T}}}(r)=\displaystyle \frac{{I}_{{\rm{T}}}}{\pi {a}^{2}}(1+\alpha ){\left[1-{\left(\displaystyle \frac{r}{a}\right)}^{2}\right]}^{\alpha },\end{eqnarray} \tag{ 6 }$$

where r is the minor radius coordinate and α > −1 is a free parameter that has a one-to-one correspondence to ℓ_i ∈ (0, ∞).

3.1. Numerical Modeling

In order to study the stability of flux ropes, we investigate the force acting on a toroidal volume element around the apex, shown in gray in Figure 1(a). The hoop force acting on this element is numerically calculated using

$$\begin{eqnarray}&&{F}_{{\rm{h}}}=\displaystyle \frac{1}{{R}_{0}{\rm{\Delta }}T}{\int }_{-{\rm{\Delta }}T/2}^{{\rm{\Delta }}T/2}{dT}{\int }_{0}^{2\pi }d\theta {\int }_{0}^{a}{dr}({{rR}}_{0}{J}_{{\bf{T}}}{B}_{{yi}}),\end{eqnarray} \tag{ 7 }$$

where R₀ is the radius of curvature, θ is the poloidal angle, T is the toroidal angle, ΔT is a small toroidal angle around the apex, and B_yi is the y-component of B_Pi, which is found via the Biot–Savart law,

$$\begin{eqnarray}&&{B}_{{yi}}=\displaystyle \frac{{\mu }_{0}}{4\pi }\int \hat{y}\cdot \left[{{\boldsymbol{J}}}_{{\boldsymbol{T}}}({\boldsymbol{x}}^{\prime} )\times \displaystyle \frac{{\boldsymbol{x}}-{\boldsymbol{x}}^{\prime} }{| {\boldsymbol{x}}-{\boldsymbol{x}}^{\prime} {| }^{3}}\right]{d}^{3}x^{\prime} .\end{eqnarray} \tag{ 8 }$$

In the space z > 0, J_T follows the apex of the rope with a magnitude given by Equation (6). The rest of the integral in Equation (8) depends on the current closure chosen. In the fixed current closures, the currents are carried by thin wires, so a line current is used. For the image current, the same current is used except with the sign of the x-component flipped. Then the normalized hoop force, F_h/F_N, can be found as a function of z_ap/x_f, ε, and ℓ_i. The hoop force has been normalized by ${F}_{{\rm{N}}}\equiv {\mu }_{0}{I}_{{\rm{T}}}^{2}/(4\pi {x}_{{\rm{f}}})$, and an example of this calculation is shown in Figure 3(a).

Zoom In Zoom Out Reset image size

In order to find the critical decay index, ${n}_{\mathrm{cr}}^{}$, for the onset of the torus instability, first consider an MFR in equilibrium with an apex height of z_ap. The total J × B force acting on the apex can be decomposed as F = F_h + F_s + F_t, where F_h, F_s, and F_t are the hoop, strapping, and tension forces as defined in Table 1. In order to find ${n}_{\mathrm{cr}}^{}$, we ask when Equation (1) is marginally satisfied,

$$\begin{eqnarray}&&{\left.\displaystyle \frac{\partial F}{\partial z}\right|}_{F=0}\mathop{=}\limits^{!}0\Rightarrow \displaystyle \frac{\partial {F}_{{\rm{s}}}}{\partial z}=-\displaystyle \frac{\partial }{\partial z}({F}_{{\rm{h}}}+{F}_{{\rm{t}}}).\end{eqnarray} \tag{ 9 }$$

The derivative of F_s can be recast in terms of decay indices,

$$\begin{eqnarray}&&-\displaystyle \frac{z}{{F}_{{\rm{s}}}}\displaystyle \frac{\partial {F}_{{\rm{s}}}}{\partial z}=-\displaystyle \frac{z}{{I}_{{\rm{T}}}{B}_{{\rm{s}}}}\displaystyle \frac{\partial }{\partial z}({I}_{{\rm{T}}}{B}_{{\rm{s}}})={n}_{I}+{n}_{{\rm{s}}},\end{eqnarray} \tag{ 10 }$$

where ${n}_{I}=-z\,\partial (\mathrm{ln}{I}_{{\rm{T}}})/\partial z$ is the decay index of the toroidal current. At marginal stability, ${n}_{{\rm{s}}}={n}_{\mathrm{cr}}^{}$; thus, the critical decay index is given by

$$\begin{eqnarray}&&{n}_{\mathrm{cr}}^{}=-{n}_{I}-{\left.\displaystyle \frac{z}{{F}_{{\rm{h}}}+{F}_{{\rm{t}}}}\displaystyle \frac{\partial }{\partial z}({F}_{{\rm{h}}}+{F}_{{\rm{t}}})\right|}_{z={z}_{\mathrm{ap}}}\end{eqnarray} \tag{ 11a }$$

$$\begin{eqnarray}&&={n}_{I}-{\left.\displaystyle \frac{{{zF}}_{{\rm{N}}}}{{F}_{{\rm{h}}}+{F}_{{\rm{t}}}}\displaystyle \frac{\partial }{\partial z}\left(\displaystyle \frac{{F}_{{\rm{h}}}+{F}_{{\rm{t}}}}{{F}_{{\rm{N}}}}\right)\right|}_{z={z}_{\mathrm{ap}}},\end{eqnarray} \tag{ 11b }$$

where we have also normalized all of the forces by F_N in the second step. This normalization of the hoop and tension forces removes their I_T dependency and thus isolates the effect of I_T (and its derivatives) to the first term. In the experimental results presented in Section 5, n_I ≈ 0 due to the large external inductance in the circuit driving the MFRs.

From the numerically calculated F(z_ap) and Equation 11, ${n}_{\mathrm{cr}}^{}$ can be evaluated, an example of which is shown in Figure 3(b). The numerical ${n}_{\mathrm{cr}}^{}$ generally increases with z_ap/x_f and ε, while it slightly decreases with ℓ_i. The new value of ${n}_{\mathrm{cr}}^{}$ is reduced from the standard ${n}_{\mathrm{cr}}^{}=1.5$ for a full torus because the general reduction of F_h with height is counteracted by the emergence of a larger fraction of the torus into the z > 0 space. Since the partial torus effects cause F_h to decay slower with height than it would for a full torus, B_s must also decay slower if stability is to be maintained. This reduces the value of ${n}_{\mathrm{cr}}^{}$. The value of ${n}_{\mathrm{cr}}^{}$ for the full torus differs from a constant because, in Equation (2), we have defined n_s with respect to z rather than R, as is often done when studying tokamaks.

Flux ropes were created inside the Magnetic Reconnection Experiment (MRX; Yamada et al. 1997). The experimental setup is shown in Figure 1(b) and described in more detail in Myers et al. (2016). The flux ropes were created by a discharge between two copper electrodes (ranging from 7.0 to 7.5 cm in radius) that were inserted into the MRX vessel with a variable foot-point half-separation of x_f = 15.5–18.0 cm. The plasma was generated by injecting a small amount of neutral hydrogen (30–40 mTorr) at the vessel wall and directly at the electrodes. Breakdown was then achieved by biasing the electrodes with a capacitor bank charged to a voltage in the range of 3–4 kV. In order to separate the plasma region above the electrodes from the space behind, a glass substrate is placed below the electrodes. Each discharge lasts about 1 ms, during which time current is quasi-statically injected into the rope with a characteristic driving time of τ_D ∼ 150 μs, which is considerably longer than the Alfvèn time of τ_A ≈ 3–8 μs. This separation of scales is what allows for MHD instabilities to be studied in MRX.

Along with the electrodes, two sets of magnetic field coils were inserted into MRX to compliment the two external coils already present. These coils determine the configuration of the vacuum fields and are arranged such that the guide and strapping fields are independently controlled by two sets of coils each. This pairing allows for control of both the strength and decay index of the fields. The magnetic field within each rope was measured by a 2D array of over 300 magnetic pickup coils inserted at the rope apex. By controlling the magnitude and direction of the currents in the four coils, both q_a and n_s can be selected for each shot.

In addition to the numerical model of the hoop force discussed in Section 3, numerical corrections have been added to F_h due to eddy currents that formed in the vessel wall well above the ropes at a height of z_w = 68 cm. Due to the large scale separation between the wall skin time (τ ∼ 3 ms) and the flux rope driving time (τ_D ∼ 150 μs), the wall can be considered a perfect conductor on the timescale of the rope's lifetime. This makes calculating the mutual inductance and driven eddy currents rather straightforward. There are also eddy currents that can be driven in the center stack that runs through MRX and under the ropes. However, the skin time of the center stack (τ ∼ 75 μs) is shorter than the driving time and therefore neglected.

In Figure 4, each shot taken during this campaign is plotted based on its value of the stability parameters for the kink and torus instabilities, q_a and n. The results of almost 2000 discharges are represented in the plot. Shots that erupted multiple times (described below) are shown by red triangles, while noneruptive shots are blue circles. The quadrants of the parameter space are labeled with their stability behavior.⁸ We are currently interested in the onset of the torus instability and so will focus on the transition from the "Failed Kink" to the "Eruptive" quadrants. To this end, we will limit our analysis to those shots with a safety factor below an empirical value of q_a < 0.67, marked in Figure 4 as a dashed line. This value was chosen based on statistical arguments that are presented in Section 5.2. The decay index that is plotted is the decay index from the total field, ${n}_{}=-z\,\partial (\mathrm{ln}{B}_{})/\partial z$, rather than just B_s, where the total vacuum field is B = B_s + B_g. This is done in part to be more consistent with how similar plots are made from solar observations (e.g., Jing et al. 2018).

Zoom In Zoom Out Reset image size

When a rope becomes unstable and erupts, its total current would normally drop. However, there is a large external inductance in the capacitor bank driving the plasma current. This inductance prevents the total current from quickly changing; so, instead of an erupting rope dissipating, a new rope is formed at a lower height. Since the vacuum fields have not significantly changed, the new rope is often unstable as well and quickly erupts. As an example, a cross section of the toroidal current distribution at the apex is shown for two times in Figure 5(a). In the first panel, only one rope exists, but since it is unstable, it quickly rises. Eventually, a new rope is formed at a lower height, causing two distinct ropes to briefly coexist. This new rope then rises and erupts, repeating the pattern. The apex is defined by nulls in the in-plane (poloidal) field and can be multivalued if multiple nulls exist. The poloidal field is represented by arrows and the nulls by black dots. When multiple ropes exist, there are usually three nulls: two O-points, one at each rope apex, and an X-point between them. An eruption was thus defined as the coexistence of two ropes, and a rope is considered "eruptive" if this occurs more than once. An example of the repeated eruptions is shown in Figure 5(b). This plot shows the evolution of the apex height versus time along with a stable example for reference.

Zoom In Zoom Out Reset image size

5.1. Prediction of ${n}_{\mathrm{cr}}^{}$

Using our more complete numerical model of the hoop force and applying Equation 11, we can predict the critical decay index, ${n}_{\mathrm{cr}}^{}$, for each shot. Each shot is fit to the model based on its experimental parameters, z_ap/x_f, ε, and ℓ_i. The ranges of these parameters that are represented by the experimental data are discussed in Section 5.3. In Figure 6, each shot is plotted based on this predicted ${n}_{\mathrm{cr}}^{}$ and the experimentally measured n. In this type of plot, unstable shots should lie above the line $n={n}_{\mathrm{cr}}^{}$ (also plotted), with stable ones below it. With some exceptions that will be discussed below, this is generally confirmed by the figure. This plot shows numerical evidence for the ${n}_{\mathrm{cr}}^{}\sim 0.8$ empirical value of the critical decay index observed previously in Myers et al. (2015).

Zoom In Zoom Out Reset image size

There are some caveats to Figure 6 that deserve explanation. In order to get good experimental agreement, the horizontal closure from Figure 2 was used instead of the semicircular closure that is more representative of the experimental setup. The real experimental closure is the physical path taken by the wires leading from the foot-point electrodes to the capacitor bank and is somewhere between these two options, which causes a slight overestimation of ${n}_{\mathrm{cr}}^{}$.

There are a number of shots in Figure 6 that do not exactly follow the predicted $n={n}_{\mathrm{cr}}^{}$ stability line. Most of these exceptions are isolated stable shots above the line but do not alter the statistical evidence of our overall prediction. One potential cause of the discrepancy comes from the value of B that is used in n. When using the value of B_y in n_s, it is assumed that the rope is exactly perpendicular to the measurement array. If, instead, the rope's axis were tilted at a small angle, θ, relative to $\hat{x}$, the nominal guide field direction, then the correct effective strapping field would be ${B}_{{\rm{s}},\mathrm{eff}}={B}_{y}\cos \theta -{B}_{x}\sin \theta$, where B_y and B_x are, respectively, the strapping and guide fields when there is no tilt. Due to experimental constraints, the guide field, B_g, is often considerably larger than the strapping field, B_s, so even small tilt angles can significantly change the decay index. Therefore, for the robustness of our results, the total field is used here to calculate ${n}_{\mathrm{cr}}^{}$, which is consistent with some previous work (Liu 2008; Myers et al. 2015). Several attempts were made to measure the tilt angle of a given rope from the available data; however, nothing convincing was found. Further tests should be preformed in future experiments to better measure the tilt angle and thus B_s,eff.

5.2. Statistical Analysis

As mentioned in Section 2.3, we have eliminated the effect of the failed torus region by only focusing on shots with an edge safety factor below a set value, ${q}_{\max }$. In order to investigate the effect of this choice of ${q}_{\max }$, the statistical validity of our ${n}_{\mathrm{cr}}^{}$ prediction was compared for different values of ${q}_{\max }$. In Figure 7, we compare two metrics, the TPR (the fraction of eruptive shots that were correctly predicted as eruptive) and the false-positive rate (FPR; the fraction of noneruptive shots that were incorrectly predicted as eruptive). As ${q}_{\max }$ is increased, the TPR quickly rises and then plateaus near 96%. However, increasing ${q}_{\max }$ also causes the FPR to increase as more shots are included that are not purely eruptive and behave more like a failed torus. For the best results, a value of ${q}_{\max }$ must be chosen to balance these two effects. By inspecting Figure 7, a value of ${q}_{\max }=0.67$ was chosen because for larger values, the FPR continues to rise for little gain in the TPR. This value of ${q}_{\max }$ results in values of TPR = 94% and FPR = 40%. It should also be noted that the total number of shots represented in Figure 7 is about half that of Figure 4. This is because a number of shots were taken with the measurement array not perpendicular to the MFR axis. The angle removes our ability to predict ${n}_{\mathrm{cr}}^{}$, and thus they cannot be included in the statistics.

Zoom In Zoom Out Reset image size

5.3. Experimental Parameter Ranges

As discussed in Section 3, the calculation of a predicted ${n}_{\mathrm{cr}}^{}$ involves three dimensionless experimental parameters: the normalized apex height, z_ap/x_f; the inverse aspect ratio, ε; and the internal inductance, ℓ_i. The shots shown in Figure 6 define ranges of these parameters based on experimental constraints. Scatter plots of these values versus the predicted ${n}_{\mathrm{cr}}^{}$ using a horizontal current closure for each of the shots are shown in Figure 8. The parameters fall within the ranges of z_ap/x_f ∈ (0.75, 1.75), ε ∈ (0.5, 1.0), and ℓ_i ∈ (0.3, 0.9). Care should be taken when extrapolating our results beyond these values. Some correlations can also been seen between the parameters in these plots and ${n}_{\mathrm{cr}}^{}$. The strongest of these correlations is between z_ap/x_f and ${n}_{\mathrm{cr}}^{}$, which is not surprising based on the predicted ${n}_{\mathrm{cr}}^{}$ curves presented in Figure 9.

Zoom In Zoom Out Reset image size

Zoom In Zoom Out Reset image size

It should be noted that some shots had a measured ε > 1. This is possible due to a varying minor radius over the length of the rope, which is not captured by our model. Since, especially at large ε, there are only minor changes in ${n}_{\mathrm{cr}}^{}$ with ε (see Figure 9(a)), these shots have had their aspect ratio adjusted to ε = 1.

The laboratory flux ropes we have created are most directly comparable to solar flux ropes with some considerations that are laid out in Table 2. There are three main differences that we will focus on: (1) the evolution of the total toroidal current, (2) the subsurface closure of the current, and (3) eddy currents in the outer vessel wall.

Table 2.  Qualitative Comparison of Experimental and Solar Conditions with Resulting ${n}_{\mathrm{cr}}^{}$ Adjustment

Parameter	Solar Value	Laboratory Value	${n}_{\mathrm{cr}}^{}$ Adjustment
Current return path	Image currents on the solar surface	Fixed wires	Decrease of ${n}_{\mathrm{cr}}^{}$ for low-lying ropes; ${n}_{\mathrm{cr}}^{\mathrm{solar}}\approx {n}_{\mathrm{cr}}^{\mathrm{lab}}-0.3$
Lower boundary condition	Conductive surface at z = 0; line-tied everywhere	Conductive, line-tied foot points; dielectric elsewhere	Affects the poloidal flux conservation but no direct effect on ${n}_{\mathrm{cr}}^{}$
Total toroidal current	Set by poloidal flux conservation; nI ∼ 0.5	External inductance forces constant current; nI ≈ 0	Increase of ${n}_{\mathrm{cr}}^{}$ due to nI; ${n}_{\mathrm{cr}}^{\mathrm{solar}}={n}_{\mathrm{cr}}^{\mathrm{lab}}+{n}_{I}\sim {n}_{\mathrm{cr}}^{\mathrm{lab}}+0.5$
Upper boundary condition	Open boundary at z → ∞ allows ropes to erupt away	Conductive wall prevents ropes from fully erupting	Minor effects except in posteruption stages; ${n}_{\mathrm{cr}}^{\mathrm{solar}}\approx {n}_{\mathrm{cr}}^{\mathrm{lab}}-0.1$

Download table as:  ASCIITypeset image

The total toroidal current as a rope rises behaves very differently on the Sun when compared to the laboratory case. In the laboratory, the total plasma current is fixed due to the large external inductance in the driving circuit. However, on the Sun, the current is determined by the conservation of magnetic flux between the rope and the solar surface. This affects ${n}_{\mathrm{cr}}^{}$ through the current decay index, ${n}_{\mathrm{cr}}^{\mathrm{solar}}={n}_{\mathrm{cr}}^{\mathrm{lab}}+{n}_{I}^{\mathrm{solar}}$. The value of n_I depends on the model being used, but for the standard large aspect ratio full torus, ${n}_{I}^{\mathrm{solar}}\approx 0.5$. This value is also similar to that used in other models. When the toroidal current is not able to decrease with height, the MFRs are less stable to the torus instability.

The laboratory return path of the current through fixed wires rather than a conductive surface as on the Sun also has an effect on ${n}_{\mathrm{cr}}^{}$. The effect is largest on low-lying ropes, and when the changes between the horizontal closure and image current models are compared near z_ap/x_f ∼ 1, we find that ${n}_{\mathrm{cr}}^{}$ is decreased by about −0.3.

The effect of eddy currents in the outer wall located at z ≈ 3.8x_f has also been numerically calculated. The wall has the largest effect on tall ropes, decreasing ${n}_{\mathrm{cr}}^{}$ by about −0.1. The effect of the wall is much larger on ropes that would fully erupt away. However, since we are focused on the onset of eruption and not the dynamics afterward, this only has a minimal effect on our results. The specific change in ${n}_{\mathrm{cr}}^{}$ due to these last two effects is shown for each shot in Figure 10.

Zoom In Zoom Out Reset image size

The net effect on ${n}_{\mathrm{cr}}^{}$ when translating from laboratory to solar conditions is therefore

$$\begin{eqnarray}&&{n}_{\mathrm{cr}}^{\mathrm{solar}}\approx {n}_{\mathrm{cr}}^{\mathrm{lab}}+{n}_{I}-0.4\sim {n}_{\mathrm{cr}}^{\mathrm{lab}}+0.1.\end{eqnarray} \tag{ 12 }$$

The modified ${n}_{\mathrm{cr}}^{\mathrm{solar}}$ is shown in Figure 9. These plots show ${n}_{\mathrm{cr}}^{\mathrm{solar}}$ as a function of z_ap/x_f, using an image current closure and n_I=0.5, for various values of ε with fixed ℓ_i = 0.5 and then for various values of ℓ_i with fixed ε = 0.5. These 1D cuts help in breaking down the 3D input space. It can be easily seen that the effects of changing either ε or ℓ_i is very minimal for tall ropes with z_ap/x_f > 1.5 but can be important for low-lying ropes. From these cuts, it is seen that ${n}_{\mathrm{cr}}^{\mathrm{solar}}$ increases with increasing ε and decreases with increasing ℓ_i. That is, MFRs with lower aspect ratios and less peaked currents are more stable to the torus instability. The nonmonotonic curves in Figure 9 are caused by two competing effects as the MFR rises: F_h decreases overall, and a larger fraction of the torus is present, increasing F_h. The latter effect becomes smaller as the MFR rises, since the relative increase is less.

Another difference between the laboratory and solar conditions is the state of the volume surrounding a flux rope. On the Sun, the external volume is the ambient plasma of the corona with frozen-in flux. However, in our experiments, the plasma is generated by the arc discharge that forms the rope, so most of the external volume is near vacuum with small neutral gas pressure. On the Sun, an erupting rope must displace this ambient plasma, while there is no such force in our experiments. It has been seen in simulations that when this background plasma contains a large guide field, otherwise unstable ropes can be stabilized (Kliem et al. 2014). However, we have not seen this effect in our experiments, possibly because of the lack of ambient plasma. Since we are currently focusing on ropes with small q_a, the ambient guide field is also small, and this effect likely would not be applicable.

The onset criteria of the torus instability in arched, line-tied flux ropes have been investigated by laboratory experiments and the application of a simplified theoretical model. Arched, line-tied flux ropes are modeled by a shifted circle with fixed foot points (see Figure 2), and the resulting forces are derived from simple Biot–Savart calculations. The effect on the torus instability of a partial torus model was previously considered in Olmedo & Zhang (2010). However, the use of numerical calculations in the hoop force of a partial torus allows for more accurate predictions of its onset and the application to low aspect ratio tori with small edge safety factors.

When applied to our experimental data, our model has predicted critical decay index values within the range of ${n}_{\mathrm{cr}}^{\mathrm{lab}}\sim (0.65,1.1)$, which is consistent with the previous empirical value of ${n}_{\mathrm{cr}}^{\mathrm{lab}}\sim 0.8$ (Myers et al. 2015). Thus, our new model provides quantitative evidence for the previous empirical description of the flux rope stability regimes. The deviation of this range from the standard value of 1.5 (Kliem & Török 2006) is caused by taking into account the full effects of the partial torus when applied to the conditions of the laboratory experiments. A large portion of the decreased ${n}_{\mathrm{cr}}^{\mathrm{lab}}$ is due to the external inductance in the flux rope–driving circuit, which prevents the total plasma current from changing. While this is a major change from solar conditions, its effect is easily quantified and separated, allowing for direct comparisons. Further changes are caused by the modification of the hoop force from partial torus effects. In particular, the emergence of more of the torus as the rope rises counteracts some of the decrease of F_h.

While the predicted ${n}_{\mathrm{cr}}^{\mathrm{lab}}$ largely agrees with our experimental results, there are still a number of exceptions that were discussed in Section 5.1. These exceptions can possibly be explained by a small tilting of the flux ropes relative to the measurement array. However, since this angle could not be found through the measurements taken, there may, of course, be other explanations that have yet to be considered. Further investigation into this portion of the parameter space may yield a better understanding of these ropes and their unexpected behavior. Even with these exceptions, our model is able to predict eruption ropes with a TPR of 94%. This rate of prediction is very substantial when considering the complexity of MFR eruptions. However, there is still a relatively high FPR of 40% that warrants further investigation.

Our analysis has also assumed small values of q_a to eliminate the failed torus regime. The twist number is difficult to measure in solar observations, but the maximum value we have chosen, ${q}_{\max }=0.67$ (or, equivalently, $| {N}_{\min }| =1/{q}_{\max }=1.49$), is relevant to many solar MFRs (Duan et al. 2019). When translated to solar conditions, we find that most of our predicted critical decay index values fall near ${n}_{\mathrm{cr}}^{\mathrm{solar}}\sim 0.9$ and within a larger range of ${n}_{\mathrm{cr}}^{\mathrm{solar}}\sim (0.7,1.2)$, which is consistent with some solar observations (Jing et al. 2018). This range provides a more accurate description of ${n}_{\mathrm{cr}}^{\mathrm{solar}}$, and our methodology can help explain the discrepancies seen in previous numerical and observational studies. This new methodology allows for the novel application of laboratory experiments directly to the conditions present on the Sun. The ability to perform direct, local measurements in the laboratory gives experiments a unique advantage over other methods of studying solar flux ropes.

We would like to thank P. Sloboda, A. Jones, J. Latham, L. Bevier, and R. Cutler for their help with the experimental setup. This research is supported by Department of Energy contract Nos. DE-SC0019049 and DE-AC0209CH11466 (Princeton) and NASA grant Nos. 80HQTR17T0005, 80NSSC17K0016, 80NSSC19K0860, and 80NSSC18K1705 and the German Research Foundation (DFG). Sandia National Laboratories is a multimission laboratory managed and operated by National Technology & Engineering Solutions of Sandia, LLC, a wholly owned subsidiary of Honeywell International Inc., for the U.S. Department of Energy's National Nuclear Security Administration under contract DE-NA0003525. This paper describes objective technical results and analysis. Any subjective views or opinions that might be expressed in the paper do not necessarily represent the views of the U.S. Department of Energy or the United States Government.