Low- and Mid-Latitude Ionospheric Response to the 2013 St. Patrick’s Day Geomagnetic Storm in the American Sector: Global Ionosphere Thermosphere Model Simulation

In this study, the low-and mid-latitude ionospheric response to the main phase of the 2013 St. Patrick’s Day geomagnetic storm in the American sector on the dayside has been investigated using the ground-based measurements and the Global Ionosphere Thermosphere Model (GITM). First, it is found that the observed ionospheric response can be well reproduced by GITM when it is driven by the electric potential and electron precipitation patterns derived from the Assimilative Mapping of Ionospheric Electrodynamics (AMIE) technique. The AMIE-driven GITM simulation also significantly improves the data-model comparison as compared with the simulation driven by the high-latitude empirical models. Second, it is found that the transport process associated with the neutral wind is largely responsible for the observed ionospheric response. Specifically, the traveling atmospheric disturbances (TADs) propagating from the opposite hemisphere play an important role in the formation of the negative storm phase at low and middle latitudes. Third, it is found that the asymmetric negative storm phases occurred at the nominal equatorial ionization anomaly (EIA) peak region in the afternoon sector are mainly attributed to the interaction of the TADs launched in different hemispheres with different phase speeds. Specifically, stronger Joule heating deposited in the northern hemisphere (NH) generates TADs with faster phase speeds than those launched in the southern hemisphere (SH). Consequently, the locations where the TADs originated from the different hemispheres interact are asymmetric about the geomagnetic equator, leading to the formation of asymmetric ionospheric negative storm phases. This study highlights the importance of accurately specifying high-latitude electrodynamic forcings in global I-T simulations and provides a new insight into the cause of the interhemispheric asymmetry phenomena during geomagnetic storms.


INTRODUCTION
Earth's ionosphere and thermosphere (I-T) system can be significantly disturbed during geomagnetic storms when the high-latitude electric field and electron precipitation undergo significant modifications. The topology and intensity of the high-latitude electric field are directly influenced by the orientation of the interplanetary magnetic field (IMF) and solar wind conditions. Sudden changes in the IMF and solar wind could cause the high-latitude electric field to penetrate towards to the equatorial region. A sudden southward turning of the IMF could produce an eastward prompt penetration electric field (PPEF) on the dayside  which reinforces the dayside fountain effect and can significantly contribute to the ionospheric positive storm effect (e.g., an increase of the F-region electron density) at low and middle latitudes (e.g., Tsurutani et al., 2004;Mannucci et al., 2005;Huang et al., 2005;Lei et al., 2008a;Lu et al., 2012). Meanwhile, a tremendous amount of electromagnetic energy from the magnetosphere is deposited into the ionospherethermosphere (I-T) system during geomagnetic storms, and the majority of this energy input is converted to Joule heating which is an important energy source of the I-T system (e.g., Lu et al., 1995;Thayer et al., 1995;Knipp et al., 2004;Deng et al., 2011;Deng et al., 2018;Knipp et al., 2021). Strong heating at high latitudes could generate equatorward-propagating disturbance neutral winds, which push the plasma upward along the magnetic field line and contribute to a positive ionospheric storm phase at low and middle latitudes (e.g., Lin et al., 2005;Lu et al., 2008). In addition, the storm-time disturbance neutral wind also generates a westward electric field on the dayside (Blanc and Richmond, 1980), which could suppress the dayside fountain effect . Joule heating not only modifies the thermospheric circulation but also alters thermospheric neutral temperature, neutral density and neutral composition, and the changes in thermospheric composition can also change the ion and electron densities through chemical processes (e.g., Fuller-Rowell et al., 1994;Burns et al., 1995;Zhang et al., 2003;Cai et al., 2021). Furthermore, impulsive Joule heating could induce large-scale traveling atmospheric disturbances (TADs) and traveling ionospheric disturbances (TIDs) (Richmond and Matsushita, 1975;Bauske and Prölss, 1997;Balthazor and Moffett, 1999;Lee et al., 2004;Shiokawa et al., 2007;Lei et al., 2008b;Bruinsma and Forbes, 2009;Guo et al., 2014;Lu et al., 2016;Ren and Lei, 2017;Lyons et al., 2019;Zhang et al., 2019;Lu et al., 2020;Sheng et al., 2021;Wang et al., 2021;Pham et al., 2022).
This study focuses on the ionospheric response to the 2013 St. Patrick's Day geomagnetic storm, which is an intense storm that has drawn substantial community interests and lots of significant processes have been made (e.g., Kalita et al., 2016;Liu et al., 2016;Lyons et al., 2016;Verkhoglyadova et al., 2016;Yue et al., 2016;Dmitriev et al., 2017;Guerrero et al., 2017;Marsal et al., 2017;Xu et al., 2017;Dang et al., 2019;Ferdousi et al., 2019;Kumar and Kumar, 2019;Lin et al., 2019;Zhai et al., 2020;Amaechi et al., 2021). However, it is still challenging for the community to accurately describe and simulate the ionospheric response, including both positive and negative phases, to the geomagnetic storms. For example, the simulation work shown in Yue et al. (2016) utilized the National Center of Atmospheric Research Thermosphere-Ionosphere-Electrodynamics General Circulation model (NCAR-TIEGCM, Qian et al., 2014) to study the ionospheric response to the storm. In their study, high-latitude empirical models were used to specify the highlatitude electrodynamic forcings in the NCAR-TIEGCM. However, the complex coexistence of the ionospheric positive and negative storm phases during the main storm phase in the American sector were not well reproduced in their simulations, preventing us from a better understanding of such interesting ionospheric response at low and middle latitudes. The datamodel discrepancies in their study may be attributed to the fact that high-latitude empirical models often fail to capture the dynamic spatial and temporal variations of high-latitude electrodynamic forcings (Heelis and Maute, 2020). In this work, whether and to what extent the more realistic highlatitude electrodynamic forcings could help improve the datamodel comparison are assessed. Meanwhile, the physical processes contributing to the ionospheric variations during the 2013 St. Patrick's Day storm are examined.
In this study, the numerical simulations using the Global Ionosphere Thermosphere Model (GITM, Ridley et al., 2006) are performed, driven by two types of high-latitude electrodynamic forcings: the electric and electron precipitation patterns derived from the Assimilative Mapping Ionospheric Electrodynamics (AMIE) procedure (Richmond and Kamide, 1988;Richmond, 1992) and the patterns provided by a recently developed empirical model, Auroral Spectrum and High-Latitude Electric field variabilitY (ASHLEY, Zhu et al., 2021). The simulation results are compared with the groundbased measurements.

MODEL AND DATA DESCRIPTIONS
Global Ionosphere Thermosphere Model GITM is a global model for the Earth's upper atmosphere which self-consistently solves for the density, velocity and temperature of neutrals, ions and electron. Unlike most global models used in the community, GITM relaxes the hydrostatic assumption and allows the propagation of acoustic waves (Deng et al., 2008;Deng and Ridley, 2014;Zhu et al., 2017;Deng et al., 2021). The version of GITM used in this study is the one coupled with the ionospheric electrodynamo solver developed by Maute and Richmond (2017), and the coupling has been described in Zhu et al. (2019). The spatial resolution of GITM is 5°in geographic longitude and 2.5°in geographic latitude and 1/3 scale height in altitude ranging from 100 to 900 km, and the time step is 2 s.

High-Latitude Electrodynamic Forcings in Global Ionosphere Thermosphere Model
In this study, the high-latitude electrodynamic forcings in GITM are specified in two different ways: by an empirical model and by the data assimilative patterns. The empirical model used in this study is the Auroral Spectrum and High-Latitude Electric field variabilitY (ASHLEY) model . Specifically, its electric potential (ASHLEY-E) and electron precipitation (ASHLEY-A) components are utilized, which are developed based on the in-situ ion drift and electron precipitation data from the Defense Meteorological Satellite Program (DMSP) satellites, respectively. Both ASHLEY-A and ASHLEY-E use interplanetary magnetic field (IMF) and solar wind data as inputs.
The AMIE ionospheric electric potential and electron precipitation patterns are also used for this study. AMIE provides an optimal estimation of high-latitude electrodynamic fields based on a variety of ground-based and space-based measurements, so it provides more realistic high-latitude electrodynamic forcing than empirical models. A more recent description of AMIE can be found in Lu (2017). The data inputs to AMIE for this event include horizontal magnetic perturbations measured by 217 ground stations (among them 173 stations were in the NH and 44 were in the SH), cross-track ion drift data measured by four DMSP satellites (F15, F16, F17, and F18), electron precipitation inferred from the Special Sensor of Ultraviolet Spectrographic Imager (SSUSI) onboard three DMSP satellites (i.e., F16, F17, and F18), line-of-sight ion drifts measured by Super Dual Auroral Radar Network high frequency radar network, and the horizontal magnetic perturbations measured by the Iridium satellite constellation and provided by the Active Magnetosphere and Planetary Electrodynamics Response Experiment (AMPERE) dataset . The temporal resolution of AMIE pattern is 5 min, and the spatial resolution of AMIE pattern is 0.67 h in magnetic local time (MLT) and 1.67°in magnetic latitude (MLAT).

Ionospheric Dataset
In this study, the vertical total electron content (TEC) data estimated from the ground-based dual-frequency global navigation satellite system (GNSS) receivers are used to validate the simulation results. First, the line-of-sight TEC (i.e., line-integrated ionospheric electron density) is estimated by analyzing the processed L1 and L2 pseudorange and phase data (Rideout and Coster, 2006). After removing the ground receiver and satellite biases, the line-of-sight TEC with an elevation angle above 30°is then converted to the vertical TEC using a mapping function (Rideout and Coster, 2006;Vierinen et al., 2016). This study focuses on the vertical TEC data derived from the GNSS receivers distributed in the American Sector are used.
In addition to the TEC data, ionosonde measurements from four stations in the American sector are analyzed. The four stations are: Wallop Island (75.5°W, 37.9°N), Eglin (86.5°W, 30.5°N), Ramey (67.1°W, 18.5°N), and Jicamarca (76.8°W, 12.0°S). In particular, the F2 peak electron density (NmF2) and height (hmF2), which are determined from the measured ionograms using the SAO-Explorer software package (Huang and Reinisch, 1996), are used. Figure 1 summarizes the geophysical conditions during the 2013 St. Patrick's Day storm, including the variations of the IMF B y and B z components, the Dst and AE indices, along with the crosspolar-cap potential (CPCP) and hemispheric power (HP) in the northern hemisphere (NH). The latter four variables are from the AMIE outputs. The AMIE Dst index is similar to the SYM-H index but using the magnetometer data from 56 stations between −40°and 40°MLAT. The AE index is derived from 98 stations between |55°| and |76°| MLAT in both hemispheres. The northern hemisphere (NH) AMIE outputs are used. The 2013 St. Patrick's Day storm was triggered by an interplanetary coronal mass ejection (ICME). As illustrated in Verkhoglyadova et al. (2016), the shock of the ICME arrived at the Earth's bow shock at around 6 UT on 03/17/2013 and was followed by the ICME sheath region. A magnetic cloud arrived at around 15:30 UT and lasted for more than 24 h. During the passage of the sheath region, the IMF B y and B z components intensified and exhibited large oscillations. Meanwhile, the Dst index dropped to−120 nT at around 12 UT and gradually recovered until the arrival of the magnetic cloud. During the passage of the magnetic cloud, the IMF B z turned southward and remained southward for about 6 h. The Dst index decreased again until 20:30 UT when it  Figure 2 shows the comparisons of the storm-quiet differential TEC in the American sector at different UTs, with the GNSS data shown in the top row and by simulation results shown in the second and third rows. The TEC data on 03/16 is used as the quiet-time reference. We also tried several other quiet-time FIGURE 2 | Comparisons of the storm-quiet differential TEC in the American sector at different UTs. At each UT, from top to bottom, the results correspond to those calculated from the GNSS data, ASHLEY-driven and AMIE-driven GITM simulation outputs, respectively. Note that the color scale is smaller for the differential TEC calculated from the GITM simulation output. references, e.g., the 5-day averaged TEC prior to 03/17, and the results do not have significant differences.

Data-Model Comparisons in the American Sector
As shown in the GNSS data, a negative storm phase occurs in the equatorial region at 11 UT and extends to higher latitudes afterwards. However, a positive storm phase occurs at around −30°magnetic latitude (MLAT) and lasts for several hours. Around 15 UT, a positive storm phase starts to occur near 30°M LAT and maintains and intensifies throughout the remaining time of the main phase. However, the negative storm phase still remains around the geomagnetic equator until 17 UT when a positive storm phase appears around the geomagnetic equator. Between 17 and 19 UT, negative storm phases occur at the nominal equatorial ionization anomaly (EIA) peak region (±15°MLAT) and are asymmetric about the geomagnetic equator. More specifically, the negative storm phase in the SH is more intense and occurs in a broader region than that in the NH at 17 and 18 UT. At 70°W, the differential TECs at 3°N (~15°MLAT) are +5 TECU and −9 TECU) at 17and 18 UT, respectively, while the differential TECs at 27°S (~15°MLAT) are −5 TECU and −17 TECU at 17 and 18 UT, respectively. However, the negative storm in the NH is more intense than its SH counterpart at 19 UT. The differential TECs are −17 TECU and −10 TECU at 15°and −15°MLAT at 70°W, respectively. At 20 UT, positive storm phases cover almost the entire American sector except for the region poleward of 50°MLAT, where a strong storm-enhanced density (SED) plume cuts through a negative storm phase which has been lasted for several hours.
The second row for each UT shows the ASHLEY-driven GITM simulation results. Note that for GITM simulation results, the color scale is half of that used for the GNSS data since the GITM TEC is roughly half of the GNSS TEC during both quiet and storm times. Several reasons may be responsible for the underestimation of TEC such as no plasmasphere contribution and inaccuracy of the top boundary conditions. Improving the quantitative representation of the electron density in GITM, especially during geomagnetic storms, is an important aspect of our future model development plan. It is clear that there are large discrepancies between the observational and ASHLEY-driven simulation results, particularly the negative storm phase shown in the GNSS data is not well captured by the ASHLEYdriven GITM simulation. We have also utilized the Weimer electric potential model (Weimer, 2005) and Fuller-Rowell and Evans electron precipitation model (Fuller-Rowell and Evans, 1987), which are the default high-latitude empirical models in GITM, to specify high-latitude electrodynamic forcings in GITM. However, that GITM simulation cannot well reproduce the observed feature either. In contrast, the major storm-related features shown in the GNSS data are well captured in the AMIE-driven simulation (third row for each UT) although some discrepancies do exist. For example, the AMIE-driven GITM simulation shows that the negative storm phase occurs in a larger area than that shown in the GNSS TEC data before 14 UT (e.g., −40°~−20°MLAT). In addition, the AMIEdriven GITM simulation shows that the positive storm phase occurs in a larger area than that exhibited in the GNSS TEC data between 15 and 17 UT (e.g., around the geomagnetic equator). Interestingly, the SED plume which cuts through the negative storm phase occurred in the region poleward of 50°MLAT (marked by the black arrow) is roughly reproduced in the AMIE-driven simulation at 19 and 20 UT. However, the SED plume in GITM also occurs at 17 and 18 UT which is not seen in the data. Meanwhile, the simulated SED plume at 19 and 20 UT is not as narrow and strong as the observed one, which may be attributed to the fact that the subauroral plasma stream (SAPS) electric field is not well specified in AMIE  and also to the fact that the GITM grid resolution used in this study is quite coarse. Improving the representation of the SAPS electric field in AMIE would be an interesting future work. The data-model comparisons shown in Figure 2 demonstrate the value of using more realistic high-latitude electrodynamic forcings to investigate the storm-time ionospheric response. In the remaining part of this paper, we will focus on the comparisons of the observations with the AMIE-driven GITM simulation results. Figure 3 shows detailed data-model comparisons at 70°W (LT~UT−5 h, LT = local time), with the GNSS differential TEC shown in the top panel and the GITM differential TEC shown in the bottom panel. The latitudinal range we focused on here is 60°S-40°N in terms of the geographic latitude (GLAT), which roughly corresponds to −50°~50°MLAT. As shown in Figure 3A, before 14 UT (~9 LT), a negative storm phase occurs at almost all latitudes except for the region around 40°S and 20°N, where a positive storm phase occurs. The negative storm phase is pretty well captured by GITM while the positive storm phase occurred around 20°N before 12 UT (~7 LT) and around 40°S between 12 and 14 UT (~7-9 LT) are not well captured in the GITM simulation ( Figure 3B). Between 14 and 16 UT (~9-11 LT), the GNSS data indicates that a positive storm phase occurs from 20°S to 45°S and from 10°N to 20°N, while a weak positive storm phase occurs near the geomagnetic equator. These positive storms are reasonably well captured by GITM except that the positive storm phase at middle latitudes extends to higher latitudes in the NH in GITM simulation than that shown in the GNSS data. After 16 UT Frontiers in Astronomy and Space Sciences | www.frontiersin.org May 2022 | Volume 9 | Article 916739 5 (~11 LT), the GNSS data shows that positive storm phases occur poleward of the 10°N and 40°S, respectively. In addition, a positive storm phase also occurs around the geomagnetic equator (11°S) after 17:30 UT, which is preceded by a negative storm phase between~16:00 and 17:30 UT. Moreover, sizable negative storm phases can be seen around 5°N as well as around 30°S after 17 UT following the positive storm phase, with the one occurred around 30°S having larger latitudinal extension and stronger magnitude in general. Again, GITM simulation captures the main TEC variations shown in the observation. Figure 4 shows the data-model comparisons of the NmF2 and HmF2 at four stations located around 70°W. The red and blue lines represent the parameters on 03/17 and 03/16, respectively. The solid lines represent the parameters from the GITM simulation and the lines with circles represent the measured parameters. In general, GITM qualitatively captures the variations of NmF2 and HmF2 shown in the observation, but there are some quantitative discrepancies between the observation and simulation, which may result from the inaccurate specifications of the O + flux at the top boundary of Frontiers in Astronomy and Space Sciences | www.frontiersin.org May 2022 | Volume 9 | Article 916739 6 the model (Lu et al., 2001). Improving the top boundary condition of GITM is an important aspect of our future model development plan as discussed earlier.

Physical Processes Responsible for the Electron Density Variations
As illustrated in the above subsection, the AMIE-driven GITM simulation reproduces the observed ionospheric response to the 2013 St. Patrick's Day geomagnetic storm reasonably well, providing us an opportunity to explore the physical processes responsible for the ionospheric response during this storm. Figure 5 shows the storm-quiet difference of the electron density, the meridional neutral wind, TEC, NmF2, and HmF2 at five selected latitudes along 70°W. One of them is located at the geomagnetic equator and the other four are located at ±15°MLAT and ±30°MLAT, respectively. The differences of the electron density and meridional neutral wind are shown as a function of UT and altitude while the other three parameters are show as a function of UT. In addition to the storm-time equatorward disturbance neutral winds (negative in NH and positive in SH), there are also poleward disturbance neutral winds at ±15°M LAT and ±30°MLAT, indicating the existence of the TADs. The TEC and ionospheric peak parameters also display quasiperiodic oscillations, which may be associated with the TADs. In particular, the variations of TEC and NmF2 are similar in general, while the correlation between the NmF2 and HmF2 oscillations is quite complicated. The oscillations of NmF2 and HmF2 are generally anti-correlated at ±15°MLAT but they are roughly positively correlated with some phase delay at other latitudes. As discussed in Lu et al. (2001), the correlation between the NmF2 and HmF2 depends on the relative importance of the vertical shear of the meridional neutral wind and the change of recombination near the F2 peak. The negative and positive correlations occur when the former and latter processes dominate, respectively. Figure 5 also shows that the reduction of the electron density roughly follows the poleward disturbance wind at low and middle latitudes, indicating a significant contribution of the meridional wind to the formation of the negative storm phase. Then a question arises: how does the neutral wind contribute to the ionospheric response? To address this question, it is necessary to investigate the evolution of the neutral wind at 70°W. Figure 6 shows the storm-quiet differences of height-integrated Joule heating, meridional neutral wind along with the vertical ion drifts associated with the neutral wind and electric field. From the meridional neutral wind and the vertical drift associated with the neutral wind ( Figures 6B,C), it is clear that several TADs are launched during this geomagnetic storm. Of particular interest are three pairs of TADs highlighted by the thick black dashed lines. The TADs launched in the NH and SH are labeled by the numbers of 1-3 and 1′-3′, respectively. The first pair of TADs (1 and 1′) launch right after the onset of the geomagnetic storm when there is a significant enhancement of Joule heating, traveling at high phase speeds (~710 m/s) and having larger magnitudes (in terms of wind perturbations) than the TADs launched afterwards. At around 10 UT, when the 70°W meridian rotates to the dawn side, they encounter each other at the geomagnetic equator. After that, they cross each other and travel to middle and high latitudes in the opposite hemisphere. The poleward disturbance meridional winds associated with the TADs push the plasma downward along the field line ( Figure 6C) where the recombination rate becomes large and thus cause a reduction of the electron density between 11and 14 UT at low and middle latitudes (Figure 3). This may also explain why the electron density reduction takes place at low latitudes at first and expands to higher latitudes as shown in Figures 2 and 3. The second pair of TADs (2 and 2′) launch around 12-13 UT, which are also able to travel to the opposite hemisphere before being dissipated. However, they are not as strong as the first pair in terms of wind perturbations so they cannot penetrate to latitudes as high as the first pair of TADs do. Moreover, the phase speed of TAD 2 (710 m/s) is~40% higher than that of TAD 2′ (510 m/s). The third pair of TADs (3 and 3′) cause poleward disturbance meridional winds in the hemispheres where they are launched and are associated with the rarefaction of TADs 2 and 2′. The phase speed of TAD 3 (540 m/s) is~13% higher than that of TAD 3′ (460 m/s). The phase speeds of those three pairs of TADs fall in the range of large-scale TADs categorized by Hunsucker (1982), and faster TAD launched in the NH has also reported in previous studies (e.g., Valladares et al., 2009;Lu et al., 2014;Pradipta et al., 2016). TAD 2 interferes with TAD 3′ at Location B (~−20°MLAT) while TAD 2′ interferes with TAD 3 at Location A (~5°MLAT) between 16 and 17 UT, causing strong poleward disturbance meridional winds at both locations. As noted above, the poleward disturbance meridional wind pushes the plasma downward along the magnetic field line, consequently reducing the plasma density later on. Although plasma transport by the neutral wind through neutral-ion coupling can explain the observed ionospheric response, its relative importance with respect to other transport and chemical processes remains unclear. For example, Figure 6D shows that the vertical drift associated with the electric field may also be able to contribute to the negative storm phase occurred between 9 and 14 UT at low and middle latitudes shown in Figure 3. To delineate the relative contributions by the different dynamical and chemical process, following the procedures described in Lei et al. (2008a), a term analysis is applied to the ion continuity equation in GITM. Since the O + is the dominant ion species in the F region and is the only species advected in GITM default setup (Ridley et al., 2006), the term analysis is carried out only for the O + . The continuity equation of O + is written as: here, N is the number density of O + , q is the production, β is the loss coefficient and V is the ion drift vector. On the right-hand side, βN and −∇ · (NV) represent the loss and transport processes, respectively, and the latter is a combination of transports by the E×B drift, by the neutral wind and by the ambipolar diffusion (Lei et al., 2008a;Liu et al., 2016;Lu et al., 2020). Thus, Eq. 1 can be rewritten as: zN zt (Production − Loss) + Trans E×B + Trans wind + Trans diff (2) Figure 7 shows the storm-time differences of the electron density and all the terms listed in Eq. 2 as a function of UT and Frontiers in Astronomy and Space Sciences | www.frontiersin.org May 2022 | Volume 9 | Article 916739 8 altitude at the five locations same as those shown in Figure 5. At the geomagnetic equator (middle column), when compared to the rate of O + density changes shown in the second row, the transport due to the E×B drift plays a more important role before 15 UT and after 19 UT. Transport due to the neutral wind is more dominant between 15 UT and 19 UT. Transport due to the neutral wind increases the electron density at first around 15 UT and reduces the electron density between 16:30 and 17:30 UT and increases the electron density afterward. Although the neutral wind cannot effectively move the plasma upward or downward near the geomagnetic equator ( Figure 6C), the convergence and divergence of meridional neutral wind ( Figure 6B) are mainly responsible for the enhancement and reduction of the neutral wind transport term at the geomagnetic equator shown in Figure 7. This is consistent with the findings by Ren and Lei (2017) who reported that both the electric field and horizontal neutral wind transport could significantly contribute to the nightside HmF2 enhancement at Sao Luis (near geomagnetic equator) during the October 2013 geomagnetic storm. At ±30°M LAT, the main electron density changes take place in the topside ionosphere, which are mainly associated with the neutral wind transport until 19 UT. At ±15°MLAT, a major electron density reduction takes place in the topside ionosphere for the negative storm phase occurred before 14 UT, which is mainly associated with the neutral wind transport. Between 16 and 19 UT, the changes of electron density occur in both topside and bottomside ionosphere. The neutral wind transport term is still the main contributor to the electron density variation in the topside ionosphere while the combination of the chemical process and neural wind transport dominates in the bottomside ionosphere. However, the variations in the chemical process are also related to the TAD (Prölss, 2011). As shown in Figure 5, while the topside electron density decreases between 16 and 17 UT due to the poleward disturbance neutral wind, the bottomside electron density increases significantly, leading to a net enhancement of TEC. The bottomside electron density enhancement is mainly due to the vertical shear of the meridional wind (Lu et al., 2001). After 17 UT, the bottomside electron density undergoes a reduction while the reduction of the topside electron density subsides and even the topside electron density increases. However, the TEC undergoes a reduction in general since it is a measure of the vertical column electron density. Altogether, the neutral wind transport plays an important role affecting the ionospheric electron density variations during this storm, and TADs play an important role in the formation of the negative storm phase.
One noteworthy feature shown in Figure 3 is that the negative storm phase is stronger in the SH than the NH between 17 and 19 UT, which is caused by the stronger downward ion drift associated with the disturbance meridional wind in the SH. As shown in Figure 6B, the disturbance neutral winds that are responsible for the negative storm phase are quite comparable in the NH and SH between 16 and 17 UT, the difference in the vertical ion drift shown in Figure 6C is mainly associated with the location where the different TADs encounter. The location where TADs 3 and 2′ encounter (Location A) is closer to the geomagnetic equator than the location where TADs 2 and 3′ encounter (Location B). Since the inclination angle at the location A is smaller than that at the location B. the effect of vertical ion drift on plasma transport is more pronounced at the location B (Rishbeth and Garriott, 1969). The reason why the location B is further from the geomagnetic equator than the location A is mainly attributed to the fact that TADs 2 and 3 (launched in NH) travel faster than TADs 2′ and 3′ (launched in the SH), and the reason why TADs 2 and 3 travel faster than TADs 2′ and 3′ may lie in the fact that the Joule heating deposited in the NH and SH is asymmetric (Richmond and Matsushita, 1975). As shown in Figure 6A, the Joule heating deposited at 70°W in the NH is generally stronger than that deposited in the SH between 12 and 13 UT (i.e.,[7][8]. The cause of the interhemispheric asymmetry of Joule heating between 12 and 13 UT will be discussed in detail in the next paragraph. If the phase speeds of TADs 2 and 3 are faster or/and the phase speeds of the TADs 2′ and 3′ are slower, the negative phase occurred between 17 and 19 UT at low latitudes may be located even southward, which is the case shown in the observation (Figure 3). If this is the case, opposite storm phases can take place in the magnetic conjugate locations. Therefore, the northward shift of the negative storm phases in the GITM simulation may be caused by the underestimation of the phase speeds of TADs 2 and 3 or/and the overestimation of the phase speeds of TADs 2′ and 3′.
To investigate the cause of the interhemispheric asymmetry of Joule heating between 12 and 13 UT, especially on the dawn side and morning side that is responsible for the TADs shown in Figure 6, the distributions of high-latitude electric potential, electron precipitation and height-integrated Joule heating in the NH and SH are examined, and Figure 8 shows an example at 12:30 UT for illustrations. Firstly, because the offset between the geomagnetic and geographic poles is different in the different hemispheres, the electric potential and electron precipitation are located more toward the dayside in the NH than in the SH. As a result, the NH receives more solar illuminations than the SH in the region where |MLAT| > 50°at 12: 30 UT even this storm occurs very close to the equinox. Secondly, the IMF B y is not negligible (about −5 nT) and is comparable with the IMF B z at 12:30 UT (Figure 1), which is responsible for the asymmetric distribution of the electric potential between the two hemispheres. As discussed in Hong et al. (2021), even if the auroral electron precipitation is symmetric and the geomagnetic field is a pure dipole field, the Joule heating displays a significant interhemispheric asymmetry at the equinox if the IMF B y is not negligible. The stronger Joule heating deposited on the dawn side in the NH also results from the stronger electron precipitation in the NH, which may be associated with the IMF B y and the substorm effects. Finally, it is worth noting that the ground-based measurements in the SH are less abundant as compared with the NH, which may partly contribute to the asymmetric distributions of high-latitude electrodynamic forcings and Joule heating in the different hemispheres shown in Figure 8. However, the satellite measurements in both hemispheres are comparable and abundant, so that the lack of ground-based measurements in the SH may have a smaller effect compared to the first two points on the inter-hemispheric asymmetry of the Joule heating. Altogether, the interhemispheric asymmetry of Joule heating between 12 and 13 UT may mainly result from the differences of the magnetosphere-ionosphere coupling in the different hemispheres and the geomagnetic field configuration.

SUMMARY
In this study, the dayside ionospheric response to the main phase of the 2013 St. Patrick's Day geomagnetic storm in the American sector has been studied by a combination of data analysis and numerical simulation. It is found that the GITM simulation driven by the AMIE electric potential and electron precipitation patterns is able to capture the observed ionospheric response reasonably well. In comparison, the GITM simulation driven by the electric potential and electron precipitation patterns from the empirical model ASHLEY does not reproduce the observations very well. Our efforts highlight the significance of properly specifying the high-latitude electrodynamic forcing in GCMs. The physical processes that are responsible for the observed ionospheric response have also been investigated utilizing the GITM simulation. It is found that the most observed storm-related features shown in the GNSS TEC are largely attributed to the transport process associated with the neutral winds. In particular, the negative storm phase is largely associated with the TADs launched in the opposite hemisphere. A noteworthy phenomenon is the interhemispheric asymmetric negative storm phases occurred around the nominal EIA peak latitudes in the local afternoon, with the one occurred in the SH is stronger than that in the NH. The negative storm phase is a consequence of the interference of the TADs launched in both hemispheres and their hemispheric asymmetry is caused by the different phase speeds of TADs due to the asymmetric Joule heating depositions. Particularly, stronger Joule heating deposited in the NH launches TADs with faster phase speeds than those launched in the SH. As a result, the locations where the TADs interfere in the NH is closer to the geomagnetic equator than that in the SH. Consequently, the TAD can more effectively affect the ionosphere in the SH and cause a stronger reduction of electron density in the SH than in the NH. This study also provides a new insight for how the high-latitude electrodynamic forcing can contribute to the interhemispheric asymmetry of the ionospheric response at low and middle latitudes during geomagnetic storms.

DATA AVAILABILITY STATEMENT
The IMF and solar wind data can be found at the NASA/GSFC's Space Physics Data Facility's OMNIWeb (https://omniweb.gsfc. nasa.gov). The GNSS TEC data is available at the CEDAR Madrigal database (http://cedar.openmadrigal.org/). The NmF2 and HmF2 data are available at the Global Ionospheric Radio Observatory (http://giro.uml.edu). The AMIE outputs along with the simulation outputs are stored at the NCAR GDEX repository, which can be accessed via the link: https://doi.org/10.5065/asf8-3d59.

AUTHOR CONTRIBUTIONS
QZ led the work, carried out the simulation and analyses and wrote the paper. GL created the AMIE electric potential and electron precipitation patterns. GL and YD contributed discussions. All co-authors contributed comments, edits to the manuscript.