© ESO 2022

1. Introduction

Solar flares and coronal mass ejections (CMEs) are violent eruptions on the Sun. Essentially they are a manifestation of the explosive release of magnetic energy in the solar corona. It is commonly believed that solar magnetic fields are generated at the base of the convection zone, often in the form of thin, intense flux tubes, which slowly emerge outwards, eventually traveling into the corona through the solar surface (i.e., the photosphere). Thereafter, the coronal magnetic fields are continuously, but rather slowly, dragged at their footpoints by the photospheric surface motions, often in a organized way, such as via shearing and rotational flows. This process is also inherent to the flux emergence process, allowing for magnetic free energy in the corona to gradually accumulate and with the magnetic field configuration often building up to a highly stressed one of an S shape (e.g., the coronal sigmoids frequently observed in X-ray and EUV images) prior to an eruption. Before the eruption onset, the coronal system is in an equilibrium state, in which the outward magnetic pressure of the low-lying, strongly stressed flux is balanced out by the inward magnetic tension of the overlying, mostly unsheared flux. At a critical point, the eruption is suddenly initiated with a catastrophic disruption of this force balance, during which the free magnetic energy is rapidly converted into impulsive heating and fast acceleration within the plasma.

The way in which solar eruptions are initiated, namely, how the force balance before eruption is suddenly destroyed and what drives the eruption, remains an open question. A number of theories have been proposed (Forbes et al. 2006; Shibata & Magara 2011; Chen 2011; Schmieder et al. 2013; Aulanier 2014; Janvier et al. 2015), often divided between two categories: one based on the ideal magnetohydrodynamic (MHD) instability and the other on magnetic reconnection. The first category generally requires the pre-existence of a magnetic flux rope (MFR), a group of twisted magnetic field lines wound tightly about a common axis. The ideal instabilities of the MFR, such as kink instability and torus instability, can initiate eruptions (Kliem & Török 2006; Török & Kliem 2005; Fan & Gibson 2007; Aulanier et al. 2010; Amari et al. 2018). In the second category, the most frequently mentioned models are the breakout model and the tether-cutting reconnection model. The breakout model requires a quadrupolar magnetic configuration in which a magnetic null point is situated above the core of sheared magnetic flux. It is proposed that the reconnection at the null point removes the overlying restraining flux to trigger an eruption (Antiochos et al. 1999; Aulanier et al. 2000; Lynch et al. 2008; Wyper et al. 2017). The tether-cutting reconnection model relies on only a single sheared arcade, namely, a bipolar magnetic field. With an increase in magnetic shear, a current sheet (CS) will be formed slowly at a low altitude above the photospheric magnetic PIL. Initially, the magnetic reconnection at that CS slowly reduces the downward magnetic tension force by “cutting the magnetic tethers”, and then the upward magnetic pressure force is unleashed and pushes up the flux, making it rise; this, in turn, enhances the tether-cutting reconnection. After a short interval, the process becomes runaway in nature, triggering the eruption, and the quickly rising flux stretches the surrounding envelope magnetic field upward, forming a new elongated CS above the PIL. The magnetic reconnection of the newly formed CS further speeds up the eruption to form a CME (Moore & Labonte 1980; Moore & Roumeliotis 1992; Moore et al. 2001). Compared with other models, the tether-cutting scenario is the simplest one in terms of magnetic topology, since it relies on a single magnetic arcade (corresponding to a pair of opposite polarities at the photosphere) without any additional special topology, such as a null point or MFR. However, unlike other models that have been extensively realized in numerical 3D MHD simulations¹, the tether-cutting model has not yet been validated in any 3D MHD simulations and thus remains a conjectural “cartoon”.

In fact, early simulations in 2D or translational invariant geometries (Mikic & Linker 1994; Amari et al. 1996; Choe & Lee 1996) show that via the continuous shearing of its footpoints, a single magnetic arcade asymptotically approaches an open state, which contains a CS. This is consistent with the Aly-Sturrock conjecture (Aly 1991; Sturrock 1991). When we take finite resistivity into account, the system experiences a global disruption once reconnection sets in at the CS, which specifically begins at the point with the largest current density in the CS. Such a simple and efficient mechanism of eruption initiation has only recently been established in fully 3D by Jiang et al. (2021), with an ultra-high accuracy MHD simulation. That simulation is initialized with a bipolar potential field. Through surface shearing motion along the PIL, a vertical CS forms quasi-statically above the PIL. Once the CS is sufficiently thin so that the ideal MHD is broken down, reconnection sets in and instantly triggers the eruption. The simulation shows that the reconnection not only cuts the magnetic tethers, but also results in a strong upward tension force; it is the latter that plays the key role in driving the eruption, that is, the slingshot effect of the reconnection is the main driver of the eruption. We note that this mechanism is different from the tether-cutting model in a twofold sense. Firstly, the tether-cutting model proposed that the reconnection only plays the role of cutting the confinement of the field lines and it is the unleashed magnetic pressure that drives the eruption. Secondly, the tether-cutting model assumes that before the onset of the eruption (i.e., the start of impulsive phase of eruption), there is a relatively long phase of slow reconnection from more than ten minutes to a few hours, which gradually cuts the tethers. It is not until a “global instability” occurs that the eruption can begin (Moore et al. 2001). Such slow reconnection is not runaway and it does not exist in the simulation from Jiang et al. (2021). In this sense, the model demonstrated by Jiang et al. (2021) stands alone with the original tether-cutting model; hereafter, we refer to this mechanism as the BASIC model, where the acronym “BASIC” refers to the key ingredients as involved in the mechanism: a Bipolar magnetic Arcade as sheared evolves quasi-Statically and forms Internally a Current sheet.

As demonstrated previously, the solar eruption can be initiated from a single bipolar field by way of the BASIC mechanism. However, a natural question arises related to how the mechanism operates differently with different flux distributions of bipolar field on the photosphere. Since Jiang et al. (2021) carried out numerical experiments for only one set of magnetic flux distribution on the bottom surface, in this paper, we investigate the BASIC mechanism with different magnetic flux distributions. In particular, we aim to find out what kind of flux distribution of the bipolar field on the photosphere is favorable for producing major eruptions and why. This study is motivated by a well-known fact based on the observation that major solar eruptions occur predominantly in source regions having a pronounced PIL with both an elongated distribution of flux along it and a strong gradient of field across it (Schrijver 2007; Toriumi & Wang 2019). Such a PIL is mostly often found in a particular type of sunspot group, namely, the δ sunspots, which are highly flare-productive, as was first found by Künzel (1959). By examining the magnetic properties of regions associated with almost 300 M- and X-class flare, Schrijver (2007) found that the magnetic flux distribution in the flare site of these regions is often of bipolar configuration with a characteristic pattern: a relatively elongated, strong-gradient PIL with strong magnetic shear. Such a pattern is even employed in developing empirical models of flare forecast by calculating the length of the high-gradient PIL (Falconer et al. 2002; Falconer 2003). Thus, we consider why a strong-gradient and long PIL, or simply a strong PIL, is more favorable for producing eruption. To the best of our knowledge, the physical reason behind the correlation of this property of the PIL and the eruption-productiveness has never been explained explicitly. This paper is devoted to answering this question on the basis of the BASIC mechanism and performing a series of 3D MHD simulations similar to those of Jiang et al. (2021), but with different photospheric flux distributions (magnetograms) for the purposes of a comparative study. Some of magnetograms exhibit a strong PIL while others have a weak PIL (i.e., a short and weak-gradient one). As we go on to demonstrate, nearly all the simulations follow the BASIC scenario of quasi-static formation of CS and triggering of eruption by reconnection, which demonstrates the robustness of the BASIC mechanism, and, strikingly, we find the magnitude of the eruption is highly dependent on the strength of the PIL.

This paper is organized as follows. In Sect. 2.1, we define the magnetograms for the bipolar fields with different magnetic flux distributions. Before showing the simulation results, in Sect. 2.2 we analyze the key parameters associated with the non-potentiality and CS of the open field configuration corresponding to these different bipolar fields, since the BASIC mechanism is closely linked to the open field. Then we compare the results of MHD simulations for the different bipolar fields in Sects. 2.3 and 2.4. We give our discussion and conclusion in Sect. 3.

2. Numerical experiments

2.1. Setting of magnetograms

Following Amari et al. (2003a) and Jiang et al. (2021), we modeled the photospheric magnetogram of a bipolar field via the composition of two Gaussian functions,

$$\begin{array}{c}\hfill B_z(x,y,0)=B_0{e}^{-x^2/{\sigma }_x^2}(e^{-(y^2-y_c^2)/{\sigma }_y^2}-e^{-(y^2+y_c^2)/{\sigma }_y^2}),\end{array}$$(1)

where σ_x and σ_y control the extents of the magnetic flux distribution in x and y directions, respectively, and y_c controls the distance between the two magnetic polarities in the y direction. By adjusting these controlling parameters, we obtain different magnetograms with different flux distributions, some of which have a long and strong-gradient PIL while other have short and weak ones. For example, by increasing σ_y but fixing σ_x and y_c, the shape of the polarities gradually changes from flat to round and the strong-gradient PIL becomes shorter, as can be seen in Fig. 1 with four values of σ_y. In another way, by increasing y_c but fixing σ_x and σ_y, the two opposite polarities become further away from each other and, thus, the field gradient across the PIL decreases, as shown in Fig. 2 for four values of y_c. For a reasonable comparison, we used the same value for the total unsigned magnetic flux of all the different magnetograms; as such, the different magnetograms only differ in the distribution of the same amount of magnetic flux. Therefore, the larger the area of magnetic polarity, the weaker the average magnetic field strength – and vice versa. For the specific values of the parameters, we chose [ − 10, 10] × [−10, 10] with a length unit of 14.4 Mm as our computational domain. The values for the parameters as shown in Figs. 1 and 2 are given in their captions.

Fig. 1. Magnetic flux distribution and surface rotation flow at the bottom end (i.e., z = 0). The background is color-coded by the vertical magnetic component Bz, and the vectors show the rotation flow. The four panelsA, B, C, and D refer to magnetogram with σy = 0.5, 1.0, 1.5, and 2.0 respectively, while σx = 2.0 and yc = 0.8 are fixed. We note that the values of σx, σy, and yc are normalized by a unit length of 14.4 Mm.

Fig. 2. Magnetic flux distribution and surface rotation flow at the bottom end. Same details as in Fig. 1 but with different values of yc = 0.1, 0.6, 1.1, and 1.6, respectively, while σx = 2.0 and σy = 1.0 are fixed.

2.2. The open field

The BASIC mechanism is closely linked to the fully open field, as by ideally shearing a magnetic arcade, it asymptotically approaches the open field state. In Fig. 3, we show an example of the open field (and compared it with the corresponding potential field) for a bipolar magnetogram. The procedure of the open field computation is as follows. First, we calculate the potential field B_pot, u(x, y, z) based on a unipolar magnetogram defined by B_z,u(x,y,0) = sign(y)B_z(x,y,0). The potential field is solved using Green’s function method. Thus, all the field lines are open, going outwards from the bottom boundary to infinity. Then, we reverse the direction of the field lines that root in the original negative polarity, and in our case, this can be simply done by defining B_open(x, y, z) = sign(y)B_pot, u(x, y, z) according to the symmetry of the field. In this way, all the field is still open and current free, except at the interface (i.e., the y = 0 plane) between the inverse-directed field lines, which forms the CS. Therefore, this comprises the open field corresponding to the bipolar flux distribution.

Fig. 3. Potential field and the open field with the same magnetic flux distribution on the bottom surface, and the distribution of CS in the open field. A: 3D prospective view of the potential field lines. The colored thick lines represent magnetic field lines and the colors denote the height. The background shows the magnetic flux distribution on the bottom boundary. B: 3D prospective view of the fully opened field lines, with same footpoints shown in panel A. The red iso-surface represents the CS, whose current density is J = 0.34 × 10−3 A m−2. C: current distribution and magnetic field lines on the central cross section (i.e., the x = 0 plane). D: profile of current density along z axis. We note that here with finite grid resolution, the CS has a finite thickness of 360 km. Thus the current density is not infinite.

The open field has important implications on the content of energy that a sheared bipolar field of the same magnetogram can store, as well as the intensity of CS that the sheared bipolar field can form. According to the Aly-Sturrock conjecture, the energy of such open field is the upper limit of the energy of all possible force-free fields with a given magnetic flux distribution on the bottom and a simply connected topology (Aly 1991; Sturrock 1991). Therefore, the upper limit of free magnetic energy that can be reached in a sheared arcade is the open field energy subtracted by the corresponding potential energy of the same magnetogram:

$$\begin{array}{c}\hfill E_{\mathrm{uf}}=E_{\mathrm{open}}-E_{\mathrm{pot}},\end{array}$$(2)

where E_uf denotes the upper limit of the free energy. If we are using the ratio of free magnetic energy to the potential energy as a measure of the non-potentiality of the field, namely, N = E_free/E_pot, the upper limit of the non-potentiality of a sheared arcade can achieve is N_max = E_uf/E_pot. The degree of non-potentiality, N, of a field has been suggested as a critical factor in producing major eruptions (Moore et al. 2012; Sun et al. 2015), and, thus, N_max is an important indicator for the eruption capability of a magnetogram. For example, Moore et al. (2012) studied a large number of active regions and found that the free energy the field can hold has a sharp upper limit (which increases with the active region’s magnetic flux content) and that most active regions approach this limit when producing a coronal mass ejection or flare eruption. In particular, using a proxy of magnetic shear for the free energy, they concluded that the non-potentiality, N, is on the order of one for the active region’s core field, that is, the field rooted around the flare PIL, when the field is close to eruption. Therefore, for a sheared bipolar field that can produce major eruption, it should have its upper limit of non-potentiality somewhat close or above one, namely, N_max ≳ 1, and by calculating this parameter for a given magnetogram, we can evaluate whether it is capable of generating a significant eruption or not.

The potential field B_pot and the open field B_open are uniquely defined by:

$$\begin{array}{c}\hfill B_{z|\mathrm{pot}}(x,y,0)=B_{z|\mathrm{open}}(x,y,0)=B_z(x,y,0),\end{array}$$(3)

and an asymptotic decay at infinity. These fields have energies given, respectively, via the standard relations (e.g., Amari et al. 2003a)

$$\begin{array}{cc}& E_{\mathrm{pot}}=\frac{1}{16{\pi }^2}{\displaystyle {\int }_{S\times S^{\prime }}}\frac{B_z(x,y,0)B_z(x^{\prime },y^{\prime },0)}{|\mathit{r}-{\mathit{r}}^{\prime }|}\mathrm{d}s\mathrm{d}{s}^{\prime },\hfill \\ \hfill & E_{\mathrm{open}}=\frac{1}{16{\pi }^2}{\displaystyle {\int }_{S\times S^{\prime }}\frac{|B_z(x,y,0)B_z(x^{\prime },y^{\prime },0)|}{|\mathit{r}-{\mathit{r}}^{\prime }|}\mathrm{d}s\mathrm{d}{s}^{\prime }.}\hfill \end{array}$$(4)

The upper limit of free magnetic energy can be obtained by Eq. (2), and then the non-potentiality, N_max, of each magnetogram can be calculated. Furthermore, the magnetic field (potential field and open field) in volume can be obtained from the bottom magnetogram, and then the magnetic energy can be obtained by integration $E=\frac{1}{8\pi }\int B^2\mathrm{d}V$. However computing the magnetic field in the full volume is very time consuming. For example, with Green’s function method, the computing time scales with the grid number as N⁵ (assuming the volume is a cube and the length of each side is N). That is why we chose to use Eq. (4), which is much faster in its calculation, since the magnetic energy can be obtained directly from the magnetogram without knowing the magnetic field in the full volume. The magnetic energy obtained by these two methods is basically the same. For the magnetogram in Fig. 1A, the potential field energy and open field energy obtained by using the first method are 8.923 × 10²⁹ erg and 2.059 × 10³⁰ erg, respectively, and the ones obtained by Eq. (4) is 8.603 × 10²⁹ erg and 2.090 × 10³⁰ erg, respectively, with relative errors of 3.586% and 1.483%, respectively.

In the open field, all the magnetic field lines have one end rooted at the bottom surface and the other end extending up to infinity, so the field lines on the two side of the PIL run antiparallel, and in between them, a CS forms (see Fig. 3B and C). The current density outside CS is zero, and all free magnetic energy is stored through the CS (but not in the CS). In the 3D prospective view of the open field (as shown in Fig. 3), the structure of the CS is denoted by the red iso-surface of current density J = 0.34 × 10⁻³ A m⁻². The central cross section of the open field and the profile of the current density along the central vertical line are also shown in Fig. 3. Since the magnetic field is discretized with a finite resolution, the current density in the CS is not infinite but rather changes with height; it first increases from nearly zero, reaches a peak value at a certain height, and then decreases towards zero again. Therefore, the maximum current density and its height can be obtained by calculating the current density in the CS. We consider that the CS of open field can be used as a proxy of the CS formed in the core field in our simulations and, in particular, the location of the maximum current density in the open field CS indicates the position where reconnection most likely starts to trigger an eruption and the maximum current density itself should be related to the reconnection rate. Therefore, the two parameters, namely: the peak value of current density of the open field CS and its height, are essential indicators of the strength of the CS, which is correlated with the strength of the eruption once being triggered. In order to calculate the current density of the open field CS, we first calculate the open field by the Green’s function method and then use the second-order central difference to get the current density in the CS.

Figure 4 shows the result for different magnetograms specified by the parameter σ_y changing from 0.5 to 2.5 with an increment of 0.1 (while σ_x = 2.0 and y_c = 0.8 are fixed) and Fig. 5 shows the results for magnetograms with a parameter y_c changing from 0.1 to 2.1, with an increment of 0.1 (with σ_x = 2.0 and σ_y = 1.0 fixed). Specifically, Fig. 4A shows the variation of different magnetic energies and the non-potentiality, N_max, with the parameter σ_y. With the increase in σ_y, all energies decrease, and N_max also shows an overall decrease from about 1.5 to below 1.0, which clearly shows that, overall, the more elongated the magnetic polarity, the higher the upper limits of free energy and non-potentiality the magnetogram can reach (except that a too small value of σ_y actually reduces N_max, which is discussed below). Figure 5A shows the variation of magnetic field energies and the non-potentiality N_max with parameter y_c. A similar decrease pattern is also seen (except that the potential energy increases mildly) with the increase of y_c, which shows that the closer the two magnetic poles, the more free magnetic energy (and non-potentiality) the magnetogram can reach.

Fig. 4. Magnetic energies, non-potentiality, Nmax, and distribution of CS varies with σy (with σx = 2.0 and yc = 0.8 fixed). A: black, blue and green curves represent potential energy, open field energy, and free magnetic energy, respectively. Red represents the non-potentiality, Nmax. B: black and red curves respectively represent the height and value of the maximum current density.

We note that in Fig. 4A, we can see that N_max has a maximum at σ_y ∼ 0.8 and a smaller σ_y actually reduces N_max. This is because as the polarity distance is fixed at y_c = 0.8, a too low value of σ_y ≤ 0.8, that is, too much concentration of the polarity in the y direction, will actually cause the strong magnetic flux to be distributed farther away from the PIL; thus the magnetic field gradient near the PIL would be reduced. This hints at the gradient on the PIL acting as a more important indicator for N_max. To confirm this, we calculate the line integral of gradient of the vertical field B_z across the PIL for each magnetogram, named L, which is defined by

$$\begin{array}{c}\hfill L={\displaystyle {\int }_{\mathrm{PIL}}\frac{\partial B_z}{\partial y}\mathrm{d}x.}\end{array}$$(5)

Figure 6A shows the diagram of L versus N_max for all the different sets of parameters shown in Figs. 4 and 5. Strikingly, N_max increases monotonically with the increase of L, which confirms our argument. Hereafter, we refer to the parameter L as simply the strength of the PIL. We note that L for σ_y = 0.5 is less than that of σ_y = 0.8 and, thus, the N_max value of the former is less than that of the latter. We also note that with the decrease in L, the value of N_max decreases more and more quickly.

Fig. 5. Magnetic energies, non-potentiality, Nmax, and distribution of CS varies with yc (with σx = 2.0 and σy = 1.0 fixed). Same details as in Fig. 4.

For the characterization of the CS, Fig. 4B shows that with the increase in the value of σ_y, the peak value of current density decreases and its height increases. This indicates that the more elongated the magnetic polarity, the stronger the CS can form. In Fig. 5B, with the increase of y_c, the maximum current value decreases and its height increases. This shows the closer the magnetic polarities, the stronger the CS can form. Again, the PIL strength, L, is crucial. As can be seen in Fig. 6, with the increase of L, the peak value of current density increases monotonically overall (except some points due the discretization errors), and the height decreases systematically, which means that the PIL strength is positively correlated with the strength of the CS.

Fig. 6. Relationship between the non-potentiality Nmax and distribution of CS and the line integral of the Bz gradient across the PIL on the magnetogram. A: diamond represents magnetic energies, and black, blue, and green represent the potential energy, open field energy, and free magnetic energy, respectively. Red solid circle represents non-potentiality, Nmax. The three dashed vertical lines represent the values of three different parameters σy (with σx = 2.0 and yc = 0.8 fixed). Magenta, cyan, and gold represent σy = 0.5, 0.7, and 0.8, respectively. B: black squares represents the maximum current density and the red solid curves represents the height of the maximum current density.

2.3. Settings of the numerical model

For each magnetogram, we carried out an MHD simulation which begins with the potential field and is driven continually by a rotational flows at each magnetic polarity at the lower boundary, which creates magnetic shear along the PIL. The rotational flow is defined as

$$\begin{array}{c}\hfill v_x=\frac{\partial \psi (B_z)}{\partial y};v_y=\frac{\partial \psi (B_z)}{\partial x},\end{array}$$(6)

with ψ given by

$$\begin{array}{c}\hfill \psi =v_0{B}_z^2{e}^{-(B_z^2-B_{z,\max }^2)/B_{z,\max }^2},\end{array}$$(7)

where B_z, max is the largest value of the photosphere B_z, and v₀ is a constant for scaling such that the maximum of the surface velocity is 4.4 km s⁻¹, close to the typical flow speed in the photosphere (∼1 km s⁻¹). The flow speed is smaller than the sound speed by two orders of magnitude and the local Alfvén speed by three orders, respectively, thus representing a quasi-static stress of the coronal magnetic field. The flow pattern is shown in Figs. 1 and 2. This is an incompressible, anti-clockwise rotational flow that does not change with time and it will not modify the flux distribution at the bottom end.

We numerically solved the full MHD equations with both coronal plasma pressure and solar gravity included, in a 3D Cartesian geometry by an advanced conservation element and solution element (CESE) method implemented on an adaptive mesh refinement (AMR) grid (Jiang et al. 2010, 2016, 2021; Feng et al. 2010). Since the controlling equations, the numerical code, and, essentially, all the setting of initial and boundary conditions are the same as used in Jiang et al. (2021), we refer to that paper (in particular, the section on methods) for details. We note that no explicit resistivity is used in the magnetic induction equation, but magnetic reconnection can still occur due to numerical diffusion when a current layer is sufficiently narrow such that its width is close to the grid resolution.

The computational volume spans a Cartesian box of approximately ( − 270, −270, 0) Mm ≤(x, y, z)≤(270, 270, 540) Mm (where z = 0 represents the solar surface). The volume is large enough such that the simulation runs can be stopped before the disturbance from the simulated eruption reaches any of these boundaries. The full volume is resolved by a block-structured grid with AMR in which the base resolution is Δx = Δy = Δz = Δ = 2.88 Mm, and the highest resolution of Δ = 360 km is used to capture the formation of CS and the subsequent reconnection.

In total, we ran nine experiments for a selected set in the parameter space as listed in Table 1. Specifically, CASE I to CASE IV represent σ_y = 0.5, 1.0, 1.5, and 2.0, respectively (with fixed σ_x = 2.0 and y_c = 0.8). CASE V to CASE VIII represent y_c = 0.1, 0.6, 1.1, and 1.6, respectively (with fixed σ_x = 2.0 and σ_y = 1.0). In order to verify the analysis of free magnetic energy maximum in Fig. 4, we set CASE IX, where σ_x = 2.0, σ_y = 0.8 and y_c = 0.8.

2.4. Simulation results

The BASIC scenario as demonstrated in Jiang et al. (2021) is that by continuously shearing a bipolar coronal field, a CS forms slowly within the arcade and once reconnection sets in, the whole arcade explodes and forms a fast-ejecting magnetic flux rope (i.e., CME). Here, we first show briefly an example of such process in CASE I and then we analyze the different runs based on the different magnetograms.

Figure 7 and online Movie 1 show evolution of magnetic field lines, current density, and velocity during the whole simulation process. The time unit is τ = 105 s (all the times mentioned in this paper are expressed with the same time unit). After a period of surface flow driving, the magnetic field structure evolves from the initial potential field to a configuration with a strong shear immediately above the PIL, where the current density is enhanced, while the envelope field is still current-free (t ≤ 76). It can be clearly seen that the entire magnetic field structure inflates during the energy injection stage, since the magnetic pressure of the core field increases gradually by the continuous shear. As a result, it stretches the envelope field outward, making the bipolar magnetic arcade tend to approach an open field configuration (Fig. 7B). During this quasi-static evolution process, the current is squeezed from the volumetric distribution into a vertical, narrow layer extending above the PIL, forming a vertical CS (Fig. 7C).

Fig. 7. Evolution of magnetic field lines, electric currents, and velocity in the whole simulation process of CASE I. A: top view of magnetic field lines. The colored thick lines represent magnetic field lines and the colors denote the value of nonlinear force-free factor defined as α = J ⋅ B/B2, which indicates the extent to which the field lines are of non-potential character. The background shows the magnetic flux distribution on the bottom boundary (i.e., plane of z = 0), and contours of Bz = ( − 48, −32, −16, 0, 16, 32, 48) G are shown. B: 3D prospective view of the same field lines shown in panel A. C: current density, J, normalized by magnetic field strength, B, in vertical cross section (i.e., the x = 0 slice). D: same as C, but with a large area. E: magnitudes of velocity. The largest velocity and Alfvénic Mach number are also denoted.

As a critical point, when the thickness of the CS decreases down to the grid resolution, magnetic reconnection sets in and triggers an eruption. This transition from the pre-eruption to eruption onset is clearly manifested in the evolution of energies, as shown in Fig. 8 (see the curves colored in magenta), which have a sharp transition at t = 78. The kinetic energy increases impulsively to nearly 7% original magnetic potential energy in a time duration of Δt = 5. Meanwhile, the magnetic energy releases quickly during the eruption. The onset of the eruption can be more clearly shown by the time profiles of the magnetic energy release rate and the kinetic energy increase rate, and both of them have a sharp increase at the beginning of the eruption (Fig. 8B).

Fig. 8. Temporal evolution of magnetic energy and kinetic energy in the simulations of different σy. The color of the lines represents five experiments, magenta, red, cyan, blue, and green represent CASE I, IX, II, III, IV, respectively. A: evolution of magnetic energy EM (solid lines) and kinetic energy EK (dotted lines). B: releasing rate of magnetic energy and increasing rate of kinetic energy of five experiments. The vertical dashed lines are shown for denoting the transition time of from pre-eruption to eruption. In the ordinates of all panels, the normalized unit is the potential field energy based on each magnetogram.

With the onset of reconnection, a plasmoid (i.e., MFR in 3D) originates from the tip of the CS and rises quickly, leaving behind a cusp structure separating the reconnected, post-flare loops from the un-reconnected field (Fig. 7C and D). The plasmoid expands quickly and meanwhile, an arc-shaped fast magnetosonic shock is formed in front of the plasmoid. All of these evolving structures are proof of the typical coronal magnetic eruption leading to CME. The shock marks the front edge of the CME and its average speed is about 603.4 km s⁻¹ (Figs. 7D and 9).

Fig. 9. Time stack map of the current distribution around x, y = 0 in CASE I, which can reveal the evolution speed of the CME.

Figures 8 and 10 show the temporal evolution of magnetic and kinetic energies (and their changing rate) of all the cases from CASE I to IX. We note that in each CASE, the energies are normalized by the corresponding potential field energy, namely, the magnetic energy at t = 0. Overall, all the different runs (except CASE VIII) show a similar evolution pattern: magnetic energy first increases monotonically for a long time, approaching the open field energy, while the kinetic energy remains very low level; at a critical point, the magnetic energy begins to decrease rapidly along with an impulsive rise of the kinetic energy and the evolutions of the two energies are closely correlated in time, which indicates that the free magnetic energy is released and accelerates the plasma. In the online Movies 2 and 3, we show the evolution of current density on the central cross-section for all the cases. As can be seen, they all follow the same BASIC scenario; that is, first a CS forms during the magnetic energy increasing phase and then reconnection sets in and triggers eruption. Therefore, these different runs demonstrate the robustness of the BASIC mechanism.

Fig. 10. Temporal evolution of magnetic energy and kinetic energy in the simulations of different yc. Same details as in Fig. 8 but with different yc. The color of the lines represents four experiments, magenta, cyan, blue, and green represent CASE V, VI, VII, and VIII, respectively. We note that CASE VIII did not produce an eruption.

Nevertheless, the magnitudes (or intensities) of the eruptions in the different cases are different, as can be seen by comparing the energy conversion rates during the impulsive phase. For example, in the five experiments with increasing σ_y (as shown in Fig. 8), the maximum release rate of magnetic energy increases first, reaching the largest at σ_y = 0.8, and then it decreases with higher σ_y, which is exactly consistent with the dependence of N_max on σ_y as shown in Fig. 4. Again, Fig. 10 shows that with the increase of y_c, the eruption intensity decreases, consistent with the dependence of N_max on y_c as shown in Fig. 5. In Fig. 11A and B, we further show how the strength of the PIL, that is, L, is related to the intensity of the eruption, as quantified by the peak values of the kinetic energy increasing rate and the magnetic energy releasing rate, as well as the speed of the leading edge of the CME (i.e., the shock). It clearly shows that the eruption intensity is correlated positively with the PIL strength. Furthermore, in Fig. 11C, we show the non-potentiality, N, at the onset time of the eruption, as compared with the corresponding N_max. The non-potentiality increases overall with increase of the PIL strength, consistent with (but a bit slower than) that of the N_max. Interestingly, we note that the non-potentiality at the eruption onset is mostly close to one, which is in striking agreement with the statistical study of observed active regions (Moore et al. 2012) and it also hints that the BASIC mechanism is responsible for those eruptions.

Fig. 11. Intensity of the eruption, as quantified by (A) the peak value of the kinetic energy increasing rate and the magnetic energy releasing rate as well as (B) the largest speed of the leading edge of the CME, depend on L. C: non-potentiality, N, at the eruption onset time. Dotted line represents the theoretical calculation value, same as in Fig. 6. D: location of CS and the maximum current density in CS at the eruption onset time.

Figure 12 shows the central vertical cross section of current density, J, (normalized by magnetic field strength B) at a time that is immediately close to the eruption onset time for the different experiments. The location of the CS and the maximum current density in the CS are also shown in Fig. 11D. These results are consistent with the calculation of the open field in Sect. 2.2. Furthermore, we find that the lower the position of the CS and the higher the current density, the greater the eruption intensity (Fig. 11).

Fig. 12. Vertical cross-section of the current density, J, (normalized by magnetic field strength, B) at the eruption onset time for the eight experiments. From left to right and top to bottom: CASE IX, II, III, IV, V, VI, VII, and VIII. The location of the CS and the maximum current density in CS are denoted on each panel. We note that there is no eruption for CASE VIII.

The only case in our experiments that did not reach an eruption is CASE VIII. This is because during the quasi-static shearing process, the field expands fast in the later phase (e.g., t > 100) and strongly presses upon the numerical boundaries before a CS is formed (or before a sufficient amount of free magnetic energy is accumulated to approach a open field); therefore, we had to stop the simulation run to avoid too much influence coming from the numerical boundaries upon the results. Ideally, with a larger computational box, a free expansion of the field driven by the surface shearing motion will create a CS, but its height is too great (and the current density is too low) to trigger an efficient eruption.

It is interesting to note that the ratio E_open/E_pot of around 1.7 is apparently an adequate threshold for starting an eruption. A similar ratio has been found in Amari’s simulations on flux rope formation and instability (Amari et al. 2000, 2003a,b). Actually, for all the different distributions of magnetic flux that we consider in this paper, the lowest ratio of E_open/E_pot is 1.7 (see Fig. 6A, where N_max = E_open/E_pot − 1). The reason why E_open/E_pot of 1.7 appears to be a threshold for an eruption to start is that in our scenario, the CS can only form (and thus trigger an eruption) when the magnetic field is sufficiently sheared, such that its energy is close to the open field energy. Furthermore, as we have deduced, the E_open/E_pot (or N_max) should be large enough to let the CS to form at a low height and with a large current density, such that the reconnection can be efficient to produce an eruption. Otherwise, if the E_open/E_pot is too small, the free energy that can be attained is small, and the CS will form at a too great a height (and with too low current density) to trigger an efficient eruption, as our experiment CASE VIII shows. The smaller the ratio of E_open/E_pot, the higher the CS is expected to form and the smaller the free energy that can be reached; thus, the less efficiently the eruption can be produced. For example, in the extreme case when E_open/E_pot reaches its lower limit of E_open/E_pot = 1, which corresponds to the flux distribution with two opposite polarities being infinitely far away from each other, it cannot produce any eruption no matter how strongly the polarities are rotated.

While using a different scenario (i.e., first assuming a flux rope is created by surface converging or cancellation (or both) and then the flux rope runs into instability to initiate an eruption), Amari’s simulations also show the same behavior: that their eruption occurs only when the magnetic energy is close to the open field energy. Thus, both of our and Amari’s simulations using a single bipolar flux distribution suggest that the ratio E_open/E_pot should be larger than 1.7 to start an eruption. Amari’s simulations consist of two important phases of energizing. Their first phase is the same as ours, namely, via the rotation of the polarities to inject free energy into the field. The key difference is that in the simulations in Amari et al. (2000, 2003a,b), the rotation is stopped before a CS is formed; then, in the second phase, they modified the flux content by opposite flux emerging (Amari et al. 2000) or surface diffusion (Amari et al. 2003a), or they modified the flux distribution by surface converging flow (Amari et al. 2003b) at the bottom boundary. In this phase, the magnetic topology will be changed from a sheared arcade to a flux rope, through the slow reconnection near the bottom boundary. More importantly, the corresponding open field energy of the evolving magnetic flux distribution will change, as it is faster than that of the total magnetic energy and can eventually lead to an eruption when the total magnetic energy is close to the open field energy.

Finally, we note that the eruption onset times of the different experiments are different. This is related to the magnetic energy injection rate, which depends on the surface flow distribution. As shown in Fig. 8, the magnetic energy injection rate decreases when σ_y increases, and, thus, the eruption onset time is systematically postponed because a longer time is needed for free magnetic energy accumulation.

3. Conclusions

It has long been known that major solar eruptions mostly occur in active regions with a strongly sheared and strong-gradient PIL. There is no doubt that a strong magnetic shear is critical for producing an eruption, since it is directly related to the degree of non-potentiality of the field. However, it lacks an explanation as to why the flux distribution with a high-gradient PIL is favorable for eruption. In this paper, we provide such a physics explanation, for the first time, based on the BASIC mechanism with different photospheric magnetic flux distributions, namely, magnetograms, by combining theoretical analysis and numerical simulations. The BASIC mechanism refers to a simple and efficient scenario in which an internal CS can form slowly in a gradually sheared bipolar field and the reconnection of the CS triggers and drives the eruption (Jiang et al. 2021).

In principle, two requirements are essential to initiate a major eruption by the BASIC mechanism and they are closely related to each other. One is to accumulate a sufficient amount of free magnetic energy to power a major eruption; specifically, the field should be sufficiently non-potential in nature. The other is to build up a strong CS in the bipolar core field, that is, a CS with a high current density and formed at a low height, such that reconnection can release magnetic energy efficiently. Focusing on these two key elements, we set up a series of magnetograms with equaling unsigned flux but different flux distributions. We first analyzed the open fields corresponding to these magnetograms. This is because the open field sets an upper limit for the energy that a sheared bipolar field can store as well as the intensity of CS that the field can form. By calculating the largest non-potentiality N_max and the peak current density based on the open field, we find that magnetogram having a stronger PIL can contain more non-potentiality, N_max, and can form stronger CS, which indicates the capability to initiate a larger eruption by the BASIC mechanism. Furthermore, we find that the strength of the PIL, denoted L, can be quantified well by a line integral of a gradient of the vertical field B_z across the PIL, which should be valuable in future studies for flare forecasts based on magnetograms.

Then we selected nine representative magnetograms to conduct MHD simulations. All of the numerical experiments exhibit the same evolution pattern; magnetic energy first increases monotonically for a long time, as driven by the boundary rotational flow, and in that period, the kinetic energy remains within a very low level, indicating that the system is undergoing a quasi-static evolution process. Then, at a critical point, when the thickness of the CS decreases down to the grid resolution, reconnection sets in and triggers an eruption, during which the magnetic energy decrease rapidly along with an fast rise of the kinetic energy. An overall comparison of the eruptions in the different cases shows a strong correlation of the eruption intensity with the strength of the PIL. Specifically, with the increase of the PIL strength, both the non-potentiality of the field and the strength of the CS at the eruption onset increase, and, consequently, the eruption intensity increases, which confirms the two key conditions in the BASIC mechanism.

In summary, through the combined study of theoretical analysis and numerical simulations, we demonstrated that the bipolar field with magnetogram of a strong PIL can hold more non-potentiality and can form stronger internal CS, which are key to the initiation of strong eruption. This is the physical reason why a magnetic field with a strong PIL is capable of producing major eruption. Our study demonstrates the robustness of the BASIC mechanism, on the one hand. On the other, it also discloses the physics explanation of why a long and strong-gradient PIL is favorable for major eruption.

Movies

Movie 1 associated with Fig. 7 Access here

Movie 2 associated with Fig. 7 Access here

Movie 3 associated with Fig. 7 Access here

Acknowledgments

This work is jointly supported by National Natural Science Foundation of China (NSFC 42174200, 41822404, 41731067), and Shenzhen Technology Project JCYJ20190806142609035. The computational work was carried out on TianHe-1(A), National Supercomputer Center in Tianjin, China.

References

Appendix A: Captions for online Movies

• Caption for online Movie 1: the evolution of magnetic field lines, current, and velocity in the whole simulation process of CASE I. Details are the same as in Fig. 7.

• Caption for online Movie 2: the evolution of current in the whole simulation process of CASE I, IX, II, III, and IV. Same details as in Fig. 7D.

• Caption for online Movie 3: the evolution of current in the whole simulation process of CASE V, VI, VII, and VIII. Same details as in Fig. 7D.

All Tables

All Figures

Fig. 1. Magnetic flux distribution and surface rotation flow at the bottom end (i.e., z = 0). The background is color-coded by the vertical magnetic component Bz, and the vectors show the rotation flow. The four panelsA, B, C, and D refer to magnetogram with σy = 0.5, 1.0, 1.5, and 2.0 respectively, while σx = 2.0 and yc = 0.8 are fixed. We note that the values of σx, σy, and yc are normalized by a unit length of 14.4 Mm.

In the text

	Fig. 2. Magnetic flux distribution and surface rotation flow at the bottom end. Same details as in Fig. 1 but with different values of yc = 0.1, 0.6, 1.1, and 1.6, respectively, while σx = 2.0 and σy = 1.0 are fixed.
In the text

Fig. 3. Potential field and the open field with the same magnetic flux distribution on the bottom surface, and the distribution of CS in the open field. A: 3D prospective view of the potential field lines. The colored thick lines represent magnetic field lines and the colors denote the height. The background shows the magnetic flux distribution on the bottom boundary. B: 3D prospective view of the fully opened field lines, with same footpoints shown in panel A. The red iso-surface represents the CS, whose current density is J = 0.34 × 10−3 A m−2. C: current distribution and magnetic field lines on the central cross section (i.e., the x = 0 plane). D: profile of current density along z axis. We note that here with finite grid resolution, the CS has a finite thickness of 360 km. Thus the current density is not infinite.

In the text

	Fig. 4. Magnetic energies, non-potentiality, Nmax, and distribution of CS varies with σy (with σx = 2.0 and yc = 0.8 fixed). A: black, blue and green curves represent potential energy, open field energy, and free magnetic energy, respectively. Red represents the non-potentiality, Nmax. B: black and red curves respectively represent the height and value of the maximum current density.
In the text

	Fig. 5. Magnetic energies, non-potentiality, Nmax, and distribution of CS varies with yc (with σx = 2.0 and σy = 1.0 fixed). Same details as in Fig. 4.
In the text

Fig. 6. Relationship between the non-potentiality Nmax and distribution of CS and the line integral of the Bz gradient across the PIL on the magnetogram. A: diamond represents magnetic energies, and black, blue, and green represent the potential energy, open field energy, and free magnetic energy, respectively. Red solid circle represents non-potentiality, Nmax. The three dashed vertical lines represent the values of three different parameters σy (with σx = 2.0 and yc = 0.8 fixed). Magenta, cyan, and gold represent σy = 0.5, 0.7, and 0.8, respectively. B: black squares represents the maximum current density and the red solid curves represents the height of the maximum current density.

In the text

Fig. 7. Evolution of magnetic field lines, electric currents, and velocity in the whole simulation process of CASE I. A: top view of magnetic field lines. The colored thick lines represent magnetic field lines and the colors denote the value of nonlinear force-free factor defined as α = J ⋅ B/B2, which indicates the extent to which the field lines are of non-potential character. The background shows the magnetic flux distribution on the bottom boundary (i.e., plane of z = 0), and contours of Bz = ( − 48, −32, −16, 0, 16, 32, 48) G are shown. B: 3D prospective view of the same field lines shown in panel A. C: current density, J, normalized by magnetic field strength, B, in vertical cross section (i.e., the x = 0 slice). D: same as C, but with a large area. E: magnitudes of velocity. The largest velocity and Alfvénic Mach number are also denoted.

In the text

Fig. 8. Temporal evolution of magnetic energy and kinetic energy in the simulations of different σy. The color of the lines represents five experiments, magenta, red, cyan, blue, and green represent CASE I, IX, II, III, IV, respectively. A: evolution of magnetic energy EM (solid lines) and kinetic energy EK (dotted lines). B: releasing rate of magnetic energy and increasing rate of kinetic energy of five experiments. The vertical dashed lines are shown for denoting the transition time of from pre-eruption to eruption. In the ordinates of all panels, the normalized unit is the potential field energy based on each magnetogram.

In the text

	Fig. 9. Time stack map of the current distribution around x, y = 0 in CASE I, which can reveal the evolution speed of the CME.
In the text

	Fig. 10. Temporal evolution of magnetic energy and kinetic energy in the simulations of different yc. Same details as in Fig. 8 but with different yc. The color of the lines represents four experiments, magenta, cyan, blue, and green represent CASE V, VI, VII, and VIII, respectively. We note that CASE VIII did not produce an eruption.
In the text

Fig. 11. Intensity of the eruption, as quantified by (A) the peak value of the kinetic energy increasing rate and the magnetic energy releasing rate as well as (B) the largest speed of the leading edge of the CME, depend on L. C: non-potentiality, N, at the eruption onset time. Dotted line represents the theoretical calculation value, same as in Fig. 6. D: location of CS and the maximum current density in CS at the eruption onset time.

In the text

	Fig. 12. Vertical cross-section of the current density, J, (normalized by magnetic field strength, B) at the eruption onset time for the eight experiments. From left to right and top to bottom: CASE IX, II, III, IV, V, VI, VII, and VIII. The location of the CS and the maximum current density in CS are denoted on each panel. We note that there is no eruption for CASE VIII.
In the text