Modeling the 10 September 2017 solar energetic particle event using the iPATH model
RResearch in Astronomy and Astrophysics manuscript no.(L A TEX: ms2020-0054.tex; printed on May 6, 2020; 0:31)
Modeling the 10 September 2017 solar energetic particle event usingthe iPATH model
Zheyi Ding , , Gang Li , ∗ , Junxiang Hu , Shuai Fu , School of Geophysics and Information Technology, China University of Geosciences (Beijing), Beijing100083, China;Department of Space Science and CSPAR, University of Alabama in Huntsville, Huntsville, AL 35899,USA; [email protected]
Lunar and Planetary Science Laboratory, Macau University of Science and Technology, Macau 519020,China
Abstract
On September 10 2017, a fast coronal mass ejection (CME) erupted from the activeregion (AR) 12673, leading to a ground level enhancement (GLE) event at Earth. Usingthe 2D improved Particle Acceleration and Transport in the Heliosphere (iPATH) model, wemodel the large solar energetic particle (SEP) event of 10 September 2017 observed at Earth,Mars and STEREO-A. Based on observational evidence,we assume that the CME-drivenshock experienced a large lateral expansion shortly after the eruption, which is modelled by adouble Gaussian velocity profile in this simulation. We use the in-situ shock arrival times andthe observed CME speeds at multiple spacecraft near Earth and Mars as constraints to adjustthe input model parameters. The modelled time intensity profiles and fluence for energeticprotons are then compared with observations. Reasonable agreements with observations atMars and STEREO-A are found. The simulated results at Earth differ from observationsof GOES-15. Instead, the simulated results at a heliocentric longitude 20 ◦ west to Earth fitreasonably well with the GOES observation. This can be explained if the pre-event solar windmagnetic field at Earth is not described by a nominal Parker field. Our results suggest that alarge lateral expansion of the CME-driven shock and a distorted interplanetary magnetic fielddue to previous events can be important in understanding this GLE event. Key words:
Sun: coronal mass ejections (CMEs) ł Sun: magnetic fields ł Sun: particle emis-sion a r X i v : . [ phy s i c s . s p ace - ph ] M a y Z.-Y. Ding et al.
Solar energetic particle events are historically classified into two broad categories: impulsive and gradualevents (Cane et al. 1986; Reames 1995, 1999). In this paradigm, impulsive SEPs are accelerated at solarflares and propagate along the interplanetary magnetic field to Earth with a rapid rise and decay phase inthe particle time intensity profiles. In comparison, gradual SEPs are accelerated via the diffusive shockacceleration mechanism at CME-driven shocks. Because shocks are of large scale and last much longerthan flares, these events are characterized by a prolonged intensity profile and often higher fluences thanimpulsive events. Large gradual SEP events are of particular concern because the accompanying high-energy protons pose the most serious radiation threats to astronauts living beyond low-Earth orbit and candamage electronics on satellites in the space (Desai & Giacalone 2016). In some largest events, acceleratedparticles can reach energies up to several GeV, leading to significant increases in particle count rates throughneutron monitors at the Earth’s surface and are known as Ground Level Enhancements. Large SEP events areusually associated with fast shocks and high ambient energetic particle intensity (referred as seed particles)prior to the shock (Kahler 1996; Kahler et al. 2000). The increasing abundance of the seed particles ispossibly from preceding flares (Mason et al. 2019) or preceding CMEs (Gopalswamy et al. 2004; Li &Zank 2005). Li et al. (2012) proposed a twin-CME scenario for GLE events in which two CMEs erupt fromsame or near active region within a period of 9 hours. The preceding CME and its driven shock disturbsthe coronal and interplanetary environment and leads to enhanced turbulence level which facilitates a moreefficient acceleration at the second CME-driven shock. The scenario was later extended to large SEP events(Ding et al. 2013).The most recent GLE event occurred on September 10 2017, classified as GLE 72 (Mishev et al.2018). A X8.2 class solar flare erupted around 15:35 UT from AR 12673 at S09W88, followed by a wideand very fast Coronal Mass Ejection. The eruption was well observed by multiple spacecraft at differentlongitudes, provided a stereo view of this extreme case. The halo CME was first observed by the Solar andHeliospheric Observatory (SOHO) Large Angle and Spectrometric Coronagraph (LASCO) at 16:00 UT,with a reported linear speed of ∼ ◦ to 105 ◦ (Luhmann et al. 2018). However, this event would not be classified as a”twin-CME” (Li et al. 2012), because the closest preceding CME erupted on 23:12 UT on 9 September2017, ∼
17 hours before the main eruption, exceeding the time interval threshold of hours (Ding et al.2014) as a twin-CME. The global EUV waves (Hu et al. 2019; Liu et al. 2018), signatures of a flux ropeand a long current sheet (Seaton & Darnel 2018; Warren et al. 2018; Yan et al. 2018) were observedby the Solar Dynamics Observatory (SDO). Liu et al. (2019) studied the geometry and kinematics ofthe CME-driven shock. They argued that the close shock arrival times at Earth and Mars is a result oflarge lateral expansion of the shock early in the eruption. They suggested that a large lateral expansion ofthe shock can affect particle acceleration and consequently the observations at different longitudes. Large odeling solar energetic particle event of 10 September 2017 3 lateral expansion of the CME-driven shock for large CME eruptions have been investigated previously.Gopalswamy et al. (2012) examined the white-light CME evolution and noted a rapid and large lateralexpansion of the CME erupted on 2010 June 13. Kwon et al. (2015) noted the apparent width of haloCMEs as seen from multiple spacecraft was related to the expanding shock with a 360 ◦ envelope. Recentworks (e.g. Liu et al. 2017, 2019; Zhu et al. 2018) have suggested that lateral expansions of shock occurfrequently in large big eruptions. However, the effects of lateral expansion on solar energetic particles havenot been considered before. In this work, we will examine this problem.Energetic particles associated with this fast CME were observed at Earth and Mars as well as SolarTerrestrial Relations Observatory A (STEREO-A). GLE was observed by several neutron monitors stationsat about 16:15 UT on 10 September 2017 (Guo et al. 2018; Mishev et al. 2018). The different energychannels from the instruments on board Geostationary Operational Environmental Satellite 15 (GOES-15)show an intensive and long-lasting enhancement of energetic particles intensity (Bruno et al. 2019; Guo etal. 2018). This event observed at Earth has a soft spectrum at high energies compared with most of GLEevents in the solar cycle 23 (Bruno et al. 2019; Cohen & Mewaldt 2018; Gopalswamy et al. 2018). On 10September at around 19:50 UT, the high energy protons were detected by Radiation Assessment Detector(RAD) (Hassler et al. 2012) from Mars Science Laboratory (MSL) (Ehresmann et al. 2018; Guo et al.2018; Lee et al. 2018; Zeitlin et al. 2018). Later, high energetic protons were observed with a slow increaseat STEREO-A at around 08:00 UT on 11 September 2017 (Guo et al. 2018; Lee et al. 2018).Attempts of modeling particle acceleration and transport of SEP events have been undertaken by manygroups (e.g. Heras et al. 1995; Kallenrode 1997; Luhmann et al. 2007, 2010, 2017). These models treatCME-driven shock as a source of energetic particles without acceleration process. A dynamic onion-shellmodel of the strong shock propagation and particle acceleration has been developed by Zank et al. (2000).This model was improved by Rice et al. (2003) for shocks with arbitrary intensities, and by Li et al.(2003) to model the transport of energetic particles using a Monte-Carlo approach. Li et al. (2005) furtherextended the model to include heavy ions. A comprehensive numerical model, developed by these authors,called Particle Acceleration and Transport in the Heliosphere(PATH) model. Modeling specific SEP eventsusing the PATH model shows their reasonable agreement (Verkhoglyadova et al. 2009, 2010). Recently, Huet al. (2017) have extended the PATH to a 2D model, named improved Particle Acceleration and Transportin the Heliosphere (iPATH). The new model has the capability to study the characters of particle accelerationand transport at a 2D shock in multiple locations of the ecliptic plane. Hu et al. (2018) modelled an examplegradual SEP event as observed at multiple locations.In the iPATH model, the shock speed profile at the inner boundary is assumed to be a Gaussian formin longitude and the propagation direction of the shock is assumed to be radial. Such an assumption isfor simplicity. As we discussed earlier, large lateral expansion of shocks can be common for large events.A lateral expansion implies that the open angle of the shock will increase over time. This means that thetreatment of the inner boundary of the shock profile in the iPATH model needs to be improved. In this Z.-Y. Ding et al. work, we modify the inner boundary of the shock profile to examine the effect of the large lateral expansionon producing solar energetic particles. We note that because the background solar wind is unlikely to behomogeneous, the lateral expansion and the shock profiles do not need to be homogeneous. Therefore thereis no need to assume a symmetric eruption.In this paper, we discuss our modeling of the 10 September 2017 GLE event using the iPATH model.In section 2, we describe the iPATH model and the model setup. In section 3, we discuss the CME-drivenshock configuration and shock parameters. Modelled time intensity profile are compared with GOES-15observation results in different energy channels. We integrate the data of proton flux from the AdvancedComposition Explorer (ACE) and GOES-15 and fit the spectrum using the Band function. The modelledspectrum is compared with observation fitted results. We also compare the modelled time intensity profileswith the observation results at Mars and STEREO-A. We then discuss our results and conclude in section4.
The iPATH model is a two-dimensional MHD code plus a particle transport code developed by Hu et al.(2017). It is a continuation of the earlier 1D PATH model (Zank et al. 2000;Li et al. 2003, 2005), addressingsolar energetic particle events in the ecliptic plane. It contains two separate modules: the first module modelsthe background solar wind and follow particle acceleration at the shock front, the second module modelsthe transport of SEPs escaped from the upstream of the CME-driven shock. Here we limit ourselves to themodifications that we introduce to the iPATH in modeling the 10 September 2017 event.We model the propagation of background solar wind and CME-driven shock limited only in the eclipticplane using a 2D MHD code. For simplicity, the background solar wind is assumed to be homogeneous. Weuse the -hour averaged in-situ solar wind observation near Earth before the CME eruption as the solar windinputs. The CME-driven shock is treated by perturbing the inner boundary (proton number density n , solarwind speed V sw and temperature T ) at 0.05 AU for a short period of time (e.g. 1 hour). Such a simplifiedtreatment is similar to the WSA-ENLIL+Cone model for the space weather prediction. The inner boundaryis set at . AU ( ∼ R s ) in this work, thus we can not model particle acceleration in the low corona.Modeling CME eruption in the low corona needs a more detailed description about the magnetic field andthe corona condition. However, CME-driven shock can be formed in the low corona (Gopalswamy et al.2013, Liu et al. 2009) and some highest energy particles are produced within several solar radii. Modelingthis part of SEPs is beyond the scope of the iPATH model. In this work, we introduce two Gaussian velocityprofiles to account for both a lateral expansion of the shock and any inhomogeneity in the background solarwind. Specifically, we consider the following eruption profile: A ( φ ) = A e − ( φ − φc σ < t < D A e − ( φ − φc σ H [ φ − φ min ] ∗ H [ φ max − φ ] 0 < t − t s < D (1) odeling solar energetic particle event of 10 September 2017 5 where A , A are the perturbed model parameters (number density, speed, temperature) at longitudes φ c and φ c . The angles φ c and φ c are the central longitude of the two Gaussian distributions. The variances σ , σ are related to the full width at half maximum (FWHM) of the perturbed width. For the second profile,we identify a range given by [ φ min , φ max ] as the perturbed range associated with the lateral expansion. D , D are the perturbing duration of two eruption profiles, and t s is the start time of the second eruption,which we take as the start time of the lateral expansion. H is the Heaviside function. For simplicity, thebackground interplanetary magnetic field is assumed to be a Parker spiral, B r = B (cid:18) R r (cid:19) ; B φ = B r (cid:18) Ω r sin θu sw (cid:19) ( r (cid:29) R ) , (2)where B r and B φ are the radial and azimuthal component of the IMF at heliocentric distance r . u sw is solarwind speed. θ = 90 ◦ corresponds to the magnetic field in the ecliptic plane. B is the radial component ofthe IMF at R .After perturbing the inner boundary, the CME-driven shock is followed in the code and the shock pa-rameters are calculated at every time step. From these shock parameters we obtain the dynamic timescale t dyn = RdR/dt . Balancing t dyn with the acceleration time scale yields the maximum particle momentum p max (Drury 1983): t dyn = (cid:90) p max p ss − κU shk p dp, (3)where p is the injection momentum, s is the shock compression ratio, κ is the diffusion coefficient ofparticle, U shk is the shock speed in the upstream frame. Once p max is obtained, the particle distributionfunction in the outmost parcel ( j r , k φ ) at each time step t k is given by, f ( j, k, p, t k ) = c ∗ (cid:15) j,k n t k ,k p − β H [ p − p j,kinj ] ∗ H [ p j,kmax − p ] exp (cid:18) − EE (cid:19) (4)where β = s j,k s j,k − , (cid:15) j,k is the injection efficiency, n t k ,k is the upstream solar wind density at the time t k infront of the parcel ( j r , k φ ) . p inj is the particle injection momentum. E is the particle energy and E is thekinetic energy corresponds to the maximum particle momentum p max obtained in the equation (3). In theabove H is the Heaviside function. c is a normalization constant, c = 1 / (cid:90) p j,kmax p j,kinj p − β (cid:110) H [ p − p j,kinj ] ∗ H [ p j,kmax − p ] (cid:111) d p (5)Note that this functional form is different from previous works (e.g. Hu et al. 2017; Li et al. 2005), wherethe exponential exp ( − E/E ) tail was not included. In Drury (1983), the dynamical time t dyn in equa-tion (3) is to be understood as the average time for accelerating particles from p inj to p max . Clearly, someportion of these particles would take longer time than t dyn . However, when t > t dyn , shock parameterscan not be regarded as constant anymore, so we expect a roll-over of f ( p ) near p max . This, of course, isthe consequence of a finite acceleration time. In an earlier paper, Forman & Drury (1983) examined the Z.-Y. Ding et al. effect of finite time analytically and showed that the particle distribution function at higher momentum is ofexponential decay of the diffusion coefficient, which is a function of particle momentum. Ellison & Ramaty(1984) adopted an exponential decay tail at high energy ∼ exp ( − E/E ) to account for the effect of finiteshock size and finite acceleration time. Some recent numerical simulations by Zuo et al. (2011) and Konget al. (2019) which examined particle acceleration at a prescribed shock showed that the time dependentparticle spectra are consistent with a power law with an exponential tail. With these considerations, weadopt equation (4).To model SEP events at multiple locations simultaneously, the iPATH model uses a 2D onion shellmodule to keep track energetic particles in the shock complex. At the j -th time step, the j -th shell (theoutermost shell) is generated and divided longitudinally into different parcels ( j r , k φ ) with an angular widthof ◦ . Particles accelerated at the shock front experience convection with the parcel and diffusion betweenparcels. Note that particles only diffuse in parcels whose p j,kmax are greater than particles’ momentum p .This is because for parcels with p j,kmax < p , there is no excited wave turbulence that can trap these particles.Particles can escape when they diffuse far enough ahead of the shock. Zank et al. (2000) assumed anescape length l = 4 λ esc , where λ esc = κ rr U shk is the scatter length scale and κ rr is the particle upstreamdiffusion coefficient in the direction of shock normal. Within the length l , the excited wave density issignificantly higher than that of the ambient solar wind. If the wave intensity has no x -dependence, then theparticle distribution function decay exponentially with the escape length in the shock front. This allows oneto obtain the escaped particle distribution function. Alternatively Li et al. (2005) calculated the escapedparticle number explicitly instead of using the particle distribution function. Note, the escape boundary isrelated to each individual parcel. In this work, we follow Li et al. (2005) and obtain the escaped particlenumber of momentum p at time t k and at longitude φ related to the outmost parcel ( j r , k φ ) by, N esc ( k, p, t k ) = J ( k,p ) (cid:88) j =1 √ π ∆ R ∗ N j,k ( t k − r ∗ j,k (cid:40) exp( − x ) + √ πr ∗ j,k ∆ R ∗ [1 − erf ( x )] (cid:41) , (6)where N esc ( k, p, t k ) is the number of particles escaped from shock complex at time t k and longitude φ within d p phase space; J ( k, p ) is the shell number for which p max = p ; (∆ R ∗ ) = 4 κ j,k ∗ ( t k − t k − ) +2([ r j +1 ,k ( t k ) − r j,k ( t k )] / decides how many particles escape from the parcel ( j, k ) . This parameter islarger for high energy particles since the diffusion coefficient κ increases with momentum. N j,k ( t k − isthe particle number in the parcel ( j, k ) before diffusion, r ∗ j = ( r j,k + r j +1 ,k ) / represents the heliocentricdistance of the center of parcel ( j, k ) , x = ( r esc − r ∗ j,k ) / ∆ R ∗ is a dimensionless parameter, and r esc = r J ( k,p ) ,k + l represents the escaping boundary for particles of momentum p .Once particles escape from the shock complex, they propagate along the IMF in the relatively undis-turbed solar wind. We follow closely Hu et al. (2017) in describing the transport of escaping particles in thesecond module of the iPATH model. In particular, we also include the cross-field diffusion and use the back-ward stochastic differential equation method. The cross-field diffusion is governed by a term ∼ (cid:79) · ( κ ⊥ · (cid:79) f ) (see equation (18) of Hu et al. (2017)), where κ ⊥ is the perpendicular diffusion coefficient of energetic par- odeling solar energetic particle event of 10 September 2017 7 ticles, Hu et al. (2017) calculated κ ⊥ from the parallel diffusion coefficient κ (cid:107) based on the Non-LinearGuiding Center theory (NLGC). In the NLGC theory, the solar wind turbulence is assumed to have a ge-ometry of ”slab+2D”, where the slab component describes Alfv´enic fluctuation along the background fieldand the 2D turbulence describes how the background field line (here the Parker field) meanders. An impor-tant parameter describing the 2D component of the turbulence is the bend-over scale l D (Shalchi et al.2010). In the work of Hu et al. (2017) the radial dependence of l slab and l D was ignored. Hu et al. (2018)considered a more general case of, l slab ∼ r α ; l D ∼ r α , (7)where α is a parameter and close to . In this work, we set α = 0 . and assume the turbulent magnetic fieldsquare δB ∼ r − , thus the radial dependence of κ (cid:107) and κ ⊥ used in Hu et al. (2018) has the followingform: κ (cid:107) ∼ v B r α +3 ; κ ⊥ ∼ v B − r α − (8)where v is particle speed, B is the background magnetic field, and r the heliocentric distance. In this work,we set κ ⊥ /κ (cid:107) to be . at 1AU for 1 MeV proton. With this choice of κ ⊥ , we examine cross fieldtransport of energetic particles and discuss the observations at Mars and STA, which, as we show below,depend strongly on κ ⊥ . For our simulation, we average hours of solar wind observation before the CME eruption as a guide tothe background solar wind inputs. This yields a solar wind speed as km s − , a magnetic field of . nT and a number density of 0.8 cm − at 1 AU. Note that there were some preceding CMEs, which mostlyaffects the number density of these three parameters. To account for the recovery of the corona after thepreceding CME, we use the 8-hours averaged number density of 3.5 cm − before shock arrival. As shown inequation (1), we use a double Gaussian profile for the velocity disturbance and the corresponding parametersdescribing these two Gaussian profiles are listed in Table 1. To model the large lateral expansion, we usea large variance for the second Gaussian distribution but limit the perturbation only within the longitudesof [ φ min , φ max ] = [10 , . The perturbation duration ( D , D ) both are set to be hour and the start time( t s ) for the second perturbation is at hour. The injection efficiency is assumed to be 1 % for θ BN = 0 (i.e.parallel shock geometry). We set the ambient turbulence level δB /B to be . at 1 AU.Figure 1 shows the location of various spacecraft when the event occurred and plots the scaled numberdensity for two cases: one for the single Gaussian perturbation and the other for a double Gaussian per-turbation. When the CME erupted, Earth was located at r = 1 . AU with φ = 0 ◦ ; and Mars was locatedat . AU and φ = 157 . ◦ , STEREO-A is located at . AU and φ = 232 ◦ . Venus and Mercury areon the propagation path of this CME. However, there were no available data at these two locations. The Z.-Y. Ding et al.
Table 1: The initial CME parameters of double Gaussian distributions
Gaussian distribution Initial speed Density amplitude Perturbation center variances i km s − φ c i σ i ◦ ◦ ◦ ◦ Notes: Density amplitude refers to the ratio of the enhanced number density to that of the background at theinner boundary.
CME-driven shock arrived at Earth and Mars at ∼ ∼
59 hours, respectively, after the eruption(Guo et al. 2018). Previous work by Guo et al. (2018) of this event has obtained a speed of km s − at . R s with the central axis at 110 ◦ . In another study by Luhmann et al. (2018) employing the Conemodel, the authors obtained a radial speed of km s − at . R s with the center axis at 108 ◦ . In oursimulation with a single Gaussian profile, shown in Figure 1 (a), we assumed the central axis is along φ c = ◦ and the variance σ of the Gaussian disturbance was set to be ◦ . At the inner boundary located at0.05 AU ( ∼ R s ) we increase the solar wind speed to km s − and the number density from thebackground by a factor of at φ c , and the disturbance lasts hour. With these parameters, our simulatedshock arrival time is . hours at Earth and . hours at φ = 155 ◦ near Mars, in reasonable agreementto observations. However, with these parameters, the eastern flank of the shock, magnetically connecting toEarth, is very weak and it does not accelerate particles to high energies. As noted by Liu et al. (2019), largelateral expansion of CME-driven shock could lead to enhanced SEPs for observers that are connected tothe shock flank. To evaluate the effect of the large lateral expansion, panel (b) of Figure 1 shows the shockprofile using a double-Gaussian profile with parameters given in Table 1. As we can see, a high densityblob now appears in the eastern flank, leading to a locally enhanced shock compression ratio. The modelledshock arrival times are also similar with observations: . hours of Earth and . hours in the longitudeof 155 ◦ near Mars. The shock speed at Earth is about km s − (Liu et al. 2019), and the modelled shockspeed at Earth is about km s − .Figure 2 plots the evolution of shock parameters in the ecliptic plane. Three observers at 1 AU withdifferent longitudes: 0 ◦ ,10 ◦ ,20 ◦ are chosen. The white dash lines in the figure are interplanetary magneticfield lines assuming a solar wind speed of km/s. From panel (a) we can see that the eastern flank of theshock maintains a quasi-parallel geometry for the entire propagation, implying that the eastern flank of theshock has a higher efficiency to accelerate particles. Panel (b) shows the evolution of the shock speed. Wenote that the shock speed is smaller than km/s most of the time when Earth connected to the shockfronts. Panel (c) shows the maximum particle energy at shock fronts as the shock propagates from R s to AU. The maximum energy occurs near the Sun at a longitude ∼ ◦ . However, even along the fieldline connecting to Earth, the maximum energy can reach hundreds MeV early on.This is due to large lateralexpansion: the presence of the second Gaussian profile leads to a larger compression ratio and a larger shockspeed in the eastern flank. Panel (d) shows the compression ratio along the shock front during the shock odeling solar energetic particle event of 10 September 2017 9 propagation. There are two high-compression-ratio lobes, corresponding to the two velocity profiles. Thehigh compression ratio in the eastern flank enhances the acceleration efficiency, indicating the importanceof the large lateral expansion.In Figure 3, the modelled time intensity profiles, shown in panels (b) to (d), at longitudes of ◦ and ◦ are compared with the observations of the Energetic Proton, Electron, and Alpha Detector (EPEAD-A)instrument on board GOES-15, shown in panel (a). Note that the EPEAD-A detector has a westward-viewing angle. From our model, we choose energy channels that are similar to the observations. Becausethe shock arrival time differs for different longitudes, the x axis is normalized by T n , the shock arrivaltime for each observer. From panel (a) we can see that the observed intensity at all energy channels showsubstantial enhancements at the beginning of the event and are followed by slow decays. Panel (b) presentsthe simulated result at Earth (i.e. φ = 0 ). The intensity shows a rapid enhancement at the beginning, similarto the observation. However, the decay is also very rapid, which differ considerably from the observation.As shown in panel (b) and (d) of Figure 2, this is due to the fact that after ∼ . AU, Earth is connected to aportion of the shock that has a weaker compression ratio and a low speed. Note that this connection is underthe assumption that the IMF is given by a Parker field. There was some preceding CMEs in this event and themagnetic connection can therefore be different from the Parker field. We will discuss this later. Panel (c) isthe simulated result at longitude of 20 ◦ . It is interesting to note that simulation results at longitudes 20 ◦ showsimilar decay behavior as the observation. This observer maintains relatively long duration connecting toregions of the shock that have high compression ratios. The presence of this long-duration high compressionratio is the result of the second Gaussian profile. Clearly, the assumption of the large lateral expansion ofCME-driven shock is crucial for our results. We also note that if the preceding CMEs perturb the IMF suchthat Earth is also magnetically connected to the second Gaussian for an extended period, then we wouldexpect observers at Earth could see results similar to panel (c). Panel (d) shows the simulated result at Earthwith five times as large as κ ⊥ . The results between panel (b) and panel (d) are only slightly different, whichillustrates that by invoking a large perpendicular diffusion alone can not explain the Earth observation. Onemay think that by aligning the center axis of the second Gaussian profile closer to the Earth direction canimprove the simulation results at φ = 0 ◦ with an ambient field that is Parker like. However, the profile of thesecond Gaussian profile is constrained by the shock arrival time at Earth. If one further tunes the center axisof the second Gaussian profile toward φ = 0 ◦ , the shock arrival time will be earlier than the observation.As a consequence, one has to decrease the amplitude of the second Gaussian, which will lead to a weakershock and does not accelerate particles to high enough energy early on. Furthermore, if the second Gaussianprofile is further towards Earth, the overall CME and shock profile will be quite distorted and appear morelike two separate eruptions, which is not supported by observations.Besides the time intensity profiles for the energies, we also consider the event-integrated spectra at φ = 20 ◦ . Reasonable agreement with the observation is also obtained. Left panel of Figure 4 shows thesimulated and observed event-integrated proton fluence. For the observation, we use the lower-energy data from the Low Energy Magnetic Spectrometer-120 (LEMS-120) of the Electron, Proton and Alpha Monitor(EPAM) on board ACE, for four energy channels ranging from 330 keV to 4.75 MeV; and the high-energydata from EPEAD-A (P2-P6) and the High Energy Proton and Alpha Detector (HEPAD, P8-P10) on boardGOES-15. The integration period is from 16:05 UT on 10 September 2017 to 18:30 UT on 12 September2017.We fit both the observed and the modelled spectra with a Band-function form (Band et al. 1993), givenby, J ( E ) = CE − γ a exp (cid:18) − EE (cid:19) f or E (cid:54) ( γ b − γ a ) E ; J ( E ) = CE − γ a [( γ b − γ a ) E ] γ b − γ a exp( γ a − γ b ) f or E (cid:62) ( γ b − γ a ) E , (9)where J ( E ) is the particle fluence, C is a normalization constant. γ a and γ b are the spectral indices in thelow energy region and high energy region respectively. E is the spectral break energy, which typicallyoccurs at energies of tens MeV. Following Bruno et al. (2019), we use the calibration schemes by Sandberget al. (2014) and Bruno (2017), below and above 80 MeV, respectively, to obtain the mean energy valuesof every channel. Using different satellite instruments, calibration schemes and integrated periods leadto different fitted parameters. Cohen & Mewaldt (2018) obtained γ a ∼ γ b ∼ E ∼ γ a is caused by using data from the Ultra-Low Energy Isotope Spectrometer (ULEIS)instrument on board ACE. Bruno et al. (2019) obtained γ a ∼ γ b ∼ . and E ∼ . MeV fora longer period, from 16:05 UT on 10 September 2017 to 00:00 UT on 17 September 2017. The fittedparameters for our period in this work are γ a ∼ γ b ∼ . and E ∼ . MeV, similar to Brunoet al. (2019). In comparison, the fitted parameters of the modelled fluence observed at longitude of 20 ◦ are: γ a ∼ . , γ b ∼ . and E ∼ . MeV. These are similar to the observations. The right panelof Figure 4 shows the simulated event-integrated proton fluence at different longitudes of 0 ◦ , 10 ◦ and 20 ◦ .Note that while the time intensity profiles at these three locations are quite different as shown in Figure 3,the differences of the event-integrated spectra of these three observers are however, small. Nonetheless,comparing the three spectral indices γ b in the right panel shows that the high energy portion of the spectrabecomes softer as the longitude decreases. This is because at lower longitudes the observer connects to partsof the shock flank with a compression ratio decreases quicker.This GLE event is one of the only two GLE events identified in the solar cycle 24. The high energyspectral index of this GLE event was softer than many GLE events observed in cycle 23 (Cohen & Mewaldt2018; Gopalswamy et al. 2018). Gopalswamy et al. (2018) suggested that the soft spectrum may becaused by poor latitudinal and longitudinal connectivities. Because the iPATH model is a 2D code and islimited to the ecliptic plane, we cannot consider the latitudinal effects. As to the longitudinal connection,our results suggests that the observation at Earth (with φ = 0 ◦ ) is very close to the simulation results atlongitude of 20 ◦ . Note that our simulation assumes the background magnetic field is given by the Parkerfield. Distortion of the background field is possible due to preceding CMEs. There were multiple CMEsprior to the September 10 event. A large and fast CME occurred on the 6th, whose propagation direction is odeling solar energetic particle event of 10 September 2017 11 ◦ west to Earth. At Earth, the magnetic cloud (MC) behind the shock was observed on the 8th and endson the beginning of 9th (Werner et al. 2019). It is very likely that the MC is a manifestation of a flux ropestructure whose two feet are anchored on the solar surface. If so, as this flux rope can further propagateout, it can affect the field that is connected to Earth on the 10th. Besides the CME on the 6th, there werealso multiple smaller CMEs occurred on the 8th and 9th. The effects of these smaller CMEs on the globalmagnetic field configuration are also unclear. In our simulation, we assume the background field is of Parker,so we can not address the connection issue. Note that the large lateral expansion of the CME-driven shock iscrucial for obtaining our results. Liu et al. (2019) have suggested that the enormous lateral expansion of theshock has an important effect on SEPs production of this event. However, to our knowledge, there has beenno simulations addressing SEPs that consider the effect of lateral expansion of a CME and its driven shock.A large lateral expansion of the shock can affect both the time intensity profile and the even-integratedspectra. This is clearly illustrated in Figure 2, from which we can see that the observer at φ = 20 ◦ caninitially connect to a weak part at the edge of first Gaussian profile (compression ratio ∼ ) of the shock,and followed by a stronger part (compression ratio ∼ ) of the shock.The high energy SEPs associated with the fast CME were also detected by Mars and STA. Zeitlin et al.(2018) reported the radiation dose rate on the surface of Mars measured by the Radiation AssessmentDetector (RAD), mounted on the Mars Science Laboratory’s (MSL’s) curiosity rover. The black curveshown in Figure 5 is the derived integral flux of proton energy greater than 113 MeV observed by RAD(Zeitlin et al. 2018). The proton onset time is 19:50 on 10 September 2017 at Mars, which was ∼ MeV. The event peaked on 11 September, from about 4:00to 14:00 UT (Zeitlin et al. 2018). The modelled result shows a similar peak in the figure. From panel (b)of Figure 1 we can see that Mars was connected to the very western flank of the shock initially and asthe shock propagates out, the connection becomes ever weaker. Many SEPs observed at Mars are particlesthat undergo cross field diffusion in the IMF. They are accelerated at the stronger part of the shock. To theobserver at Mars, the event is a typical eastern event. Therefore we see a clear decay phase after the peak.The shock arrival time is ∼ ∼ . - MeV proton (red curve) from the High Energy Telescope (HET) instruments onboard STEREO-A.Both of them show a gradual rising phase after the CME eruption. The absolute value of the flux are alsosimilar. The STA was at longitude of 232 ◦ which was connected to the back side of the sun when the CMEerupted. Therefore the SEPs detected at STA must have undergone cross field diffusion as they propagateout. Our results suggest that the perpendicular diffusion coefficients and its radial dependence used in thiswork are reasonable. In this paper, we model the 10 September 2017 SEP event using the iPATH model. The most important as-sumption we made in this work is a large lateral expansion of CME-driven shock, modelled by two Gaussianvelocity profiles. In previous practices of iPATH modeling (Fu et al. 2019; Hu et al. 2017, 2018), a sin-gle symmetric Gaussian-shape velocity profile is adopted to simulate the inner boundary condition of thethe CME-driven shock. To illustrate the lateral expansion, we use a velocity profile at the inner boundaryconsisting of two Gaussian profiles with different propagation directions and different start times. We thenexamine the time intensity profiles and the event-integrated spectra at multiple observation locations corre-sponding to observers at Earth, Mars and the STEREO-A spacecraft. We adjust the initial CME perturbationparameters such that the modelled shock arrival times at Earth and Mars are close to the observations. Theintroduction of a second velocity profile to mimic the lateral expansion leads to stronger particle accelera-tion at the eastern flank of the shock, which is necessary to understand the observation at Earth. Our modelresults at φ = 0 ◦ , which is the Earth’s location, show similar event-integrated spectra as the ACE and GOESobservations, but the modelled time intensity profile shows a clear differences from that of observation. Incomparison, the model results at φ = 20 ◦ are comparable to ACE and GOES observations for both thetime intensity profiles and the event-integrated spectra. We interpret this by a possible field line distortiondue to preceding events. There were multiple smaller CMEs before this event, which can cause the IMFto be deviated from the nominal Parker configuration. We do point out that in the acceleration module ofour model a Parker-like background IMF is assumed. When the unperturbed field line is non-Parker, it willaffect the shock obliquity, and therefore the acceleration process. Modelling a CME with a non-Parker fieldis out of the scope of this work and will be pursued in a future work.We also obtain time intensity profiles at longitudes that correspond to Mars and STEREO-A. The modelresults at these two locations are in good agreements with the observations. Because both Mars and STA arenot well connected to the CME-driven shock, SEPs observed in these locations have to undergo cross-fielddiffusion. The agreements between the modeling results and the observations therefore provide a confirma-tion of our choice of the perpendicular diffusion coefficient, and its radial dependence.Our study suggests that modeling a realistic SEP event needs 1) an in-depth understanding of the in-fluence of preceding CMEs, in particular the effect of a non-Parker field, and 2) taking into account ofpossible lateral expansion. Indeed, observations have shown that large lateral expansion of CME-drivenshock is quite frequent in large eruptions, so considering the lateral expansion and capturing a realisticCME-driven shock profile, such as the shock speed and obliquity along the shock surface is also crucial forany SEP modelings. We point out that lateral expansion are intrinsically inhomogeneous because the plasmaconditions at different part of the shock flank are different. In this work, we use a second Gaussian profileto mimic the lateral expansion. Finally, fitting observations at multiple longitudes can provide very strongconstraints on the perpendicular diffusion coefficients of energetic particles. Although our simulation is forthe 10 September 2017 event, it forms a nice basis for understanding other large SEP events as well. odeling solar energetic particle event of 10 September 2017 13 Acknowledgements
References
Band, D., Matteson, J., Ford, L., et al. 1993, ApJ, 413, 281Bruno, A. 2017, Space Weather, 15, 1191Bruno, A., Christian, E. R., de Nolfo, G. A., Richardson, I. G., & Ryan, J. M. 2019, Space Weather, 17, 419Cane, H. V., McGuire, R. E., & von Rosenvinge, T. T. 1986, ApJ, 301, 448Cohen, C. M. S., & Mewaldt, R. A. 2018, Space Weather, 16, 1616Desai, M., & Giacalone, J. 2016, Living Reviews in Solar Physics, 13, 3Ding, L.-G., Li, G., Dong, L.-H., et al. 2014, Journal of Geophysical Research (Space Physics), 119, 1463Ding, L., Jiang, Y., Zhao, L., & Li, G. 2013, ApJ, 763, 30Drury, L. O. 1983, Reports on Progress in Physics, 46, 973Ehresmann, B., Hassler, D. M., Zeitlin, C., et al. 2018, Geophys. Res. Lett., 45, 5305Ellison, D. C., & Ramaty, R. 1984, Advances in Space Research, 4, 137Forman, M. A., & Drury, L. O. 1983, in International Cosmic Ray Conference, Vol. 2, International CosmicRay Conference, 267Fu, S., Jiang, Y., Airapetian, V., et al. 2019, ApJ, 878, L36Gopalswamy, N., Nitta, N., Akiyama, S., Makela, P., & Yashiro, S. 2012, ApJ, 744, 72Gopalswamy, N., Yashiro, S., Krucker, S., Stenborg, G., & Howard, R. A. 2004, Journal of GeophysicalResearch (Space Physics), 109, A12105Gopalswamy, N., Yashiro, S., M¨akel¨a, P., et al. 2018, ApJ, 863, L39Gopalswamy, N., Xie, H., M¨akel¨a, P., et al. 2013, Advances in Space Research, 51, 1981Guo, J., Dumbovi´c, M., Wimmer-Schweingruber, R. F., et al. 2018, Space Weath er, 16, 1156Hassler, D. M., Zeitlin, C., Wimmer-Schweingruber, R. F., et al. 2012, Space Sci. Rev., 170, 503Heras, A. M., Sanahuja, B., Lario, D., et al. 1995, ApJ, 445, 497Hu, H., Liu, Y. D., Zhu, B., et al. 2019, ApJ, 878, 106Hu, J., Li, G., Ao, X., Zank, G. P., & Verkhoglyadova, O. 2017, Journal of Geophysical Research (SpacePhysics), 122, 10,938Hu, J., Li, G., Fu, S., Zank, G., & Ao, X. 2018, ApJ, 854, L19Kahler, S. W. 1996, in American Institute of Physics Conference Series, Vol. 374, American Institute ofPhysics Conference Series, ed. R. Ramaty, N. Mandzhavidze, & X.-M. Hua, 61Kahler, S. W., Reames, D. V., & Burkepile, J. T. 2000, in Astronomical Society of the Pacific ConferenceSeries, Vol. 206, High Energy Solar Physics Workshop - Anticipating HESSI, ASP Conference Series,
Vol. 206., 468Kallenrode, M.-B. 1997, J. Geophys. Res., 102, 22347Kong, F. J., Qin, G., Wu, S. S., et al. 2019, ApJ, 877, 97Kwon, R.-Y., Zhang, J., & Vourlidas, A. 2015, ApJ, 799, L29Lee, C. O., Jakosky, B. M., Luhmann, J. G., et al. 2018, Geophys. Res. Lett., 45, 8871Li, G., Moore, R., Mewaldt, R. A., Zhao, L., & Labrador, A. W. 2012, Space Sci. Rev., 171, 141Li, G., & Zank, G. P. 2005, in International Cosmic Ray Conference, Vol. 1, 29th International Cosmic RayConference (ICRC29), Volume 1, 173Li, G., Zank, G. P., & Rice, W. K. M. 2003, Journal of Geophysical Research (Space Physics), 108, 1082Li, G., Zank, G. P., & Rice, W. K. M. 2005, Journal of Geophysical Research (Space Physics), 110, A06104Liu, W., Jin, M., Downs, C., et al. 2018, ApJ, 864, L24Liu, Y. D., Hu, H., Zhu, B., Luhmann, J. G., & Vourlidas, A. 2017, ApJ, 834, 158Liu, Y. D., Zhu, B., & Zhao, X. 2019, ApJ, 871, 8Liu, Y., Luhmann, J. G., Bale, S. D., & Lin, R. P. 2009, ApJ, 691, L151Luhmann, J. G., Ledvina, S. A., Krauss-Varban, D., Odstrcil, D., & Riley, P. 2007, Advances in SpaceResearch, 40, 295Luhmann, J. G., Ledvina, S. A., Odstrcil, D., et al. 2010, Advances in Space Research, 46, 1Luhmann, J. G., Mays, M. L., Odstrcil, D., et al. 2017, Space Weather, 15, 934Luhmann, J. G., Mays, M. L., Li, Y., et al. 2018, Space Weather, 16, 557Mason, G. M., Mazur, J. E., & Dwyer, J. R. 1999, ApJ, 525, L133Mishev, A., Usoskin, I., Raukunen, O., et al. 2018, Sol. Phys., 293, 136Reames, D. V. 1995, Reviews of Geophysics, 33, 585Reames, D. V. 1999, Space Sci. Rev., 90, 413Rice, W. K. M., Zank, G. P., & Li, G. 2003, Journal of Geophysical Research (Space Physics), 108, 1369Sandberg, I., Jiggens, P., Heynderickx, D., & Daglis, I. A. 2014, Geophys. Res. Lett., 41, 4435Seaton, D. B., & Darnel, J. M. 2018, ApJ, 852, L9Shalchi, A., Li, G., & Zank, G. P. 2010, Ap&SS, 325, 99Verkhoglyadova, O. P., Li, G., Zank, G. P., Hu, Q., & Mewaldt, R. A. 2009, ApJ, 693, 894Verkhoglyadova, O. P., Li, G., Zank, G. P., et al. 2010, Journal of Geophysical Research (Space Physics),115, A12103Warren, H. P., Brooks, D. H., Ugarte-Urra, I., et al. 2018, ApJ, 854, 122Werner, A. L. E., Yordanova, E., Dimmock, A. P., & Temmer, M. 2019, Space Weather, 17, 357Yan, X. L., Yang, L. H., Xue, Z. K., et al. 2018, ApJ, 853, L18Zank, G. P., Rice, W. K. M., & Wu, C. C. 2000, J. Geophys. Res., 105, 25079Zeitlin, C., Hassler, D. M., Guo, J., et al. 2018, Geophys. Res. Lett., 45, 5845Zhu, B., Liu, Y. D., Kwon, R.-Y., & Wang, R. 2018, ApJ, 865, 138 odeling solar energetic particle event of 10 September 2017 15
Zuo, P., Zhang, M., Gamayunov, K., Rassoul, H., & Luo, X. 2011, ApJ, 738, 168 a b
Fig. 1: Snapshots of the CME-driven shock configuration. The color scheme is for the normalized density nr . The simulation domain is from . to AU. The bold black circle is r = 1 AU. Planets Mercury,Venus, Earth and Mars and the spacecraft STEREO-A are marked as red dots, with the correspondingdashed lines the field lines passing them. The left panel shows a symmetric eruption with a single Gaussianvelocity profile and the right panel shows the lateral expansion of CME-driven shock with a two-Gaussianprofiles. ba dc
Fig. 2: The time evolution of the shock parameters along the shock front. The black solid curves indicateshock fronts at different times (labelled on the left part of the x-axis). The white dash lines represent theunperturbed Parker magnetic field lines connecting observers at longitudes of 0 ◦ , 10 ◦ , 20 ◦ . Earth is locatedat φ = 0 ◦ . The color schemes for panel (a) is the shock obliquity angle θ BN (the angle between shocknormal and the upstream magnetic field); for panel (b) the shock speed; for panel (c) the maximum particleenergy in the shock front; and for panel (d) the shock compression ratio. odeling solar energetic particle event of 10 September 2017 17 (a) (b) (d) (c) Fig. 3: panel (a) shows the proton time intensity profile measured by the GOES-15 from onset to shockarrival. Other panels show the simulation results before shock arrival at different longitudes ( 0 ◦ , 20 ◦ ). Thecurves of different colors correspond to different particle energies. The x-axis of all panels are the timesince the eruption, normalized to T n , where T n is the shock arrival time for different observers. Fig. 4: Time interval-integrated proton fluence. Left panel: The yellow solid dots represent the fluences fromthe instrument EPAM-LEMS120 on board the ACE. The red solid triangle dots and black solid triangle dotsrepresent the fluences with the instrument EPEAD-A and HEPAD on board the GOES-15. The blue curveis the fit to the observation using the band function. The fitting parameters are shown in blue. The redcurve is the modelled fluence at longitude of 20 ◦ with the red parameters from the band function fitting.Right panel: The modelled fluences at three different longitudes (0 ◦ , 10 ◦ , 20 ◦ ) and the fitting parametersare labelled with corresponding colors. odeling solar energetic particle event of 10 September 2017 19 Fig. 5: The proton time intensity profiles from the observation and simulation at Mars. The black curveshows the integral flux in the energy bands >
113 MeV of the instrument RAD on the MSL (Adapted fromZeitlin et al. 2018). The red curve shows the modelled result in the energy bands >
119 MeV.
Fig. 6: The proton time intensity profiles of the observation and simulation observed by STEREO-A. Thered line represents the flux in the energy bands . − MeV of the instrument HET on board theSTEREO-A. The blue line is the modelled result in the energy bands . − .3