Open Access
Issue
J. Space Weather Space Clim.
Volume 7, 2017
Article Number A17
Number of page(s) 25
DOI https://doi.org/10.1051/swsc/2017015
Published online 01 August 2017

© C. Guennou et al., Published by EDP Sciences 2017

Licence Creative Commons
This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

1. Introduction

Solar flares are sudden brightenings of the solar atmosphere occurring over the whole spectral range, from radio to X-rays, resulting from the release of a huge amount of magnetic energy. Flares can be classified in terms of the strength of the X-ray flux measured close to Earth. They can be divided into two general subcategories: eruptive flares, accompanied by a coronal mass ejection (CME), large ejection of magnetized plasma released into the interplanetary space, and confined flares, for which only electromagnetic radiation is emitted. Both CMEs and solar flares can impact the near-Earth environment and in particular human technologies. The most energetic solar flares can heat and ionize the upper Earth atmosphere, engendering disruption in radio communication and Global Positioning System (GPS) inaccuracies. CMEs can be accompanied by the burst of energetic particles and can be associated with strong geomagnetic storms, inducing ground-level electric fields likely to produce power grid damages. Astronauts, pilots, and satellites can also be affected by radiations originating from flares or CMEs. The ability to accurately predict both solar flares and CMEs is therefore fundamental to protect both space and ground-based technologies from strong space-weather events.

The European H2020 research project FLARECAST (Flare Likelihood and Region Eruption Forecasting – http://flarecast.eu/) aims to develop a fully automated solar flare forecasting system. FLARECAST will automatically extract various magnetic field parameters of solar active regions (ARs) from vector solar magnetograms to produce accurate predictions using the state-of-the-art forecasting techniques based on data-mining and machine learning. Various systems of prediction invoking different categories of models have been developed in the past decades, as e.g., the “Theophrastus” tool (McIntosh 1990), the linear-prediction system of Gallagher et al. (2002) used by SolarMonitor, the discriminant analysis of Barnes et al. (2007), the Automated Solar Activity Prediction (ASAP) of Colak & Qahwaji (2009), based on machine learning, and more recently the statistical learning technique of Yuan et al. (2010). A recent comparison between the current forecasting tools using line-of-sight (LOS) magnetograms has been performed by Barnes et al. (2016), showing that none of them substantially outperformed all others. All these forecasting techniques, including the future FLARECAST forecasting tool, require scalar quantities derived from photospheric magnetic field to be able to make flare and/or eruption predictions. In this paper, we focus on the predictability capabilities of these scalar quantities, in order to determine the most reliable for flare and eruption forecasting use.

It is now widely accepted that the energy source of solar flares is stored in highly non-potential magnetic fields. The storage of energy, typically estimated in the 1028–1032 erg range (Schrijver et al. 2012), arises from a long phase of magnetic stress and free magnetic energy build-up before sudden flare or/and eruption. A large variety of models have been recently developed to explain this storage and released mechanism (see e.g., the recent reviews of Aulanier 2014; Janvier et al. 2015; Lin et al. 2015; Schmieder et al. 2015), but it is still unclear why some ARs trigger a flare, while others remain quiet. In this context of absence of a clear physical scenario for flare triggering, efforts have been concentrated on empirical flare prediction methods, investigating the behavioral patterns of ARs.

Mapping the 3D coronal magnetic field is not systematically practicable, while the photospheric magnetic field is more easily measurable using spectro-polarimeters such as the Helioseismic and Magnetic Imager (HMI; Schou et al. 2012) on board the Solar Dynamics Observatory (SDO). Consequently, flare forecasting tools are mainly based on quantities derived directly from solar magnetograms. As flare occurrence is correlated with size, variability, and complexity of the ARs, i.e., the associated non-potentiality degree, many attempts have been made in order to find reliable photospheric eruptivity indicators. Global photospheric features, such as magnetic field gradient, currents, magnetic geometry, and magnetic free energy, have been proposed to establish a link between photospheric observations and coronal activity.

However, it has been shown that observing clear pre-signature of flare and eruption in observations is challenging. Many studies investigated the link between changes in some photospheric parameters and the coronal eruptive and flaring activity, by investigating prior temporal changes, performing superposed epoch analysis, or forecasting the likelihood of the flaring events (see e.g., Bao et al. 1999; Leka & Barnes, 2003a; Schrijver, 2007; Jing et al. 2010; Mason & Hoeksema, 2010; Falconer et al. 2011; Ahmed et al. 2013; Bobra & Ilonidis, 2016, and references therein), with moderate results. More recently, Al-Ghraibah et al. (2015) studied the magnetic parameters of about 2000 ARs to search for the best eruptive predictor. They used a wide range of parameters, including a wavelet analysis to resolve multiple-scale changes, and found that that magnetic field properties alone are not sufficient for powerful flare forecasting, although the magnetic-gradient related features appear to be one of the most reliable indicators for flare predictions.

Despite all these previous studies, none of the current photospheric indicators, or combinations thereof, have yet demonstrated an unambiguous eruptive or flaring criterion. However, controlled cases (e.g., originating from numerical datasets) have barely been used to test the predictability of these parameters. Kusano et al. (2012) presented a parametric analysis of eruption onset, by varying two parameters associated with the magnetic structure. They showed that these two parameters are able to discriminate between eruptive and non-eruptive ARs, although no comprehensive scalar quantity directly measurable is provided. In this work, we use MHD numerical simulations of the formation of stable and unstable magnetic flux ropes (Leake et al. 2013, 2014) in order to systematically investigate the pre-eruptive signature included in different magnetic parameters. This series of numerical experiments is based on the emergence of a convection zone magnetic flux tube into a solar atmosphere. The interaction of the emerging magnetic flux with the pre-existing magnetic field plays a key role for triggering impulsive events. This class of simulations is thus able to explain and reproduce multiple active solar phenomena (see e.g., the review of Cheung & Isobe 2014). Time series of magnetograms from parametric simulations of stable and unstable flux emergence, i.e., corresponding to respectively quiet AR versus eruptive flare configurations, are used to compute a large range of parameters. This list includes parameters previously used for operational forecasting, as well as parameters used for the first time for this purpose, such as helicity (Pariat et al. 2017), the current-weighted magnetic inversion line, and the length of the strong-current magnetic inversion line portions.

The remainder of this paper is as follows. The set of eruptive and non-eruptive simulations is summarized in Section 2. Section 3 describes the magnetograms and the magnetic parameters extracted from our simulation sets. The results are presented in Section 4, and a discussion on the impact of data masking on the analysis is proposed in Section 5. A basic analysis of the noise impact on our results is presented in Section 6. Section 7 provides a parametric study of the inversion line related parameters and the influence of the different thresholds on eruption predictability. Finally, Section 8 summarizes our work and presents our main conclusions.

2. Numerical datasets: parametric simulations of eruptive and non-eruptive active regions

2.1. MHD simulations

In order to test the reliability of photospheric eruptive predictors, we used the three-dimensional visco-resistive magnetohydrodynamic (MHD) simulations of Leake et al. (2013, 2014). These simulations describe the partial emergence process of a twisted magnetic flux tube into a stratified solar atmosphere, where a coronal arcade field is present. Both stable and non-stable flux ropes are formed as a result of the flux emergence, depending on the choice of certain parameters. These simulations are thus analogous to respectively quiet ARs and eruptive flare productive ARs, where the newly formed flux ropes are ejected higher in the solar corona.

The evolution of the system is described by the visco-resistive MHD equations (see Eqs. (1)–(4) from Leake et al. 2013), and the plasma is assumed to be fully ionized. The MHD equations are solved using the Lagrangian-remap code Lare3D (Arber et al. 2001), using an irregular Cartesian grid. The initial conditions consist of a hydrostatic background atmosphere, stratified such as the solar convection zone, the photosphere/chromosphere, the transition region, and the corona. An arcade field covering the entire simulation domain is imposed on this background atmosphere (cyan solid lines in Fig. 1), and a right-hand twisted flux tube is inserted in the solar convection zone (red solid lines in Fig. 1), aligned along the y axis. The arcade field is transitionally invariant along the y axis, generated by a source much deeper in the solar convection zone than the initial horizontal flux tube. This background coronal field is designed to reproduce the magnetic field of an old decaying active region. In these simulations, the initially buoyant twisted flux tube, which is line-tied at the side boundaries, partially emerges from the convection zone into the corona due to a perturbation in the flux tube’s pressure and density (see Eqs. 18 and 19 from Leake et al. 2013) at the center of the flux tube. The later evolution of the emerging field involves surface shearing and rotation, and the formation of a new coronal flux rope above the new active region. The ultimate state of this new flux rope depends on the choice of the overlying field parameters.

thumbnail Fig. 1.

Emergence of the convection magnetic flux tube into the solar atmosphere at t = 110 t 0 for both the stable SD (A) and the unstable SD E (B) simulations. The cyan lines trace the coronal arcade magnetic field, originating from the lower boundaries, whereas the red lines belong to the magnetic flux tube and originate from the y = ±maxy side boundaries. The gray scale indicates the magnetic field strength at the surface.

The emergence process of the flux tube is examined through a range of initial coronal arcade strengths and orientations, in order to investigate the effects of the external magnetic field on the formation mechanisms of both stable and unstable flux ropes. The flux rope, self-consistently formed by the emergence of the flux tube, is either ejected, i.e., corresponding to an eruptive flare, or confined in the corona, i.e., equivalent to a quiet stable AR, depending on the orientation of the coronal arcade. For the non-eruptive simulations, the direction of the arcade field is the same as that of the top of the flux tube, maximizing the confinement of the flux tube by the arcade field, as shown in Figure 1A. Conversely, the eruptive simulations are driven by the magnetic reconnection between the arcade and the top of the flux tube’s axis, due to the opposite orientation of the magnetic field in the two structures. In this case, the horizontal magnetic field Bx changes sign at the interface between the flux tube and the arcade, making this separatrix a favorable location for magnetic reconnection. For each of the two arcade orientations, the arcade magnetic field strength is varied, leading to three different simulations: strong (SD), medium (MD), and weak (WD), corresponding to the initial surface field strength of 26, 19.5, and 13 G. An additional simulation where no coronal field is present (ND) is also performed, resulting in the formation of a stable flux rope.

For the eruptive simulations, the coronal magnetic field strength varies in the same way, namely {WD E, MD E, SD E}, where E corresponds to “Eruptive” simulations. The continued reconnection process between the emerging field and the arcade starts as the flux rope reaches the corona, changing the connectivity of the system and forming magnetic field lobes on both sides of the emerging flux rope. Figure 1B shows the magnetic field lines configuration of the flux tube emergence process for the SD E simulation, at t = 110 t0 (with t0 = 55.7 s, see Sect. 2.2 for details), slightly before the ejection of the flux rope. The formation of a shearing quadrupole configuration above the photosphere is clearly visible, with a central arcade of the emerging structure (red lines) now allowed to expand higher in the corona, whereas horizontal expansion is limited by the lobes on either sides (in red). Once the central arcade reaches a sufficient height, internal reconnection takes place at about t ~ 120 t0 (depending on the overlying arcade strength, see Leake et al. 2014, for details) beneath the flux rope axis, and the newly formed flux rope accelerates its vertical expansion, and the flux rope is immediately ejected. The objective of this work is to investigate if some pre-eruptive flare variations in one or more physical parameters are detectable using only photospheric magnetograms, at a reasonable stage before t ~ 120 t0, providing thus a reliable eruptive indicator.

This set of seven simulations treats the cases of both eruptive and non-eruptive flux rope formation. The presence and the orientation of the coronal arcade are critical parameters for the eruption onset. The ratio of the arcade flux to emerging flux is also a fundamental parameter, controlling the flux rope vertical acceleration, its size, and the amount of reconnection allowed to occur.

2.2. Scaling

The visco-resistive MHD equations are non-dimensionalized, using a normalization constant for each variable. In this work, we choose a slightly different normalization than the original paper, in order to increase the representative size and the magnetic flux of the active regions, to obtain more representative ARs. For the length, we rescale the simulations using the normalizing value L0 = 8.5 × 105 m, corresponding to five times the original value adopted by Leake et al. (2013, 2014; see Sect. 2.2). The two other normalizing values, the magnetic field (B0 = 1300 G) and gravitational acceleration (g0 = gsun = 274 ms−2), remain unchanged. The other normalizing variables that are affected by this rescaling are, namely the density (ρ0 = 5.77 × 10−5 kg m−3), the velocity (v0 = 15.26 × 103 ms−1), the time (t0 = 55.7 s), the temperature (T0 = 28.2 × 103 K), the current density (j0 = 0.122 Am−2), the viscosity (ν0 = 7.49 × 105 kg m−1s−1), and the resistivity (η0 = 16.3 × 103 Ωm). The Reynold’s number and the magnetic Reynolds number are kept in the same order as in the initial simulations, both evaluated to 100.

Using a larger spatial scaling than that of the initial analysis allows us to obtain larger active regions and larger magnetic flux, therefore more representative of the observed eruptive active region. Still, we are limited to the study of a relatively small active region with a characteristic flux of 1–2 × 1021 Mx during the early emergence process, and a characteristic size of about 30 Mm, roughly giving an area of about 296 μhs (i.e., microhemisphere – 1 μhs = 3.04 × 106 km2). According to Sammis et al. (2000), our simulations belong to the smallest flaring active region. For comparison, the authors have shown that during the 1989–1997 time frame, all the X4 flares only originated from active regions with an area greater than 1000 μhs. It should be noted that not all active regions in this size range are flaring: size seems to be a necessary condition to X4 flares, but not sufficient. Given our framework, our study is thus limited in terms of size and complexity of active regions that we cannot explore: the reliability of eruptive indicators is tested only for a given class of active region size and complexity. Nonetheless, it is worth noting that a recent study by Toriumi et al. (2017) concludes that a complex δ-sunspot is not a necessary condition for flaring ARs, even for the X-class flares. Thus, this kind of controlled study is important in the sense that the size and the complexity of an active region are not discriminant parameters for eruptivity: small active regions can still produce flares, while two active regions in the same size and complexity range do not necessarily have the same flaring likelihood.

However, if the rescaling of the simulation set enables us to have greater size, it also influences other quantities used in the simulations. Using our scaling, the transition region is 8.5 Mm thick (1.7 Mm in the original simulations), which is significantly higher than its characteristic size of about 0.1 Mm. However, given the spatial resolution imposed by the total domain size and the computational cost, such a thickness is required to resolve the large temperature and density gradient occurring in this part of the solar atmosphere.

3. Analysis: extracting the physical properties of ARs

In the following, the methodology that we follow to investigate the reliability of the eruptive flare indicators is detailed. From the series of 2D-plane vector magnetograms of the three eruptive and four stable simulations, we compute a series of parameters in order to detect some particular behavior associated only with the eruptive simulations. A predictor can be considered as reliable if it shows significant change(s) in a reasonable timescale prior to the eruption. Since flares and eruptions are driven by a store-and-release mechanism of magnetic energy, eruptive indicators must provide indications of such process occurring. Such signature includes threshold beyond which the system becomes unstable, most likely generating an eruption or a flare. Therefore a reliable predictor should present a different behavior respective to the eruptive or the non-eruptive nature of the simulations, but a similar behavior for simulations of the same nature. This signature should also be significant enough to be measured, and detectable in a sufficient time prior to the eruption to be able to then perform operational forecasting. Here, we focus only on determining which quantities could potentially be proficient eruptive flare predictors, i.e., parameters associated with higher values for the three eruptive simulations {WD E, MD E, SD E} than that of the stable {ND, WD, MD, SD} prior to the eruption. The associated thresholds needed for performing operational predictions require the use of (1) real photospheric magnetogram observations and (2) substantially more AR samples to be statistically significant.

3.1. Magnetograms

From the time series of the 3D MHD simulations, we first extract photospheric-like 2D-plane vector magnetograms, by interpolating the cube domain at z = 0 (as defined by Leake et al. 2013), onto a regular grid using a resolution of 0.86 L0 (730 km pixel−1) in both directions, which is equivalent to the instrument resolution of about 1′′. This resolution is the same as that of the SDO/HMI magnetograms. In our simulations, the surface (z = 0) is defined as the beginning of the temperature minimum region, i.e., corresponding to the region of the lowest pressure scale height. Since the eruptive flare occurs around t ~ 120 t0 for the three eruptive simulations {WD E, MD E, SD E}, we restrict the time windows of our analysis from t = 0 t0 to t = 150 t0, using a time sampling of Δt = 5 t0, where t0 = 55.7 s. Therefore, we finally obtain 31 magnetogram series for each vector magnetic field component Bx, By, Bz, and for each of the seven parametric numerical simulations {ND, WD, WD E, MD, MD E, SD, SD E}.

We first apply a mask to the entire magnetogram series, excluding the pixels for which , where Bmask = 30 G. This Bmask threshold has been chosen in order to be greater than the background magnetic field due to the coronal arcade field, respectively, about 26, 19.5, and 13 G for the SD, MD, and WD simulations (see Sect. 2.1). In observational analysis, a threshold is usually applied to exclude noisy data and only retain strong field areas. Our data mask threshold is relatively low compared to the typical uncertainties measured in observed magnetograms and so far the usual threshold applied to data. This is due to the relatively small typical values of magnetic field in our simulations, with Bz values in the [−800 G, 800 G] range, due to the small size and flux of our ARs (see Sect. 2.2). Many different masking methods have been used in previous analyses, using different threshold values applied on either B or Bz magnetograms, and will be briefly introduced in Section 5. The impact of applying a data mask on the detection of photospheric eruptive signatures will also be discussed in more detail in Section 5. The effects of noise, using a random perturbation of the magnetograms, will also be investigated in Section 6.

Figure 2 shows an example of such magnetograms, for both eruptive and non-eruptive simulations. Figure 2A displays Bz magnetograms for the MD simulation, either masking (right) or not (left) of the data, while Figure 2B exhibits the same photospheric maps for the MD E simulations. Masking data has a great impact on the area of the region of interest, used then for the computation of the different eruptive indicators. Eruptive and non-eruptive magnetograms are very similar, except near the external edge of the polarities. This is due to the quadrupolar nature of the eruptive simulations: as mentioned before, the coronal arcade has an opposite orientation with respect to the flux tube, i.e., the arcade Bx component has the opposite sign, to favor reconnection above the flux tube axis. From such time series of masked magnetograms, we compute a set of scalar parameters, most of them typically used in standard solar flare operational forecasting methods.

thumbnail Fig. 2.

Magnetograms of the Bz component at t = 100 t0 for both the MD (A) and MD E (B) simulations. Both non-masking (left) and masking (right) magnetograms are displayed, in order to show the reduction of the region of interest for further analysis. Eruptive and non-eruptive simulations present very similar magnetograms, except near the external edge of the polarities, due to the quadrupolar geometry above the surface.

3.2. The parameters

All the quantities considered in this work are scalar and can be derived using exclusively vector magnetic field data. They described the overall physical condition of an active region and the ongoing evolution of its magnetic field. Most of them have been chosen based on their previous identification as potential flare predictors, but new quantities such as helicity or WLsc, the current integral along the magnetic polarity inversion line, are also tested. Each derived quantity is parametrized such as one single number is able to characterize the state of the whole active region at a given time. For quantities that are spatially distributed, such as Bz (x, y), we compute the four first moments as described in Leka & Barnes (2003a): the mean, the standard deviation, the skewness, and the excess kurtosis (i.e., the kurtosis −3, for comparison to normal distribution). The skewness describes the asymmetry of the spatial distribution, while the kurtosis accounts for small changes in extremal values. As pointed out by Leka & Barnes (2003a), computing the four first moments of the spatial distribution allows us to quantify the subtle evolution of the different parameters.

3.2.1. Magnetic field

This first class of parameters is directly related to magnetic field evolution. All these parameters, listed in Table 1, come from the work of Leka & Barnes (2003a, hereafter LB03), and the readers are referred to Section 3 of their paper for a full detailed description. The spatial/scalar column specifies if the computed parameters are spatially distributed or not. If so, the four first moments are computed, easily identifiable by the notation. The temporal evolution of the total, vertical, and horizontal magnetic fields, respectively, noted B, Bz, Bh, and their associated moments provide insight about how the emerging field changes the distribution of B, and the overall polarity imbalance of the active region. The associated horizontal spatial gradients , and the corresponding moments indicate how much the magnetic field is sheared and distorted, an important indicator of the AR non-potentiality, especially along the neutral line. We also compute four quantities related to the AR flux: the total unsigned flux ϕtot, the net flux ϕnet, and the positive and negative flux Φ±, providing insights into the overall magnetic flux and the flux imbalance during the partial emergence process. Finally, we compute the photospheric free magnetic energy, which refers to the energy directly available for eruptive activity, i.e., measured with respect to the potential field energy (see the corresponding equation in Table 1). It is worth noting that in our case, the free magnetic energy is computed from surface integrals, not from the whole AR volume, which would require 3D magnetic field extrapolations.

Table 1.

Parameters relying on magnetic field tested in the present study.

3.2.2. Magnetic field geometry

This second class of parameters, listed in Table 2, quantifies the morphology of the AR magnetic field and its deviation from the potential configuration. The inclination angle γ measures how much a magnetic field is inclined relative to the vertical axis. It reflects the magnetic field orientation, providing hints about the flux emergence evolution. In the present work, we adopt the following definition(1)

Table 2.

Parameters relying on magnetic field geometry tested in the present study.

As defined, the inclination angle tends toward small values when the magnetic field is approaching the vertical, while great values correspond to an almost horizontal magnetic field. This quantity indicates the morphology evolution of the magnetic field during the emerging process.

In a force-free field framework, the Lorentz force is null and the magnetic field follows(2)where α is the local twist parameter, referring to the torsion of each individual field line. Using magnetogram data, α can be approximated by (see e.g., Pevtsov et al. 1994)(3)where μ0 is the permeability in free space. Using a set of 133 of flaring ARs, Nindos & Andrews (2004) found that the pre-flare value of the twist parameter α was, in general, higher for M-class eruptive flares, suggesting the twist parameter as a reliable indicator for CMEs.

The magnetic shear angle has been extensively used in previous studies to measure the non-potentiality of active regions. Initially introduced by Lu et al. (1993) and Hagyard et al. (1984), the three-dimensional shear angle Ψ measures the angle between the observed magnetic field and its potential component, while its horizontal projection, namely the “horizontal” or “planar” shear angle , quantifies the azimuthal difference between the observed and potential magnetic field. The potential magnetic field is computed following the method of (Valori et al. 2012). We also tested the additional shear angle-related parameters A[Ψ > 80°] and A[ > 80°] proposed by Leka & Barnes (2003a), corresponding to the AR area where respectively the 3D and the horizontal shear angle exceeds 80°.

Magnetic helicity measures to what extent the magnetic field lines are wrapped around each other, and how much the individual magnetic field lines are twisted and writhed, relatively to their lowest energy state. This parameter provides a quantitative estimation of the geometric properties of the magnetic field lines. Because helicity is a conserved MHD quantity, even in resistive MHD where the dissipation is very small (Pariat et al. 2015), im portant efforts have been carried out concerning its estimation and its relation with solar flares (see e.g., Nindos & Andrews 2004; Démoulin 2007; Démoulin & Pariat 2009; Park et al. 2010). In the present study, helicity has been computed using the Gθ proxy for the helicity flux density from Pariat et al. (2005), using the series of photospheric magnetograms. In this framework, the time variation of the relative magnetic helicity can be written as(4)where is the photospheric surface, Bz is the magnetic field normal to this surface, r = x – x′ is the vector between the two photospheric positions x and x′, with u and u′ the associated flux transport velocities of the photospheric field line footpoints. While the flux transport velocity could be directly extracted from the simulations, the velocity u has been estimated using the differential affine velocity estimator for vector magnetograms (DAVE4VM) from Schuck (2008) in order to analyze the simulation data as observations. Using simulations, Welsch et al. (2007) have tested the robustness of the DAVE algorithm, the line-of-sight (LOS) magnetogram equivalent of the DAVE4VM algorithm (Schuck 2008), and found a good agreement between the simulated and the DAVE-derived velocities.

3.2.3. Current properties

Current-carrying magnetic fields are understood to be the building block for understanding the flares and CMEs drivers. The current properties computed here (see Table 3) allow us to quantify the energy stored by the magnetic field relative to its lowest energy state given its boundary distributions, i.e., potential. The vertical current density component Jz and the associated moments can be estimated through the classical Ampère’s law. The total, net, positive, and negative currents, respectively, noted as Itot, Inet, and I±, are estimated using various reckonings through the AR area. The heterogeneity and chirality current density components are also derived, according to Zhang (2001) (5)with B = Bb the first term of equation (5) refers to the current of chirality and the second relates to the heterogeneity, perpendicular by construction to the magnetic field. Thus, in the case of the heterogeneity current dominates over the chirality component, the AR magnetic field is far from the force-free field hypothesis. The associated four moments for each component are also evaluated, as well as the total and net heterogeneity and chirality currents, , . We also computed the net current originating from each polarity, , as described in Leka & Barnes (2003a) (6)as well as the vertical contribution to the current helicity density hc along the vertical axis. Indeed, the current helicity density is defined as B·J but since only the vertical current component Jz can be deduced from the observations, the current helicity computed here is thus limited to the vertical contribution alone. From this partial estimation of the helicity density, we then derived the total and net partial helicity, and , characterizing the current helicity imbalance over the ARs. Bao et al. (1999) examined the relation between flare activity and the vertical contribution of the current helicity, and found that the time variation of partial current helicity is higher in flaring in ARs, suggesting partial and most probably total current helicity as a valuable eruptive indicator.

Table 3.

Current properties tested in the present study.

Determining whether or not the net electric current is neutralized over individual polarities of ARs is crucial for some theoretical flare and CMEs models (e.g., Melrose 1991; Parker 1996; Titov & Démoulin 1999; Forbes 2010). If the AR currents are fully neutralized, the net current integrated over one photospheric polarity is set to zero. Recent observations and simulations tend to confirm the existence of non-neutralized ARs (e.g., Török et al. 2014, Dalmasse et al. 2015, and references therein). Forbes (2010) argued that neutralized ARs may inhibit the eruption process. Hence, the measurement of direct and return current within ARs could provide insight into the flaring activity. The computations of both the direct current Id and the return current Ir have been carried out following Török et al. (2014), who used the same simulations as this present work. They carried on the computation at the base of the corona, while we only concentrated on the surface, i.e., photospheric measurements. Several combinations of Id and Ir have also been tested: Id + Ir, |Id/Ir| as well as the same quantities normalized by the total flux Φtot of the AR (see Sect. 3.2.1).

3.2.4. Lorentz forces

The Lorentz forces are thought to be an important parameter for solar eruptions. In particular, Fisher et al. (2012) argued that in eruptive flare ARs, the impulse coming from the Lorentz forces is dominant over all other forces, suggesting that an observed change in the Lorentz force before the eruption could be a pre-eruptive flare signature. To test this proposition in our simulations, we used the expressions derived by Fisher et al. (2012) in order to compute the horizontal and vertical components Fx, Fy, Fz of the Lorentz forces. We also estimated the total Lorentz force F and the normalized horizontal and vertical components , , (Table 4).

Table 4.

Lorentz force parameters tested in the present study.

3.2.5. Magnetic polarity inversion line (MPIL) properties

Within an AR, the magnetic field is sheared, stressed, and twisted, thus deviating from the potential state of minimum energy. In particular, near the magnetic polarity inversion line (MPIL), the magnetic field can be almost perpendicular to the potential field (see e.g., Hagyard et al. 1990; Falconer et al. 2003). This region of strong magnetic shear is particularly well observed in X-rays, corresponding to the sigmoidal structures with an overall shape of “S” or an inverse “S” for bipolar ARs. Evaluating the properties of the magnetic field and currents near the MPIL provides a quantitative measure of the AR’s overall non-potentiality, even though some studies, such as Li et al. (2005), did not find obvious strongly sheared magnetic neutral line or strong overlying arcade associated with a X3 flare event. Since the degree of non-potentiality seems correlated with eruptive flare productivity (Falconer et al. 2002), the MPIL and near-MPIL areas’ properties are potentially proficient predictors.

We computed a series of parameters relying on MPIL lengths, most of them introduced by Falconer et al. (2003, 2006, 2008) (Table 5). The total length of the strong field neutral line Ls is defined as the length of the MPIL where is greater than a given threshold . The quantity Lss defined the length of the MPIL for which the observed horizontal magnetic field Bh is greater than the given threshold and the shear angle is greater than . The total strong-gradient length of the MPIL Lsg is specified for MPIL regions where is greater than the previously defined threshold and where the horizontal gradient of the vertical magnetic field .We introduce here a new measure of the AR non-potentiality through the total strong-current length Lsc, measured along the MPILs, where and the current density Jz is greater than a given threshold .

Table 5.

Magnetic Polarity Inversion Line (MPIL) parameters tested in the present study.

The quantities WLss, WLsg were first introduced by Falconer et al. (2008) while we define the new parameter WLsc, measuring the current along the MPIL. WLss,WLsg represent the integral of, respectively, the shear angle and the horizontal magnetic gradient , along the MPIL strong horizontal magnetic field portions. The WLsc parameter corresponds to the current density Jz-integral along the MPIL portions of strong Bh. For the WLss and WLsc parameters, only the MPIL portions where are taken into account, while WLsg considers only the regions where . The two parameters L[Ψ > 80°] and L[ψ > 80°] from Leka & Barnes (2003a) are also tested, corresponding to the portions of the MPILs where, respectively, the 3D and the horizontal shear angles Ψ and ψ are greater than 80°. Finally, the unsigned flux near the MPIL is computed according to Schrijver (2007), given characteristic values around log R ~ 5.

Falconer et al. (2008) chose the following thresholds of G,  = 45° and G Mm−1 while we initially adopt lower thresholds for G and G Mm−1 given the small size and flux of our simulated ARs (see Sect. 3.1). For the Lsc parameter, we imposed a threshold  = 12 m Am−2. It is worth noting that the choice of these thresholds can have substantial impact on the detection of eruptive flare signatures as we will show in Section 7, where a parametric study is proposed to estimate their influence.

4. Eruptive indicators evolution with B mask = 30 G

For each output magnetogram of each simulation, we have computed the 99 parameters described in Section 3.2. For clarity, in this section, we only present the evolution of the most representative and interesting parameters among them. The parameters shown in Figures 35 have all been mentioned in previous studies as potential good indicators for eruptive and non-eruptive flaring activity. Even though the analysis method used in all these aforementioned works may differ from ours, i.e., using superposed epoch or forecasting methods, we used these studies as guidance to present our results and avoid quantitative comparisons. The displayed time window starts at t = 50 t 0, once the flux tube has significantly emerged, the emergence starting at t = 35 t 0. Below t = 50 t 0, none of the parameters exhibits a significant variation, and for clarity, we thus restrained the displayed time window from t = 50 t 0 to t = 150 t 0, allowing us to analyze both prior and post eruptive flare (occurring at t ~ 120 t 0) physical conditions. In the following, we do not attempt to provide physical interpretation of the full variation of the parameters, which is beyond the scope of this paper; rather we focus on the detection of significant changes in one or more parameters, according to the eruptive nature of the simulations only and occurring at a reasonable stage prior to the eruption.

thumbnail Fig. 3.

Some parameter evolution for the seven simulations, as a function of time in t 0 units (t 0 = 55.5 s). From (A–F), the following predictor evolution is respectively displayed: the total current helicity , the total current I tot, the mean value of the free energy , the total flux Φ tot, the total Lorentz force F and the Schrijver R value. The warm colors indicate the evolution for the three eruptive simulations {WD E, MD E, SD E}, respectively, corresponding to the magenta, orange, and red solid lines. Conversely, the evolution of the non-eruptive simulations {ND, WD, MD, SD} is displayed using the cold colors, respectively, as follows: purple, deep blue, light blue, and green solid lines. The eruption time is indicated by the vertical dashed gray line.

Figure 3 displays the evolution of some physical properties (see caption for details). The total partial current helicity (Fig. 3A), the total current I tot (Fig. 3B), the free magnetic energy (Fig. 3C), and the total Lorentz force F (Fig. 3E) are in the top-ranked eruptive indicators by Bobra & Couvidat (2015), whereas the total unsigned flux ϕ tot (Fig. 3B) and have been ranked in the best predictors for 48 hours flare forecasting by Bobra & Ilonidis (2016). The R value (Fig. 3F) was empirically established as a powerful eruptive indicator by Schrijver (2007). The parameters , I tot, ϕ tot, F, and the R value exhibit a very similar evolution, possibly due to an existing correlation between each parameter (Barnes et al. 2016). There is no evidence of physical changes prior to the eruptive flare, at about t ~ 120 t 0 (vertical dashed gray line) for the eruptive simulations.

On the other hand, the mean free magnetic energy displays different behavior for the set of simulations, with a different peak value at t = 70 t 0. The magnetic free energy is highly dependent on the region of interest selected using the B-mask described in Section 3.1. This selected area increases more rapidly than the free magnetic energy per surface unit after t = 70 t 0, due to the emergence process. Therefore, the mean excess magnetic energy is observed to slowly decrease, but this is an artifact engendered by the size of the AR area. Using a constant mask as a test, the free magnetic energy is increasing as expected during the formation of a flux rope. However, the pattern does not correlate with the flare-eruptive nature of the simulations (warm vs. cold colors solid lines), but rather with the coronal arcade magnetic field strength. Higher peak values correspond to weaker coronal arcade, although that for a given coronal arcade magnetic field strength, the non-eruptive simulation shows systematically slightly higher values. Thus, the behavior of the indicator is not only dependent on the eruptive nature of the simulations, but is also coronal-field-strength dependent. This kind of behavior is not useful for flare forecasting, as we look for clear changes in some physical properties between stable and non-stable simulations only, in a way that some thresholds can, for example, be imposed. The property changes do not allow such an approach, since the same free magnetic energy peak could correspond to both eruptive and non-eruptive ARs, with different coronal arcade magnetic strengths.

Figure 4 presents the evolution of some magnetic field properties, such as geometry and current (see caption for details). The time evolution of the kurtosis of the twist parameter κ(α) is displayed in Figure 4A. The kurtosis of the horizontal magnetic field κ(B h ) (Fig. 4C) is among the more efficient predictors in the four-variable discriminant analysis of Leka & Barnes (2003b), applied on flaring/non-flaring ARs. The mean gradient of the horizontal magnetic field (Fig. 4F) has been found by Bobra & Ilonidis (2016) to be the best predictive variable for 24 h flare forecasting, while the mean shear angle (Fig. 4D) is also in the top-ten best variables for both the 24 and 48 h predictions of the same study. We also examine the predictive capabilities of the relative magnetic helicity time variation (Fig. 4E) and the direct current I d (Fig. 4B) since both have been suggested to play an important role in eruptive flare mechanisms (see Sect. 3.2 and e.g., Nindos & Andrews 2004; Dalmasse et al. 2015, and references therein).

thumbnail Fig. 4.

Same as Figure 3, but for different parameters. From (A–F), the following indicator evolutions are respectively displayed: the kurtosis of the twist parameter κ(α), the direct current , the kurtosis of the horizontal magnetic field κ(B h ), the mean value of the shear angle , the time variation of the relative magnetic helicity dH m/dt and the mean value of the vertical magnetic field B z .

As before, none of these series of parameters exhibits a clear eruptive signature. The direct current I d and the relative helicity variation show very similar evolutions, even after the eruption at t ~ 120 t 0. The kurtosis of the twist parameter κ(α) shows slightly higher values for the SD simulations, but apart from that difference, the evolution is almost the same for the whole simulation set. The κ(B h ) presents a strong peak at t = 60 t 0, with distinct values as a function of the simulations, but as for the , this behavior does not depend on the eruptive nature. Instead, the peak intensity depends on the strength of the coronal arcade, with stronger coronal magnetic field associated with stronger κ(B h ) peak. The horizontal gradient also presents a strong peak, at about t ~ 66 t 0, whose magnitude depends once again on the overlying field strength. As the coronal field gets smaller, the horizontal gradient increases, no matter the stable/unstable nature of the simulation. Notwithstanding, the mean value of the shear angle displays a slightly distinct behavior between eruptive and non-eruptive simulations. The parameter exhibits somewhat greater values for the eruptive numerical experiments, notably for the MD E and SD E simulations compared to that of the stable simulations. Both the eruptive and the non-eruptive remain stable once the flux tube has significantly emerged, i.e., after t = 70 t0. However, there is no behavioral change before and after the eruption at t = 120 t0, a required feature for a parameter to be a reliable predictor. As the eruption cannot be detected observing only the parameter, the mean shear angle is not a fully discriminant parameter.

Figure 5 shows the evolution for six of the MPIL properties extracted from the magnetogram series. All of them have been established as a powerful measure of the non-potentiality of ARs by Falconer et al. (2003, 2006, 2008), and are therefore potentially good predictors for both eruptive and non-eruptive flare activity. For this class of parameters, clear eruptive flare signatures are observable, and the evolution of eruptive and non-eruptive simulations is significantly different, except perhaps for the length of the strong-shear MPIL Lss (Fig. 5A). For each of these parameters, a rise is undeniably observed during the early stage of the emergence for the three eruptive simulations, i.e., between t = 60 t0 and t ~ 95 t0. During the whole emergence process, the parameter evolutions of the SD E and MD E simulations have a significantly different comportment relative to the others, while those of the WD E simulation become similar to the non-eruptive simulations once the eruption is imminent, after t ~ 95 t0. For the strong-shear length Lss, both the SD E and MD E Lss lengths are significantly longer, by respectively of factors 2 and 3, than that of the stable simulations, whereas the WD E Lss length is barely 1.5 longer and only at the very early stage of the emergence, over a short time range. As a result, the Lss is not a very strong predictor, since the eruptive signature quickly disappears, at least for weak coronal magnetic field. The parameter WLsg (Fig. 5D) shows a similar trend, with a peak at t ~ 80 t0 only about 1.2 times that of the non-eruptive simulations, making it a poorly robust predictive parameter.

thumbnail Fig. 5.

Same as Figure 3, but for different parameters. From (A–F), the following indicator evolutions are respectively displayed: the length of the strong-shear MPIL Lss, the integral of shear angle along the MPIL WLss, the length of the strong-gradient MPIL Lsg, the integral of the Bz horizonal gradient along the MPIL WLsg, the length of the strong-current MPIL Lsc and the integral of the current along the MPIL WLsc.

However, the other parameters, WLss (Fig. 5B), Lsg (Fig. 5C), Lsc (Fig. 5E), and WLsc (Fig. 5F), are clearly more robust predictors, showing a distinct enhancement prior to the eruption, even for the WD E simulations. The WLss parameter is 3–8.5 times greater than that of the stable simulations, while the Lsg length is 2–3 times higher, but on a longer time range, between t = 75 t0 and t = 95 t0. The Lsc parameter presents a sharp peak at t = 80 t0, between 2.8 and 5.4 times greater than the non-eruptive values. The WLsc MPIL property is the most efficient eruptive indicator, since the three eruptive simulations reach similar values at the beginning of the emergence process, being about 8 times greater than the non-eruptive WLsc measurements. For this latter parameter, the influence of the coronal arcade is reduced, at least during the early emergence stage.

All these MPIL properties provide comprehensible eruptive flare signatures in these models, prior to the eruption. Apart from the Lss and WLsg parameters, for which the signature exists but is weak, the Lsg, WLss, Lsc, and WLsc parameters provide robust measurements of pre-flare eruptive conditions. As the MPIL properties are correlated with the magnetic configuration of the flux emergence, these parameters provide measurements of the magnetic complexity of the ARs. Since the eruptive flare formation and ejection are related to the quadrupolar versus bipolar nature of the emergence for this simulation set, a clear eruptive flare signature can be detected. The newly defined Lsc and WLsc MPIL-current properties, as well as the WLss parameter appear to be the better predictor, with significant changes between eruptive and non-eruptive simulations. However, none of the other 93 parameters, related to current, magnetic field geometry, or properties, provides unambiguous eruptive flare signatures, since no significant changes have been detected in their evolution.

5. Influence of data masking

The eruptive indicators are highly sensitive to the AR area selected to perform the computation, as seen in Section 4, where the flux emergence process makes the AR area increase very rapidly. Bobra & Couvidat (2015) already highlighted that the AR parameters are highly sensitive to the data masking, and pointed out that a study to optimize such parameters is necessary. In this section, we present, as an example, the same results as in Section 4, but this time using a Bmask = 100 G. This new threshold leads to a reduced AR area, and because most of the parameters are AR-area dependent, most of the results are affected.

A variety of data masking methods can be found in previous studies. The automated system Solar Monitor Active Region Tracker (SMART) for detecting ARs and their associated properties using the Solar and Heliospheric Observatory (SOHO)/Michelson Doppler Imager (MDI; Scherrer et al. 1995) uses a static Bz threshold of 70 G to remove the background (Higgins et al. 2011). The global features provided by the Space-Weather HMI Active Region Patches (SHARPs; Bobra et al. 2014) pipeline select only the pixels for which the magnetic field strength is above the disambiguation threshold of about 150 G (see also Hoeksema et al. 2014). Using the ground-based data from the University of Hawai’i Vector Magnetograph (Mickey et al. 1996), Leka & Barnes (2003a) used only pixels above 3σ detections. Falconer et al. (2008) imposed the Bz component to be greater than 150 G, while Falconer et al. (2011) used either a threshold of 25 or 35 G, both studies using the MDI data. Al-Ghraibah et al. (2015) also used a 3σ binary mask to exclude noisy MDI data, while Mason & Hoeksema (2010) imposed a static threshold of 100 G.

Figure 6 displays the same parameter evolutions as Figure 5, but excluding pixels for which the total magnetic field B is lower than 100 G. The seven simulations exhibit the same evolution over time, reaching similar values close to the eruption starting time. In comparison with the previous results using a lower mask threshold, the rise observed for the six MPIL-related parameters of the eruptive simulations is lost, and the characteristic eruptivity signature is no longer discernible. This is due to the higher threshold masking imposed to the initial data, excluding initially from the dataset a fraction of the MPIL, and therefore reducing the computation domain of the MPIL properties.

thumbnail Fig. 6.

Same as Figure 5, but using a mask threshold Bmask = 100 G.

In order to visualize the data masking impact on the detection of the MPIL, Figure 7 displays the Bz magnetograms for the SD (Figs. 7A7C) and the SD E (Figs. 7D7F) simulations, at t = 100 t0. The color scale is saturated between −400 and 400 G to highlight the two polarities. Since the orientation of the arcade is opposite in SD versus SD E simulations (see Sect. 2 for details), the magnetic configuration is purely bipolar for the non-eruptive simulations, whereas the eruptive magnetic topology is quadrupolar. This can be clearly observed in Figures 7A and 7D, where the MPIL is only located between the two polarities for the stable SD simulation (Fig. 7A), whereas the MPIL possesses additional portions at the polarity external edges for the SD E simulation. However, as the masking threshold increases, the external MPIL is reduced, and for the case where Bmask = 100 G, both MPILs, either in the eruptive or non-eruptive simulation, are almost exactly the same. Consequently, the deviation between stable and unstable simulations observed in the MPIL properties for Bmask = 30 G (see Fig. 5) is significantly reduced with Bmask = 100 G (see Fig. 6). The eruptive and non-eruptive simulations cannot be distinguished in this case, and thus the MPIL indicators become inefficient. This shows how much the initial mask threshold Bmask should be carefully chosen, since it may be crucial to be able to detect significant eruptive signatures and therefore make reliable flare forecasting.

thumbnail Fig. 7.

Magnetograms of the SD (A–C) and SD E (D–F) simulations, including the MPIL (red lines). The MPIL width has been dilated by a factor 4 in order to increase its visibility. The masking threshold increases from left to right, with respective values of 30, 50, and 100 G.

6. Impact of the noise on the detection of pre-eruptive signatures

Noise is an important issue for observational analysis. There are many different sources of error affecting the photospheric magnetic field measurements, e.g., photon noise, detector noises, spacecraft radial velocity inducing periodic systematic uncertainties (Hoeksema et al. 2014), or uncertainties associated with the inversion of the Stokes parameters. In the present work, these different sources are not individually treated, which is beyond the scope of this paper. Instead, we simply include a random Gaussian noise to our data in order to evaluate noise impact on the detection of pre-eruptive flare signatures.

We used a Monte-Carlo scheme, randomizing the magnetograms using Gaussian perturbations, with a standard deviation of 3.5 G. The error associated with HMI magnetic field data is typically about 10 G for the line-of-sight (LOS) component (V. Bommier, private communication) and can be generally one order of magnitude larger for the horizontal components (Hoeksema et al. 2014). However, our characteristic magnetic field values are about 3–4 times smaller than those that are typically observed (see Sect. 2.2). Therefore, the standard value of 10 G is not representative for our time series magnetograms. Accordingly, we adopted the characteristic noise standard deviation of 3.5 G, in order to remain consistent between our simulated observations and the relative noise scale. For this preliminary study, we assume the same error for the three components. Noise levels of 1 G and 5 G have also been tested with no sensitive differences with the results presented here. All the parameters investigated in this work are then derived from these time series of noisy magnetograms, following the same approach. First, for each simulation and for each time step, 50 different noisy magnetograms are computed by randomizing the initial one. The indicators are computed 50 times from these 50 noisy magnetograms and then averaged together. The error associated with each parameter for each time step is assumed to be the standard deviation corresponding to the 50 computed noisy indicator series, while the mean value provides the simulated magnetic field measurements.

Figure 8 displays the same parameters as Figure 5 but including the Monte-Carlo estimation of the measurement errors. The overall behavior remains unchanged, although some specific changes occur for the Lsc and the WLsc (Figs. 8E and 8F) indicators. The Lss and the WLss predictors for the MD E and SD E simulations are still significantly higher than that of the non-eruptive simulations. However, for these two parameters, in the WD E simulation, the pre-eruptive signature is lost. The sharp increase of the Lsg and WLsg indicators, observed between t = 60 and t = 100 t0, is still observable, making them robust to noise perturbation. However, for the two current-related predictors Lsc and WLsc, the impact of noise is more important, and for instance the pre-eruptive signature provided by the WLsc indicator is completely dominated by the noise. For the Lsc quantity, the trend is different, and a peak at t ~ 75 t0 is still detectable for the eruptive simulations. However, we observe a very different behavior for the WD and ND simulations compared to the clean data (see Fig. 4). This is due to complex effect of the noise, which generates additional pseudo-MPIL in strong-current regions.

thumbnail Fig. 8.

Same as Figure 5, but including the noise perturbation.

To summarize, the presence of noise may affect the detection of the potential pre-eruptive signature, with impact depending on the parameter considered. The Lsg and WLsg parameters are mostly not affected, and their pre-eruptive sharp increases are still detectable. Lss and WLss and Lsc are slightly affected, showing weaker pre-eruption peaks, and moderate effects for the WD E simulations. However, WLsc is strongly affected and the pre-eruptive signature is no longer observable.

7. Parametric study of MPIL properties

As described in Section 4, given our controlled-case study, only the MPIL features are able to provide a clear eruptive signature. Seven of the parameters characterizing the MPILs, namely Ls, Lss, Lsg, Lsc, WLss, WLsg, and WLsc depend on the four different thresholds , , and (see Sect. 3.2.5 for details). These thresholds all refer to the portion of the MPILs taken into account in the computation of the various MPIL properties. As for the data masking process used to isolate the AR core area, the calculation of the eruptive indicators is highly sensitive to these values. In order to optimize the choice of these criterion and quantify their influence, we present in this section a parametric study of the MPIL properties. The threshold values explored for each parameter have been chosen based on typical values commonly used in observational data. Since our simulations correspond to small AR, we explore the parameter space using smaller and characteristic values used in previous studies.

To illustrate how much the choice of the thresholds potentially affects the detection of the MPILs, Figures 911 represent the MPIL portions (displayed as white line) detected as a function of the threshold , , and for the SD (top rows) and SD E (bottom rows) simulations, in a similar fashion as Figure 7. The length of the MPIL (white line) depends on these different parameters. Figure 9 displays Bh maps on the background, and the parameter varies from left to right, while Figure 10 represents maps as a function of the threshold. Figure 11 displays Jz maps and the threshold changes from 0 (left) to 12 m Am−2 (right). In this section, the objective is to start to investigate the effects of each thresholding parameter on the MPIL parameters that have been demonstrated to be efficient eruptive predictors, under the condition that the initial mask thresholding was low enough (see Sect. 4). In this section, we do not initially mask the data (see Sect. 3.1, i.e., Bmask = 0 G) as we focus on threshold impact. Obviously, it is worth noting that the conjugate effects of high mask threshold (see Sects. 3.1 and 5) and individual parameter thresholdings are even stronger on the detection of the eruptive signature and therefore highly worsen the results.

thumbnail Fig. 9.

Maps of the horizontal magnetic field Bh, for the SD (A–C) and the SD E simulation (D–F) at t = 100 t0. The white line represents the portion of the MPIL for which . The threshold increases from left to right, with values of  = 0, 50, and 100 G. No initial data masking is applied here, in order to highlight only the impact of the threshold on the computation of Lss, Lsc, WLss, and WLsc.

thumbnail Fig. 10.

Maps of the quantity, using the same layout as Figure 9. The white line now represents the portion of the MPIL for which the horizontal gradient of the magnetic field is higher than . From left to right the threshold increases from 0–100 G Mm−1.

thumbnail Fig. 11.

Maps of the vertical current Jz, using the same layout as Figure 9. The white line now represents the portion of the MPIL for which the current density . From left to right, the threshold increases, with respective values 0, 3, and 12 mA m−2.

The threshold has a strong effect on the detection of the whole MPIL: as seen in Figure 9, the quadrupolar configuration, corresponding to outer MPIL on the edge of the polarities (see Fig. 9D), is no longer detected as threshold is greater than 50 G, therefore corresponding to similar MPIL (white line in Figs. 9C9F). As such, any MPIL properties computed using such a threshold will not provide reliable eruptivity signature. The threshold is less restrictive, and a small fraction of the MPIL outer portions (white line in Fig. 9F) is still detected for the SD E simulations. However, it is worth noting that the MPIL represented in Figure 10 is dilated to increase the visibility; actually the MPIL is thinner and only a few additional pixels are detected. On the other hand, the thresholding does not affect that much the larger measurements of eruptive MPIL: as seen in Figure 11, the external portions of the MPIL are still observed, even using a threshold of 12 mA m−2. From these series of maps as a preliminary coarse evaluation of thresholding effects on the MPIL properties, we can conclude that the current thresholding is the most flexible parameter, since the quadrupolar configuration of the eruptive simulations can still be detected even when increasing threshold.

To further investigate the influence of the physical thresholds on the MPIL predictors, we compute the variations of the Lss, Lsg, Lsc, WLss, WLsg, and WLsc indicators as a function of their associated threshold. The three first quantities depend on two different thresholds, making therefore their evolution even more sensitive to the chosen cutting values. Figure 12 displays the Lss variations for both eruptive SD E (solid lines) and non-eruptive SD (dashed lines) simulations, as a function of the and thresholds. The is increasing from panel to panel, while the is fixed for a given color line. If the is below 25 G (Figs. 12A and 12C), we can still observe longer Lss lengths for the SD E simulation, whatever the threshold. In these three cases, the SD E Lss parameter is respectively about 4, 3, and 2 times longer than that of the SD simulation, a parameter difference still measurable. However, for higher , this difference is reduced and to discriminate between eruptive and non-eruptive simulations becomes difficult. On the other hand, the influence of the threshold is weak and whatever the chosen value, SD E and SD simulations can still be distinguished.

thumbnail Fig. 12.

Parametric evolution of the total length of the strong-shear MPIL Lss, as a function of the two thresholds and . The varies for each figure, with respective values of 0, 10, 25, 50, 75, and 150 G, from (A–F). For each panel, the Lss quantity is represented for the strong arcade eruptive (SD E; solid lines) and strong arcade (SD; dashed lines) simulations. The threshold is changed for each curve, using the values of 30° (blue lines), 45° (red lines), and 60° (magenta lines). The eruption time is denoted by the vertical gray dot dashed line.

Figure 13 is the same as Figure 12, but exploring the variations of the Lsg predictor as a function of and thresholds. As before, each panel corresponds to a given threshold, and each color line corresponds to a given . For this parameter, the influence of the threshold is lower, while the impact of the is significantly stronger. For example, for the given G Mm−1 (red lines on each panel), whatever the imposed, the SD E Lsg length is still longer than that of the SD simulation. However, as the increases, the divergence between the SD E and SD curves rapidly disappears, making the eruptive signature undetectable. For Figures 12A and 12B, it is worth noting that the length Lsg (dashed blue line) of the SD simulation remains constant over the whole numerical experiment. This is due to the low imposed, below the photospheric arcade magnetic field mean value, conjugated to the G Mm−1 threshold, allowing therefore to measure the whole MPIL on the entire numerical domain. Since SD is a non-eruptive simulation, the magnetic field configuration is bipolar, and the MPIL is almost kept constant.

thumbnail Fig. 13.

Same as Figure 12 but for the parameter Lsg, as a function of the two thresholds and . The varies in the same way as for Figure 12. The threshold is changed for each curve, using the values of 0 (blue lines), 10 (red lines), 25 (magenta lines), 50 (green lines), 75 (yellow lines), and 100 G Mm−1 (black lines).

Figure 14 shows in the same way as Figures 13 and 12 the variations of the newly introduced Lsc parameter, as a function of its two associated thresholds and . For this MPIL property, the influence of both thresholds is important. As for the Lss parameter, detecting the eruptive signature, i.e., a longer Lsc for the eruptive simulations, becomes harder as exceeds 25 G. For G, a longer Lsc for the SD E simulation is however still detectable, but the difference between eruptive and non-eruptive simulations is less obvious. For higher (Figs. 14E and 14F), discriminating the simulations is impossible, whatever the current . On the other hand, increasing the also has a strong effect on the detection of the eruptive signature, even though small differences between eruptive and non-eruptive simulations persist. For instance, if G (Figs. 14B) and  = 3 mA m−2 (blue lines), Lsc is almost four times longer before the eruption for the SD E simulation than that of the SD. If we increase to 24 mA m−2 (yellow lines), the Lsc SD E to SD ratio decreases to 1.7.

thumbnail Fig. 14.

Same as Figure 12 but for the parameter Lsc, as a function of the two thresholds and . The threshold is changed for each curve, using the values of 3 (blue lines), 6 (red lines), 12 (magenta lines), 18 (green lines), and 24 mA m−2 (yellow lines).

From these three plot analyses, general dependence trends can be deduced. The Lss parameter is strongly dependent on the threshold, and for greater than 50 G, the eruptivity signature is no longer measurable. However, the influence of the shear-angle threshold is rather weak, and whatever the threshold assumed, Lss parameters still allow us to discriminate between the SD E and SD simulations. Conversely, the Lsg threshold is weakly dependent on the threshold, while the gradient threshold has a strong impact on the detection of eruptive ARs. For G Mm−1, the distinction between both AR types is hardly observable. Finally, the Lsc quantity is impacted by the choice of both its associated thresholds, even if small remnants of eruptive signature are still detectable for high thresholding values.

The WLss, WLsg, and WLsc parameters only depend on the threshold (see Table 5). Figure 15 displays the variations of these indicators for both the SD E (solid lines) and the SD (dashed lines) simulations, varying the threshold, in order to optimize the choice of this parameter. Figure 15A displays the WLss parameter evolutions, and as the threshold increases, the WLss parameter decreases for the SD E simulation while it remains stable for that of the SD non-eruptive. Using a low excluding value allows us to detect the external portion of the MPIL (blue, red, and magenta lines) and therefore to discriminate between the eruptive and non-eruptive simulations, while thresholds larger than 50 G (green, yellow, and black lines) remove this eruptive signature, characterized by a larger WLss for eruptive ARs.

thumbnail Fig. 15.

Parametric evolution of the WLsc (A), WLsg (B), and WLss (C) parameters, as a function of the threshold . Solid lines represent the parameter evolutions for the SD E simulation, while dashed lines are computed for the SD stable simulation. The varies from 0 (blue lines) to 150 G (black lines), with the following intermediate values 10 (red lines), 25 (magenta lines), 50 (green lines), and 75 G (yellow lines).

Figures 15A and 15B present respectively the WLsc and WLsg and parametric analysis. For these two parameters, the thresholding impact is weak, and even using high threshold allows to distinguish between eruptive and non-eruptive simulations. The WLsg property (see Fig. 5D) is about 2–2.5 times higher for the SD E simulation than that of the SD, whatever the adopted threshold imposed on . The newly defined WLsc parameter shows sharp variations for the eruptive simulation, increasing quickly by a factor about 4, whereas the non-eruptive evolution remains smooth over the whole emergence process. These sharp evolutions are observed whatever the Bh threshold imposed on the horizontally observed magnetic field, demonstrating that the noise-free WLsc is a robust predictor for detecting flare producing ARs.

From these parametric studies of the MPIL properties WLss, WLsg, and WLsc, clear conclusions can be drawn. The WLss predictor is highly dependent on the threshold process, and the eruptive signature, corresponding to higher WLss in eruptive ARs, is rapidly lost as increases. Consequently, if the WLss is able to detect imminent eruption in an AR under certain conditions (i.e., low masking threshold Bmask, see Sect. 5), this may not be the most reliable eruptive indicator to be used for flare forecasting, unless different thresholds are simultaneously tested. On the other hand, both noise-free WLsg and WLsc parameters appear as very robust eruptive predictors since the choice of the threshold does not affect the detection of the eruptive signature.

8. Summary and conclusions

The predictability of magnetic properties has been barely estimated using numerical simulations (see Sect. 1). In this work, we have studied the reliability of some eruptive flare indicators using the 3D parametric MHD simulations of flux emergence in Leake et al. (2013, 2014), leading to the self-consistent formation of either stable or unstable flux ropes. A set of four stable and three eruptive simulations have been used, corresponding to the partial emergence of a twisted magnetic flux tube into a stratified atmosphere, where a coronal overlying field is present. The eruption is triggered by a combination of external magnetic reconnection between the dipole and the overlying magnetic field, followed by internal reconnection allowing the vertical expansion and thus the ejection of the flux rope. The behavior of the emergence has been investigated through a range of initial coronal arcade strengths and orientations. Depending on the coronal arcade orientation, the magnetic configuration is either dipolar or quadrupolar, leading to respectively stable or unstable coronal flux rope (see Sect. 2.1).

Parameters have been examined in order to detect whether some physical changes specific to eruptive flaring ARs could be detected and measured, thus establishing a certain relation to flare occurrence. We rigorously analyzed the simulations as observations, including data masking and magnetic flux transport velocity derivation from the DAVE4VM code (Schuck 2008), using time series of 2D-plane magnetograms extracted from the 3D numerical datasets. A list of 99 parameters has been tested, including indicators currently used for both eruptive and non-eruptive flare operational forecasting as well as new quantities, such as helicity or current-weighted neutral line lengths Lsc and WLsc.

From the 99 predictors’ list used in this work, only the parameters relying on MPIL properties (see Table 5) demonstrated predicting eruptive flare capabilities. The other parameters showed a very similar evolution for both simulations, and no eruptive signature was detectable (see Sect. 4). The strong-shear MPIL length Lss, the strong-gradient MPIL length Lsg and the strong-current MPIL length Lsc are significantly longer for the unstable simulations due to the quadrupolar magnetic configuration generating additional external MPIL (see e.g., Fig. 7). The length increase rate depends on the coronal magnetic field strength: stronger coronal dipole is associated to longer property-weighted MPIL lengths. However, a strong rise of these three lengths is systematically observed during the early stage of the emergence, whatever the arcade field strength. The WLss, WLsg, and the WLsc parameters present very similar trends, making them also promising eruptivity predictors.

The detection of eruptive flare signatures, i.e., higher MPIL properties for eruptive ARs, is dependent on the initial masking process. Indeed, to isolate the ARs from the background magnetic field, a thresholding mask is usually applied when observations are analyzed. Notwithstanding, we demonstrated that the choice of the initial masking threshold is crucial for detecting MPIL property variations. The same analysis has been presented for both Bmask = 30 and 100 G (see Sect. 5), showing that the eruptive flare signatures disappear using the higher mask threshold. Noise in the measurements can also scramble the eruptivity signatures: through a Monte-Carlo scheme, we estimated the noise influence by including random perturbation to our time series magnetograms. Our results show that the four Lss, Lsg, WLss, and WLsg parameters are not strongly impacted by noise and their associated peaks prior to the eruption are still detectable. However, the Lsc and WLsc are more strongly impacted: the increase of Lsc is weaker, although still measurable, but the WLsc parameter no longer allows to discriminate between eruptive and stable simulations.

In addition, the MPIL properties, associated with valuable eruptive flare predictabilities, depend not only on the initial masking process and the noise, but also on additional thresholds imposed on physical properties, such as current, magnetic field, or gradients. We also investigated the impact of these physical thresholds through a parametric study measuring the detectability of eruptive signature as a function of the thresholding process (see Sect. 7). Results show that the three MPIL lengths Lss, Lsg, and Lsc and the WLss parameter are strongly sensitive to the choice of the physical threshold, whereas WLsg and WLsc are robust relative to threshold changes.

Our study is limited in terms of AR size and complexity that we cannot explore given the computational cost of such simulations. Given the small scale of our MHD simulations, we do not catch the same distribution of magnetic field strengths as observations. Because of the small flux (in the 1021 Mx range), and the absence of interaction between granulation and magnetic field, the intermediate strength fields are not caught. Hence, by applying the B-mask, which is done at these intermediate strengths, some features from the simulations may be lost. Besides, eruptive flare indicators have only been tested for a given AR flux and size, corresponding to the smallest observed flaring AR class (see Sect. 2.2 for details). Still, this type of analysis is yet relevant to connect flare physical models and observations, and provide a comprehensive perspective of what can be done using actual observations. Future work will be extended to AR simulations with different sizes and fluxes.

We therefore conclude that from our 99 predictors’ list, WLsg and Lsg are the best eruptive flare indicators for these model tests. The other parameters relying on MPIL properties tested in this study, namely Lss and Lsc, WLsc and WLss should be tested using various physical thresholds. The current-related parameters seem to be more sensitive to noise, even if a more detailed analysis of uncertainty sources is needed. Apart from the physical thresholds, we also recommend the testing of various masking processes for actual observational analysis, using both low and high values in order to detect complex MPILs, potentially indicating an imminent eruptive flare activity.

Acknowledgments

This work has received partial support by the European Unions Horizon 2020 Research and Innovation Programme under Grant Agreement No. 640216 (Flare Likelihood and Region Eruption Forecasting [FLARECAST]). JEL is supported by NASA’s Heliophysics program. Simulations used in this study were carried out using a grant of computational time from the DoD’s High Performance Computing Program. The editor and the authors thank two anonymous referees for their assistance in evaluating this paper and improving its quality.

References

Cite this article as: Guennou C, Pariat E, Leake JE & Vilmer N. Testing predictors of eruptivity using parametric flux emergence simulations. J. Space Weather Space Clim., 7, A17, 2017, DOI: 10.1051/swsc/2017015.

All Tables

Table 1.

Parameters relying on magnetic field tested in the present study.

Table 2.

Parameters relying on magnetic field geometry tested in the present study.

Table 3.

Current properties tested in the present study.

Table 4.

Lorentz force parameters tested in the present study.

Table 5.

Magnetic Polarity Inversion Line (MPIL) parameters tested in the present study.

All Figures

thumbnail Fig. 1.

Emergence of the convection magnetic flux tube into the solar atmosphere at t = 110 t 0 for both the stable SD (A) and the unstable SD E (B) simulations. The cyan lines trace the coronal arcade magnetic field, originating from the lower boundaries, whereas the red lines belong to the magnetic flux tube and originate from the y = ±maxy side boundaries. The gray scale indicates the magnetic field strength at the surface.

In the text
thumbnail Fig. 2.

Magnetograms of the Bz component at t = 100 t0 for both the MD (A) and MD E (B) simulations. Both non-masking (left) and masking (right) magnetograms are displayed, in order to show the reduction of the region of interest for further analysis. Eruptive and non-eruptive simulations present very similar magnetograms, except near the external edge of the polarities, due to the quadrupolar geometry above the surface.

In the text
thumbnail Fig. 3.

Some parameter evolution for the seven simulations, as a function of time in t 0 units (t 0 = 55.5 s). From (A–F), the following predictor evolution is respectively displayed: the total current helicity , the total current I tot, the mean value of the free energy , the total flux Φ tot, the total Lorentz force F and the Schrijver R value. The warm colors indicate the evolution for the three eruptive simulations {WD E, MD E, SD E}, respectively, corresponding to the magenta, orange, and red solid lines. Conversely, the evolution of the non-eruptive simulations {ND, WD, MD, SD} is displayed using the cold colors, respectively, as follows: purple, deep blue, light blue, and green solid lines. The eruption time is indicated by the vertical dashed gray line.

In the text
thumbnail Fig. 4.

Same as Figure 3, but for different parameters. From (A–F), the following indicator evolutions are respectively displayed: the kurtosis of the twist parameter κ(α), the direct current , the kurtosis of the horizontal magnetic field κ(B h ), the mean value of the shear angle , the time variation of the relative magnetic helicity dH m/dt and the mean value of the vertical magnetic field B z .

In the text
thumbnail Fig. 5.

Same as Figure 3, but for different parameters. From (A–F), the following indicator evolutions are respectively displayed: the length of the strong-shear MPIL Lss, the integral of shear angle along the MPIL WLss, the length of the strong-gradient MPIL Lsg, the integral of the Bz horizonal gradient along the MPIL WLsg, the length of the strong-current MPIL Lsc and the integral of the current along the MPIL WLsc.

In the text
thumbnail Fig. 6.

Same as Figure 5, but using a mask threshold Bmask = 100 G.

In the text
thumbnail Fig. 7.

Magnetograms of the SD (A–C) and SD E (D–F) simulations, including the MPIL (red lines). The MPIL width has been dilated by a factor 4 in order to increase its visibility. The masking threshold increases from left to right, with respective values of 30, 50, and 100 G.

In the text
thumbnail Fig. 8.

Same as Figure 5, but including the noise perturbation.

In the text
thumbnail Fig. 9.

Maps of the horizontal magnetic field Bh, for the SD (A–C) and the SD E simulation (D–F) at t = 100 t0. The white line represents the portion of the MPIL for which . The threshold increases from left to right, with values of  = 0, 50, and 100 G. No initial data masking is applied here, in order to highlight only the impact of the threshold on the computation of Lss, Lsc, WLss, and WLsc.

In the text
thumbnail Fig. 10.

Maps of the quantity, using the same layout as Figure 9. The white line now represents the portion of the MPIL for which the horizontal gradient of the magnetic field is higher than . From left to right the threshold increases from 0–100 G Mm−1.

In the text
thumbnail Fig. 11.

Maps of the vertical current Jz, using the same layout as Figure 9. The white line now represents the portion of the MPIL for which the current density . From left to right, the threshold increases, with respective values 0, 3, and 12 mA m−2.

In the text
thumbnail Fig. 12.

Parametric evolution of the total length of the strong-shear MPIL Lss, as a function of the two thresholds and . The varies for each figure, with respective values of 0, 10, 25, 50, 75, and 150 G, from (A–F). For each panel, the Lss quantity is represented for the strong arcade eruptive (SD E; solid lines) and strong arcade (SD; dashed lines) simulations. The threshold is changed for each curve, using the values of 30° (blue lines), 45° (red lines), and 60° (magenta lines). The eruption time is denoted by the vertical gray dot dashed line.

In the text
thumbnail Fig. 13.

Same as Figure 12 but for the parameter Lsg, as a function of the two thresholds and . The varies in the same way as for Figure 12. The threshold is changed for each curve, using the values of 0 (blue lines), 10 (red lines), 25 (magenta lines), 50 (green lines), 75 (yellow lines), and 100 G Mm−1 (black lines).

In the text
thumbnail Fig. 14.

Same as Figure 12 but for the parameter Lsc, as a function of the two thresholds and . The threshold is changed for each curve, using the values of 3 (blue lines), 6 (red lines), 12 (magenta lines), 18 (green lines), and 24 mA m−2 (yellow lines).

In the text
thumbnail Fig. 15.

Parametric evolution of the WLsc (A), WLsg (B), and WLss (C) parameters, as a function of the threshold . Solid lines represent the parameter evolutions for the SD E simulation, while dashed lines are computed for the SD stable simulation. The varies from 0 (blue lines) to 150 G (black lines), with the following intermediate values 10 (red lines), 25 (magenta lines), 50 (green lines), and 75 G (yellow lines).

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.