EDP Sciences logo
Subscriber Authentication Point
Search
Menu
Home All issues Volume 626 (June 2019) A&A, 626 (2019) A91 Full HTML
Free Access
Issue
A&A
Volume 626, June 2019
Article Number A91
Number of page(s) 13
Section The Sun
DOI https://doi.org/10.1051/0004-6361/201935225
Published online 19 June 2019

© ESO 2019

1. Introduction

Coronal mass ejections (CMEs) are the manifestations of the evolution of the magnetic field on the Sun. The magnetic field of the Sun is produced by the dynamo process in the deep layers of the convection zone and rises to the surface due to magnetic buoyancy (Parker 1997; Worden & Harvey 2000). Because the density of the magnetic energy is greater than the thermal, kinetic, and gravitational energy density, the plasma behavior is dominated by the magnetic field, especially in the active regions. The measurement and reconstruction of the magnetic field structures in the solar atmosphere are very helpful to improve our understanding of physical processes in the solar atmosphere. Currently, the photosphere magnetic field can be regularly measured by many ground-based and space-based observatories, including Synoptic Optical Long-term Investigations of the Sun (SOLIS), Global Oscillation Network Group (GONG), and Solar Dynamics Observatory (SDO). However, there are still many difficulties associated with the coronal magnetic field measurements due to the fact that the plasma is optically thin and the recorded radiation is integrated over the line of sight. To acquire the magnetic field in the corona, solar physicists have developed lots of theoretical models by using the photospheric data as boundary values and defining some reasonable assumptions about the state of corona.

A magnetic flux rope is defined as a locally cylindrical bundle of field lines that have a helical structure except for the axial field line (Gopalswamy 2016). The field strength on the axis of flux rope has the maximum value and decreases towards the edges. According to the flux rope ejection models, flux rope can be formed in two ways (i) dynamically formed from the emergence below the photosphere, (ii) formed from the topological changes in the corona. In the former category, a twisted flux rope is assumed in the convection zone, which can emerge due to magnetic buoyancy. However, magnetohydrodynamic (MHD) models show that it is difficult for a flux rope to emerge as a whole from the underneath of the photosphere because the dense material is trapped in the dips of the field line (Schmieder et al. 2015). In the second category, the rotation of the two polarities and the flux cancellations lead to the formation of a flux rope because of the shearing motions. The flux cancellation is defined observationally as the mutual disappearance of magnetic fields of opposite polarities at the neutral line separating them (Linker et al. 2003). Martens & Zwaan (2001) explain the formation of a flux rope through flux convergence and cancellation.

There are many mechanisms to explain the flux rope ejections. The simulations of Wang & Sheeley (1999) showed that the eruption of the flux rope is due to the emergence and the shearing motions along the polarity inversion line. The kink instability of a flux rope can be triggered when the twist (winding of magnetic field lines around the rope axis) exceeds a threshold, for instance, one full turn (Liu et al. 2016; Wang et al. 2017). Xue et al. (2017) use the observations of New Vacuum Solar Telescope (NVST) and SDO to study the formation of the flux rope in NOAA AR 11967. They find the features of this active region are in agreement with the situations produced by the tether-cutting reconnection, such as the changes of the magnetic connection, the brightening in multiple wavelengths, and the outflows from the reconnection site. Simulation results from Syntelis et al. (2017) reveal that the flux rope became torus unstable when the strength of the magnetic envelope field above the flux rope dropped fast. Wyper et al. (2017) simulate the coronal jet driven by filament ejection. They find the energy release came from the magnetic breakout process.

Based on the observed magnetograms, there have been many attempts to simulate flux rope ejections (Inoue 2016; Green et al. 2018). Roussev et al. (2007) use the radial magnetic field provided by SOHO/MDI as the boundary conditions to investigate the dynamics of the CMEs from ARs 9906 and 10069. Zuccarello et al. (2012) apply twisting and shearing motions on the boundary to study the evolution of the flux rope. Jiang et al. (2013) used the non-linear force free (NLFF) model to extraploate the magnetic field of AR 11283, then the results were put into the MHD model to mimic the eruption of a flux rope. Amari et al. (2014) model the solar eruption from AR 10930 by using the vector magnetograms taken from the Hinode/SOT. Rodkin et al. (2017) coupled the NLFF and MHD model to produce the CME on 2 August 2011. They use the NLFF to obtain a flux rope. The evolution of the ejection of the flux rope was described in MHD simulation. Pagano et al. (2018) use the global non-potential magnetic fields deduced from the observed magnetograms as the initial condition for MHD simulation to study the dynamic evolution of erupting flux ropes.

In this paper, we present a three dimensional (3D) MHD simulation model to test the eruption mechanism of the coronal mass ejection from NOAA AR 12371. Many studies about driving a CME are based on the artificial flux rope (Borovikov et al. 2017; Fan 2017; Fan et al. 2018; Jin et al. 2017a,b; Singh et al. 2018; Török et al. 2018). Here we have applied the data-driven MHD method in a real active region to simulate the CME’s eruption process. The results are compared with the eruption models to explain the trigger mechanism of this event.

The content of this paper is as follows. In Sect. 2, we briefly describe the active region and the data used in our model. In Sect. 3, we describe the detailed setting-up of the numerical model. In Sect. 4, we describe the simulation results and compare them with the observations. Finally, a summary of conclusions is presented in Sect. 5.

2. Observations

The eruption of NOAA AR 12371 (8W, 12N) on 2015 June 22 led to a halo CME (appearing to surround the occulting disk of the observing coronagraph) with the projected speed (on the plane of the sky) of about 1000 km s−1. A two-ribbon flare (showing two bright bands of emission in Hα along the PIL) was observed which started at 17:39 universal time (UT), peaked at 18:23 UT, and ended at 18:51 UT according to the Geostationary Operational Environmental Satellite (GOES) 1.6−2.4 keV soft-X ray (SXR) flux.

We follow the evolution of the active region from 11:00 UT to 20:00 UT on 2015 June 22. During this time, the AR was almost at the center of the solar disk. The SDO spacecraft provides continuous observations and makes AR 12371 as a suitable candidate for this study. Figure 1 shows the basic magnetic configurations of this AR, which consists of two main negative polarities and two positive polarities. The leading negative polarity is set as N1 and the following one set as N2. The leading positive polarity is set as P1 and the following one set as P2.

thumbnail Fig. 1.

Magnetic fields of AR 12371 at 11:00:08 UT. The background shows the radial magnetic fields with white and black regions representing positive and negative polarities, respectively. The red line represents the polarity inversion line (PIL).

The HMI vector magnetogram provides photosphere field measurement projected and remapped onto the heliosphere coordinates with 12 min cadence and 0.5 arcsec per pixel. The vector magnetic field data are obtained using the Very Fast Inversion of the Stokes Vector (VFISV) algorithm, which assumes the Milne-Eddington model of the solar atmosphere. The 180° ambiguity in the azimuth angle of the magnetic field is resolved by the minimum energy method. The vector field data are then transformed to heliocentric spherical coordinate system.

3. The numerical model

In this study, we have used the CESE-MHD method to simulate the evolution of NOAA AR 12371. The equations are as follows

∂ ρ ∂ t + ∇ · ( ρ v ) = 0 $$ \begin{aligned}&\frac{\partial \rho }{\partial t} + \nabla \cdot (\rho {\boldsymbol{v}})=0 \end{aligned} $$(1)

∂ ( ρ v ) ∂ t + ∇ · ( ρ v v + ( p + 1 2 B 2 ) I − B B ) = − ( ∇ · B ) B + ρ g $$ \begin{aligned}&\frac{\partial (\rho {\boldsymbol{v}})}{\partial t}+\nabla \cdot (\rho {\boldsymbol{v}} {\boldsymbol{v}} + \left(p + \frac{1}{2}{\boldsymbol{B}}^2\right){\boldsymbol{I}}-{\boldsymbol{B}}{\boldsymbol{B}})=-(\nabla \cdot {\boldsymbol{B}}){\boldsymbol{B}}+ \rho {\boldsymbol{g}} \end{aligned} $$(2)

∂ p ∂ t + ∇ · ( γ p v ) = ( γ − 1 ) v · ∇ p $$ \begin{aligned}&\frac{\partial p}{\partial t} + \nabla \cdot (\gamma p {\boldsymbol{v}}) = (\gamma -1){\boldsymbol{v}}\cdot \nabla p \end{aligned} $$(3)

∂ B ∂ t + ∇ · ( v B − B v + ψ I ) = − ( ∇ · B ) v $$ \begin{aligned}&\frac{\partial {\boldsymbol{B}}}{\partial t} + \nabla \cdot ({\boldsymbol{v}}{\boldsymbol{B}}-{\boldsymbol{B}}{\boldsymbol{v}}+\psi {\boldsymbol{I}})=-(\nabla \cdot {\boldsymbol{B}}){\boldsymbol{v}} \end{aligned} $$(4)

∂ ψ ∂ t + c h 2 ∇ · B = − v · ( ∇ ψ ) − c h 2 c p 2 ψ $$ \begin{aligned}&\frac{\partial \psi }{\partial t} + c_{\rm h}^2\nabla \cdot {\boldsymbol{B}} = -{\boldsymbol{v}}\cdot (\nabla \psi )-\frac{c_{\rm h}^2}{c_{\rm p}^2}\psi \end{aligned} $$(5)

where B is the magnetic field, v is the velocity, ρ is the plasma density, p is the pressure, and ψ is the scalar potential. The last equation is used to reduce deviation from ∇ ⋅ B = 0, which is introduced by Dedner et al. (2002). The parameter ch is chosen to be the largest speed of all domain, cp set as 0.18 c h $ \sqrt{\mathrm{0.18} c_{\mathrm{h}}} $. Additionally, the Powell’s source terms −∇⋅B[(0, B, 0, v, 0)]T are added to the equations to deal with the divergence of the magnetic field. The magnetic monopoles can be advected with the velocity of the plasma (Powell et al. 1999). γ is the ratio of the specific heats. According to the previous studies (Mikić et al. 1999; Riley et al. 2006; Hu et al. 2008; Jiang et al. 2011), we set γ = 1.05 to describe the thermodynamics for corona plasma.

Quantities are normalized in the following units: Rs = 6.96 × 108 m, B0 = 20 G, ρ0 = 1.67 × 10−12 kg m−3, va = B0/ μ 0 ρ $ \sqrt{\mu_0\rho} $ = 1379 km s−1, and τ0 = Rs/va = 504 s. The simulation is done in a spherical domain (see Fig. 2). In the radial direction, r ranges from Rs to 2 Rs (Rs is the radius of the Sun). The colatitude, θ, spans from θ = 61.8° to θ = 91.8°, and the longitude, ϕ, spans from 16.2° to 61.2°. This domain is resolved by the grid of 256 × 128 × 128, with 256 grid points in r direction, 128 grid points in θ direction, and 128 grid points in ϕ direction.

thumbnail Fig. 2.

Full-disk of SDO AIA 94 Å image of the Sun observed at 2015 June 22. The bottom of computational domain is colored by red curves.

3.1. Initial condition

Following the approach of Tóth et al. (2011), we apply a potential field extrapolation to achieve the initial magnetic field, in which the synoptic map updated on June 22, 2015 from SDO/HMI1 is adopted as the bottom condition. The initial atmosphere of the domain is assumed to be in the state of hydrostatic equilibrium with density and pressure profiles given by

ρ = ρ 0 exp [ − g R s R T 0 ( 1 − R s r ) ] $$ \begin{aligned} \rho =\rho _0 \mathrm{exp} \left[-\frac{g R_{\rm s}}{R T_0}\left(1-\frac{R_{\rm s}}{r}\right)\right] \end{aligned} $$(6)

p = ρ R T 0 , $$ \begin{aligned} p=\rho R T_0 , \end{aligned} $$(7)

where R is the gas constant of the plasma, g = 274 m s−2 is the gravitational acceleration of the solar surface, the density at the base ρ0 = 1.67 × 10−12 kg m−3, Rs = 6.96 × 108 m is the radius of the Sun, and the initial temperature T0 is 2 MK. The velocity field is initially imposed to zero (Fan 2010; Pagano et al. 2015).

3.2. Boundary conditions

Because the CME expands laterally during its outward propagation, the setting of the computational domain should be a compromise between the limited computational resources available and the requirement to reduce the effect of the boundaries on the evolution of the CME. The size of vector magnetograms of AR 12371 provided by the SDO/HMI is about 12°–26°. The bottom of computational domain is 30°–45°, which is significantly larger than the size of the vector map. In addition, we aligned the computational domain with the vector map and keep their centers coincident with each other.

On the bottom boundary, we applied two methods to prescribe the conditions according to whether the region is covered with the vector magnetograms. For the area covered with the magnetograms, the three components of the magnetic fields on the bottom boundary are updated with the vector magnetograms after some treatments. The strength of magnetic field in the original observed magnetograms, as high as thousands of Gauss in the active region, is too strong to be directly employed in the MHD simulation. We first smoothed the magnetograms with Gaussian window of σ = 2 arcsec. To remove the temporal oscillations, we also smoothed the data in time with Gaussian window of σ = 4 × 12 min (Yamamoto & Kusano 2012; Jiang et al. 2016). Then, following the approach proposed by Wu et al. (2006), we reduced the strength of the magnetic field by dividing the smoothed magnetograms by 40 to get the final magnetic field maps that are input into the numerical model. This can guarantee that the plasma beta over an active region ranges from 10−3 to 10−1 at the lower corona (Gary 2001). Figure 3 shows the three components of the vector magnetogram observed by SDO/HMI and the map preprocessed by the above techniques. The three-order Lagrange interpolation is employed to achieve the magnetic fields on the computational grids at the observational instants from the preprocessed vector magnetograms due to the consideration that the resolution of the vector magnetograms is higher than that of the computational grids. We call these interpolated magnetograms as the “spatially-interpolated magnetic maps”. The magnetic fields on the bottom boundary are updated at each time step by utilizing the linear interpolation between two sequential spatially-interpolated magnetic maps. The velocity field on the bottom boundary can be updated by applying the algorithm of Differential Affine Velocity Estimator for vector magnetograms (DAVE4VM: Schuck 2005, 2006, 2008). This method is an approach to infer the magnetic footpoint velocities from a sequence of magnetograms by exploiting a variational principle to minimize deviations in the magnitude of the magnetic induction equation constrained by an affine velocity profile. Figure 4 shows the vector velocity field derived from the DAVE4VM algorithm. We have employed the compatible equations associated with the eigenvalues of vr − vfr and vr (vr is the radial component of the plasma flow velocity, vfr is the fast magnetosonic speed in the r-direction, vfr = 1/ 2 ( γ p + B · B ) / ρ + [ ( γ p + B · B ) / ρ ] 2 − 4 γ p B r 2 / ρ 2 $ \sqrt{2}\sqrt{(\gamma p\,{+}\,{\boldsymbol{B}}\cdot{\boldsymbol{B}})/\rho\,{+}\,\sqrt{[(\gamma p\,{+}\,{\boldsymbol{B}}\cdot{\boldsymbol{B}})/\rho]^2 - 4\gamma p B_{\mathrm{r}}^2/\rho^2}} $) in the projected normal characteristic (PNC) approach to determine the density and thermal pressure. The variable ψ is set to zero (Susanto et al. 2013).

For the area without vector magnetograms, we applied the PNC method to determine the boundary conditions. The PNC method, first proposed by Nakagawa (1980, 1981a,b) and Nakagawa et al. (1987), was widely employed in many numerical studies (Wang et al. 1982, 2011; Wu & Wang 1987; Wu et al. 2001, 2006; Hayashi 2005; Hayashi et al. 2006, 2018; Feng et al. 2010, 2012; Yang et al. 2012; Li et al. 2018; Li & Feng 2018). The method can keep the physical consistency and reduce the unphysical vibrations near the sub-Alfvénic solar surface boundary to an acceptable level (Yang et al. 2012).

In practice, the variable ψ is set to zero, and the other eight variables are solved by the PNC method. The approach finds the Eigenvalues of the Jacobian matrix with each of them representing one basic wave and the associated equations. For the outgoing waves that propagate from the computational domain to the outside, the associated compatible relations can be used to prescribe the boundary conditions. For typical radial speeds on the bottom boundary, there are three outgoing waves. We used the first-order upwind differencing method to calculate the spatial derivatives of the solution variables in the r-direction, and employed the central difference method in the θ and ϕ directions. For the incoming waves that propagate from outside of the computational domain to the inside, we used non-reflecting boundary conditions to fully determine all physical quantities (Wu et al. 2001; Li et al. 2018). For details of the Eigenvalues of the Jacobian matrix, their associated compatible equations, and the additional equations for the incoming waves, please refer to Hayashi (2005) and Yang et al. (2012).

For the four lateral planes and the top boundary, we also used the PNC method to update the variables at each time step. We note that the outward normal directions for these boundaries are different from each other. Therefore, the Eigenvalues and the related compatible relations should also be modified according to the corresponding Jacobian matrices and outward normal directions.

4. Results

4.1. Observation results

To examine the evolution of the magnetic flux, we divided the magnetic flux into positive and negative parts according to the sign of the radial component of the magnetograms. The magnetic flux is acquired by calculating the product of the area (A) and the radial magnetic field (Br). We used the equation F+ = ∑Br >  0BrA and F− = ∑Br <  0BrA to calculate the positive and negative flux, respectively. The positive and negative magnetic flux as a function of time are presented in Fig. 5. We find the magnetic flux of the active region undergoes two stages. In the first stage (from 11:00 to 13:00 UT), both the positive and negative flux decrease. This is probably due to the shear and convergence motions of the sunspots that result in cancellation between opposite polarities. In the second stage (from 14:00 to 20:00 UT), both the positive and negative flux increase significantly. It’s worth noting that the positive magnetic flux undergoes quick rise before the flare (about 16:00 UT). This process lasts for about an hour. We performed a linear fit to the time profiles of magnetic flux between 16:20 UT and 17:39 UT and calculate the slopes. The estimated positive and negative flux increase rates are 3.879 × 1020 Mx h−1 and 9.697 × 1019 Mx h−1, respectively. During the flare, both the positive and negative flux increase sightly. On the other hand, the increase rate of the negative flux is higher at the peak of GOES 1−8 Å soft X-ray flux.

thumbnail Fig. 3.

Three components of the vector magnetogram observed by SDO/HMI (left column) and the preprocessed maps (right column).

thumbnail Fig. 4.

Vector velocity field derived from the DAVE4VM algorithm. The background represents the radial velocity and the arrows represent the horizontal velocity.

thumbnail Fig. 5.

Positive (red) and negative (blue) magnetic flux evolution of active region 12371. The vertical black dashed line represents the flare peak time at 18:23 UT, and the gray bar represents the flare span time from 17:39 UT to 18:51 UT.

We selected three areas in the active region (see Fig. 1). The Area1 mainly includes sunspot N1, the Area2 contains sunspots P1, P2, and N2, the Area3 contains sunspots P1 and N2. Figure 6 shows the magnetic flux evolution of these three areas. The positive flux decreases slowly in the Area1 from 11:00 to 13:00 UT, while the negative magnetic flux decreases quickly during this time. We can see the positive flux decreases from 11:00 to 13:00 UT in the Area2. There is no obvious change of the positive and negative flux in the Area3 from 11:00 to 13:00 UT. We can find many small magnetic polarities around the four sunspots N1, N2, P1, and P2. The magnetic flux cancellation may occur there. From 14:00 to 17:30 UT, the increase of positive magnetic flux of the active region comes from the Area1 and Area2, while the Area2 mainly contribute to the increase of negative flux. During the flare, the magnetic flux in the Area1, Area2, and Area3 all change significantly. This is possibly due to the magnetic flux injection and cancellation.

thumbnail Fig. 6.

Positive (red) and negative (blue) magnetic flux evolution in different areas of active region 12371. The vertical black dashed line represents the flare peak time at 18:23 UT, and the gray bar represents the flare span time from 17:39 UT to 18:51 UT.

We obtain the velocity flow from the sequence of the HMI vector magnetograms on the photosphere by using the DAVE4VM algorithm (Schuck 2005, 2006, 2008). In this analysis, we used a FWHM of the resolving window of 6″.

Figure 7 displays the selected velocity field maps. From the maps, we can clearly distinguish the persistent shearing motion between the positive polarity sunspot P1 and the negative polarity sunspot N2. However, both the northern parts of the two sunspots N2 and P1 are subjected to the northeast direction shearing motions. The southern part of the N2 is subjected to the shearing motions directed to the southwest and the southern part of polarity sunspot P1 moves to the southeast directions. We can also find the converge motions towards the polarity inversion line from the center of the negative sunspot N2 (the yellow circle area). In the northeast part of the polarity, the flow arrow suggests that the sunspot rotates in the anticlockwise direction. The negative polarity N2 also shows flows towards the west. In the west, the sunspot N1 expands to the surrounding. Some fragment opposite polarities are distributed around the sunspot N1, N2, and P1. The motion of sunspots can interact with these opposite polarities.

thumbnail Fig. 7.

Horizontal velocity derived from the DAVE4VM technique overlaid on Br maps. The arrow shows transverse velocity of flux motions, and the white and black regions show positive and negative radial component of magnetic polarities, respectively. North is at the top and west to the right.

The active region shows complex interaction phenomena. These persistent converge shearing motions, flux cancelation, rotation and emergence are likely responsible for the intense flare and violent eruption.

From the observations of SDO/AIA 94 Å (see Fig. 8), two J-shaped structure can be recognized (tagged as L1 and L2, respectively) at 13:46 UT. These two J-shaped loops last for several hours until 17:32 UT. The sudden brightening arises near the neutral line. Subsequently, the area of the brightening part increases gradually. An inverse S-shaped structure can be clearly distinguished and the former two J-shaped structures disappear at 18:54 UT.

thumbnail Fig. 8.

Snapshots of SDO/AIA 94 Å maps at four different times that show the development of the inverse S-shaped flux rope in active region.

We added two RHESSI images to investigate the magnetic reconnection sites during the flare. Figure 9 shows the AIA 94 Å and HMI maps with contours of the RHESSI X-ray 6−12 keV and 12−25 keV using the CLEAN algorithm (Hurford et al. 2002) with 30 s integration time. The contour levels are at 50%, 70% and 90% of the peak intensity. Interestingly, the X-ray sources are co-spatial with sunspots N2 and P1, where the brightening in the AIA channels is detected (see Fig. 8). These suggest that the magnetic reconnection most likely occurred at the location between the sunspots of opposite polarities (N2 and P1).

thumbnail Fig. 9.

RHESSI images at 6−12 keV (blue) and 12−25 keV (yellow) energy range using the CLEAN reconstruction algorithm overlaid on AIA 94 Å and HMI images near the peak phase of the flare (integrated for 30 s). The contour levels are 50, 70 and 90% of the peak intensity.

4.2. Simulation results

Figure 10 shows the radial velocity of the plasma. We can see that the maximum velocity is located in the vicinity of the PIL and between P1 and N2. At the time of t = 6.6 h, the maximum radial velocity is about 700 km s−1. We note in particular that in the front there is a fan-shaped structure with higher velocity driven by the propagation of the coronal mass ejection. This fan-shaped structure in our simulation resembles shock waves caused by CMEs.

thumbnail Fig. 10.

Radial component of the velocity on the r − ϕ plane through the center of the PIL at t = 6.6 h and t = 7.7 h.

To study the evolution of the energy in the whole domain, we calculate the magnetic energy ( E Mag = ∫ V B 2 8 π d V $ E_{\mathrm{Mag}} = \int_{\mathrm{V}} \frac{{\boldsymbol{B}}^2}{8\pi}\ \mathrm{d}V $) and kinetic energy ( E Kin = ∫ V ρ v 2 2 d V $ E_{\mathrm{Kin}} = \int_{\mathrm{V}} \frac{\rho {\boldsymbol{v}}^2}{2}\ \mathrm{d}V $), where B is the magnetic field, V denotes the volume of the computational domain, ρ is the plasma density, and v is the plasma velocity. Figure 11 shows the magnetic and kinetic energy as a function of time. Initially, the magnetic energy in the domain increases gradually after we input the magnetic maps. About 5.5 h later, the growth curve of the magnetic energy changes much steeper as a result of the flux emergence keeping a rapid rate of growth during this period (see Fig. 5). The changes of kinetic energy show a similar trend with magnetic energy. The kinetic energy keeps in a lower level at the beginning stage and then gradually increases with time. We analyze the evolution of magnetic energy above the sunspots P1 and N2 (see Fig. 12), the region is as follows: 1 Rs <  r <  1.3 Rs, 75.4 ° < θ <  81.4°, 34.5 ° < ϕ <  40.5°. The magnetic energy begins to decrease at about 5.8 h. This may be caused by the magnetic reconnection or cancellation.

thumbnail Fig. 11.

Evolution of the magnetic and kinetic energy in the whole domain.

thumbnail Fig. 12.

Evolution of the magnetic energy above sunspots P1 and N2.

To follow the propagation of the plasma, Fig. 13 shows the slice of density distribution along r − ϕ plane through the middle of the PIL at different times. From this sequence, we can clearly see a higher dense bow front, a subsequent low density cavity, and a densest core in the center. These distributions of plasma correspond to the three parts of standard CMEs. In the early phase, the plasma gradually converges near the base of the corona. At the time of t = 5.4 h, a loop structure composed of plasma is formed. We find a low density region under the loop, which is followed by a more dense core. Besides the forward movements, the CME also expands laterally.

thumbnail Fig. 13.

Snapshots of the plasma density distribution on the r − ϕ plane through the center of the PIL. The black dot denotes the front of CME, where the distance from the Sun’s surface and the longitude have been given.

Figure 14 shows the magnetic field lines at different times. From the map we can identify two bundles of the wrested J-shaped loops at the time of t = 2.8 h, which are tagged as L1 (P2 to N2) and L2 (P1 to N1), respectively. The PIL is located between P1 and N2. In our simulation, we also discover the similar phenomenon with 17:32 UT corresponding to the time t = 6.0 h in the simulation. The primary loop L1 (connect P2 and N2) and loop L2 (connect P1 and N1) disappear and turn into a longer inverse S-shaped loop L3 (connect P2 and N1). This situation is consistent with the observation which is similar to the tether-cutting reconnection. However, from the map of Fig. 6, the flux near the PIL grows quickly from 14:00 UT, especially during the flare. The increasing rate of the magnetic flux reaches its maximum. This shows a strong link between the eruptions and the magnetic emergence. Thus, we confirm that the occurrence of the coronal mass ejection and flare are caused by the combined effects of the tether-cutting reconnection and the flux emergence.

thumbnail Fig. 14.

Two perspectives of the selected field lines (black) of the active region at different times. The bottom boundary is colored with radial magnetic field.

Figure 15 shows snapshots of the three-dimensional magnetic field evolution. The active region is covered by different magnetic arcades. Here we use different colors to distinguish different magnetic flux systems. The green magnetic field lines denote the magnetic field lines that connect P1 and N2. The two J-shaped magnetic field lines are colored pink and baby blue, respectively. The magnetic field lines that connect P2 and N1 formed in the later phase are denoted as black field lines. In the early phase (t = 2.8 h) of the dynamics, the pink, green, and baby blue magnetic systems can be seen. The pink and baby blue field lines show the J-shaped morphology. Overlying the PIL we can see the green field lines, which is also on the top of one foot of the pink field lines and baby blue field lines. Between the times of t = 2.8 h and t = 6.0 h, part of the pink and baby blue field lines are reconnected. These two sets of field lines change their connectivity. The newly formed field lines are shown in black, and connect the polarity of P2 and N1. After reconnection, the blue field lines across the PIL are raised by the rising black field lines. The pink and baby blue field lines become more perpendicular towards the bottom. At the time of t = 7.7 h, the configuration of the ejected flux ropes become significantly more complex and parts of the field lines change their connectivity. As the flux rope continues to move up, they interact with each other. The magnetic reconnection keeps happening. Then at t = 8.8 h the ejection reaches the outer boundary and the dynamics of the field lines are no longer followed.

thumbnail Fig. 15.

Snapshots of three dimensional magnetic field lines of AR 12371 at different times. Different flux systems are shown with different colors. Left column: magnetic field lines viewed from the top. Right column: magnetic field lines viewed from the side. The bottom boundary is colored with radial magnetic field. The domain is outlined by the red curves.

We calculated the decay rate of the corresponding radial magnetic field with height along the radial direction through the center of the AR by using the equation n = − h B r ∂ B r ∂ h $ n = -\frac{h}{B_{\mathrm{r}}}\frac{\partial{B_{\mathrm{r}}}}{\partial h} $ (Jiang et al. 2014). Where Br is the radial component of the magnetic field, h is the distance to the solar surface. The critical value for the torus instability is about ncr = 1.5 (Kliem & Török 2006; Démoulin & Aulanier 2010). The flux rope L3 (marked in Fig. 14) gradually reached the height (about 0.11 Rs) of the critical value for torus instability at about 5.1 h (see Fig. 16). We infer that the tether-cutting reconnection pushes the newly formed flux rope to the height where torus instability can happen.

thumbnail Fig. 16.

Height of the critical value n = 1.5 (red line). The green line shows the height of the flux rope L3.

5. Discussion

Many explosive phenomena are associated with flux cancellation. In the mean time, magnetic reconnection takes place in the corona. Wang & Shi (1993) observed flux cancellation in the photosphere, and they pointed out that the magnetic energy converted by reconnection is important for the following explosive solar activities. Zhang et al. (2001) analyzed the process of magnetic evolution of NOAA 9077 by using observational data, and they find that the magnetic flux cancellation in the vicinity of the neutral line is responsible for the filament eruption and the flare’s onset.

Emerging process which transfers the free energy from the convection zone to the surface is very important for the buildup of the magnetic energy in the corona (Forbes 2000). Magara (2006) study the activities on the photosphere and solar atmosphere by using 3D MHD simulation. Their results showed that a twisted flux rope is created by the flux emergence. Kumar et al. (2013) analyzed the multiwavelength data of the NOAA AR 11112. They find that the emerging process is essential for producing solar flares and might result in sunspots rotation.

Tether-cutting reconnection has close relationship with the flux rope ejection (DeForest et al. 2013; Chen et al. 2016; Moore et al. 2018). Liu et al. (2010) investigated a sigmoid eruption, and they found a twisted flux rope formed via tether cutting. The simulation performed by Aulanier et al. (2010) showed that the tether-cutting reconnection was important for flux rope formation, but the eruption was triggered by torus instability. The AIA and RHESSI observations (see Figs. 8 and 9) suggest the occurrence of the tether-cutting reconnection. Our analysis also reveals that the flux reached a region that can trigger torus instability.

6. Summary

In this paper, we present a 3D MHD data-driven simulation for the CME eruption from the NOAA AR 12371. Initially a current-free potential filed is made to represent the basic magnetic structure in the corona. Then, the time sequences of the vector magnetograms and the vector velocity maps are input into the bottom to simulate the evolutions of coronal magnetic fields in responds to the photospheric variations. Our aim is to investigate the physical mechanisms for the flux rope formation and the cause of the CME eruption under the real background.

The active region 12371 contains four large sunspots, one negative (N1) in the west side, and the following one also has negative polarity (N2). Besides the negative N2, there is a positive sunspot (P1). The negative sunspot P2 is in the east. The distribution of the sunspots is shown in Fig. 1. With the motion in the photosphere, such as shear and convergence, two J-shaped flux ropes were formed with one connecting the polarity N1 and P1, the other one connecting the polarity N2 and P2. The configurations of the two J-shaped flux ropes are consistent with the description of the tether-cutting reconnection proposed by Moore et al. (2012, 2018). We analyzed the magnetic flux of the active region. The results show that both the negative and positive magnetic flux increased before the flare. These indicate that the energy is injected into the active region which can provide the fuels for the eruption. The emergence of the magnetic flux contributed to the change of the magnetic configuration and induced the magnetic field to reconnect with the pre-existing ones. These could reduce the confinement of the ambient field above the flux rope and prompted the flux rope to rise to a higher place. During the phase of the flare, the magnetic flux of the four sunspots showed great variations, which indicates that the flux injection and cancellation occurred there. The horizontal velocity maps (see Fig. 7) show the complicated movements of the sunspots in this active region with the shearing and converging motions near the PIL. From the observation of the AIA 94 Å, two J-shaped flux ropes can be identified before the flare onset. During the time of the flare, we can find a brighter area above the PIL. The RHESSI X-ray contour maps (see Fig. 9) show the high energy release site above the PIL. These phenomena can be evidence of the magnetic reconnection. After the flare, two J-shaped flux ropes were reformed into an S-shaped flux rope. This scenario is similar to the tether-cutting reconnection. Our simulation results reproduce these processes as seen from Figs. 14 and 15. The results of our simulation are consistent with the AIA and RHESSI observations. These reveal that our model is successful in simulating this eruption event. We also analyzed the decay index of the magnetic field, and found that the flux rope gradually rose to the height where the decay index of the surrounding magnetic field close to 1.5. Therefore, we concluded that the tether-cutting reconnecting caused the lifting of the flux rope when the flux ascend to the height where the torus instability occurred. The flux rope finally erupted due to the multiple factors mentioned above.

1

Available from http://jsoc.stanford.edu/

Acknowledgments

This work is jointly supported by the National Natural Science Foundation of China (41474153, 41531073, 41374176, 41731067, 41274179, and 41674170) and the Youth Innovation Promotion Association CAS under grant No. 2016133. We acknowledge the use of HMI and AIA data of SDO. We thank Dr. C. Q. Xiang for his assistance in improving this paper.

References

  1. Amari, T., Canou, A., & Aly, J.-J. 2014, Nature, 514, 465 [NASA ADS] [CrossRef] [Google Scholar]
  2. Aulanier, G., Török, T., Démoulin, P., & DeLuca, E. E. 2010, ApJ, 708, 314 [Google Scholar]
  3. Borovikov, D., Sokolov, I. V., Manchester, W. B., Jin, M., & Gombosi, T. I. 2017, J. Geophys. Res.: Space Phys., 122, 7979 [NASA ADS] [CrossRef] [Google Scholar]
  4. Chen, H., Zhang, J., Li, L., & Ma, S. 2016, ApJ, 818, L27 [NASA ADS] [CrossRef] [Google Scholar]
  5. Dedner, A., Kemm, F., Kröner, D., et al. 2002, J. Comput. Phys., 175, 645 [Google Scholar]
  6. DeForest, C. E., Howard, T. A., & McComas, D. J. 2013, ApJ, 769, 43 [NASA ADS] [CrossRef] [Google Scholar]
  7. Démoulin, P., & Aulanier, G. 2010, ApJ, 718, 1388 [NASA ADS] [CrossRef] [Google Scholar]
  8. Fan, Y. 2010, ApJ, 719, 728 [NASA ADS] [CrossRef] [Google Scholar]
  9. Fan, Y. 2017, ApJ, 844, 26 [NASA ADS] [CrossRef] [Google Scholar]
  10. Fan, Y., Gibson, S., & Tomczyk, S. 2018, ApJ, 866, 57 [NASA ADS] [CrossRef] [Google Scholar]
  11. Feng, X., Yang, L., Xiang, C., et al. 2010, ApJ, 723, 300 [Google Scholar]
  12. Feng, X., Jiang, C., Xiang, C., Zhao, X., & Wu, S. T. 2012, ApJ, 758, 62 [Google Scholar]
  13. Forbes, T. G. 2000, J. Geophys. Res.: Space Phys., 105, 23153 [NASA ADS] [CrossRef] [Google Scholar]
  14. Gary, G. A. 2001, Sol. Phys., 203, 71 [NASA ADS] [CrossRef] [Google Scholar]
  15. Gopalswamy, N. 2016, Geosci. Lett., 3, 8 [NASA ADS] [CrossRef] [Google Scholar]
  16. Green, L. M., Török, T., Vršnak, B., Manchester, W., & Veronig, A. 2018, Space Sci. Rev., 214, 46 [Google Scholar]
  17. Hayashi, K. 2005, ApJS, 161, 480 [NASA ADS] [CrossRef] [Google Scholar]
  18. Hayashi, K., Benevolenskaya, E., Hoeksema, T., Liu, Y., & Zhao, X. P. 2006, ApJ, 636, L165 [NASA ADS] [CrossRef] [Google Scholar]
  19. Hayashi, K., Feng, X., Xiong, M., & Jiang, C. 2018, ApJ, 855, 11 [Google Scholar]
  20. Hu, Y., Feng, X., Wu, S., & Song, W. 2008, J. Geophys. Res.: Space Phys., 113, A3 [Google Scholar]
  21. Hurford, G., Schmahl, E., Schwartz, R., et al. 2002, Sol. Phys., 210, 61 [NASA ADS] [CrossRef] [Google Scholar]
  22. Inoue, S. 2016, Prog. Earth Planet. Sci., 3, 19 [NASA ADS] [CrossRef] [Google Scholar]
  23. Jiang, C., Feng, X., Fan, Y., & Xiang, C. 2011, ApJ, 727, 101 [NASA ADS] [CrossRef] [Google Scholar]
  24. Jiang, C., Feng, X., Wu, S. T., & Hu, Q. 2013, ApJ, 771, L30 [NASA ADS] [CrossRef] [Google Scholar]
  25. Jiang, C., Wu, S. T., Feng, X., & Hu, Q. 2014, ApJ, 780, 55 [NASA ADS] [CrossRef] [Google Scholar]
  26. Jiang, C., Wu, S. T., Yurchyshyn, V., et al. 2016, ApJ, 828, 62 [NASA ADS] [CrossRef] [Google Scholar]
  27. Jin, M., Manchester, W. B., van der Holst, B., et al. 2017a, ApJ, 834, 173 [CrossRef] [Google Scholar]
  28. Jin, M., Manchester, W. B., van der Holst, B., et al. 2017b, ApJ, 834, 172 [NASA ADS] [CrossRef] [Google Scholar]
  29. Kliem, B., & Török, T. 2006, Phys. Rev. Lett., 96, 255002 [NASA ADS] [CrossRef] [PubMed] [Google Scholar]
  30. Kumar, P., Park, S.-H., Cho, K.-S., & Bong, S.-C. 2013, Sol. Phys., 282, 503 [NASA ADS] [CrossRef] [Google Scholar]
  31. Li, H., & Feng, X. 2018, J. Geophys. Res.: Space Phys., 123, 4488 [Google Scholar]
  32. Li, C., Feng, X., Xiang, C., et al. 2018, ApJ, 867, 42 [NASA ADS] [CrossRef] [Google Scholar]
  33. Linker, J., Mikić, Z., Lionello, R., Amari, P. R., & Odstrcil, D. 2003, Phys. Plasmas, 10, A5 [Google Scholar]
  34. Liu, R., Liu, C., Wang, S., Deng, N., & Wang, H. 2010, ApJ, 725, L84 [NASA ADS] [CrossRef] [Google Scholar]
  35. Liu, R., Kliem, B., Titov, V. S., et al. 2016, ApJ, 818, 148 [NASA ADS] [CrossRef] [Google Scholar]
  36. Magara, T. 2006, ApJ, 653, 1499 [NASA ADS] [CrossRef] [Google Scholar]
  37. Martens, P. C., & Zwaan, C. 2001, ApJ, 558, 872 [Google Scholar]
  38. Mikić, Z., Linker, J. A., Schnack, D. D., Lionello, R., & Tarditi, A. 1999, Phys. Plasmas, 6, 2217 [Google Scholar]
  39. Moore, R. L., Falconer, D. A., & Sterling, A. C. 2012, ApJ, 750, 24 [NASA ADS] [CrossRef] [Google Scholar]
  40. Moore, R. L., Sterling, A. C., & Panesar, N. K. 2018, ApJ, 859, 3 [NASA ADS] [CrossRef] [Google Scholar]
  41. Nakagawa, Y. 1980, ApJ, 240, 275 [NASA ADS] [CrossRef] [Google Scholar]
  42. Nakagawa, Y. 1981a, ApJ, 247, 707 [NASA ADS] [CrossRef] [Google Scholar]
  43. Nakagawa, Y. 1981b, ApJ, 247, 719 [NASA ADS] [CrossRef] [Google Scholar]
  44. Nakagawa, Y., Hu, Y. Q., & Wu, S. T. 1987, A&A, 179, 354 [NASA ADS] [Google Scholar]
  45. Pagano, P., Mackay, D. H., & Poedts, S. 2015, J. Astrophys. Astron., 36, 123 [Google Scholar]
  46. Pagano, P., Mackay, D. H., & Yeates, A. R. 2018, J. Space Weather Space Clim., 8, A26 [Google Scholar]
  47. Parker, E. N. 1997, Sol. Phys., 176, 219 [NASA ADS] [CrossRef] [Google Scholar]
  48. Powell, K. G., Roe, P. L., Linde, T. J., Gombosi, T. I., & Zeeuw, D. L. D. 1999, J. Comput. Phys., 154, 284 [NASA ADS] [CrossRef] [Google Scholar]
  49. Riley, P., Linker, J. A., Mikić, Z., et al. 2006, ApJ, 653, 1510 [Google Scholar]
  50. Rodkin, D., Goryaev, F., Pagano, P., et al. 2017, Sol. Phys., 292, 90 [NASA ADS] [CrossRef] [Google Scholar]
  51. Roussev, I. I., Lugaz, N., & Sokolov, I. V. 2007, ApJ, 668, L87 [NASA ADS] [CrossRef] [Google Scholar]
  52. Schmieder, B., Aulanier, G., & Vršnak, B. 2015, Sol. Phys., 290, 3457 [Google Scholar]
  53. Schuck, P. W. 2005, ApJ, 632, L53 [NASA ADS] [CrossRef] [Google Scholar]
  54. Schuck, P. W. 2006, ApJ, 646, 1358 [NASA ADS] [CrossRef] [Google Scholar]
  55. Schuck, P. W. 2008, ApJ, 683, 1134 [NASA ADS] [CrossRef] [Google Scholar]
  56. Singh, T., Yalim, M. S., & Pogorelov, N. V. 2018, ApJ, 864, 18 [NASA ADS] [CrossRef] [Google Scholar]
  57. Susanto, A., Ivan, L., Sterck, H. D., & Groth, C. 2013, J. Comput. Phys., 250, 141 [Google Scholar]
  58. Syntelis, P., Archontis, V., & Tsinganos, K. 2017, ApJ, 850, 95 [NASA ADS] [CrossRef] [Google Scholar]
  59. Török, T., Downs, C., Linker, J. A., et al. 2018, ApJ, 856, 75 [NASA ADS] [CrossRef] [Google Scholar]
  60. Tóth, G., Van der Holst, B., & Huang, Z. 2011, ApJ, 732, 102 [Google Scholar]
  61. Wang, J., & Shi, Z. 1993, Sol. Phys., 143, 119 [NASA ADS] [CrossRef] [Google Scholar]
  62. Wang, Y.-M., & Sheeley, Jr., N. R. 1999, ApJ, 510, L157 [NASA ADS] [CrossRef] [Google Scholar]
  63. Wang, S., Hu, Y., & Wu, S. T. 1982, Scientia Sinica Series Mathematical Physical Technical Sciences, 25, 1305 [NASA ADS] [Google Scholar]
  64. Wang, A. H., Wu, S. T., Tandberg-Hanssen, E., & Hill, F. 2011, ApJ, 732, 19 [NASA ADS] [CrossRef] [Google Scholar]
  65. Wang, W., Liu, R., & Wang, Y. 2017, ApJ, 834, 38 [NASA ADS] [CrossRef] [Google Scholar]
  66. Worden, J., & Harvey, J. 2000, Sol. Phys., 195, 247 [Google Scholar]
  67. Wu, S., & Wang, J. 1987, Comput. Methods Appl. Mech. Eng., 64, 267 [NASA ADS] [CrossRef] [Google Scholar]
  68. Wu, S. T., Zheng, H., Wang, S., et al. 2001, J. Geophys. Res.: Space Phys., 106, 25089 [Google Scholar]
  69. Wu, S. T., Wang, A. H., Liu, Y., & Hoeksema, J. T. 2006, ApJ, 652, 800 [NASA ADS] [CrossRef] [Google Scholar]
  70. Wyper, P. F., Antiochos, S. K., & DeVore, C. R. 2017, Nature, 544, 452 [NASA ADS] [CrossRef] [Google Scholar]
  71. Xue, Z., Yan, X., Yang, L., Wang, J., & Zhao, L. 2017, ApJ, 840, L23 [NASA ADS] [CrossRef] [Google Scholar]
  72. Yamamoto, T. T., & Kusano, K. 2012, ApJ, 752, 126 [NASA ADS] [CrossRef] [Google Scholar]
  73. Yang, L. P., Feng, X. S., Xiang, C. Q., et al. 2012, J. Geophys. Res.: Space Phys., 117, A8 [Google Scholar]
  74. Zhang, J., Wang, J., Deng, Y., & Wu, D. 2001, ApJ, 548, L99 [NASA ADS] [CrossRef] [Google Scholar]
  75. Zuccarello, F. P., Meliani, Z., & Poedts, S. 2012, ApJ, 758, 117 [NASA ADS] [CrossRef] [Google Scholar]

All Figures

thumbnail Fig. 1.

Magnetic fields of AR 12371 at 11:00:08 UT. The background shows the radial magnetic fields with white and black regions representing positive and negative polarities, respectively. The red line represents the polarity inversion line (PIL).
In the text
thumbnail Fig. 2.

Full-disk of SDO AIA 94 Å image of the Sun observed at 2015 June 22. The bottom of computational domain is colored by red curves.
In the text
thumbnail Fig. 3.

Three components of the vector magnetogram observed by SDO/HMI (left column) and the preprocessed maps (right column).
In the text
thumbnail Fig. 4.

Vector velocity field derived from the DAVE4VM algorithm. The background represents the radial velocity and the arrows represent the horizontal velocity.
In the text
thumbnail Fig. 5.

Positive (red) and negative (blue) magnetic flux evolution of active region 12371. The vertical black dashed line represents the flare peak time at 18:23 UT, and the gray bar represents the flare span time from 17:39 UT to 18:51 UT.
In the text
thumbnail Fig. 6.

Positive (red) and negative (blue) magnetic flux evolution in different areas of active region 12371. The vertical black dashed line represents the flare peak time at 18:23 UT, and the gray bar represents the flare span time from 17:39 UT to 18:51 UT.
In the text
thumbnail Fig. 7.

Horizontal velocity derived from the DAVE4VM technique overlaid on Br maps. The arrow shows transverse velocity of flux motions, and the white and black regions show positive and negative radial component of magnetic polarities, respectively. North is at the top and west to the right.
In the text
thumbnail Fig. 8.

Snapshots of SDO/AIA 94 Å maps at four different times that show the development of the inverse S-shaped flux rope in active region.
In the text
thumbnail Fig. 9.

RHESSI images at 6−12 keV (blue) and 12−25 keV (yellow) energy range using the CLEAN reconstruction algorithm overlaid on AIA 94 Å and HMI images near the peak phase of the flare (integrated for 30 s). The contour levels are 50, 70 and 90% of the peak intensity.
In the text
thumbnail Fig. 10.

Radial component of the velocity on the r − ϕ plane through the center of the PIL at t = 6.6 h and t = 7.7 h.
In the text
thumbnail Fig. 11.

Evolution of the magnetic and kinetic energy in the whole domain.
In the text
thumbnail Fig. 12.

Evolution of the magnetic energy above sunspots P1 and N2.
In the text
thumbnail Fig. 13.

Snapshots of the plasma density distribution on the r − ϕ plane through the center of the PIL. The black dot denotes the front of CME, where the distance from the Sun’s surface and the longitude have been given.
In the text
thumbnail Fig. 14.

Two perspectives of the selected field lines (black) of the active region at different times. The bottom boundary is colored with radial magnetic field.
In the text
thumbnail Fig. 15.

Snapshots of three dimensional magnetic field lines of AR 12371 at different times. Different flux systems are shown with different colors. Left column: magnetic field lines viewed from the top. Right column: magnetic field lines viewed from the side. The bottom boundary is colored with radial magnetic field. The domain is outlined by the red curves.
In the text
thumbnail Fig. 16.

Height of the critical value n = 1.5 (red line). The green line shows the height of the flux rope L3.
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.