EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 675 (July 2023) A&A, 675 (2023) A80 Full HTML
Open Access
Issue
A&A
Volume 675, July 2023
Article Number A80
Number of page(s) 12
Section The Sun and the Heliosphere
DOI https://doi.org/10.1051/0004-6361/202244650
Published online 03 July 2023

© The Authors 2023

Licence Creative CommonsOpen Access article, published by EDP Sciences, under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

This article is published in open access under the Subscribe to Open model. Subscribe to A&A to support open access publication.

1. Introduction

Coronal mass ejections (CMEs; e.g. Webb et al. 2000) are of primary interest for space weather as they cause the largest disturbances in the near-Earth space environment (e.g. Kilpua et al. 2017; Zhang et al. 2007). Forecasting their impacts on our planet requires modelling their properties and evolution from the Sun to the Earth. Both semi-empirical and physics-based space weather models require realistic input parameters to constrain their magnetised CME models (e.g. Kilpua et al. 2019). In particular, recent advancement includes runs with spheromaks and flux rope CMEs with 3D magnetohydrodynamic heliospheric models such as the EUropean Heliospheric FORecasting Information Asset (EUHFORIA; e.g. Pomoell & Poedts 2018; Scolini et al. 2019, 2020; Verbeke et al. 2019; Asvestari et al. 2021) and Space-weather-forecast-Usable System Anchored by Numerical Operations and Observations (SUSANOO; e.g. Shiota & Kataoka 2016). The goal is to predict a time series of magnetic field vectors impinging the Earth with a good level of accuracy.

Currently, the magnetic field in the corona is difficult to routinely measure reliably and with global coverage due to the relatively weak field strength and the plasma being hot and tenuous, resulting in the broadening of spectral lines (Lin et al. 2000). The information of intrinsic CME flux rope parameters must thus be obtained by other means, using either indirect observational proxies or data-driven coronal modelling (e.g. Kilpua et al. 2019; Palmerio et al. 2017, and references therein). The knowledge of the magnetic field of a CME is not only important for space weather but also for fundamental investigations of the formation, destabilisation, and eruption of solar magnetic flux ropes (e.g. Green et al. 2018; Welsch 2018; Kilpua et al. 2019).

From a space weather forecasting point of view, a data-driven model for eruptive coronal fields has to be computationally fast. Time-dependent magnetohydrodynamic (MHD) simulations are in principle well suited for investigating flux rope dynamics in the corona (e.g. Jiang et al. 2016), but they are computationally expensive, and not all necessary boundary conditions are routinely available from observations. An alternative approach is to simplify or even neglect the thermodynamics, which is justified by the dominance of magnetic forces in the low corona. One such approach, the non-linear force free field (NLFFF) method, uses photospheric magnetograms to extrapolate a static snapshot of the coronal magnetic field at a particular instance in time (Wiegelmann & Sakurai 2012).

The magnetofrictional method (MFM; Yang et al. 1986) is a particular method for computing NLFFF configurations. The basic principle of the method is the addition of a frictional term −νv, where ν is the frictional coefficient and v the plasma velocity, to the MHD momentum equation and assuming a low-beta (pressure gradient ignored) and quasi-static evolution. With these assumptions, the MHD momentum equation reduces such that an explicit relation for the magnetofrictional velocity is obtained and given as , which is further used to evolve the magnetic field according to Faraday’s law via an electric field obtained by invoking Ohm’s law. The velocity, which is proportional to the Lorentz force, evolves the system towards a force-free minimum-energy state in the absence of a net Poynting flux.

The MFM can be extended to a time-dependent version using time-varying photospheric boundary conditions (van Ballegooijen et al. 2000). Several studies have demonstrated that data-driven time-dependent magnetofrictional modelling (TMFM) can successfully describe the formation and initial rise of solar flux ropes when using a time-sequence of boundary conditions derived from photospheric observations (Cheung & DeRosa 2012; Gibb et al. 2014; Fisher et al. 2015; Yardley et al. 2018; Pomoell et al. 2019; Price et al. 2019, 2020; Kilpua et al. 2021; Lumme et al. 2022).

The sole driving boundary condition required by TMFM is the photospheric electric field. The so-called inductive electric field component (EI) is obtained directly by inverting Faraday’s law using a time series of photospheric magnetograms as input, while to specify the remaining non-inductive (curl-free) component one has to incorporate additional information, such as using Dopplergrams and optical flow methods (Kazachenko et al. 2015) or, alternatively, ad hoc prescriptions (Lumme et al. 2017; Welsch 2018; Yeates & Bhowmik 2022). Previous studies have shown that the non-inductive electric field component has a significant contribution to the total electric field, and ignoring it estimates only a certain percentage of the Poynting flux (Fisher et al. 2010; Kazachenko et al. 2014; Price et al. 2019). The non-inductive component has also been found to be critical for the formation and rise of flux rope configurations in magnetofrictional simulations (e.g. Pomoell et al. 2019; Cheung & DeRosa 2012).

The optimisation method described in Lumme et al. (2017) estimates the non-inductive electric field by employing ad hoc assumptions on the functional form of the sources responsible for generating the electric field. In the process, free parameters in the non-inductive electric field are introduced and are constrained by finding the best match with the energy injection computed by Differential Affine Velocity Estimator for Vector Magnetograms (Schuck 2008). While well defined, this approach yields an electric field based on a single metric while ignoring other possible metrics, such as the injection of helicity or the behaviour of the coronal magnetic field when used to drive a coronal model.

It is an open question how sensitive the formation and evolution of coronal flux ropes and their key parameters are to the choices made in estimating the photospheric non-inductive electric field. The purpose of this study is to explore how the choice of the non-inductive electric field in the TMFM simulation affects the flux rope formation and early evolution, as well as their key properties. We performed a detailed study of the effects of the driving electric field on AR 12473 from December 22, 2015 to January 2, 2016. Price et al. (2020) performed a detailed study of the magnetic evolution of AR 12473 using TMFM with the energy-optimised non-inductive driving electric field. In the present study, we show how varying the optimisation parameters affects the evolution and eruption of the flux rope. The paper is arranged as follows. In Sect. 2, we describe the active region (AR) we studied and the preparation of the data for TMFM. We give a detailed description of the method used to perform the TMFM. We describe the analysis of the simulated data in Sect. 3 and provide results obtained. Finally, we conclude this study in Sect. 4 with the summary.

2. Data and methods

2.1. Event overview

We studied the evolution of AR 12473 from December 22, 2015 to January 2, 2016. The location of this AR on the aforementioned days was at S22E661 and S21W862 (heliographic coordinates) as seen with the Helioseismic and Magnetic Imager (HMI; Scherrer et al. 2012) on board the Solar Dynamics Observatory (SDO; Pesnell et al. 2012). The sunspot associated with the AR was of a bipolar (β) nature. There were several C- and M-class flares associated with this AR (e.g. Mulay et al. 2021).

On December 28, 2015 at approximately 11:30 UT, an eruption was visible in the multi-wavelength observations from the Atmospheric Imaging Assembly (AIA; Lemen et al. 2012) on board the SDO. On the day of this eruption, the AR was located close to the disc centre (S23W11). An M1.9 class flare was associated with this eruption. The start, peak, and end times of this flare were 11:20 UT, 12:45 UT, and 15:00 UT, respectively. There was a halo CME associated with the eruption as viewed from Earth3, which was first seen at 01:12 UT by the C2 coronagraph of the Large Angle Spectroscopic Coronagraph (Brueckner et al. 1995) on board on the Solar and Heliospheric Observatory (Domingo et al. 1995). The linear and second-order speeds of the CME were ≈1212 km s−1 and ≈1471 km s−1, respectively4. An interplanetary CME (ICME) was recorded on December 31, 2015 with the shock and leading edge at ≈00:50 UT and ≈17:00 UT5, respectively, with the Advanced Composition Explorer (ACE; Chiu 1998).

2.2. Data

For the data-driven coronal simulation (Pomoell et al. 2019), we used a time sequence of disambiguated vector magnetograms obtained by HMI. The vector magnetograms have a spatial and temporal resolution of ≈0.5″ in the plane of sky and ≈700 s, respectively. The HMI data were downloaded from the Joint Science Operations Center (JSOC)6 using the ELECTRIC field Inversion Toolkit (ELECTRICIT) developed by Lumme et al. (2017).

2.3. Methods

2.3.1. Vector magnetogram processing

The data were processed using ELECTRICIT to make them suitable as input for the simulation. For this, firstly, the very strong and weak magnetic field vector pixels (B ≥ 750 Mx cm−2 and B ≤ 250 Mx cm−2, respectively) were disambiguated using the minimum energy method (Metcalf 1994) and the potential field acute angle method (Liu et al. 2017), respectively. The spurious pixels in the vector magnetograms were then removed by checking the formal error of total magnetic field strength by the Stokes inversion module (Hoeksema et al. 2014), which was fixed as σB = 750 Mx cm−2. These spurious pixels were then replaced with the median of their surrounding good pixel values.

This disambiguated dataset was used to produce a time series of re-projected magnetogram cutouts of the AR from the start to the end date of the observation. Using ELECTRICIT, we tracked AR 12473 across the disc and produced a time series of cutouts in a coordinate system with the AR remaining in the centre of the cutout. The tracking is similar to the Space-weather HMI Active Region Patches (SHARPs; Bobra et al. 2014; Hoeksema et al. 2014), which tracks each HMI Active Region Patch (HARP) using a fixed rotation rate. We also masked the noise-dominated weak-field pixels. The data was re-projected using a Mercator map projection centred at the region of interest to a local Cartesian system.

The re-projected time series cutouts were produced from 18:00 UT December 22, 2015 to 22:00 UT January 02, 2016, which included the heliographic longitude from ≈ − 50° to ≈ + 50°. The previous and later data were discarded due to the poor quality of magnetograms closer to the limb (Sun & Norton 2017). The dataset contained the emergence of AR 12473, the evolution of the AR, and the eruption that took place on December 28, 2015 at ≈11:30 UT. This time series was used for optimising the free parameters of the electric field inversion.

2.3.2. Electric field inversion methodology

The photospheric electric field E is decomposed into two components: (i) inductive EI; and (ii) non-inductive −∇ψ, as

(1)

The inductive component is determined from Faraday’s law directly using the set of photospheric magnetic field measurements, .

In order to determine the non-inductive component which is not constrained by Faraday’s law, we employed ad hoc assumptions on the form of the divergence of the non-inductive electric field. They are given by Cheung & DeRosa (2012), Cheung et al. (2015), Lumme et al. (2017), Pomoell et al. (2019):

(2)

(3)

(4)

where, Ω (units rad s−1) and U (units in Mm s−1) are free parameters in the non-inductive electric field, which are assumed to be spatially and temporally constant (see their physical interpretation in Lumme et al. 2017; Pomoell et al. 2019). These assumptions are henceforth referred to as the (i) zero- (Eq. (2)); (ii) Ω- (Eq. (3)); and (iii) U- (Eq. (4)) assumption, respectively.

2.3.3. Choice of U and Ω

To obtain the values of the constant U and Ω, we optimised them with respect to an independent metric, namely the injection of photospheric magnetic energy as a function of time, EM(t). As the target, we used the magnetic energy injection produced by the photospheric electric field computed using the ideal Ohm’s law:

(5)

where VDAVE4VM is the photospheric flow velocity vector determined by DAVE4VM (Schuck 2008) on the basis of optical flow methods.

The photospheric magnetic energy injection to the upper solar atmosphere is computed by integrating the vertical Poynting flux (Sz) over time and over the area (A) of the magnetogram (Kazachenko et al. 2015):

(6)

where, μ0 is the magnetic constant. For a given choice of ad hoc assumption and value for the associated U and Ω, the photospheric electric field can be inverted and the magnetic energy injection determined. Thus, a given reference magnetic energy injection is sufficient to determine the non-inductive component of the driving electric field. The optimised values of the U and Ω were found with an automated procedure. The procedure consists of a bracketing-like procedure in which the parameter space is discretised and for each value comparing the resulting energy injection to the reference DAVE4VM estimate by computing the root mean square (rms) difference until a minimum deviation is found. This was done for the data where the magnetic field was masked with a threshold of < 250 Mx cm−2 (hereafter called masked data). This procedure (i.e. Eqs. (5) and (6)) was also done for the data without any magnetic field masking and is hereafter called unmasked data. Using this method, the value of optimised U and Ω was found to be 120 Mm s−1 and 0.06 rad s−1, for unmasked data, respectively, and 200 Mm s−1 and 0.13 rad s−1, for masked data. The accuracy in the estimates of U and Ω is the step of increment and decrements used for finding the optimised values, which were ±5 ms−1 and ±0.01 × 2π day−1, respectively, for U and Ω runs.

The aim is to study the effects that the choice of U and Ω values have on the data-driven coronal simulation. By varying U and Ω values, we can explore the changes in flux rope properties and evolution. In order to achieve this, we prepared additional electric field datasets using modified values for U and Ω. Specifically, we constructed eight additional electric field sets by scaling the optimised U and Ω values by the factors 0.25, 0.5, 2, and 4.

For the last stage of preparing the data, following our previous work, we smoothed the magnetograms and rebinned. We used 4 × 4 pixel rebinning, which resulted in a pixel size of ≈1.46 Mm. The electric fields were re-inverted using the smoothed data and with the different values of Ω and U obtained during optimisation.

2.4. TMFM of AR 12473

The TMFM simulations were carried out from 23:36 UT on December 23, 2015, when the centre of AR 12473 was approximately at E50 on the disc until 15:36 UT on January 2, 2016, at which time the AR had reached ≈W70. The HMI data were prepared as mentioned in the previous sections. The simulation input magnetograms were padded with 25 pixels of zeroes at the boundary to reduce boundary effects due to the finite size of the computational domain. The simulation box height was chosen to be 200 Mm (≈0.29 R⊙). The boundary conditions of the simulation box were chosen to be open (Pomoell et al. 2019).

The simulations were run using the input magnetograms prepared for various different ad hoc U and Ω values. The optimised values of U and Ω were obtained using the procedure explained in Sect. 2.3.3. In addition, inputs were also prepared for the electric fields inverted using scaled values for the inversion parameters. Table 1 shows the optimised values for all investigated cases in bold. The errors shown for the optimised parameters in Table 1 represent the step of the increment and decrements used in the optimisation procedure. The simulations were also run using the magnetograms prepared by setting the non-inductive electric field to zero.

Table 1.

Summary of the key results for different U and Ω runs.

For each value of U and Ω, as well as for the case of zero non-inductive electric field, the simulation was run with both unmasked and masked magnetograms. In the latter case, we used the masking threshold B = 250 Mx cm−2, that is all pixels for which B ≤ 250 Mx cm−2 were assigned to 0. As an example, Fig. 1 shows the excerpts of vertical component (Bz) of raw, re-projected magnetograms cutouts from HMI during the eruption (11:36 UT on December 28, 2015). The left panel contains both unmasked and masked Bz components of magnetograms (threshold: B ≤ 250 Mx cm−2). The right panel shows the Bz component of simulation-ready cutouts (re-binned and padded) produced with ELECTRICIT at the same time for one of the simulations for both the masked and unmasked cases.

thumbnail Fig. 1.

Bz components of magnetogram cutout used in the analysis and simulation taken at 11:36 UT on December 28, 2015 (the time when the eruption was observed). These are (a) the re-projected cutout where no masking is applied, (c) the re-projected cutout where masking is applied with a threshold of B ≤ 250 Mx cm−2, (b) the simulation-ready cutouts produced with ELECTRICIT where no masking is applied, and (d) the simulation-ready cutouts produced with ELECTRICIT where masking is applied with a threshold of B ≤ 250 Mx cm−2.

3. Results

In this section, we present the simulation results for AR 12473 using different values of U and Ω used in the photospheric electric field inversion (Table 1) and investigate their effects on the emergence, evolution and eruption of the flux rope.

3.1. Photospheric energy and relative helicity injection

First, we investigated the photospheric energy injection, computed using Eq. (6), for various values of the parameters U and Ω. The results were compared with the energy injection obtained from the DAVE4VM optical flow-based electric field estimate that was also used to find the optimal values for U and Ω (see Sect. 2.3.3 and Table 1).

Figure 2 shows the energy injection for all the U and Ω values as summarised in Table 1. As noted in Sect. 2.3.1, during the optimisation process the full-resolution and unsmoothed vector magnetograms were used. The energy injection closely matches with the DAVE4VM reference curve (black) for the optimised values throughout the displayed time for all four optimised datasets. At the time of the observed eruption (i.e. 28 December, 2015 at 11:36 UT, vertical dashed line), for 0.25× the optimised value the energy injection is about ≈80 ± 2% less than the DAVE4VM value, and for 0.5× the optimised value is about 54% less. The percentages do not vary significantly between the U and Ω assumptions, nor with the masked and unmasked cases. We note that the optimisation procedure minimises the rms difference for the entire considered time-interval and not for a specific instance of time. It is also important to note that the value of the optimised U and Ω values differ for the unmasked and masked cases. For 2× of the optimised values, the energy injection is about 94% larger, and for 4× of the optimised values it is about 290% larger, that is, the increase is approximately linear with the increase in the parameter value. The overestimation for 2× and 4× optimised values increases at later times. The dot and star markers in the figures represent for each simulation the time of the formation of the flux rope and the time at which the flux rope reaches 100 Mm, respectively, and is detailed in Sect. 3.2.

thumbnail Fig. 2.

Temporal evolution of total photospheric magnetic energy injection for AR 12473 from 23:36 UT on December 23, 2015 to 08:36 UT on January 2, 2016. The black dot and black star represent the flux rope formation time and the time when the flux rope reaches ≈100 Mx for various values of U and Ω. These are (a) the unmasked U-assumption, (b) the masked U-assumption, (c) the unmasked Ω-assumption, and (d) the masked Ω-assumption. The free parameters in the non-inductive electric field values vary from 0.25× to 4× of the optimised values. The DAVE4VM reference is also shown in solid black. The vertical dashed line is the M1.9 class flare peak time on Dec. 28, 2015.

Figure 3 shows the temporal evolution of the photospheric injection of relative helicity for the same electric field inversions as those provided in Fig. 2, together with the corresponding DAVE4VM result. The total photospheric helicity injection is larger than the estimate from DAVE4VM for all Ω runs and masked U runs (Eq. (18) in Pomoell et al. 2019). Regarding the injected helicities for 0.5× and in particular for 0.25×, the optimised values are, however, very close (±20%) to those obtained from DAVE4VM (see Fig. 3c). In addition, for the unmasked U case for 0.25× and 0.5×, the optimised values were lower than the DAVE4VM estimate, and the energy injection for the optimised U also matched relatively closely with DAVE4VM. With 2× and 4×, the optimised values led to a drastic overestimation of the helicity injection when compared to DAVE4VM. We also noticed that for the optimised datasets, the helicity estimates were 2 − 3× higher in Ω runs than U runs. This effect was seen during the flux rope evolution (discussed in detail in Sect. 3.2). Similar to what is shown in Fig. 2, we did not notice much relative variation in the helicity injection between masked and unmasked Ω datasets.

thumbnail Fig. 3.

Temporal evolution of total photospheric helicity injection for AR 12473 from 23:36 UT on December 23, 2015 to 08:36 UT on January 2, 2016. The black dot and black star represent the flux rope formation time and the time when the flux rope reaches ≈100 Mm for various values of U and Ω. These include (a) the unmasked U-assumption, (b) the masked U-assumption, (c) the unmasked Ω-assumption, and (d) the masked Ω-assumption. The free parameters in the non-inductive electric field values vary from 0.25× to 4× of the optimised values. The DAVE4VM reference is also shown in solid black. The vertical dashed line is the M1.9 class flare peak time on December 28, 2015.

3.2. Formation and evolution of the flux rope

3.2.1. Flux rope identification

For the identification of the flux rope in the simulation, we used the twist number Tw, which quantifies the winding of two infinitesimally close field lines around each other. It is defined as (Berger 2003), where J∥ is the current density parallel to the magnetic field and ds is the arc length increment along the magnetic field line. To facilitate locating the flux rope using the twist number, we compute the value of Tw in a plane that is approximately perpendicular to the local axis of the flux rope. Thus, we placed the plane at y = 0, that is vertically to the XZ plane. For all simulation runs using a non-zero, non-inductive electric field component, we were able to identify a coherent region of Tw ≥ 1.5 that rose in height through the simulation domain with increasing time. We considered the volume of space filled by the field lines passing through the coherent Tw region as the flux rope (FR). With the FR defined, we tracked the evolution of the structure until the height of 200 Mm (upper edge of the simulation domain) and determined various FR parameters, such as accumulated helicity, axial flux, and axial magnetic field. For the cases where the non-inductive component of the driving photospheric electric field was set to zero (last two rows of Table 1), a flux rope as defined above did not form during the simulation time.

Snapshots of the Tw planes at the times when the apex of the identified FR reached ≈100 Mm are given in Fig. 4, while a visualisation of selected magnetic field lines through the Tw ≥ 1.5 region for each corresponding time is shown in Fig. 5. Figure 4 clearly demonstrates that all simulations, irrespective of optimisation value, generated a consistently positively twisted flux rope. The figure also shows that the twist values are higher for higher values of Ω and U. Although the overall appearance of the FRs are very similar between U and Ω runs, the twist structure and visualised magnetic field lines of the FRs appear more coherent for the U runs. In particular for some of the Ω runs, there are some field lines that exit the simulation domain through the lateral boundaries that are associated with the FR. The twist maps also reveal a significant negative region above the apex of the positively twisted FR.

thumbnail Fig. 4.

Snapshots of the twist number Tw maps (green-purple) for different TMFM runs.

thumbnail Fig. 5.

Excerpts from different simulations showing the time when the flux rope reaches the mid height in the simulation domain.

3.2.2. Magnetic energy and relative helicity in the corona

Next, we quantified the effect of the different driving photospheric electric fields on the evolution of the coronal magnetic field in the TMFM simulation. We computed the total magnetic energy (εM) and free magnetic energy (εfree) in the entire coronal simulation volume as

(7)

(8)

where Bp is the associated potential field.

Figure 6 shows that the evolution of total energy varies significantly depending on the driving electric field. For unmasked and masked 0.25×, 0.5× and 1× optimised U runs, the total energy first decreases and reaches a minimum around the time of flux rope formation (Sect. 3.2), and then it slowly increases. For 2× and 4× optimised U, the energy initially shows a sharp increase before reaching a plateau. For Ω runs, the changes in total energy for 0.25× and 0.5× optimised value runs are subtler, while the 4× optimised run does not reach the plateau during the tracked flux rope evolution (see Sect. 3.2 for details).

thumbnail Fig. 6.

Temporal evolution of total magnetic energy for AR 12473 from 23:36 UT on December 23, 2015 to the time the flux rope reached ≈100 Mm for various values of U and Ω using the processed simulation data (top panel). These are (a) the unmasked U-assumption, (b) the masked U-assumption, (c) the unmasked Ω-assumption, and (d) the masked Ω-assumption. The black dot and black star represent the flux rope formation time and the time when the flux rope reaches ≈100 Mm for various values of U and Ω. The free parameters in the non-inductive electric field values vary from 0.25× to 4× of the optimised values. Temporal evolution of the ratio of the free to the total magnetic energy (bottom panel). The vertical dashed line is the M1.9 class flare peak time on December 28, 2015.

The magnetic energy followed a very similar pattern for all of the runs with lower values of optimised Ω and U. For U 0.25× and 0.5× runs, the magnetic energy reaches a minimum during the same time frame, which indicates the flux ropes would have formed at very close times. These volume metrics show similar pattern for all the U values. However, for Ω runs, the total magnetic energy for higher values becomes very high after the flux rope reaches ≈100 Mm. The 4× runs never reached a minimum, which indicates that the flux rope was formed at the very beginning of the simulation (see Table 1 for details). The εM peaks at ≈23:36 UT on 25 Dec., for all U runs (except 4×), for both masked and unmasked values. It reaches the minima at ≈18:36 UT on 27 Dec. (except 4×). Contrary to this, for different values of Ω, εM always increased between flux rope formations until it reached the middle of the simulation box.

The ratio of the free magnetic energy to total magnetic energy (εfree/εM) for different values of U and Ω are shown in Fig. 6 (lower panel). For higher values of U, εfree/εM increased during the FR emergence until it reached the height ≈100 Mm. For lower values of U, the εfree/εM value remained the same during the time of FR emergence until it reached the height ≈100 Mm. Contrary to this, εfree/εM increased for all the values of Ω, irrespective of the mask effect.

The relative helicity (HR) and the current-carrying helicity (Hj) were calculated in the simulation as in Berger (2003):

(9)

(10)

(11)

where A and Ap are the vector potential of magnetic field B and potential magnetic field Bp, respectively, and Hpj is the mutual helicity between Hj and HR.

Figure 7 (upper panel) shows the temporal evolution of the relative helicity for various values of U and Ω. We found that HR always increased for all the runs between the FR formation time and the time it reached the middle of the simulation domain, irrespective of the optimisation parameter values. The ratio of the current-carrying helicity to the relative helicity (Hj/HR) is shown in Fig. 7 (bottom panel). It can be seen that Hj/HR shows a similar rise and then decline to εfree/εM (see Fig. 6, bottom panel). However, for all the U and Ω values, the Hj/HR ratio had increased between the time of FR formation and until it reached the height ≈100 Mm. It was previously reported by Zuccarello et al. (2018) that for their torus unstable flux ropes, the Hj/HR has a threshold of ≈0.29 ± 0.01, but that it likely is not universal. We found the Hj/HR value close to this threshold only when the optimisation parameters were close to the optimised values. Otherwise, the Hj/HR value varied from ≈10 − 40% from the threshold value found by Zuccarello et al. (2018). This means that in the TMFM runs the flux rope can rise through the simulation domain for a relatively large range of Hj/HR values. We note that Rice & Yeates (2022) found, using two-dimensional magnetofrictional simulations, that this eruptivity index has only a weak predictive ability and depends on the orientation of the overlying magnetic field with that of the flux rope.

thumbnail Fig. 7.

Temporal evolution of relative helicity for AR 12473 from 23:36 UT on December 23, 2015 to the time the flux rope reached ≈100 Mm for various values of U and Ω using the processed simulation data (top panel). These are (a) the unmasked U-assumption, (b) the masked U-assumption, (c) the unmasked Ω-assumption, and (d) the masked Ω-assumption. The black dot and black star represent the flux rope formation time and the time when the flux rope reaches ≈100 Mm for various values of U and Ω. The free parameters in the non-inductive electric field values vary from 0.25× to 4× of the optimised values. Temporal evolution of the ratio of the current-carrying to the relative helicity (bottom panel). The vertical dashed line is the M1.9 class flare peak time on December 28, 2015.

3.2.3. Flux rope parameters

The following properties were investigated for the identified FRs and they are listed in Table 1: (i) flux rope formation time, (ii) the time when the FR apex reached ≈100 Mm height, (iii) the rise time of the FR to ≈100 Mm, (iv) the accumulated helicity, (v) the total axial flux, (vi) the axial magnetic field, and (vii) the flux rope radius. In addition, Figs. 8 and 9 show the flux rope appearance and various flux rope parameters for different values of the optimisation parameters, respectively. The total accumulated relative helicity (see Eq. (6)) was calculated for the complete simulation box. The axial flux was computed as ΦA = ∫AB ⋅  dA, where A is the area of integration in a plane normal to the flux rope axis. The flux rope radius was used to determine the extent of the flux rope (see Fig. 5). The flux rope radius was determined at the apex by fitting a vertical line to its extent in the Z direction.

thumbnail Fig. 8.

Time (in hours) for the flux rope to (a) first appear in the simulation domain, (b) reach the middle height of the simulation domain, and (c) from first appearance to reach the middle of the simulation domain.

thumbnail Fig. 9.

Various flux rope properties when the flux rope apex had reached the altitude of ≈100 Mm; these are (a) the helicity accumulated from the time the flux rope formed to the observed eruption time, (b) the axial flux, (c) the axial magnetic field magnitude, (d) the flux rope radius, and (e) the peak current-carrying helicity.

Table 1 and Fig. 9 show that the FR forms earlier for Ω runs than for U runs. This is the case for all investigated values of Ω and U, although the difference decreases with the increasing value. The formation time also clearly decreases with an increasing value of the optimised parameter. For example, in the simulation runs performed using 0.25×, the optimised value of the FR forms 15−24 h later than in the corresponding optimised value run. The difference in the formation time is smaller for the Ω runs than for the U runs and for the masked runs than for the un-masked runs. For the values larger than the optimised values, the FRs form earlier than for the optimised value runs. In the runs conducted with 4×, the optimised value of the FR forms 10−39 h earlier than in the corresponding optimised value run, the largest difference being for the unmasked U run.

Approximately similar trends are also found for the times when the FR reaches ≈100 Mm in altitude and the rise time of the FR. This height is reached earlier for the Ω runs than for the U runs. For all the runs made with 0.25× the optimised value, the evolution of the FR was so slow that its apex did not reach the ≈100 Mm altitude during the simulation time. It is also noteworthy that the rise of the FR is very slow for 0.5× the optimised value, but the difference in the rise time is not that large between the runs made with the optimised value of 2× and 4× the optimised value.

For the U runs (both masked and unmasked), the total current-carrying helicity, total axial flux, axial magnetic field, and current-carrying helicity increase monotonically with the increasing value of U. For Ω runs in turn some variations are detected. This could be due to a less coherent FR structure in Ω runs, as discussed previously. The FR radius increases monotonically only for the unmasked U run. For masked runs, it is lower for 2× and 4× optimised U runs than for the optimised U run.

4. Summary and discussion

We performed data-driven time-dependent magnetofrictional modelling (TMFM) of the active region AR 12473 and studied the effects of the values of U and Ω used in determining the non-inductive photospheric electric fields on the formation and early evolution of the flux ropes and on some of their derived key parameters. First the optimised values of U and Ω parameters were found by comparing photospheric energy injection with the energy injection from DAVE4VM. The simulations were performed with [0.25, 0.5, 1, 2, 4]× the optimised values and both for the masked (threshold B = 250 G) and unmasked magnetograms.

The key findings can be summarised as follows:

  • Flux ropes formed in every simulation where the non-inductive component was not zero. They also had an overall similar appearance and evolution, with the largest difference being in the formation and evolution times. For the runs made by setting the non-inductive component to zero, the flux rope did not form at all (both for masked and unmasked magnetograms) during the simulation time. This is consistent with previous studies (Cheung & DeRosa 2012; Pomoell et al. 2019) and emphasises the importance of including the non-inductive electric field in the TMFM simulations.

  • Flux ropes formed earlier and evolved considerably faster in the Ω runs than in the U runs using the same scaling value of the free parameters with respect to their optimised values. This could be related to a significantly higher photospheric helicity injection (an order of magnitude) for Ω runs than for the U runs, presumably due to the strong twisting motion it presents. In contrast, the photospheric energy injection for the corresponding electric fields were similar.

  • Flux ropes formed earlier and rose faster in the masked simulations than in unmasked simulations. The difference was, however, smaller than between the Ω and U simulation runs. Photospheric helicity injection was found to be larger for masked cases than for unmasked cases, while no difference was found between the energy injections.

  • The flux rope formed earlier and rose faster with the increasing values of free parameters in the non-inductive electric field. It is not straightforward to match the flux rope evolution in the simulation with the real eruption time. For the studied event, the flux ropes formed at the bottom of the simulation domain considerably earlier than the real eruption time (from 1 to 4 days). The height of ≈100 Mm was also reached earlier than the real eruption time for most of the runs.

  • The current-carrying helicity increases and reaches its maximum just before the eruption or when it reaches the height ≈100 Mm, and then it decreases with time.

  • The accumulated helicity, magnetic fields, and axial flux were always higher for Ω simulation than U. It was also noted that increasing the optimisation values increases all these parameters for Ω and U runs.

  • The size of the flux ropes (in terms of radius) were also increased with a higher value of Ω and U in most of the runs.

  • The twist was found to be consistently positive for all flux ropes. The derived parameters (accumulated helicity, axial flux, and current-carrying helicity) varied between the simulations but were of the same order of magnitude.

The results of this study indicate that flux rope is formed and has overall similar evolution and properties with a large range of ad hoc free parameters needed to determine the non-inductive electric field component that is critical for energising and introducing twist to the coronal magnetic field. As discussed in the introduction, the flux rope parameters derived from the TMFM simulation can be used to constrain magnetised CME flux rope models inserted in heliospheric simulations and semi-empirical CME models. However, this study shows that, irrespective of the values, flux ropes are formed and erupted. Therefore, data-driven TMFM can be used to estimate flux rope properties early in their evolution without the need to employ a lengthy optimisation process.

Acknowledgments

The authors acknowledge the European Research Council (ERC) under the European Union’s Horizon 2020 Research and Innovation Programme Project SolMAG 724391 and the Academy of Finland Project 310445. This research used version 0.0.7 (software citation) of the SunPy open source software package. A.K.’s research was supported by an appointment to the NASA Postdoctoral Program at the NASA Goddard Space Flight Center (GSFC).

References

  1. Asvestari, E., Pomoell, J., Kilpua, E., et al. 2021, A&A, 652, A27 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  2. Berger, M. A. 2003, in Topological Quantities in Magnetohydrodynamics, eds. A. Ferriz-Mas, & M. Núñez, 345 [Google Scholar]
  3. Bobra, M. G., Sun, X., Hoeksema, J. T., et al. 2014, Sol. Phys., 289, 3549 [Google Scholar]
  4. Brueckner, G. E., Howard, R. A., Koomen, M. J., et al. 1995, Sol. Phys., 162, 357 [NASA ADS] [CrossRef] [Google Scholar]
  5. Cheung, M. C. M., & DeRosa, M. L. 2012, ApJ, 757, 147 [NASA ADS] [CrossRef] [Google Scholar]
  6. Cheung, M. C. M., De Pontieu, B., Tarbell, T. D., et al. 2015, ApJ, 801, 83 [Google Scholar]
  7. Chiu, M. C., von-Mehlem, U. I., Willey, C. E., et al. 1998, Space Sci. Rev., 86, 257 [NASA ADS] [CrossRef] [Google Scholar]
  8. Domingo, V., Fleck, B., & Poland, A. I. 1995, Sol. Phys., 162, 1 [Google Scholar]
  9. Fisher, G. H., Welsch, B. T., Abbett, W. P., & Bercik, D. J. 2010, ApJ, 715, 242 [NASA ADS] [CrossRef] [Google Scholar]
  10. Fisher, G. H., Abbett, W. P., Bercik, D. J., et al. 2015, Space Weather, 13, 369 [NASA ADS] [CrossRef] [Google Scholar]
  11. Gibb, G. P. S., Mackay, D. H., Green, L. M., & Meyer, K. A. 2014, ApJ, 782, 71 [NASA ADS] [CrossRef] [Google Scholar]
  12. Green, L. M., Török, T., Vršnak, B., Manchester, W., & Veronig, A. 2018, Space Sci. Rev., 214, 46 [Google Scholar]
  13. Hoeksema, J. T., Liu, Y., Hayashi, K., et al. 2014, Sol. Phys., 289, 3483 [Google Scholar]
  14. Jiang, C., Wu, S. T., Feng, X., & Hu, Q. 2016, Nat. Commun., 7, 11522 [NASA ADS] [CrossRef] [Google Scholar]
  15. Kazachenko, M. D., Fisher, G. H., & Welsch, B. T. 2014, ApJ, 795, 17 [NASA ADS] [CrossRef] [Google Scholar]
  16. Kazachenko, M. D., Fisher, G. H., Welsch, B. T., Liu, Y., & Sun, X. 2015, ApJ, 811, 16 [NASA ADS] [CrossRef] [Google Scholar]
  17. Kilpua, E. K. J., Balogh, A., von Steiger, R., & Liu, Y. D. 2017, Space Sci. Rev., 212, 1271 [CrossRef] [Google Scholar]
  18. Kilpua, E. K. J., Lugaz, N., Mays, M. L., & Temmer, M. 2019, Space Weather, 17, 498 [NASA ADS] [CrossRef] [Google Scholar]
  19. Kilpua, E. K. J., Pomoell, J., Price, D., Sarkar, R., & Asvestari, E. 2021, Front. Astron. Space Sci., 8, 35 [NASA ADS] [CrossRef] [Google Scholar]
  20. Lemen, J. R., Title, A. M., Akin, D. J., et al. 2012, Sol. Phys., 275, 17 [Google Scholar]
  21. Lin, H., Penn, M. J., & Tomczyk, S. 2000, ApJ, 541, L83 [NASA ADS] [CrossRef] [Google Scholar]
  22. Liu, Y., Hoeksema, J. T., Sun, X., & Hayashi, K. 2017, Sol. Phys., 292, 29 [NASA ADS] [CrossRef] [Google Scholar]
  23. Lumme, E., Pomoell, J., & Kilpua, E. K. J. 2017, Sol. Phys., 292, 191 [CrossRef] [Google Scholar]
  24. Lumme, E., Pomoell, J., Price, D. J., et al. 2022, A&A, 658, A200 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  25. Metcalf, T. R. 1994, Sol. Phys., 155, 235 [Google Scholar]
  26. Mulay, S. M., Tripathi, D., & Mason, H. 2021, MNRAS, 504, 1201 [CrossRef] [Google Scholar]
  27. Palmerio, E., Kilpua, E. K. J., James, A. W., et al. 2017, Sol. Phys., 292, 39 [Google Scholar]
  28. Pesnell, W. D., Thompson, B. J., & Chamberlin, P. C. 2012, Sol. Phys., 275, 3 [Google Scholar]
  29. Pomoell, J., & Poedts, S. 2018, J. Space Weather Space Clim., 8, A35 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  30. Pomoell, J., Lumme, E., & Kilpua, E. 2019, Sol. Phys., 294, 41 [NASA ADS] [CrossRef] [Google Scholar]
  31. Price, D. J., Pomoell, J., Lumme, E., & Kilpua, E. K. J. 2019, A&A, 628, A114 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  32. Price, D. J., Pomoell, J., & Kilpua, E. K. J. 2020, A&A, 644, A28 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  33. Rice, O. E. K., & Yeates, A. R. 2022, Front. Astron. Space Sci., 9, 849135 [NASA ADS] [CrossRef] [Google Scholar]
  34. Scherrer, P. H., Schou, J., Bush, R. I., et al. 2012, Sol. Phys., 275, 207 [Google Scholar]
  35. Schuck, P. W. 2008, ApJ, 683, 1134 [NASA ADS] [CrossRef] [Google Scholar]
  36. Scolini, C., Rodriguez, L., Mierla, M., Pomoell, J., & Poedts, S. 2019, A&A, 626, A122 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  37. Scolini, C., Chané, E., Temmer, M., et al. 2020, ApJS, 247, 21 [Google Scholar]
  38. Shiota, D., & Kataoka, R. 2016, Space Weather, 14, 56 [CrossRef] [Google Scholar]
  39. Sun, X., & Norton, A. A. 2017, Res. Notes Am. Astron. Soc., 1, 24 [Google Scholar]
  40. van Ballegooijen, A. A., Priest, E. R., & Mackay, D. H. 2000, ApJ, 539, 983 [NASA ADS] [CrossRef] [Google Scholar]
  41. Verbeke, C., Pomoell, J., & Poedts, S. 2019, A&A, 627, A111 [NASA ADS] [CrossRef] [EDP Sciences] [Google Scholar]
  42. Webb, D. F., Lepping, R. P., Burlaga, L. F., et al. 2000, J. Geophys. Res., 105, 27251 [NASA ADS] [CrossRef] [Google Scholar]
  43. Welsch, B. T. 2018, Sol. Phys., 293, 113 [Google Scholar]
  44. Wiegelmann, T., & Sakurai, T. 2012, Liv. Rev. Sol. Phys., 9, 5 [Google Scholar]
  45. Yang, W. H., Sturrock, P. A., & Antiochos, S. K. 1986, ApJ, 309, 383 [NASA ADS] [CrossRef] [Google Scholar]
  46. Yardley, S. L., Mackay, D. H., & Green, L. M. 2018, ApJ, 852, 82 [NASA ADS] [CrossRef] [Google Scholar]
  47. Yeates, A. R., & Bhowmik, P. 2022, ApJ, 935, 13 [NASA ADS] [CrossRef] [Google Scholar]
  48. Zhang, J., Richardson, I. G., Webb, D. F., et al. 2007, J. Geophys. Res., 112, A10102 [NASA ADS] [CrossRef] [Google Scholar]
  49. Zuccarello, F. P., Pariat, E., Valori, G., & Linan, L. 2018, ApJ, 863, 41 [Google Scholar]

All Tables

Table 1.

Summary of the key results for different U and Ω runs.

In the text

All Figures

thumbnail Fig. 1.

Bz components of magnetogram cutout used in the analysis and simulation taken at 11:36 UT on December 28, 2015 (the time when the eruption was observed). These are (a) the re-projected cutout where no masking is applied, (c) the re-projected cutout where masking is applied with a threshold of B ≤ 250 Mx cm−2, (b) the simulation-ready cutouts produced with ELECTRICIT where no masking is applied, and (d) the simulation-ready cutouts produced with ELECTRICIT where masking is applied with a threshold of B ≤ 250 Mx cm−2.
In the text
thumbnail Fig. 2.

Temporal evolution of total photospheric magnetic energy injection for AR 12473 from 23:36 UT on December 23, 2015 to 08:36 UT on January 2, 2016. The black dot and black star represent the flux rope formation time and the time when the flux rope reaches ≈100 Mx for various values of U and Ω. These are (a) the unmasked U-assumption, (b) the masked U-assumption, (c) the unmasked Ω-assumption, and (d) the masked Ω-assumption. The free parameters in the non-inductive electric field values vary from 0.25× to 4× of the optimised values. The DAVE4VM reference is also shown in solid black. The vertical dashed line is the M1.9 class flare peak time on Dec. 28, 2015.
In the text
thumbnail Fig. 3.

Temporal evolution of total photospheric helicity injection for AR 12473 from 23:36 UT on December 23, 2015 to 08:36 UT on January 2, 2016. The black dot and black star represent the flux rope formation time and the time when the flux rope reaches ≈100 Mm for various values of U and Ω. These include (a) the unmasked U-assumption, (b) the masked U-assumption, (c) the unmasked Ω-assumption, and (d) the masked Ω-assumption. The free parameters in the non-inductive electric field values vary from 0.25× to 4× of the optimised values. The DAVE4VM reference is also shown in solid black. The vertical dashed line is the M1.9 class flare peak time on December 28, 2015.
In the text
thumbnail Fig. 4.

Snapshots of the twist number Tw maps (green-purple) for different TMFM runs.
In the text
thumbnail Fig. 5.

Excerpts from different simulations showing the time when the flux rope reaches the mid height in the simulation domain.
In the text
thumbnail Fig. 6.

Temporal evolution of total magnetic energy for AR 12473 from 23:36 UT on December 23, 2015 to the time the flux rope reached ≈100 Mm for various values of U and Ω using the processed simulation data (top panel). These are (a) the unmasked U-assumption, (b) the masked U-assumption, (c) the unmasked Ω-assumption, and (d) the masked Ω-assumption. The black dot and black star represent the flux rope formation time and the time when the flux rope reaches ≈100 Mm for various values of U and Ω. The free parameters in the non-inductive electric field values vary from 0.25× to 4× of the optimised values. Temporal evolution of the ratio of the free to the total magnetic energy (bottom panel). The vertical dashed line is the M1.9 class flare peak time on December 28, 2015.
In the text
thumbnail Fig. 7.

Temporal evolution of relative helicity for AR 12473 from 23:36 UT on December 23, 2015 to the time the flux rope reached ≈100 Mm for various values of U and Ω using the processed simulation data (top panel). These are (a) the unmasked U-assumption, (b) the masked U-assumption, (c) the unmasked Ω-assumption, and (d) the masked Ω-assumption. The black dot and black star represent the flux rope formation time and the time when the flux rope reaches ≈100 Mm for various values of U and Ω. The free parameters in the non-inductive electric field values vary from 0.25× to 4× of the optimised values. Temporal evolution of the ratio of the current-carrying to the relative helicity (bottom panel). The vertical dashed line is the M1.9 class flare peak time on December 28, 2015.
In the text
thumbnail Fig. 8.

Time (in hours) for the flux rope to (a) first appear in the simulation domain, (b) reach the middle height of the simulation domain, and (c) from first appearance to reach the middle of the simulation domain.
In the text
thumbnail Fig. 9.

Various flux rope properties when the flux rope apex had reached the altitude of ≈100 Mm; these are (a) the helicity accumulated from the time the flux rope formed to the observed eruption time, (b) the axial flux, (c) the axial magnetic field magnitude, (d) the flux rope radius, and (e) the peak current-carrying helicity.
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.