© ESO 2020

1. Introduction

Coronal mass ejections (CMEs) are some of the most energetic solar events that expel magnetic field-loaded plasma into the interplanetary space. These events release energies of the order of 10³² ergs (e.g. Forbes 2000) and they leave the solar environment with speeds between 20 km s⁻¹ and 2500 km s⁻¹ (Webb & Howard 2012). If headed towards Earth, CMEs can provoke serious geomagnetic disturbances, induce electrical currents in power grids, disrupt the global positioning system, and endanger the life of astronauts. One way of mitigating these effects is to forecast the arrival of the associated interplanetary coronal mass ejection (ICME) at our planet, by identifying its source, direction, and speed. If source signatures such as flares, filament eruptions, coronal dimmings, or extreme ultraviolet (EUV) waves are not clearly visible on the side of the Sun directed towards Earth, it is usually assumed that the event is backsided. However, there have been reported cases of eruptions lacking clear source locations, which were still Earthward directed. These situations are extremely puzzling for space weather forecasters, since geomagnetic storms caused by such CMEs can occur unpredicted and can potentially cause disruptions. Ejections without clear low coronal signatures are referred to as ‘stealth’ CMEs and were first investigated by Robbrecht et al. (2009a). The authors reported a slow CME eruption caused by a streamer blowout on 1 and 2 June 2008, which spanned only 54° in angular width and lacked the typical low-coronal signatures. Nevertheless, using the multi-viewpoint capability provided by the twin Solar Terrestrial Relations Observatory (STEREO) spacecraft (Kaiser et al. 2008) and the onboard coronagraphs COR1 and COR2 of the Sun Earth Connection Coronal and Heliospheric Investigation (SECCHI) imaging suite (Howard et al. 2008), these authors reconstructed the source of the CME and located it in a quiet-Sun front-sided region, contrary to the intuitive assumptions. The authors attributed the absence of traces to the large lift-off height and concluded that a CME originating from above 1.4 solar radii (R_⊙) should not leave any associated dimming, hence creating a stealth eruption. Other studies reporting CMEs without low coronal signatures include those of Ma et al. (2010), D’Huys et al. (2014), Kilpua et al. (2014), and Nitta & Mulligan (2017). Most of the studies on such events agree on two of their prevalent properties: a small angular width, of mostly <60°, and low initiation and propagation speeds, which are usually below 300 km s⁻¹. Another study conducted by Alzate & Morgan (2017) used advanced image processing techniques and located sources for the eruptions presented in D’Huys et al. (2014), confirming the low kinetic energies observed so far. These results emphasise the issues faced by space weather forecasters, which are the need for better imaging instruments, continuous observations from different vantage points, and a better understanding of the causes of such problematic eruptions.

Whether CMEs have a clear source region is closely related to the eruption mechanism, the height of reconnection, the overlying magnetic field, the limitation of current instrument capabilities (Howard & Harrison 2013), and other solar atmospheric conditions. It is generally accepted within the solar community that photospheric shearing motions present one of the driving mechanisms for CMEs, as they increase the magnetic helicity and consequently, the free magnetic energy. Linker & Mikic (1995), Van der Holst et al. (2006, 2007), DeVore & Antiochos (2000), and Lynch et al. (2016) are just some examples of studies using numerical simulations with imposed shearing motions to obtain CMEs and investigate their dynamical properties, ranging from the simplest 2.5D bipolar magnetic configurations, to full Sun 3D multi-polar structures.

Proposed mechanisms for stealth CME initiation include streamer blowout (e.g. Robbrecht et al. 2009a), empty filament channel eruption (e.g. Pevtsov et al. 2012), and coronal magnetic field reconfiguration, usually in the trailing current sheet of a first CME (e.g. Bemporad et al. 2012). We modelled and analysed the latter type of event by means of shearing motions imposed at the lower boundary (low corona), and also focussed on the effect of variation of the shear amplitude on the resulting eruptions.

It could be argued that this particular type of stealth events fits in the ‘plasma blobs’ category, as they originate in the current sheet created behind a first CME. These blobs were observed and analysed by several authors, including Webb & Cliver (1995), Ko et al. (2003), and Webb & Vourlidas (2016), who reported such narrow structures flowing along the current sheets that form after fast CMEs; these blobs were also numerically simulated by Riley et al. (2007). Whether these blobs can be called CMEs or not is still a topic of debate, but their eruption mechanism is very similar, and therefore, for this study the distinction is not important.

This work is similar to that of Zuccarello et al. (2012) and Bemporad et al. (2012) owing to the initial magnetic configuration and shearing profile imposed, but the main difference is that we did not require a change in magnetic field strength to obtain the second eruption. Instead, small changes applied to the shearing speed resulted in the coronal magnetic field reconfiguration, creating the second stealthy CME. Furthermore, as opposed to having a fixed grid, we use an adaptive mesh refinement (AMR) routine that increases the resolution only in regions of interest, which are areas where the electrical current is enhanced.

2. Observations

The observational basis of this work consists of a multiple coronal mass ejections event (MCMEs; Bemporad et al. 2012) that occurred on 21 and 22 September 2009. The first CME (hereafter CME1) was associated with a western limb prominence eruption and originated from an approximate latitude of 37° south, as seen by SECCHI-Extreme UV Imager (EUVI; Wuelser et al. 2004) on board the STEREO-B spacecraft. This CME was correlated with a small active region located at 38°S and 15° W, which is not recorded in the NOAA catalogue. The ejected material left the EUVI field of view on 21 September, at around 19:37 UT, as seen in Fig. 1. The second CME (hereafter CME2) entered the COR1-B (Thompson et al. 2003; Howard et al. 2008) coronagraph field of view on 22 September 2009 at ∼04:05 UT, almost eight hours after CME1. We investigated running difference images in different wavelengths and with different time steps, but we could not identify any plasma motions that could be distinguished clearly from the background noise and that would indicate the origin of the second CME. Because of the difficulty of locating the source of CME2, the lack of obvious low coronal signatures, and the fact that it was visible only when reaching higher distances from the solar surface, we consider it a stealth CME. Even though it was fainter and narrower than CME1, CME2 propagated very similarly to the first event, originating from approximately the same coronal region and presenting a strong equatorial deflection. This behaviour was explained and numerically simulated by Zuccarello et al. (2012). Their paper also contains a more detailed description of the eruptions and positions of STEREO-A and B spacecraft during the event.

Fig. 1. Image from STEREO-B EUVI (304) of the erupting prominence, taken on 21 September 2009, at 16:36 UT. The image was scaled to enhance the prominence and a mask was applied from 1.35 R⊙ outwards to remove the large noise at the edges.

Bemporad et al. (2012) calculated the velocities of both CMEs within 4 R_⊙ using height versus time maps (‘ht-plots’), and obtained values below 164 km s⁻¹. Since we were more interested in the direction of propagation of the two CMEs after the deflection, we performed 3D reconstructions using the forward modelling technique higher in altitude (above 6 R_⊙), in the COR2 field of view. The geometrical shape chosen for the CME fitting was the graduated cylindrical shell (GCS) model (Thernisien et al. 2006, 2009; Thernisien 2011) available as a SolarSoftWare (SSW) package¹, which approximates the CME structure to a bent cylinder connected to the Sun by two cones, representing the CME front and its footpoints, respectively. The shape can be changed through six free parameters, which we adjusted such that the predefined configuration matches the eruption in white-light running difference images (Fig. 2). For a more accurate computation, the reconstruction was performed when the CMEs were best seen in the observations, i.e. at two different times for CME1, and at three times for CME2. Table 1 shows the resulting longitude, latitude, and height in the Stonyhurst system of coordinates (Thompson 2006). All five measurements were made for the same day, that is 22 September 2009.

Fig. 2. Top: running difference images from COR2 used in the GCS reconstruction model. Bottom: reconstructions (green wire) using GCS superposed onto the upper corresponding panels in the case of CME1 (left) and CME2 (right).

The deprojected speed of a CME is mainly influenced by its longitude, which in this case had average values of 5.82° W for CME1, and 6.7°W for CME2. The smallest longitude value found for CME2 (5.6°, see Table 1) is probably due to errors from the reconstruction method and not a real variation. We previously measured the projected heights of the fronts of both CMEs in COR2-B plane of sky, when the angle between STEREO-B and the Sun-Earth line was −55.6°. We thus computed the deprojected distances of both CME fronts from the Sun centre. We plotted the deprojected height of the CME1 leading edge as a function of time from 6.3 R_⊙ to 16.7 R_⊙ (see Fig. 3) and fitted the data with a first order polynomial, resulting in a deprojected velocity of 257 ± 69 km s⁻¹. This speed is higher than 164 km s⁻¹ measured in the low corona within 4 R_⊙, and this acceleration above 6 R_⊙ is also seen in the automated CACTus CME list (Robbrecht et al. 2009b) based on COR2-B observations, where the velocity is listed as 277 km s⁻¹. As seen in the ht-plot, CME2 was faster and propagated in the same field of view with a deprojected velocity of 349 ± 70 km s⁻¹, mainly because it was expelled into a depleted solar wind following the passage of CME1. Given these slow speeds, a good assumption would be that the CMEs eventually reached the slow solar wind speed (as seen in Fig. 4) due to the drag force exerted on them. The error bars in Fig. 3 were calculated by selecting the CME front several times and extracting the largest variation between measurements. We computed the standard deviation in speed using the derivsig.pro function available in IDL, based on the height and time measurements and their error bars. The function computes the speed uncertainty values from the errors in height using an error propagation formula.

Fig. 3. Height-time plots of the two studied CMEs, erupting between 21 and 22 September 2009, based on STEREO-B coronagraph observations.

Fig. 4. Solar wind speed recorded at Earth, taken from the OMNI database. The two red lines indicate the expected arrival time interval of the two CMEs.

In order to obtain the speed of the wind into which the CMEs were ejected, we used the in-situ measured solar wind speed as extracted from NASA/GSFC’s OMNI data set through OMNIWeb² between 3 and 8 days after the eruption. This parameter is needed as input for the simulations described in Sect. 3. In this time interval, we carried out a back-calculation for each OMNI speed measurement for the date and time at which the embedded CMEs should have originated from the Sun, in the assumption they had approximately the same velocity as the solar wind. We then compared this with the real eruption time and concluded the arrival time of the CMEs should be 27 and 28 September, as indicated by the red vertical lines in Fig. 4. In this period, the solar wind speed did not present any major jumps, as expected owing to the low speeds of the CMEs, and had an average value of ≈330 km s⁻¹. The small jump in speed on 26 September at 13:00 is not correlated to any other in situ ICME signature and this jump also arrived too early, therefore it was not considered to be the beginning of the ICME. We performed a more detailed study of other solar wind parameters as well, such as magnetic field components, to determine the arrival of discussed CMEs at Earth, and this topic shall be the focus of a follow-up paper.

3. Simulations set-up and method

Our goal was to investigate the effect of variations in the amplitude of shearing motions on the dynamics of CMEs, starting from the observed event discussed in Sect. 2. For this purpose, we performed 2.5D spherical axisymmetric numerical simulations using the Message Passing Interface – Adaptive Mesh Refinement – Versatile Advection Code (MPI-AMRVAC; Keppens et al. 2012; Porth et al. 2014; Xia et al. 2018), which resolves the system of partial differential equations for magnetohydrodynamics (MHD) on a 2D logarithmically stretched grid. The mesh has an initial resolution of 480 × 240 cells and the domain spans between the north and south poles of the Sun, and the Earth, that is (r, θ) ∈ [1, 215] R_⊙ × [0°, 180°], where r is the radial distance from the centre of the Sun and θ is the heliographic colatitude. Even though the grid extends until 1 AU, in this study we only focussed on the near-Sun region, below 20 R_⊙. The propagation of the simulated CMEs to Earth will be discussed in a follow-up paper. In order to achieve much better detail of the dynamics of the eruption and of the structures that form, without increasing the resolution of the entire domain, an AMR protocol was used. The grid was refined to a maximum of two additional levels and only in the regions of interest, where the electrical current presented an enhancement, that is at potential magnetic reconnection sites. This increase is quantified through a parameter taken from Karpen et al. (2012) and also used by Hosteaux et al. (2018) as follows:

$$\begin{array}{c}\hfill c\equiv \frac{|\int {\int }_S\mathrm{\nabla }\times \mathit{B}·\mathrm{d}\mathit{a}|}{{\oint }_C|\mathit{B}·\mathrm{d}\mathit{l}|}=\frac{|{\oint }_C\mathit{B}·\mathrm{d}\mathit{l}|}{{\oint }_C|\mathit{B}·\mathrm{d}\mathit{l}|}=\frac{|{\sum }_{n=1}^4{B}_{t,n}{l}_n|}{{\sum }_{n=1}^4|B_{t,n}{l}_n|},\end{array}$$(1)

where B_t, n represents the tangential component of the magnetic field along the segment l_n. The line integral at the numerator is obtained by applying the Stokes’ theorem on the surface integral in the left hand-side term, and in this case the curve (contour) and the surface are those of a grid cell. In their discrete form, the line integrals are rewritten as the sum of the product between the tangential magnetic field and the length of the edges, along all four sides of a cell. Therefore, c is simply the ratio of the magnitude of the electrical current passing through the surface S spanning loop C to the sum of the absolute value of all of its components. In this definition, the parameter c varies from 0 to 1, depending on the non-potentiality of the magnetic field. The blocks of the grid are refined if and where c crosses a certain threshold, which we chose as c >  0.02. On the other hand, if c <  0.01, the grid is coarsened to a lower level, since there are no strong current-carrying structures in that region. In between these two values, the AMR routine does not impose any constraints, so the grid retains its resolution from previous time steps. The refinement is also fixed to a maximum number of two levels (on top of the base grid) for distances below 2.5 R_⊙ and for θ ∈ [27°, 153°], such that there is no change in diffusivity close to the boundary, which would have affected the eruption dynamics and introduced artificial phenomena. Because of the stretching of the grid, the scale of cells is kept constant, and the size ratio of the widths (or heights) of the farthest to closest cells to the inner boundary is ≈215 for the same mesh refinement level.

The MHD equations are temporally discretised using a two-step predictor corrector, which is suitable for a finite-volume discretisation method called the total variation diminishing Lax-Friedrichs (TVDLF) scheme. We used the most diffuse and stable type of slope limiter, minmod, and a CFL number of 0.3. The magnetic field solution is kept divergence free using the generalised Lagrange multiplier (GLM; Dedner et al. 2002) method, which introduces a new variable into the system, that damps and transports the unphysical magnetic monopoles outside the computational domain.

The initial conditions (at the start of shear) consist of a realistic bimodal solar wind created by introducing extra source terms to the momentum and energy equations to account for gravity and heating mechanisms. This solar wind model is similar to that used in previous studies by Jacobs et al. (2005), Chané et al. (2006, 2008) and Hosteaux et al. (2019). The volumetric heating function was initially used by Groth et al. (2000) and Manchester et al. (2004), and takes the following empirical form:

$$\begin{array}{c}\hfill Q=\rho q_0(T_0-T)exp[-\frac{{(r-1R_{\odot })}^2}{{\sigma }_0^2}],\end{array}$$(2)

where ρ is the mass density, the amplitude of the volumetric heating q₀ = 10⁶ ergs g⁻¹ s⁻¹ K⁻¹, r (R_⊙) is the distance from the centre of the Sun, and T (K) is the temperature. The values T₀ (K) and σ₀ (R_⊙) represent the target temperature and the heating scale height, respectively, and these are both dependent on the value of a critical angle, θ₀ (measured from the North pole), as follows: from θ₀ equator-ward, T₀ = 1.32 × 10⁶ K and σ₀ = 4.5 R_⊙, and from θ₀ pole-ward, T₀ = 2.3 × 10⁶ K and σ₀ = 4.5[2 − sin²(θ)/sin²(θ₀)] R_⊙. To better reproduce the slow/fast solar wind features, θ₀ is also defined with a distance dependence as follows:

• sin²(θ₀) = sin²(17.5°) + cos²(17.5°)(r − 1 R_⊙)/8 R_⊙, for r <  7 R_⊙,

• sin²(θ₀) = sin²(61.5°) + cos²(61.5°)(r − 7 R_⊙)/40 R_⊙, for 7 R_⊙ ≤ r <  47 R_⊙,

• sin²(θ₀) = 1, for 47 R_⊙ ≤ r.

The temperature and mass density are fixed at the inner boundary to 1.32 × 10⁶ K and 1.66 × 10⁻¹⁶ g cm⁻³, respectively. The temperature value at the boundary is chosen based on the impact it has on the solar wind speed through the momentum equation, provided that the observed speed at Earth was ≈330 km s⁻¹ (as discussed in Sect. 2). Therefore, the temperature is adjusted to provide a realistic slow solar wind speed at Earth of ≈330 km s⁻¹, and a fast solar wind speed at the poles of ≈735 km s⁻¹. The radial component of the momentum is extrapolated in the ghost cells, while the latitudinal component (v_θ) is set to 0 at the boundary. The azimuthal component (v_ϕ) is set such that it resembles the differential rotation of the Sun. In order to introduce a dipole field, the term r²B_r is fixed at the inner boundary, while r⁵B_θ and B_ϕ are extrapolated from the first inner cell. At the outer supersonic boundary, the variables r²ρ, r²ρv_r, ρv_θ, rv_ϕ, r²B_r, B_θ, rB_ϕ, and T are continuous as well.

A commonly encountered magnetic field structure on the Sun depicts a triple arcade system embedded in a helmet streamer (see Fig. 5), which was also simulated by Bemporad et al. (2012) and Zuccarello et al. (2012) to study the deflection and dynamics of the chosen multiple event, and by Karpen et al. (2012) to analyse a breakout event. We applied the same configuration in our model, by taking the curl of the following vector potential:

$$\begin{array}{c}\hfill A_{\varphi }=\frac{A_0}{r^{4}sin\theta }{cos}^2\left[\frac{180°(\lambda +11.5°)}{2\mathrm{\Delta }a}\right],\end{array}$$(3)

through which we obtained the magnetic field components

$$\begin{array}{cc}& B_r=\frac{A_0}{r^{5}sin\theta }\frac{180°}{\mathrm{\Delta }a}cos\left[\frac{180°(\lambda +11.5°)}{2\mathrm{\Delta }a}\right]sin\left[\frac{180°(\lambda +11.5°)}{2\mathrm{\Delta }a}\right],\hfill \end{array}$$(4)

$$\begin{array}{cc}& B_{\theta }=\frac{3A_0}{r^{5}sin\theta }{cos}^2\left[\frac{180°(\lambda +11.5°)}{2\mathrm{\Delta }a}\right],\hfill \end{array}$$(5)

where Δa = 37.2° represents half the width of the arcade system, λ = 90° −θ the solar latitude, $A_0=0.73\mathrm{G}·R_{\odot }^5$, and the whole arcade system is shifted by 11.5° to the south. These components were only added to the dipole magnetic field in the region determined by the latitude λ ∈ [− 48.7°,25.8°]. This configuration provides a polar magnetic field strength of 1.8 G (or 1.8 × 10⁵ nT), and a maximum arcade strength of 1.57 G (or 1.57 × 10⁵ nT), measured at the first cell of the domain.

Fig. 5. Initial magnetic field configuration – Br (colour scale) and selected magnetic field lines in the meridional plane.

Once the solar wind has reached a steady-state solution after ≈200 h, CMEs are obtained by applying time-dependent shearing motions at the inner boundary, thus magnetically stressing the corona above the selected region. The shear was introduced through the following additional azimuthal component of the speed:

$$\begin{array}{c}\hfill v_{\varphi }=v_0{({\alpha }^2-\mathrm{\Delta }{\theta }^2)}^{2}sin\alpha sin[180°(t-t_0)/\mathrm{\Delta }t],\end{array}$$(6)

where α = λ − λ₀. This flow is almost symmetric with respect to the latitude of the southernmost polarity inversion line, which is approx. λ₀ = −41°, and spans over 2Δθ = 17.2°. It is applied for Δt = 16 h starting from t₀ = 0 h and it has a slow increase and decrease, with a maximum at half the time interval. Throughout the simulations, the scaling factor v₀ is given such that v_ϕ does not exceed 10% of the local Alfvén speed. A representation of the shearing profile can be seen in Fig. 6.

Fig. 6. Normalised shearing profile (vϕ/vmax) as a function of latitude, at different simulation times from the start of shear.

4. Results

In order to make the system erupt, we applied shearing motions to the southernmost polarity inversion line and modified the amplitude of the shear by adjusting the scaling factor v₀. Thus, a first erupting flux rope (CME1) is created approximately 8 h after the start of the shear (see Figs. 7a–c), and its evolution and propagation are mostly independent of the small variations we impose on v₀, as seen in Figs. 8c,e,g. On the other hand, for the second and third flux ropes, changes within 1% of the initial lowest shearing speed (measured at the first cell of the domain) result in three different dynamical scenarios:

Fig. 7. Simulated relative density (grey scale) and selected magnetic field lines during the formation phase of flux ropes (FRs), in the case of: first row – CME1; second row – CME2 (double eruption case); third row – second FR (FR2) which falls back to the Sun; fourth row – stealth eruption; fifth row – failed eruption. The green circles indicate reconnection sites. The relative density is ${\rho }_{\mathrm{rel}}=\frac{\rho (t)-{\rho }_0}{{\rho }_0}$, where ρ0 is the density of the initial relaxed state before the shear.

Fig. 8. First column: snapshots of observations and MPI-AMRVAC simulations taken when CME1 front is located at 10 R⊙. Second column: snapshots of observations and simulations taken when CME2 front is located at 6 R⊙ (except (h)). Panels a and b: running difference images taken from EUVI, COR1, and COR2 instruments on board STEREO-B and visualised in JHelioviewer software. Snapshots of the simulated relative density and selected magnetic field lines: panels c and d indicate the double eruption scenario, panels e and f show the stealth eruption scenario, and panels g and h indicate the failed eruption scenario. The parameter Δt represents the time interval from the eruption of CME1, chosen as 21 September 2009 21:00 UT, or 9.6 h simulation time after the start of shearing motions, when the centre of the flux rope crosses 1.5 R⊙ and enters COR1-B field of view.

1. The double eruption case where $\left|v_{\varphi }^{\max }\right|=37.4\text{km}{\text{s}}^{-1}$. After the first CME, a second flux rope also erupts (CME2) from the shear applied at the boundary (see Figs. 7d–f and 8c,d).

2. The stealth eruption case, where $\left|v_{\varphi }^{\max }\right|=37\text{km}{\text{s}}^{-1}$. After the first CME, the second flux rope (FR2) created from the shearing motions falls back to the Sun (Figs. 7g–i); a third flux rope emerges from the reconfiguration of the coronal magnetic field (Figs. 7j–l) and erupts in the trail of the preceding eruption (see Figs. 8e,f).

3. The failed eruption case, where $\left|v_{\varphi }^{\max }\right|=\text{37}.2\text{km}{\text{s}}^{-1}$. After the first CME, the second flux rope (FR2) created from the shearing motions falls back to the Sun (Figs. 7g–i); a third flux rope emerges from the reconfiguration of the coronal magnetic field (Figs. 7m–o) and reconnects with the northern closed arcade, thus not erupting (Figs. 8g,h).

Values of $\left|v_{\varphi }^{\max }\right|$ below the smallest or above the largest values of those mentioned result in only one eruption, or multiple eruptions, respectively, with the stealth occurring solely in the threshold specified. The amplitude of the applied shear is larger than the observed photospheric motions, which have typical speeds of around 5 km s⁻¹ (Manchester 2007; Malherbe et al. 1983) and usually lead to flux rope formation, system destabilisation, and eventually eruption. Nevertheless, our simulated inner boundary describes the low coronal environment, and the observed magnitude of the shearing speed increases with height in the solar atmosphere, as shown by Athay et al. (1982, 1985) (transition region: 20 km s⁻¹) and by Chae et al. (2000) (low corona: 20–50 km s⁻¹). Therefore, the amplitude of our imposed shear flow is in agreement with the computed velocities from coronal observations. Initially, Zuccarello et al. (2012) and Bemporad et al. (2012) applied the shear for 36 h in their simulations; this duration was calculated based on observed photospheric motions. In order to save computation time and resources, in this work we decrease the time interval during which the additional v_ϕ is applied from 36 h to 16 h, the latter being the value used in the discussed cases. This leads to an increase in the shear magnitude, but with almost no effect on the simulated eruptions.

In all three cases listed above, the formation of the first flux rope and its eruption are essentially the same. Firstly, the applied azimuthal flow increases the magnetic pressure inside the southern arcade, making it rise and expand. Secondly, the imbalance between the magnetic pressure gradient and tension compresses the magnetic field locally. This triggers reconnection between the sides of the arcade (green circle, Fig. 7b) and creates the flux rope, which is deflected towards the equator owing to the southern polar magnetic pressure (Figs. 7a–c), and propagates radially inside the current sheet of the northern helmet streamer approximately after 5 R_⊙ (Fig. 8). For a more detailed view of the evolution, supporting videos for all three simulations and of the observed eruptions are available online. In the simulation movies, extra flux ropes can be seen in the wake of the CMEs, and these flux ropes are somewhat comparable to those appearing in the coronagraph observations, but the occurrence of these is highly influenced by the initial coronal magnetic field configuration and the background solar wind. The magnetic configuration simulated is surely more simplified than the real configuration, and the solar wind is only adjusted such that the speed matches that measured at 1AU; therefore we do not compare these small eruptions any further since they would require simulation from a data-driven MHD model. The time in the simulation videos is not normalised to a real unit system and the conversion equation is t_real = t_sim × 1.9325 h, where t_sim is the time shown in movies. Also, the distances on axes are expressed in solar radii. The movie of the observations was visualised in the JHelioviewer software, using running difference images from the EUVI, COR1, and COR2 instruments on board STEREO-B.

The fastest shear we applied (37.4 km s⁻¹) results in the formation of another flux rope in the wake of the first CME, much in the same manner as the previous flux rope. The arcade expands again, and the sides pinch (green circle, Fig. 7d) and create the flux rope, which then erupts into the equatorial current sheet to produce CME2 (Figs. 7d–f). The evolution and propagation can be visualised in Figs. 8c,d. This second flux rope also forms in the next two scenarios, but it fails to erupt and falls back to the Sun, dissipating at the boundary (Figs. 7g–i).

The lowest speed we used (37 km s⁻¹) produced CME1 in the same way, but the magnetic flux introduced was insufficient to inflate the arcade long enough to make it detach once more. Instead, the unstable current sheet formed in the trail of CME1 reconnects (green circle, Fig. 7j) and creates a third flux rope after ≈17.5 h of shearing, which is then slowly deflected towards equator (Figs. 7j–l). Given the high altitude at which the reconnection sets in (1.4 R_⊙) and the low speed of the process, this second eruption most likely does not leave any clear signatures, placing it in the stealth CME category (Robbrecht et al. 2009a). Compared to CME2 from the previous case, this one is narrower, slower, and occurs later, as seen in Figs. 8e,f.

Applying an intermediate shearing speed (37.2 km s⁻¹) again gives rise to CME1, to a second flux rope that falls back to the Sun, and to a third flux rope emerging through the same mechanism as the stealth, from coronal magnetic field reconnection (green circle, Fig. 7m) ≈19 h after the start of the shearing. Its evolution, however, is different since the small amount of additional flux translates into a slightly more energetic reconfiguration, which increases the speed of the flux rope. As a consequence of its higher acceleration, slightly higher speed, and lower formation height (1.55 R_⊙, taken at the centre of the flux rope, as opposed to 1.7 R_⊙ for the stealth), the CME cannot be deflected as easily towards the equator by the magnetic pressure gradient, and eventually reconnects with the large northern arcade, creating a failed (confined) eruption (Figs. 7m–o). We analysed and compared the dynamics of the three scenarios in more detail, and these results are presented in the next section.

5. Analysis of the simulations

In order to compare the dynamics of the simulated and observed CMEs, we plotted the deprojected heights of the eruption fronts and the locations of simulated leading edges in Fig. 9. The real distances from Sun to the observed CMEs were required in this step due to the 2.5D topology of the simulated eruptions, which would exclude any longitudinal deflection. The CMEs on 21 and 22 September 2009 were tracked only in the COR2-B field of view, from which they began to propagate radially. Below 6 R_⊙, the strong deflection and diffuse fronts of the CMEs prevented us from extracting accurate information.

Fig. 9. Location of the CMEs fronts in the meridional plane as a function of time. The black dots indicate the observed CME1, black squares indicate the observed CME2, blue triangles represent the simulated CME1, red crosses indicate the simulated CME2 (double eruption scenario), and red diagonal crosses represent the simulated stealthy CME2.

In the simulation data shown in Fig. 9 there is only one representation of the first CME because it evolved similarly for all three cases, as described in Sect. 4. We note that the first simulated eruption was slightly faster than the observed eruption, with a speed of 279 km s⁻¹ corresponding to the last measured height in COR2-B, where the speed was 257 km s⁻¹. Given the fact that the background magnetic configuration and solar wind were not perfectly simulated to match remote sensing and in situ observations, this value is a good approximation. The second CME of the double eruption case correlates better with its observed counterpart, which has a simulated speed of 342 km s⁻¹ at the last data point at a height of 18 R_⊙, as compared to the observed 349 km s⁻¹. The stealth CME front was starting to merge with the trail of the first eruption at that distance, and therefore we were unable to detect it, but the speed of the centre of the flux rope was 339 km s⁻¹. All three values closely resemble the observed values, therefore managing to realistically reproduce two slow CMEs, with two different erupting scenarios.

Next, we focussed more on the simulations as an independent result, rather than comparing them to the observations, and we analysed the differences between the three cases. As shown in Fig. 10, all three secondary eruptions present a stronger acceleration compared to CME1 when close to the Sun, owing to the depletion of background wind, and to the already open field lines created from the passage of the first CME. We note that CME2 from the double eruption case presents a slow initial rising phase. All three of these secondary eruptions emerge in the current sheet created by the first eruption, whereas CME1 needs to overcome the magnetic tension of the overlying field, hence the smaller acceleration and slow rise. The higher initial acceleration makes it harder for these narrower flux ropes to be deflected towards the equator, so they reach higher latitudes than CME1, as shown in Fig. 11. In the case of the failed eruption, the slightly larger momentum of the flux rope increases the latitudinal deflection, making it unable to escape the reconnection with the northern arcade. The timeline of its evolution can be seen in Figs. 7m–o, and is described as follows: 20.5 h after the initiation of the shearing motions, the flux rope enters the northern streamer; 20 min later, it gets trapped by the closed arcades, and after 23.5 h, it reconnects entirely and disappears. In the double eruption scenario, the second flux rope begins to form around 13.5 h after the start of the shearing motions, and is accelerated away from the Sun 3 h later. In the case of the stealth CME, the flux rope becomes embedded into the streamer after 19 h (Fig. 7k), a step which induces deceleration for all the scenarios, at their corresponding times.

Fig. 10. Total speed of the simulated CMEs, calculated at the centre of the flux rope. The blue triangles indicate CME1, the red crosses indicate CME2 (double eruption scenario), the red dots represent CME2 (failed eruption), and the red diagonal crosses indicate stealthy CME2.

Fig. 11. Latitude of the simulated CMEs with respect to the Sun’s equator as a function of time, taken at the centre of the flux rope. The legend is the same as in Fig. 10.

Upon reaching the highest latitude, the polar magnetic pressure together with the northern streamer force all CMEs to deflect equatorward, thus also contributing to the deceleration of the CMEs. Once the propagation becomes radial and the CMEs enter the solar wind, they accelerate again and reach speeds of up to 350 km s⁻¹. Even though the CME formation and cause of eruption are different, it is interesting to point out the similarities in the overall speed profile of CME2 for the stealth and double eruption scenarios.

We also performed Poynting flux analysis by calculating the radial component of its associated vector close to the inner boundary through a sphere located at 1.14 R_⊙ as follows:

$$\begin{array}{c}\hfill S_r=\frac{{[-(\mathit{v}\times \mathit{B})\times \mathit{B}]}_r}{{\mu }_0}=\frac{B^2{v}_r-(\mathit{v}·\mathit{B})B_r}{{\mu }_0},\end{array}$$(7)

and multiplying every value with the corresponding area of the surface obtained by rotating the bottom of the grid cell around the Sun at 1.14 R_⊙. The final values consist in the summation over the entire surface of the shell around the Sun, and their variation in time is shown in Fig. 12, where time = 0 indicates the start of the shear. As expected, the induced shearing motions increase the electromagnetic energy flow through the inner boundary, as seen in the first 6 h of simulation, during which the southern arcade starts expanding. Afterwards, the lateral magnetic pressure gradient compresses the field locally to create the first flux rope, and after ≈8 h this flux rope detaches from the Sun as it is accelerated a first time, which corresponds to the ‘shoulder’ in the radial Poynting flux, indicated by the arrow in Fig. 12. Following the eruption of CME1, a second flux rope begins to form as a consequence of ongoing shearing motions and lateral pinching commences after the peak at 12 h. This evolution is the same for all three simulation cases within the first 14 h. In the next part, the dynamics start to differentiate because this recently created flux rope erupts and forms CME2 in the double eruption case. For the other two scenarios, the remaining shear does not impose sufficient energy onto the flux rope to overcome the magnetic tension of the overlying closed magnetic field, so the flux rope falls back to the Sun at ≈15 h and disperses at the inner boundary, as described in Sect. 4.

Fig. 12. Evolution in time of the radial component of the Poynting flux, summed over the whole Sun.

A particularly interesting aspect presented in Fig. 12 is that the Poynting flux profile is practically the same until near the end of the shearing motions (which were applied for 16 h), yet the coronal environment behaves completely different from that point onwards. Furthermore, the radial electromagnetic energy flux is very similar in the cases of stealth and failed eruptions, emphasising the argument stated before that these flux ropes are not caused by the shear, instead they result from the coronal magnetic field reconfiguration. This can also be seen in Fig. 10, in the difference between initiation phases, which presents a slow rise for CME2 in the double eruption case, whereas for the other two, which form higher in the corona, the moment they are created also coincides to the initial acceleration stage.

6. Summary

In this paper, the initial observational basis consisted in a multiple CME event on 21 and 22 September 2009, with an interval between the eruptions of approximately eight hours. In order to obtain the real propagation directions, we performed GCS analysis high in the corona on both CMEs, and derived an average longitude of 5.82°W for CME1 and of 6.7°W for CME2, respectively.

Subsequently, we used the MHD package of the code MPI-AMRVAC to numerically simulate consecutive CMEs by applying shearing motions at the inner boundary, starting from an initial triple arcade structure embedded into a bimodal solar wind. We changed the amplitude of the shear to understand the effects of its variation on the induced eruptions with the goal of obtaining a stealth eruption in the trailing current sheet of a preceding CME. These new results indicate that within just 1% variation of the lowest shearing speed (37 km s⁻¹), three different eruption scenarios for the second CME were created: a stealth, a failed, and a double eruption. We emphasise in particular the initiation, for the first time, of a stealth from the change in applied shear, and not only from the overlying magnetic field constraints, as previously studied. All three cases are comprised of a first eruption generated from the imposed stress at the inner boundary, and a second eruption occurring either from the reconfiguration of coronal magnetic field or from the additional v_ϕ component. We compared the height-time evolution of these eruptions with the observations described above, and achieved a good slope correlation for CME1 and CME2 from the double eruption scenario. The simulated deflection of both ejections in all cases was also studied by tracking the latitude and total speed of the flux rope centres. We would like to emphasise the physical interpretation of these results, in the sense that the solar corona can react and produce eruptions differently, given almost the same initial triggering factors, making the predictions of such eruptive processes even more difficult.

We also investigated the evolution of the radial Poynting flux close to the inner boundary, to acquire a global perspective of the electromagnetic energy introduced into the domain. The steps in dynamical evolution of the simulated CMEs could also be correlated with certain features, such as maximum, minimum or fluctuations, in the Poynting flux profile. This highlights the surprising similarities between the stealth and failed eruptions and the sensitive coronal response to stress factors. In a follow-up paper we will discuss the propagation of these simulations in the interplanetary space and compare their signatures to in situ measurements.

Movies

Movie 1 of Fig. 8 Access here

Movie 2 of Fig. 8 (stealth_eruption) Access here

Movie 3 of Fig. 8 (failed_eruption) Access here

Movie 4 of Fig. 8 (double_eruption) Access here

Acknowledgments

We thank the referee for their constructive comments and suggestions, which led to the improvement of this article. D. C. T. was funded by the Ph.D. fellowship of the Research Foundation - Flanders (FWO), contract number T1 1118918N. D. C. T., E. D. and M. M. acknowledge support from the Belgian Federal Science Policy Office (BELSPO) in the framework of the ESA-PRODEX program, Grant No. 4000120800. This work was partially supported by a PhD Grant awarded by the Royal Observatory of Belgium. E. C. was funded by the Research Foundation – Flanders (grant FWO 12M0119N). These results were also obtained in the framework of the projects C14/19/089 (C1 project Internal Funds KU Leuven), G.0A23.16N (FWO-Vlaanderen), and C 90347 (ESA Prodex). For the computations we used the infrastructure of the VSC – Flemish Supercomputer Center, funded by the Hercules foundation and the Flemish Government - department EWI. We acknowledge use of NASA/GSFC’s Space Physics Data Facility’s OMNIWeb (or CDAWeb or ftp) service, and OMNI data. We acknowledge the usage of JHelioviewer software (Müller et al. 2017) and the use of STEREO/SECCHI data (http://stereo-ssc.nascom.nasa.gov/).

References

All Tables

All Figures

	Fig. 1. Image from STEREO-B EUVI (304) of the erupting prominence, taken on 21 September 2009, at 16:36 UT. The image was scaled to enhance the prominence and a mask was applied from 1.35 R⊙ outwards to remove the large noise at the edges.
In the text

	Fig. 2. Top: running difference images from COR2 used in the GCS reconstruction model. Bottom: reconstructions (green wire) using GCS superposed onto the upper corresponding panels in the case of CME1 (left) and CME2 (right).
In the text

	Fig. 3. Height-time plots of the two studied CMEs, erupting between 21 and 22 September 2009, based on STEREO-B coronagraph observations.
In the text

	Fig. 4. Solar wind speed recorded at Earth, taken from the OMNI database. The two red lines indicate the expected arrival time interval of the two CMEs.
In the text

	Fig. 5. Initial magnetic field configuration – Br (colour scale) and selected magnetic field lines in the meridional plane.
In the text

	Fig. 6. Normalised shearing profile (vϕ/vmax) as a function of latitude, at different simulation times from the start of shear.
In the text

Fig. 7. Simulated relative density (grey scale) and selected magnetic field lines during the formation phase of flux ropes (FRs), in the case of: first row – CME1; second row – CME2 (double eruption case); third row – second FR (FR2) which falls back to the Sun; fourth row – stealth eruption; fifth row – failed eruption. The green circles indicate reconnection sites. The relative density is ${\rho }_{\mathrm{rel}}=\frac{\rho (t)-{\rho }_0}{{\rho }_0}$, where ρ0 is the density of the initial relaxed state before the shear.

In the text

Fig. 8. First column: snapshots of observations and MPI-AMRVAC simulations taken when CME1 front is located at 10 R⊙. Second column: snapshots of observations and simulations taken when CME2 front is located at 6 R⊙ (except (h)). Panels a and b: running difference images taken from EUVI, COR1, and COR2 instruments on board STEREO-B and visualised in JHelioviewer software. Snapshots of the simulated relative density and selected magnetic field lines: panels c and d indicate the double eruption scenario, panels e and f show the stealth eruption scenario, and panels g and h indicate the failed eruption scenario. The parameter Δt represents the time interval from the eruption of CME1, chosen as 21 September 2009 21:00 UT, or 9.6 h simulation time after the start of shearing motions, when the centre of the flux rope crosses 1.5 R⊙ and enters COR1-B field of view.

In the text

	Fig. 9. Location of the CMEs fronts in the meridional plane as a function of time. The black dots indicate the observed CME1, black squares indicate the observed CME2, blue triangles represent the simulated CME1, red crosses indicate the simulated CME2 (double eruption scenario), and red diagonal crosses represent the simulated stealthy CME2.
In the text

	Fig. 10. Total speed of the simulated CMEs, calculated at the centre of the flux rope. The blue triangles indicate CME1, the red crosses indicate CME2 (double eruption scenario), the red dots represent CME2 (failed eruption), and the red diagonal crosses indicate stealthy CME2.
In the text

	Fig. 11. Latitude of the simulated CMEs with respect to the Sun’s equator as a function of time, taken at the centre of the flux rope. The legend is the same as in Fig. 10.
In the text

	Fig. 12. Evolution in time of the radial component of the Poynting flux, summed over the whole Sun.
In the text