Main text

Introduction

The radiation from the sun is the most important ionization and heating energy source for the Earth’s upper atmosphere. When solar flares occur, the intensity of multi-wavelength electromagnetic emissions suddenly increases (e.g., Kane 1974). Among these, the intensification of X-rays (< 10 nm) to extreme ultraviolet (EUV; 10–121 nm) emissions accelerates the ionization and molecular dissociation of atmospheric components in the ionosphere and thermosphere, and it may cause a sharp increase in electron density (e.g., Qian et al. 2011). It is generally considered that ultraviolet (UV) emissions with a long wavelength of 100 nm or more affect the lower atmosphere, whereas short-wavelength EUVs and X-rays may also affect the mesosphere (e.g., Woods et al. 2000). Moreover, these short-wavelength UV radiation reaches the ionosphere D region at 60–90 km, thereby causing an increase in electron density in this region. A rapid increase in the total electron content (TEC) may degrade the accuracy of single-frequency global navigation satellite system (GNSS) positioning, which has been widely used in recent years. In an air navigation system, GNSS positioning is critical for navigation and safety. As augmentation systems, the satellite-based augmentation system, which covers a wide area by geostationary satellites, and ground-based augmentation system have been developed and are in operation. Temporal and spatial variations of TEC can be a threat to these systems (International Civil Aviation Organization 2010).

The communication failure caused by the absorption of the short-wave due to the variations in electron density is known as the Dellinger phenomenon (Dellinger 1937). Generally, the Dellinger phenomenon is considered to be caused by the occurrence of M-class (~ 10–5 W/m2 at 1–8 Å band) or higher solar flares and predicted by using the magnitude of the solar flares. However, reportedly, the Dellinger phenomenon sometimes occurred even in during C-class (~ 10–6 W/m2 at 1–8 Å band) flares; however, it did not occur even in during X-class (~ 10–4 W/m2 at 1–8 Å band) flares (e.g., Tao et al. 2020). From these results, it is considered that the flare emissions contributing to the occurrence of the Dellinger phenomenon may also be affected by emissions that are not proportional to the X-ray intensity. It has also been reported that fluctuations of the total electron content (TEC) depends not only on the flare class, but also the occurrence position of solar flare on the solar surface (Qian et al. 2010; Le et al. 2012). The X-ray intensity from solar corona does not change depending on the location of flare. On the other hand, EUV emissions from the lower layer of the solar atmosphere, such as the chromosphere, attenuates as the flare location goes from the center of the solar surface to the limb because the thickness of the passing solar atmosphere changes. From this, it is thought that EUV emission is the main source on TEC fluctuations. In previous studies, it is suggested that 0–450 Å EUV band (Strickland 1995), and/or 0–170 Å bands (Woods et al. 2011) is an important indicator for the EUV irradiance input to the Earth’s upper atmosphere.

In recent years, the understanding of flare triggers and their physical processes has been advanced, and physics-based flare prediction research is rapidly progressing (e.g., Kusano et al. 2020). By applying these physics-based flare prediction methods to the prediction model of space weather, we would like to make one prediction model of the flare emission effects on the space weather in the future. The combination of physics-based predictive research for flare emissions and physics-based predictive research for Earth's atmosphere is considered to be compatible because the parameters to be used are common in the physical process between the two models. On the other hand, the mainstream of research for predicting the impact to the space weather, is basically non-physical-based research. For example, the Flare Irradiance Spectral Model (FISM; Chamberlin et al. 2006, 2007, 2008) is most famous prediction model of solar flare X-ray and EUV emissions, but FISM was constructed empirically using observational X-ray (1–8 Å) data of the Geostationary Operational Environmental Satellites (GOES; Bornmann 1996), EUV data of the Thermosphere Ionosphere Mesosphere Energetics and Dynamics Solar EUV Experiment (TIMED/SEE) (Woods et al. 2005) and so on. Usually, the empirical model well reproduces the observation results, but it is impossible to investigate the cause when a special phenomenon occurs. The advantage of physics-based model research is that we can understand where the problem lies in the model. In this study, we attempted to predict and evaluate the effects from the solar flare on the Earth's atmosphere using physics-based model studies.

First, to reproduce the observed flare emissions, we constructed a new flare emission model based on physical processes (Imada et al. 2015; Kawai et al. 2020). In our proposed model, the physical process of the plasma in the flare loop is reproduced by combining the 1D hydrodynamic calculation using Coordinated Astronomical Numerical Software (CANS; http://www-space.eps.s.u-tokyo.ac.jp/~yokoyama/etc/cans/) 1D package with the CHIANTI atomic database (Dere et al. 1997, 2019). Using this model, we reproduced EUV dynamic spectra for some flare events and compared them with the observed EUV spectra (Nishimoto et al. submitted to Earth, Planets and Space) using the Extreme Ultraviolet Variability (EVE: Woods et al. 2012) onboard the Solar Dynamics Observatory (SDO). We performed our calculations with the estimated loop length from the observed flare ribbon distance (Nishimoto et al. 2020). When we compared the observed SDO/EVE Multiple EUV Grating Spectrograph (MEGS)-A (60–370 Å) spectra with our calculation results, we found that our results clearly reproduced most of the EUV lines during flare.

Then, to examine the effects of the constructed flare irradiance on the Earth’s ionosphere, we use an Earth’s global whole atmosphere and ionosphere model, GAIA (ground-to-topside atmosphere and ionosphere model for aeronomy; Jin et al. 2011). The output of the solar irradiance spectrum from some models is applied to GAIA as an input source for ionospheric photoionization and thermospheric heating, using an ionization cross-section model by Solomon and Qian (2005). The vertical TEC was calculated using the results from GAIA runs and compared with GNSS-TEC observations.

In this study, we attempted to evaluate the flare emissions following a flare event and its effect on the Earth's atmosphere by a physics-based model study. We focused on the X9.3 flare of September 6, 2017, which is the largest solar flare in the solar cycle 24. In Sect. 2, we introduce observational data and the models for the flare spectrum and TEC used in this study. Section 3 presents the results of the GAIA runs, and Sect. 4 discusses and summarizes the results.

Observational data and models

Two X-class flares occurred on AR102673 on September 6, 2017. Among them, we focused on the X9.3 flare that occurred at 11:53 UT. In this section, we present an overview of the observed solar flare emission spectra by SDO/EVE and the model calculation for its reproduction. Then we also present the observation and model calculation of the TEC fluctuation associated with the X9.3 flare.

Flare emission spectra and their prediction model

The flare emissions from the X9.3 flare have been observed in X-rays (0.5–4 and 1–8 Å) by the GOES, and SDO/EVE MEGS-B observed the EUV emissions (360–1060 Å). However, emissions within the range of 60–370 Å, which have a large effect on the Earth's atmosphere, have not been observed due to a failure of the SDO/EVE MEGS-A. Therefore, we estimated the flare emission spectrum using the method of Kawai et al. (2020).

In our model (Kawai et al. 2020), we assumed a semicircular flare loop with a diameter that is the distance between two flare ribbons observed in this flare, for simulated plasma motion in the flare loop with the 1D-hydrodynamic model. Main input parameters for CANS 1D are the half loop length and heating rate. For the X9.3 flare, the distance between the two ribbons can be estimated as 32.5 Mm from SDO/AIA 1600 Å data (http://jsoc.stanford.edu/ ajax/lookdata.html), and then the loop length can be obtained as 51 Mm. The heating rate is 120 erg cm−3 s−1, which is estimated by comparing observed and estimated two channel of GOES X-ray peak intensities. The input heating duration and width along the loop of the flare are fixed at 60 s and 6000 km, respectively. These parameters are same as those used for the X9.3 flare described in Kawai et al. (2020). Then, we obtained the EUV emissions strength and its time variation using CHIANTI. We calculated only for the X9.3 class flare in our model, and the result is shown in Fig. 1 depicted by the blue line during 11:53–12:24 UT. Hereafter, this model will be called the “1D Hydro flare model” in this paper.

Fig. 1
figure 1

Light curve of soft X-ray (1–10 Å) emission. “FISM flare” (red line) is proportional to the light curve of GOES soft X-ray (1–8 Å). “bgd” is the daily component of FISM which is a constant value. “1D Hydro flare” reproduces only the light curve of the X9.3 flare during 11:53–12:24 UT. “1D Hydro + bgd” (blue line) is the sum of the light curves of “bgd” and “1D Hydro flare”. Gray regions have no 1D Hydro calculation data

Another prediction model for flare emissions is the FISM (Chamberlin et al. 2006, 2007, 2008). The FISM is an empirical model that provides a spectrum of solar flare at a wavelength range of 1–1900 Å every minute with a wavelength resolution of 10 Å. The calculation results for the FISM are available on the Web site (https://lasp.colorado.edu/lisird/). Recently (October 2020), FISM2 with improved wavelength resolution (0.1–1900 Å at 1 Å spectral bins) was released (Chamberlin et al. 2020), however this research is verifying using FISM. FISM assumes that flare emission time evolution is proportional to GOES soft X-ray, and thus, the FISM time evolution X-ray profile of 1–10 Å X-ray, shown in Fig. 1, is similar to that of the GOES X-ray profile of 1–8 Å. This is also the case for EUV line emissions, FISM has not been able to adequately reproduce the observed delay in EUV line emissions (Nishimoto et al. 2020, submitted to Earth, Planets and Space).

Figure 2 shows the flare spectrum at 12 UT on September 6, 2017, for each emission model. The red broken line in Fig. 2 is the daily component of FISM, and this spectrum can be considered as background emissions. The red line is the FISM flare spectrum, and the green line is the 1D Hydro flare model spectrum. As the 1D Hydro flare model does not consider the continuum component, it can be observed that there is no continuum component at approximately 75–95 nm. The blue line is the sum of the 1D Hydro flare model spectrum and the FISM background spectrum.

Fig. 2
figure 2

Emission spectra of X9.3 flare at 12 UT. The red solid line shows the FISM spectrum and the green solid line shows the 1D Hydro flare model spectrum. The blue solid line is the sum of the 1D Hydro flare model spectrum (green solid line) and the background (red dashed line). Because the 1D Hydro flare model does not consider the continuum component, there is no continuum component around 70–90 nm

The flare spectra variations, which comprise Fig. 2 spectra and Fig. 1 light curves, were input to the GAIA model, which is described in Sect. 2.2. There are four types of input models: (1) FISM background (“bgd”); (2) FISM flare (“FISM”); (3) 1D Hydro flare or background (“1D Hydro/bgd”), and (4) 1D Hydro flare and background (“1D Hydro + bgd”).

TEC observation and GAIA

The TEC in the Earth’s ionosphere can be obtained from the radio propagation delay between GNSS satellites and receivers (Otsuka et al. 2002; Tsugawa et al. 2018; Shinbori et al. 2020). Figure 3 shows the global distribution of verticalized TEC processed, and provided by National Institute of Information and Communications Technology (NICT) and Nagoya University. This is the TEC distribution observed at 12 UT on September 6, 2017, when the X9.3 flare occurred. As the flare occurred around 12 UT, we can observe the variations in TEC distribution due to the flare emissions using the TEC data around 0°E longitude. In the European region, sufficient data were obtained around 0°E longitude, and therefore, we used 0°E, 25–75°N data for comparative studies.

Fig. 3
figure 3

Distribution of TEC observed by global TEC (Drawing-TEC) at X9.3 class flare (12 UT). Distribution of TEC observed by global TEC (Drawing-TEC) at X9.3 class flare (12 UT). The European region was in the daytime at the time of the flare, and has sufficient data

Figure 4 shows the time variation of the TEC distribution at 0°E and 25–75°N. The upper left part of Fig. 4 shows the average TEC during the time period from September 1 to 10, 2017, which can be used as the background TEC. The upper right of Fig. 4 shows the TEC variation observed with the X9.3 flare, and evidently, the TEC sharply increases at any observed latitude at approximately 12UT. The background data in the upper left panel in Fig. 4 show that the amount of TEC variation increases in the afternoon, whereas the flare effect on TEC is significantly stronger than the diurnal change and is approximately 1.3 times stronger at the flare peak (Fig. 4 bottom).

Fig. 4
figure 4

TEC observations at 0°E 25–75°N. (Top left) Averaged TEC from September 1 to 10, 2017. This is the background. (Top right) TEC on September 6, 2017. (Bottom) Ratio of TEC from flare (top-right figure) and the background (top-left figure). TEC significantly increases at any latitude when flare occurs

Although the actual variations in TEC due to flare emissions are displayed in Fig. 4, it is not evident which wavelength of the flare emissions causes the TEC variation. Therefore, to decipher this, we used the GAIA model which can calculate the ionospheric variation by inputting the flare emissions (discussed in Sect. 2.1).

Flare emission ionization rates in the upper atmospheres depend on the ionization cross sections and neutral species abundance. Figure 5 shows the cross sections of photoionization and photo-dissociative ionization of oxygen atoms, oxygen molecules, and nitrogen molecules, which are the main constituents in the Earth’s thermosphere. These cross sections are based on the model of Solomon and Qian (2005). Notably, the X-ray cross sections are two orders or more smaller than the cross sections of > 15 nm EUVs. In addition, evidently, the cross sections of > 60 nm EUVs are considerably smaller depending on the ionization reaction.

Fig. 5
figure 5

Photoionization and photo-dissociative ionization cross sections of major solar atmosphere atoms and molecules. Photoionization and photo-dissociative ionization cross sections for each wavelength. (Left) Photoionization cross section of oxygen atom, (middle) cross section of photoionization and photo-dissociative ionization of oxygen molecule, and (right) photoionization and photodissociation of nitrogen molecule. Cross-section data are based on Solomon and Qian (2005)

Figure 6 shows the change in ionization rate at each altitude as a result of inputting the flare spectrum of the FISM model (Fig. 2) into GAIA. From this result, the wavelength ionized at each altitude in the upper atmosphere can be identified. X-rays mainly ionize oxygen and nitrogen atoms at 100–150 km; approximately 10 nm EUVs mainly ionize nitrogen atoms at similar altitude, and ~ 30 nm EUVs mainly ionize nitrogen molecules at a slightly higher altitude. Approximately 100 nm EUVs also ionizes oxygen molecules at an altitude of ~ 100 km.

Fig. 6
figure 6

Difference in ionization rate due to X9.3 flare that occurred in September 6, 2017. The difference in ionization rate of each ion at 0°E and 30°N due to X9.3 flare on September 6, 2017, including ionization by solar irradiance and secondary photoelectrons. The FISM flare spectrum is used for the input spectrum to GAIA. Differences in ionization rate before and after flare for a O+, b O2+, c N+, d N2+, and e total ion as a function of the wavelength of solar irradiance (unit: Log10(W m−2 nm−1)). f Differences in the altitude distribution of the ionization rate of each ion. Each solid line starts immediately after the flare (12:00 UT), and each dotted line starts before the flare (11:30 UT)

Results

TEC obtained from GAIA simulations

As shown in Fig. 7, the four types of solar flare spectra presented in Sect. 2.1 were input to the GAIA model, and the TEC time variation was obtained. Run1 inputs FISM background spectrum, which is referred to as background; Run1, 2, 3 and 4 input the (1) “bgd”, (2) “FISM flare”, (3) “1D Hydro/bgd”, and (4) “1D Hydro + bgd” spectral models, respectively. Evidently, TEC significantly increased in Run2, 3, and 4 during the X9.3 flare.

Fig. 7
figure 7

TEC distribution at 0°E 25–75°N calculated by GAIA. TEC distribution derived by inputting (upper left: Run1) background, (lower left: Run2) “FISM”, (upper right: Run3) “1D Hydro/bgd”, and (lower right: Run4) “1D Hydro + bgd” spectral changes

Figure 8 shows the calculation results ratio between Run1 and Run2, 3, and 4. The increase in TEC due to flare emissions can be confirmed in all model calculations, and the increase rate was the largest in Run4, whereas the smallest in Run3. As the FISM flare model was input to Run2, a second TEC increase was observed after 12:30 UT.

Fig. 8
figure 8

Ratio of TEC at 0°E 25–75°N calculated by GAIA and the background. Upper left: background TEC distribution. Same as the top-left figure in Fig. 7. Lower left: ratio of TEC for the “FISM” input (Run2/Run1). Upper right: ratio of TEC for the “1D Hydro/bgd” input (Run3/Run1). Lower right: ratio of TEC for the “1D Hydro + bgd” input (Run4/Run1)

Comparison with observed TEC data

The time variation of TEC obtained from GAIA (Sect. 3.1) was compared with the observed TEC, which is shown in Fig. 4. Figure 9 shows the TEC time variation at 0°E and 44°N. These variations do not change much even if the latitude and longitude are changed slightly, and the increase in TEC decreases as the latitude increases due to the difference in zenith angle. In the observed data, there is another prominent TEC increase from 12:30 UT after the peak at 12 UT. This is the effect of another flare event after 12:30 UT, which is seen in the X-ray light curve in Fig. 1 and was captured by the Run2 GAIA calculation. In the case of the 1D Hydro flare model input (Run3 and 4), the background emission model was input after 12:30 UT, and therefore, the second increase is not observed. In the GAIA model calculation, a TEC increase was also observed at 15–20 UT; however, this was not noticed in the observation data and therefore an artifact of the GAIA model.

Fig. 9
figure 9

Variation in TEC at 0°E and 44°N. The solid black line is the observed value (top-right panel in Fig. 4), and the dashed black line is the background of the observation (top-left panel in Fig. 4). The red dashed line, red solid line, green line, and blue line represent the calculation results of GAIA with background input (Fig. 7: Run1), “FISM” input (Fig. 7: Run2), “1D Hydro/bgd” input (Fig. 7: Run3), and “1D Hydro + bgd” input (Fig. 7: Run4), respectively. Gray regions have no 1D Hydro model data

As the absolute TEC value calculated by GAIA does not match the observed value, we decided to compare each increment from the TEC background. The difference in the observed TEC from the background, and the difference in calculated TEC by GAIA from the FISM background is shown in Fig. 10. The error of the observed TEC value is about 0.01 TECU. The absolute value of the observation is about 1–3 TECU (Otsuka et al. 2002), but these effects can be eliminated by subtracting background value. The FISM model accurately reproduces the observed TEC, although it is slightly larger than the observed TEC due to the effects of previous X-class flare. The “1D Hydro/bgd” model and the “1D Hydro + bgd” model was 0.8 times and 2.4 times larger than the observed TEC, respectively.

Fig. 10
figure 10

Difference in TEC from background at 0°E 44°N. The black line is the difference between the observed TEC and the background. The red, green, and blue lines represent the difference in the calculation results of GAIA with “FISM” input, ‘‘1D Hydro/bgd” input, and “1D Hydro + bgd” input from the background, respectively. Gray regions have no 1D Hydro model data

Discussion and summary

We attempted to calculate TEC variations by inputting the temporal variation of the solar flare spectrum for the X9.3 flare into the GAIA model and compared it with the observed TEC. In particular, we analyzed the increments of TEC from background TEC at 0°E and 44°N, and the FISM model almost reproduced the observed values; the “1D Hydro/bgd” and “1D Hydro + bgd” models were 0.8 and 2.4 times larger than the observed value, respectively. Since the amount of TEC fluctuation depends on the zenith angle, it increases as the latitude decreases and decreases as the latitude increases, however, there was no latitude-dependent variation ratio of the difference between the flare spectrum models. Therefore, this discrepancy is considered to represent the intensity ratio of the EUV radiation that affects the TEC increment. By comparing this with the solar flare spectra in Fig. 2, we can infer the primary wavelength causing TEC enhancement.

From Fig. 2, X-ray is the flare emission with the highest intensity, which has an intensity of one order or more compared to EUV emissions. However, the fluctuation in TEC is also caused by a combination of the emission spectrum, the photoionization cross section (Fig. 5), and the distribution of atmospheric densities. Furthermore, as shown in Fig. 6e, atmosphere ionization at approximately 100–150 km is primarily attributable to X-rays. However, from Fig. 1, the background intensity of X-rays is two orders of magnitude smaller than the flare peak intensity of X-rays, indicating there is no significant difference between the “1D Hydro/bgd” and “1D Hydro + bgd” models. Figure 2 also shows that the background spectrum below 15 nm is smaller than the flare spectrum by one order or more, thereby indicating that there is no difference in the X-ray intensity between the “1D Hydro/bgd” and “1D Hydro + bgd” models. In fact, the total emission energy of < 10 nm was 1.3 times that of FISM in both models, with no difference between the presence and absence of background. This relationship implies that X-ray emissions are not the main wavelength causing TEC enhancement.

However, large spectral differences depending on the type of flare emission model were observed for long wavelengths of > 75 nm. In particular, because the 1D Hydro flare model does not consider the continuum component, the EUV flux of the “1D Hydro/bgd” model is two orders of magnitude smaller than that of the FISM model. In this wavelength band, the intensity of the “1D Hydro + bgd” model is approximately 1/3 times that of the FISM model. However, the TEC in the “1D Hydro + bgd” model at the flare peak at 0°E and 44°N is > 2 times that of the FISM model, and therefore, the intensity ratio is not consistent with that of the EUV emissions. This indicates that long-wavelength EUV radiation (> 75 nm) does not affect the TEC.

The wavelengths that become stronger in the order of “1D Hydro/bgd”, “FISM”, and “1D Hydro + bgd” are approximately 15–55 nm emission. Among this wavelength range, the EUV intensity at > 40 nm is about an order of magnitude smaller than < 40 nm (Fig. 2), and at > 35 nm, cross section of photodissociation of nitrogen molecule is an order of magnitude smaller than that at < 35 nm (Fig. 5), there is not enough ion production rate to fluctuate the TEC (Fig. 6e). Therefore, 15–35 nm is a candidate for radiation that is thought to be primarily effective in TEC fluctuation. Especially, at 15–20 nm wavelength band, the total emission energy of the “1D Hydro/bgd” and “1D Hydro + bgd” models were 0.8 and 2.1 times larger than the FISM, respectively, which is similar to the increase in each TEC. The EUV intensity in this wavelength band is ≥ one order of magnitude smaller than that of X-rays, whereas approximately one order of magnitude larger than that of long-wavelength EUVs, which indicates that this EUV wavelength mainly affects the TEC of the Earth’s ionosphere. In addition, EUVs of ~ 30 nm wavelength band have a high ionization rate at a wide altitude of 120–300 km. Ionization is neutralized immediately at lower altitudes; however, 120–300 km is a relatively high altitude, and thus, it is likely that ionization is maintained for a considerable amount of time. Therefore, it is considered that EUVs of ~ 15–35 nm are the main cause of TEC increase.

When reproducing the flare spectrum with our model’s 1 Å wavelength resolution, we can identify the dominant EUV line emission in the ~ 15–35 nm wavelength range. The strongest line emission in this wavelength range is the 304 Å He II line. In addition, lines such as Fe XII at 195 Å, He II at 256 Å, Fe XV at 284 Å, and Fe XVI at 335 Å are also stronger lines. Among these lines, the Fe lines come from several mega Kelvin (MK) plasmas, whereas He II is mainly from the chromosphere. This indicates that the flare emission model can also reproduce the TEC variations correctly if it can reproduce the flare spectrum of several MK plasmas and emissions from the chromosphere.

We discussed the TEC variation due to the X9.3 flare that occurred on September 6, 2017. As this analysis was limited to one case and quantified the specific feature of only one flare, we will conduct statistical analysis of multiple flares to obtain global relationship between flare emission spectrum and TEC in future.