Magnetic Helicity Flux Oscillations in the Atmospheres of Flaring and Nonflaring Active Regions

Analyzing the evolution of magnetic helicity flux at different atmospheric heights is key for identifying its role in the dynamics of active regions (ARs). The three-dimensional (3D) magnetic field of both flaring and nonflaring ARs is constructed using potential field extrapolations, enabling the derivation of emergence, shearing, and total magnetic helicity components at a range of atmospheric heights. An analysis of temporal oscillations of the derived components shows that the largest significant period of the three helicity fluxes are common (within ±2 hr) from the photosphere up to at least 1 Mm for flaring ARs—being consistent with the presence of a coupled oscillatory behavior that is absent in the nonflaring ARs. We suggest that large, energetic solar eruptions may have been produced in ARs when the vertical and horizontal helicity flux components became a coupled oscillatory system in the low solar atmosphere.


Introduction
The most common source of large flares and coronal mass ejections (CMEs) is the magnetically complex and highly twisted δ-type active regions (ARs) (e.g., Georgoulis et al. 2019;Toriumi & Wang 2019, and references therein). Threedimensional (3D) models of the coronal magnetic field above ARs, constrained by the underlying photospheric observations, enable a search for promising observational-based signatures of flare and CME productivity. These include the stored free magnetic energy of an AR, which is related to the frequency and intensity of an upcoming flare or a series of flares (Jing et al. 2010;Su et al. 2014). While the AR free magnetic energy is a requirement for flaring, it does not directly dictate whether the flare(s) will occur, nor directly predict a flare's intensity, since not all of the free magnetic energy is necessarily released at once (see, e.g., Gilchrist et al. 2012, and references therein). The free magnetic energy of an AR is not directly related to the complexity of the magnetic field. Another derived value, the magnetic helicity, does depend on the strength and the structure of the field.
Helicity, itself, can be visualized as the net linking or winding of individual magnetic field lines around each other (Moffatt 1969). This arrangement of the magnetic field lines may correspond to the lowest magnetic energy state because they cannot pass through one another or reconnect (Moffatt 1990;Freedman & He 1991;Berger 1993). The two primary sources of magnetic helicity are (i) the shearing of magnetic field lines by footpoint motions, and (ii) the emergence of twisted fields, which originate from the solar interior. Magnetic helicity is invariant in the highly conducting plasma of ideal magnetohydrodynamics, and is almost conserved even for a nonideal finite resistivity (Berger 1984;Pariat et al. 2015).
Using photospheric magnetic field observations is still essential to analyze the magnetic helicity evolution (Elsasser 1956) before flare and/or CME events. Park et al. (2008) discovered that helicity flux can slowly increase and then remain constant prior to some flares. Smyrli et al. (2010) and Park et al. (2012) reported that the injected helicity flux changed its sign before some very impulsive eruptive events. Furthermore, a strong correlation was found between the free magnetic energy and the magnetic helicity in flaring ARs by Tziotziou et al. (2012). They suggest that magnetic helicity is an essential ingredient, alongside free energy, for the production of intense flare and CME eruptions.
Studying the evolution of magnetic helicity in 3D can help reveal and understand the processes of an AR leading to a flare and/or a CME. Based on a 3D numerical analysis, Pariat et al. (2017) argues that the magnetic energy and the total relative helicity are not effective diagnostics for flare prediction. However, the relative magnetic helicity (Berger 2003) can be decomposed into the current-carrying component and its counterpart. Thalmann et al. (2019) had a similar conclusion to Pariat et al. (2017) but with observation-based 3D magnetic field extrapolation data analyses. Thalmann et al. (2019) found that the ratio of current-carrying to total helicity is useful to identify an eruptive AR. Furthermore, Vemareddy (2021) found that, based on the evolution of helicity in 3D, successive injection of opposite helicity could invoke CME(s) in an AR. Recently, Gupta et al. (2021) performed a comprehensive analysis of the coronal magnetic energy and helicity evolution around the time of large flares. They found that if the values of helicity and magnetic energy exceeded a certain threshold, then an AR is likely to produce large, CME-associated flare(s). However, when the helicity and magnetic energy remained below the threshold then an AR only would be a host for a confined flare.
A recent theoretical paper, based on multiresolution wavelet decomposition, by Prior et al. (2020) reported that spatial scales of magnetic field helicity is consistently additive. This theory is supported by, or at least is consistent with, Korsós et al. (2020b) and Soós et al. (2022) in their observation-based helicity analyses. They found a relationship between the flaring activities and the unique oscillatory behavior pattern of the emergence, shearing, and total magnetic helicity flux components in the photosphere. Their results motivate us to extend and apply the method of Korsós et al. (2020b) and Soós et al. (2022) to the 3D magnetic field of ARs. In particular, here, we will analyze the evolution of the emergence, shearing, and total helicity flux components at a range of atmospheric heights, in the case of 14 flaring and 10 nonflaring δ-type ARs. We investigate the evolution of these three helicity flux components and test their ability to indicate the eruptive potential of a δ-type AR in 3D.
The work is organized as follows: Section 2 describes the adopted tools for the magnetic helicity flux calculations and the selection criteria of the investigated ARs. Section 3 introduces and describes the application of our analyses for one typical flaring and one nonflaring AR. Section 4 discusses our findings in the case of a further 14 flaring and 10 nonflaring δ-type ARs. Key findings and conclusions, along with plans for future work, are given in Section 5.

Data Selection
Magnetic helicity is one of the typical scalar parameters that contribute to the description of the dynamic evolution of an AR at any moment of time in a 3D volume. Therefore, we propose to study and compare the derived evolution of the emergence, shearing, and total magnetic helicity flux components at different lower atmospheric heights, in the case of 14 flaring and 10 nonflaring ARs, which satisfy the following selection criteria: 1. The AR has a δ-spot(s), because the helicity flux injection is strong in this type of configuration. 2. The AR shows one of two distinct behaviors in terms of flare activity. Namely, the AR is either host to intense X-class flares (the "flaring" group), or the AR only produces B-and C-class flares (the "nonflaring" group). 3. The angular distance of an AR from the solar disk center is within ±60°, to optimize the quality of the data (Bobra et al. 2014).

Magnetic Field Extrapolation for 3D Analyses
Direct observations of the 3D magnetic field in the lower solar atmosphere are currently only possible for exceptional observing conditions and under certain special circumstances (e.g., Kuridze et al. 2019). Instead, the measured line-of-sight component or the full vector magnetic fields observed at the photosphere can be used as a lower boundary to model the atmospheric magnetic field using the potential free (PF), or more complex nonlinear force-free field (NLFFF) model extrapolations. Other than advanced methods that couple magnetic modeling with extreme ultraviolet observations of the distribution of AR coronal loops (Aschwanden 2016), such extrapolations lack direct observational constraints above the photosphere-thus it is a challenging task to quantify the accuracy of the model field (Wiegelmann & Sakurai 2012). Furthermore, the majority of models give a static approximation to a dynamic system, and this is a particular concern for ARs. However, without improved observational constraints , the static extrapolation models currently remain our best practical approach for approximating AR magnetic fields. Dalmasse et al. (2018) studied the 2D and 3D transfer of magnetic helicity in ARs with a connectivity-based helicity flux density method. They found that the helicity flux density calculations derived from PF and NLFFF extrapolations were highly correlated (the Spearman rank correlation coefficient is 0.72) and consistent with each other. They showed that the local sign of helicity flux given by the two different extrapolations are in good agreement (∼85%). Korsós et al. (2020b) reported that their flare prediction method performed well using either PF or NLFFF extrapolated data; however, the prediction could be done more often earlier with employing PF data overall.
Motivated by the findings of Dalmasse et al. (2018) and Korsós et al. (2020b) the 3D solar magnetic structures of 14 flaring and 10 nonflaring ARs are calculated by PF. The PF approximation serves as a tool for initial insight, but has obvious limitations for this application, e.g., the quasi-static evolution may not reflect the dynamics of the actual magnetic field in important ways. Also, there are regions in the lower atmosphere where the magnetic field is likely to be highly nonpotential; thus, care is needed to not overinterpret dependent results. One way forward is to apply the PF approach to a much larger sample. For practical reasons, we reserve this for a future study and limit this work to a smaller sample, thus showing a proof of concept. To determine the magnetic field above the photosphere with PF extrapolation, we used the hmi. sharp_cea_720s photospheric vector magnetic field measurements at every moment of time, namely the B r , B t , and B p components of the Spaceweather Helioseismic Magnetic Imager Active Region Patches (SHARPs; 6 Bobra et al. 2014).
In this work, we study the extrapolated magnetogram data from the z = 0 level (as the photosphere) up to 3.24 Mm with a 60 minute cadence. Therefore the three components of the vector magnetograms are smoothed from the original components by PF at z = 0. This allows us to work with a homogeneous data set from z = 0 up to 3.24 Mm with a step size of z = 0.36 Mm, which is the same as the SHARP pixel size. The extrapolated magnetic field strength of each voxel is recorded for each relevant frame and saved in dedicated data files. The results of two example extrapolations are shown in Figures 1(a)-(b) for ARs 11158 and 11775.

Magnetic Helicity Flux Calculation
Following Korsós et al. (2020a) and Soós et al. (2022), Equation (1) by Berger (1984) is used to calculate the time evolution of the magnetic helicity flux (i.e., the helicity injection rate) through the photosphere and over the atmospheric domain: ere, A p is determined by the surface magnetic field and the Coulomb gauge (Berger 1997;Berger & Ruzmaikin 2000), which is the vector potential of the potential magnetic field B p . B h and B z are the tangential and normal components, respectively. The v ⊥h is the tangential velocity component, while v ⊥z is the normal one. On the right-hand side of the equation, the first term is the emergence term associated with emerging twisted flux tubes. The second term is generated by the shearing of the field lines by tangential motions in the photosphere and is called the shearing term. The total magnetic helicity flux is calculated as the sum of the emergence and shearing components.
Analog to Korsós et al. (2020a) and Soós et al. (2022), the plasma velocity is calculated by applying the differential affine velocity estimator for vector magnetograms (DAVE4VM 7 ) algorithm (Schuck 2008). In this study, the physical quantities are determined at several atmospheric heights above an AR, including the photosphere; therefore the DAVE4VM algorithm is applied for each discrete height within the PF model domain. We use a window size of 19 pixels for the velocity calculations, which is derived based on the nonparametric Spearman rankorder correlation coefficients, Pearson correlation coefficients, and slopes between Δ h · (v z B h − v h B z ) and δB z /δt (Schuck 2008). The vector potential A p is calculated by MUDPACK, a multigrid software for solving elliptic partial differential equations (Adams 1993). The velocity was determined from the time sequence of the calculated magnetogram maps at a certain height. In this case, the magnetic maps themselves are changing as a function of time, where there is an inherent photospheric velocity accounted for in the observation.
Before the emergence, shearing, and total magnetic helicity flux components are calculated by Equation (1), we define lower and upper magnetic field thresholds, as was set up in Soós et al. (2022). A 200 | | G value is chosen as the lower threshold to avoid regions with low signal-to-noise (Deshmukh et al. 2020). The voxels with field magnitude larger than 2000 | | G was taken as upper boundary to damp the artificial 12 and 24 hr periods of the Solar Dynamics Observatory (SDO)/ Helioseismic and Magnetic Imager (HMI) magnetogram measurements, based on Table 1 of Smirnova et al. (2013). The artificial 12 and 24 hr periods correspond to Doppler shifts that move the spectral lines back and forth by about one tuning step every 12 h . Furthermore, this is combined with the solar rotation, which causes temporal and spatial variations of the inverted magnetic field at every 24 hr. The time series of the three magnetic helicity flux components (emergence (EM), shearing (SH), and total (T)) are further normalized by their respective largest absolute value in order to facilitate comparison on similar scales. The normalized time series are then smoothed with a 24 hr sliding window average. The resulting smoothed series is subtracted from the original normalized data series, which could further mitigate the 24 hr SDO/HMI period artifact, if it has remained in the normalized EM, SH, and T time series after the applied 2000 | | G threshold. Also, the resulting smoothed series is subtracted from the original normalized data series that helps to filter out the redundant long periods to enhance the shorter periods (McAteer et al. 2002).

Data Analysis
After data preparation, the time series of the smoothednormalized EM (dotted line), SH (solid line), and T (dasheddotted line) helicity flux components are analyzed, as is demonstrated in the cases of AR 11158 and 11775 in the top panels of Figures 2-3.
For further analyses, we use the wavelet analyses software 8 developed by Torrence & Compo (1998) to generate the wavelet power spectrum (WPS) of the time series of EM, SH, and T, with (i) a red-noise background (α = 0.72), (ii) zero padding, and (iii) a default Morlet wavelet profile. The associated global power spectrum (GPS) is also calculated as the time-averaged WPS. For each studied AR, to identify significant periodicities at every investigated height, the confidence level is one σ, as was set in Korsós et al. (2020a) and Soós et al. (2022).
To reveal how the oscillatory behavior of the three helicity flux components develop in 3D, we identify the local maxima of the identified average periods by GPS within the one-σ contours of the EM/SH/T WPS plots (red dots in Figures 2-3) at all atmospheric heights. Following Soós et al. 2022, we adopt the Python package 9 implementation of the zeroth dimensional persistent homology (Huber 2021) to find these local maxima within the one-σ contours of the EM/SH/T WPS plots. After the local maxima of the contours in the EM/SH/T WPS are successfully identified, we investigate their evolution as a function of height, see, e.g., in Sections 3.1-3.2.
Some of the local maxima are in the area of the cone of influence (i.e., the hatched area in Figures 2-3), but only those are taken into account where the corresponding period is also within the cone of influence (i.e., not in the hatched area in Figures 2-3). Here, we also note that the actual value of the local maxima may differ slightly from the value of the corresponding average period obtained by GPS. Therefore, the value of the average period by GPS is taken instead of the actual value of the local maxima, while we study the evolution of the local maxima as a function of height.

Example of Flaring AR 11158
We first demonstrate the 3D analyses of a flaring AR using the case of AR 11158. In   Similar to the examples above, the local maxima of one-σ significant periodicities are recorded for EM/SH/T at height increments of 0.36 Mm within the PF model domain. Figure 4 plots the local maxima of the identified periodicity as a function of height to reveal how the oscillatory behavior of EM (stars), SH (crosses), and T (dots) behave in the lower solar atmosphere before the X2.2 flare onset. In Figure 4, the color of each plotting symbol represents the value of the identified periodicity by GPS, as indicated by the color bar.
From Figure 4(a), we could conclude that the measured ∼18 hr period of oscillation of the EM/SH/T (in Figures 2(a)-(b)) is the largest and common period of the three helicity flux components from the photosphere up to 1.44 Mm prior to the flare. Above 1.44 Mm the largest period of the EM/SH/T is ∼26/23/23 hr prior to the flare.

Example of Nonflaring AR 11775
We now focus on AR 11775 in detail, as a nonflaring AR example. This AR was in fact host to a few low-energy C-class flares; however, such low-energy class flares we do not consider to be significant. Figure 3 shows the results of the WPS and GPS analysis for this AR. We find the following: 1. Photosphere. In Figure 3( Figure 3(c), interestingly, the EM component has a similar periodicity to the SH and T helicity flux components. In the EM time series, a 25.9 hr periodicity can be identified for a time span of ∼4.5 cycles. Meanwhile, the SH and T helicity flux components show a 24.2 hr periodicity from the beginning of the time series. This 24.2 hr period declines after ∼5.3lifetime cycles. The corresponding significance ratio intervals of the local maxima can be seen on the color bars.
The local maxima of the periods are recorded in the WPS of EM/SH/T, respectively, from the photosphere up to 3.24 Mm at every 0.36 Mm height, and are plotted in Figure 4(b). We can see that a ∼25 hr period becomes the largest and most common period of the three helicity flux components from 1.08 to 3.24 Mm. However, the ∼25 hr period is absent from the time series of EM helicity flux component, from the photosphere up to 1 Mm.

Analyses of More Active Regions
The analysis of Section 3 is extended for all the selected 14 flaring AR and 10 nonflaring ARs. Similar to AR 11158 and 11775, the evolution of the three helicity flux components vary as a function of height in all the AR cases. In the flaring AR cases, based on Figure 5(a) and Table 1, the actual value of the largest periods (the value taken by GPS) of each helicity flux component are very close to each other at z = 0 level, approximately only by a maximum of 2 hr difference. These largest periodicities remain contemporaneously the largest and most common periodicity of EM, SH, and T from z = 0 up to a certain height before an X-class flare.
In the nonflaring AR cases, the largest periods of SH and T are similar (see Figure 5(b) and Table 2). However, the largest period of EM is different compared to SH and T from the photosphere up to at least 1 Mm. Similar to AR 11775, in the case of a few nonflaring ARs, the largest periods of each three helicity flux components are approximately similar above heights of 1 Mm, e.g., ARs 11512, 11635, and 11793 in Table 2.
We also investigated the corresponding significance ratio of the largest period of the EM, SH, and T components as a function of height. In Table 3, we summarized the corresponding significance ratio intervals of the identified largest periods at z = 0 and the optimum height, in the flaring AR cases. While, for the nonflaring ARs, Table 3 includes the corresponding significance ratio of the largest period of the three helicity flux components only at z = 0. Based on the significance ratio analyses, we could not find any distinctive behavior pattern between the flaring and the nonflaring AR cases.
The appearance of a common period in certain height ranges is an important difference between the flaring and nonflaring AR cases. Based on our result, when the horizontal and the vertical helicity flux components become a coupled oscillatory system from the photosphere up to the chromosphere by photospheric footpoint motions then the stored free magnetic energy could be released as an intense flare event.
Furthermore, we also tested the conjecture of Korsós et al. (2020b) of the existence of an optimal height, where some flare precursors manifest earlier than in the photosphere. Here, we consider an optimal height where the local maxima of a common peak periodicity can be detected earlier compared to the photosphere. Based on Korsós et al. (2020a) and Soós et al. (2022), we now define the optimal height via EM because the periodic behavior of the EM helicity flux component plays such an important role in the flare activity of an AR.
In Figure 6, we can see that AR 11515 behaves differently to the other flaring ARs because the optimal height of EM is at the photosphere. However, in the other 13 flaring cases, the EM local maxima of the common period can be observed 2-8 hr earlier between 0.36 and 1.5 Mm compared to the photosphere. In three flare cases, the EM local maxima of the common period can be observed tens of hours earlier above 1.5 Mm. Here, we also note that these defined optimal heights for EM are not the dedicated optimal heights of the SH and T components, as can be seen in Figure 6, especially when the T Gain is negative (e.g., AR 11166, 11890, etc.).
Interestingly, the 2-8 hr gain time found between in 0.36-1.5 Mm height range is similar to the result of Korsós et al. (2020b), who stated that "1 and 1.8 Mm above the solar surface would improve the prediction of the flare onset time by around 2-8 hr." The conclusion of the identified optimal height range with the EM helicity flux component is not as strong as that of Korsós et al. (2020b), but still offers support to improving flare prediction methods through focus on the lower solar atmosphere.

Summary and Discussion
Recently, magnetic helicity has been found to be an invaluable tool for understanding underlying physical processes of solar activity. Therefore some helicity-based quantities could be an efficient approach for the purpose of flare and CME prediction (see, e.g., Pariat et al. 2017;Thalmann et al. 2019;Gupta et al. 2021;MacTaggart et al. 2021;Soós et al. 2022, and references therein), but it still remains a challenging task to employ this information as a practical tool in the context of flare and CME forecast.  Table 1.
Comparing the evolution of the helicity fluxes, or its components (EM, SH, and T), between the flaring and nonflaring δ-type ARs in 3D is important since it reflects the dynamic evolution of an δ-type AR before solar eruptions in the lower solar atmosphere. We focus on δ-type ARs because these types of ARs tend to produce large solar eruptions (e.g., Georgoulis et al. 2019;Toriumi & Wang 2019, and references therein). In this work, we studied the periodic patterns in the evolution of emergence, shearing, and total helicity components in the case of 14 flaring and 10 nonflaring δ-type ARs at different solar atmospheric heights from the photosphere to the lower corona.
Our aims were realized by (i) implementing a potential (PF) exploration technique and (ii) creating a sample of 3D magnetic maps of ARs using the SHARP 0.36 Mm pixel size in the zdirection, above the photosphere up to 3.24 Mm. We calculated the time series of the emergence, shearing, and total magnetic helicity flux components of the investigated 24 ARs by using the DAVE4VM algorithm (Schuck 2008) from the photosphere up to 3.24 Mm, at height increments of 0.36 Mm. After normalizing, the time series are smoothed with a 24 hr moving average window that was subtracted from the original normalized data series to reduce the problematic 24 hr SDO/HMI artifact (Higgins et al. 2022). Finally, wavelet analysis is performed on the time series. The following results are revealed by wavelet analyses: 1. In the flaring AR cases, the largest photospheric periods of the EM, SH, and T helicity fluxes are also a common period of the three helicity flux components. This common periodicity tends to appear before the occurrence of a large flare. This photospheric common period remains the largest and common period in the time series of the EM, SH, and T helicity fluxes at least up to 1 Mm or even higher, before the flare event(s) (see Figure 5(a)). 2. The 3D analyses of the helicity fluxes corroborate the finding of Korsós et al. (2020b). They define an optimum height range (1-1.8 Mm) where the potential lead-time improvement is 2-8 hr for predicting an upcoming intense flare. In the case of 13 out of 17 flare cases, the local maxima of the dedicated common period in the oneσ contours of the WPS of EMs can be measured 2-8 hr earlier between 0.36 and 1.5 Mm compared to the photospheric counterpart. This is consistent with Korsós et al.'s (2020b) findings, namely the 2-8 hr lead-time improvement within the 1-1.8 Mm optimum height range. During this work, the optimum height of each studied case was defined by the EM component because the dynamic evolution of the EM helicity flux component plays an important role in the production of intense solar Note. GOES X-flare class intensity before the largest and common period of EM, SH, and T appear; the detected largest period of the EM, SH, and T helicity flux components in the photosphere; the atmospheric height interval where the largest photospheric period remain the largest and the common period of the three helicity flux components before the flare(s); defined optimum height by EM; time difference between the moment of the local maxima at the photosphere and the moment of the local maxima at the EM optimum height, which also refer to how much is the gain time in the case of EM helicity flux component. These defined optimal heights by EM are not the dedicated optimal height of SH and T components, as you can see via negative values in the case of SH and T. Note. The columns are the NOAA number of AR; the detected largest period of the EM, SH, and T helicity flux components in the photosphere; when the largest period of the EM, SH, and T became the common period of them somewhere in the solar atmosphere; the atmospheric height interval where the largest periods are also a common period of the three helicity flux components. eruptions in δ-type AR (Korsós et al. 2020b;Soós et al. 2022). 3. In the case of nonflaring ARs, the largest period of the EM is different when compared to SH and T helicities at the photosphere. However, in some cases, it was found that the largest period of each helicity flux becomes a common period, just above 1 Mm (see Table 2).
Based on our results, we conclude that to release the stored free magnetic energy of an AR then one of the conditions could be when the horizontal and vertical helicity flux components become a coupled oscillatory system from the photosphere up to the chromosphere due to photospheric footpoint motions.
In this paper, we analyzed 14 flaring and 10 nonflaring events that gave us a first insight. One has to be careful not to overinterpret the results, and analyses on a larger number of samples may indeed be needed to confirm these initial findings. In particular, a limitation of this study is that the actual value of the identified periods (see Tables 1-2) can be affected by factors such as (i) how long the studied time series is and (ii) the implementation of the wavelet analysis (e.g., the type of mother wavelet function).
However, in Soós et al. (2022), it was shown that next to the identified fundamental periods some of the harmonics were also detectable in the time series of the three helicity flux components. Therefore in this present work, the different Note. The columns are the NOAA number of AR; the detected largest period of the EM, SH and T helicity flux components at z = 0; the significance ratio interval of the EM, SH and T periods at defined optimum height by EM in case of flaring ARs. Figure 6. Summary of the optimal heights and lead times based on the evolution of EM helicity flux component for 14 AR cases under the PF extrapolation. The three columns of each AR (line-crossed columns EM, filled columns SH, dotted columns T) show the time difference between the moments of the local maxima at the optimum height and the photosphere. The color bar represents the atmospheric height, which was identified as an optimum height by the EM helicity flux component (not the optimal heights defined by the SH and T components). For the actual values, see Table 1. global oscillatory behavior could be identified between the flaring and nonflaring cases with the wavelet analysis of a longenough data series (i.e., satisfying the requirement set by the Nyquist frequency), but further conclusions cannot be made based on the exact value of the periods. We also note that there is some measured significant power at ∼24 hr (see in Table 3) but it relates to real periods because (i) we used an upper magnetic field threshold value to damp out the artifact (Smirnova et al. 2013), and (ii) independent studies also found these types of long periods, e.g., see Table 1 of Griñón-Marín et al. (2020) for a summary.