Coronal disturbances in extreme ultraviolet (EUV), soft X-ray (SXR), and other wavebands during the solar flare and coronal mass ejection (CME) eruptive events have been observed and simulated extensively during past decades (Liu and Ofman, 2014; Warmuth, 2015). To understand their observed characteristics, the researchers proposed three kind of models based on numerical simulation studies, i. e., magnetohydrodynamic (MHD) wave/shock models (Uchida, 1970; Vršnak and Cliver, 2008; Selwa et al., 2012; Wang et al., 2021), non-wave models (Chen et al., 2002; Attrill et al., 2007; Delannée et al., 2008) and hybrid models (Chen et al., 2005; Cohen et al., 2009; Downs et al., 2012; Mei et al., 2020a; Mei et al., 2020b; Downs et al., 2021). The wave models explain the EUV disturbances as fast MHD shock during the flare/CME events. This wave interpretation is supported by lots of observations, including reflection, refraction and transmission across magnetic structures (Thompson and Myers, 2009; Shen and Liu, 2012; Kienreich et al., 2013; Shen et al., 2013; Muhr et al., 2014), broadening of shock front and its decreased amplitude (Wills-Davey et al., 2007; Muhr et al., 2011; Long et al., 2017) and quasi-periodic wave trains (Liu et al., 2011; Liu et al., 2012; Nisticò et al., 2014; Zheng et al., 2018; Shen et al., 2019). In non-wave models, the EUV disturbances were explained as adjusting of the magnetic field due to the expanding CME bubble. In the high-cadence AIA observations, researchers report lots of events with both wave and non-wave disturbances (Liu et al., 2010; Chen and Wu, 2011; Asai et al., 2012; Cheng et al., 2012; Liu et al., 2013; Cunha-Silva et al., 2018; Fulara et al., 2019). The bimodality characteristics of coronal disturbances support the hybrid models. In a typical physical scene of the hybrid model, a fast shock for the wave component of the disturbance appears in front of the upward erupting CME, the CME bubble (Downs et al., 2012) or other accompanied structures, such as helical current boundary/current shell (Delannée et al., 2008; Mei et al., 2020a) are responsible for the non-wave component.

Although the hybrid models have become widely accepted by researchers, there exist some remain problems. For example, in the numerical simulation of EUV disturbances, some physical processes, such as vortices, slow MHD shock wave (SS), and velocity separatrix (VS), has been noticed and should exist as ubiquitous as the fast MHD wave/shock in realistic observations. However, they have not been confirmed by observational studies. Forbes (1990) had shown the vortices on both sides of the CME in a 2D MHD simulation, which has also been confirmed by Wang et al. (2009) and Mei et al. (2012). Wang et al. (2009) had performed a 2D simulation and proposed the SS as one of the physical mechanisms behind coronal disturbances. Furthermore, Mei et al. (2020b) had performed 3D MHD simulation for EUV disturbance and find that the SS is associated with a VS, which separates plasma moving inward to the center of the eruptive source region and plasma moving after the fast shock (FS). These numerical simulations reflect the existence of a velocity distribution with complex structure in the eruptive source region, which results from the interaction among the CME, the FS and other structures. In this work, we utilize the Titov and Démoulin (1999) model (TD99 hereafter) to perform a 3D high-resolution MHD numerical simulation of the eruptive MFR, emphasizing the complex velocity field and corresponding EUV manifestations. Here, we use the forward modeling method (Van Doorsselaere et al., 2016) to create the synthetic SDO/AIA images to directly compare our numerical results with actual observations. In the studies of EUV disturbances and other phenomena during the eruption, the synthetic images method has been widely used to compare the numerical models with the actual observations, such as X-ray sigmoids (Roussev et al., 2012), prominence formation (Xia and Keppens, 2016), the global EUV disturbances (Downs et al., 2012, 2021; Mei et al., 2020a; Mei et al., 2020b), the CME and reconnecting current sheet in EUV emission (Lugaz et al., 2011; Pagano et al., 2014; Zhao et al., 2019; Ye et al., 2020) and white-light (Lugaz et al., 2007; Manchester et al., 2008; Jin et al., 2017). Furthermore, the Fourier local correlation tracking code (Fisher and Welsch, 2008; Fisher and Welsch, 2020) has been applied to the synthetic images to deduce the velocity field, in which the vortices, the SS, and the VS may appear. In Section 2, the setup of this simulation is given; in Section 3, the main results are presented; in the last section, we summarize this work.

Utilizing the MPI-parallelized adaptive mesh refinement code (MPI-AMRVAC) (Keppens et al., 2012; Porth et al., 2014; Xia et al., 2018; Keppens et al., 2020), this simulation’s governing resistive MHD equations are solved numerically by a three-order accurate finite volume scheme, which consists of a Harten-Lax-van Leer approximate Riemann solver (Harten, 1983), a third-order slope limiter (Čada and Torrilhon, 2009) and a three-step Runge–Kutta time-marching method. The simulation domain is a box of size −6 ≤ x ≤ 6, −6 ≤ y ≤ 6, and 0 ≤ z ≤ 12 in the Cartesian coordinate system, which resolves by 360³ uniformly distributed grid points. The dimensionless units for length, velocity, pressure, and magnetic field are 5 × 10⁹ cm, 1.2 × 10⁷ cm s⁻¹, 3.2 Pa and 6.3 G respectively.

The initial magnetic structure comes from the well-known TD99 model, which consists of an MFR to model filament/prominence, a background field to confine the MFR, and a dipole to control the twist feature of the MFR (see Figure 1 in Mei et al., 2020b). The MFR has been described by parameters major radius R, minor radius a, and total toroidal current I along the MFR axis. The background field comes from a pair of magnetic sources ± q separated by a distance L, lying on the MFR symmetry axis. The dipole introduces an extra toroidal component around the MFR to control the twist profile of the MFR. Although the analytical model is much simpler than a realistic magnetic structure, this simple model allows us to study the fundamental physical process during MFR eruption without being disturbed by other features. On the other hand, we have adopted gravity stratification atmosphere and thermal conduction, two crucial components for obtaining reliable density and temperature distributions during the MFR eruption and creating synthetic EUV image. The heat conduction aligns the magnetic field line, hence κ = κ ∥ T 2 / 5 e ̂ B e ̂ B , with e ̂ B = B ⃗ / | B | and κ _∥ = 10^–6 erg s⁻¹ cm⁻¹ K^−3.5. The super-time-stepping scheme has been used to handle the heat conduction term to significantly reduce iterations and improve stability when the plasma temperature exceeds10⁷ K (Alexiades et al., 1996; Meyer et al., 2012; Zhou et al., 2021).

We use a two-layer gravitationally stratified atmosphere, with z ≤ z _p and z > z _p representing the photosphere and the corona, respectively. Here, the photosphere provides a high-beta environment to realize a line-tied bottom boundary, where the magnetic field line foot-points are anchored into photosphere (Wang et al., 2009; Mei et al., 2012; Wang et al., 2015; Xie et al., 2019). At the height z = z _p , the plasma pressure equals 0.2 Pa, the strength of the magnetic field nearby the MFR can reach 100 G, and so the plasma beta value approximately equals 10^–4, which is close to the realistic coronal environment. In addition, for the other five boundaries, we adopt the simplest open boundary conditions, i.e., all physical quantities are deduced via the extrapolation of internal grid points. Detail formulae of components of magnetic structure and two-layer stratified atmosphere had been already given in Mei et al. (2020b). Involved parameters for initial magnetic configuration and atmosphere are the same with Mei et al. (2020a), and one can refer to this work for more detail.

The MFR starts to erupt immediately after the simulation begins because of the un-equilibrium initial magnetic structure and so a net upward Lorenz force acting on the initial MFR. In the meanwhile, the MFR experiences kink instabilities due to the high twist turn of the MFR. Figure 1 shows the evolution snapshots of the eruptive MFR. The golden twisted curves are magnetic field lines inside the MFR, which show significant expansion during the eruption. The electric current distribution on a plane y = 0 shows that several structures appear as a result of the upward rise of the MFR, including the FS, the helical current boundary (HCB) and other features. The FS is a piston-driven shock invoked by the MFR (Wang et al., 2009), which means that its outward moving speed can be much faster than the MFR, as shown in Figure 1. The HCB results from the interaction between the background magnetic field and the upward eruptive magnetic structure. Its helical shape comes from the rotation movement of the MFR because of the kink instability. After a short acceleration process in the very early stage, the MFR, the HCB, and the FS expand outward with almost constant speeds of 360 km s⁻¹, 470 km s⁻¹ and 500 km s⁻¹ respectively. These kinetic features of the eruptive magnetic structure have been significantly affected by the magnetic reconnection rate and fine structure inside a 3D current sheet (CS) (Mei et al., 2017), which grows continuously under the MFR.

The distribution of density ρ, temperature T and velocity curl ∇ × v ⃗ on cut x = 0 and z = z _p (low-left panel) at t = 2 are shown in Figure 2. The grey arrows show the velocity distribution, and its length has been scaled by the norm of velocity vector | v ⃗ | . For the cut x = 0, the velocity distribution shows that parts of plasma move toward the reconnecting CS and the CME bubble, and part of plasma moves outward, following the FS. Between two streams of plasma, there exists a VS marked by a pink curve. The top-tip of this curve connects to a vortex region, which is already noticed by preceding 2D numerical simulations (Forbes, 1990; Wang et al., 2009). For the cut z = z _p , the grey isosurface | v ⃗ | = 0.07 approximately illustrates the position of the 3D VS. Inside the VS, plasma moves toward the center of the eruptive source region. Outside the VS, plasma follows the expanding FS front. In addition, around the VS, there exists a SS, as indicated by colorful shading of velocity curl.

The structure of the velocity field at t = 3 is also shown by eight groups of streamlines in Figure 3. In all panels, the golden and red curves are magnetic field lines of the MFR and the outer boundary of the CME bubble, respectively. In the upper panels, the streamlines (light-green and light-purple curves) illustrate the velocity field structures in quadrants II and IV of the x-y coordinate system. The light-green streamlines show plasma flow related to the reconnection inflow and the CME bubble. They originate from regions R1 marked in the upper-left panel, located in the lower atmosphere near the foot-points of the MFR. The plasma in R1 comes to the side of the 3D CS through a vortex channel surrounding the MFR. Later, it enters into the CS as reconnection inflow and finally becomes parts of the CME bubble or flare loop system. Unlike light-green curves, the light-purple curves illustrate another kind of plasma stream, only related to the FS. It originates from region R2, located at the bottom of the simulation box. It indicates a stream of plasma moves upward and follows the expanding FS. In the lower panels, the streamlines show plasma flow structures in the I and III quadrants. Like the curves in upper panels, the green streamlines relate to the reconnection inflow and CME bubble, and the purple ones are associated with the FS. However, unlike the curves in upper panels, the green and purple streamlines originate from the same region R3, marked on the lower-left panel.

Furthermore, the relationship among streamlines, the VS, and the vortexes are presented in panel (a) of Figure 4. The streamlines in quadrants I and II are the same with the curves in upper panel of Figure 3. The grey iso-surfaces with v = 0.07 illustrate the location of the 3D VS in the lower atmosphere of the eruptive source region. These streamlines show that the plasma inside the VS moves inward to the CS and the CME bubble, and the plasma outside moves outward. The regions R1, R2 and R3 are nearby the VS. Two grey slices inside two boxes have been chosen, as marked in panel (a), to illustrate the detailed information of velocity distribution nearby these regions. Panels (b) and (c) have given velocity fields on these two slices.

In panel (b), we can see a vortex center labeled as C1. The surrounding velocity field arrows show that ambient plasma moves around, converges toward C1, and finally forms a vortex channel, as already demonstrated by the colorful streamlines. This vortex channel transports the plasma near the region R1 and R2 to the ambient region around the reconnection CS and the CME bubble. Part of the transported plasma moves toward the FS, and others move toward the CS and the CME bubble, contributing to the plasma composition of the CME bubble. Although we have not considered the composition and ionization state of the plasma in this work, a more realistic numerical experiment likely also exists in the vortex channel. This transport process suggests that the plasma composition at the lower atmosphere can change the plasma composition inside the CME bubble.

In panel (c), another type of vortex center is labeled as C2, that all velocity arrows rotate around and spread away from it. Also, part of the plasma moves toward the FS, and others move toward the CS and the CME bubble. Unlike C1, no vortex channels associated with C2 can be seen. Because of the slight streamline distortion at the vortex center and the low plasma beta, a noticeable vortex channel can not form, and so there is no continual plasma that has been transported into the vortex center to compensate for scattered plasma away from the center.

In Figure 5, the synthetic SDO/AIA 193Å, 171Å and 131Å images are shown to give observational features of our numerical results. These images are created by utilizing the forward modeling code FOMO (Van Doorsselaere et al., 2016), which translate the plasma density and temperature of the optically thin coronal atmosphere in numerical simulation into the EUV emission and then integrated along the line-of-sight (LOS). The coordination system of our numerical simulation connects to the rotated frame-of-reference of the observer by two angles L and B . Here, L is between the LOS and the z-axis, − B is between the LOS and the y-axis. The synthetic images plane (i.e., the plane-of-the-sky, POS) is perpendicular to the LOS. In the beginning, the MFR appears as a dark filament because of the low temperature of its internal plasma. At t = 1, the MFR has been brightened due to its internal untwisting process and resultant magnetic reconnection. The FS front can be seen clearly as a 3D bright dome. Later, the MFR expands outward adiabatic, and so its brightness decays continuously. Due to the different contribution functions for 171Å and 131Å (Lemen et al., 2012), the MFR and the FS in the lower panels are slightly different from the upper panels. At t = 2, the HCB can been seen in the 171Å and 131Å images. Different from the FS, the HCB should be classified as non-wave components of the EUV disturbance.

The velocity field is very important to diagnose the nature of the EUV disturbance. It is valuable to determine whether we can deduce a reliable plasma velocity distribution from the realistic EUV observations, and whether these deduced velocity field are similar to these shown in Figure 2. Therefore, we utilize the Fourier local correlation tracking code FLCT (Fisher and Welsch, 2008; Fisher and Welsch, 2020) to deduce velocity field. The FLCT estimates a 2D velocity field from two successive images. The first image evolves into the second image over a small time step, which usually depends on observations’ time resolution. Instead to apply the FLCT to realistic observations, here we apply them to our synthetic images with different view angles at t = 2, 2 + Δt, four and 4 + Δt. Here, Δt = 0.0023 times the dimensionless unit of time equals one second. Figure 6 show the synthetic images of AIA 193Å with log-scale. In all panels, black arrows show the deduced velocity fields. In the top panels of Figure 6, the synthetic image shows a situation in which the EUV disturbance has been observed in the limb. The velocity field inside the two boxes marked on the top panels is shown in the middle panels. Unfortunately, the deduced velocity does not consist of inward plasma flow and outward flow, significantly different from the velocity distribution on cut x = 0 in Figure 2. In the bottom panels of Figure 6, the EUV disturbance has been observed on the solar disk. The black arrows show that almost all plasma move outward to the expanding FS front, which is also different from the velocity field on cut z = z _p , shown in Figure 2. The apparent difference between the deduced velocity field and the velocity field shown in Figure 2 comes from the fact that the deduced velocity is based on the EUV emission and involves an integration along the LOS. In the deduced images, we can not see the VS, the SS and the vortex, so that explains why they have not been reported usually in realistic observations, although the numerical simulation indicates that they should present similar to the ubiquitous FS. In other words, applying FLCT on realistic AIA images, the deduced velocity field can not provide essential evidence for the existence of the VS, the SS and the vortex.

In this work, we have performed a 3D MHD simulation for the EUV disturbance during the eruptions, emphasizing the complex velocity distribution in the lower atmosphere around the eruptive source region. The TD99 model (Titov and Démoulin, 1999) has been used as an initial un-equilibrium magnetic structure, in which a magnetic fluxrope (MFR) models prominence or filament in the corona. An isothermal gravitationally stratified atmosphere has been used to model the background corona. During the MFR eruption, the current sheet (CS) under the MFR grows continuously. The magnetic reconnection inside the CS generates new magnetic field lines to attach to the expanding CME bubble. In front of the CME bubble, the fast shock (FS) and the following helical current boundary (HCB) appear. The HCB comes from the interaction between the CME bubble and the background field. To directly compare with realistic EUV observations, we created synthetic SDO/AIA images for different wavelengths. In synthetic images, the FS moves outward as a 3D dome, followed by the HCS and the MFR, which form the typical three-components CME and may also correspond to the non-wave components of EUV disturbances.

At flanks of the CME bubble, the velocity distribution of the lower atmosphere develops complex structures. Two streams of plasma exist in the lower atmosphere of the eruptive source region, divided by a 3D velocity separatrix (VS). Outside the VS, plasma moves outward to the expanding FS front. Inside the VS, the plasma moves toward the center of the source region. The interaction of two streams of plasma flows has invoked two types of vortexes and the slow shocks (SS) near the VS. The plasma around the first kind of vortex converges to the vortex center. It forms a vortex channel, which lifts plasma nearby the bottom of the simulation box and provides reconnection inflow for the CS and continual plasma for the outermost boundary of the CME bubble. For the second type of vortex, the plasma spreads out from the vortex center, located at higher position than the first kind center, and no associated vortex channel has been observed. In addition, we use the local correlation tracking method to deduce plasma velocity field based on the successive synthetic images. Regrettably, the deduced velocity distributions based on the synthetic images are significantly different from the velocity distribution on the cuts, such as x = 0 and z = z _p shown in Figure 2. For the cases of EUV disturbances observed on the solar disk or the limb, the distribution of deduced velocity shows that almost all plasma moves outward and almost no plasma moves inward, so there is no VS and evidence of the SS. The synthetic image involves the integration of EUV emission along the line-of-sight, so that the deduced velocity can not represent the complex 3D velocity field of plasma during the eruptive events. Thus, it is not easy to find the evidence for the VS, the SS, and the vortex, except we have an in-situ measure datum of eruptive events.