Global MHD simulation of magnetospheric response of preliminary impulse to large and sudden enhancement of the solar wind dynamic pressure

Abstract

A sudden increase in the dynamic pressure of solar wind generates a prominent and transient change in ground-based magnetometer records worldwide, which is called a sudden commencement (SC). The magnetic field variation due to an SC at high latitudes shows a bipolar change, which consists of a preliminary impulse (PI) and main impulse (MI). The largest recorded SC had an amplitude of more than 200 nT with a spiky waveform at low latitudes, and the mechanism causing this super SC is unknown. Here, we investigate the cause of the super SC using a newly developed magnetosphere-ionosphere coupling simulation, which enables us to investigate the magnetospheric response to a large increase in the solar wind dynamic pressure. To simulate SCs, the dynamic pressure of the solar wind is increased to 2, 5, 10, and 16 larger than that under the stationary condition, and two different types of dynamic pressure increase are adopted by changing the solar wind density only or the solar wind speed only. It was found that the magnetic field variations of the PI and MI are several times larger and faster for a jump in the speed than for a jump in the density. It is inferred that a solar wind velocity of more than 2500 km/s in the downstream shock, which cannot be directly simulated in this study, would be consistent with the super SC.

Findings

Introduction

Sudden commencement (SC) is defined by a distinct variation in the horizontal component of the magnetic field on the ground, which is caused by a rapid increase in the solar wind dynamic pressure. SC shows a prominent bipolar change at auroral latitudes, in which the first variation is called the preliminary impulse (PI) and the second variation is called the main impulse (MI) (Araki 1994). On the other hand, SC shows a prominent stepwise increase at low and middle latitudes, which is called a disturbance dominant at low latitudes (DL) (Araki 1994). The PI/MI variations are caused by a field-aligned current (FAC) generated in the dayside magnetosphere (Fujita et al. 2003a, b). The DL variation is caused by compression of the magnetosphere. It is well known that the amplitude of a DL is determined by the solar wind dynamic pressure (Siscoe et al. 1968).

These rapid magnetic field variations could give rise to large power blackouts because the current induced by these magnetic field variations, which is called a geomagnetically induced current (GIC), may damage electric power stations. On 13 March 1989, for example, the power blackout of the Hydro-Quebec system in Canada was attributed to an intense geomagnetic storm (Kappenman 2001). On 30 October 2003, a power blackout of a high-voltage transmission system occurred in southern Sweden during a space storm (Pulkkinen et al. 2005). A crucial task in the study of space weather is to predict and estimate the risk from GICs.

Araki et al. (1997) reported an SC with an anomalously large amplitude on 24 March 1991. An amplitude of more than 200 nT and a rise time of 30 s were observed at Kakioka Magnetic Observatory (geomagnetic latitude = 26.6° N). A geosynchronous magnetopause crossing (GMC) event occurred at the same time, where the magnetopause entered inside the geostationary orbit. In addition, Araki (2014) investigated the other SC with an anomalously large amplitude observed at Kakioka for the period 1924 and 2013 and Colaba-Alibag for 1868 to 1967. The anomalous SC had different features from ordinary SC events. The waveform of the low-latitude ground-based magnetometer signal was not stepwise like that of a DL, but it was an unusually spiky waveform similar to the PI/MI variations typically observed at auroral latitudes. The rise time of the anomalous SC of only 30 s was several times shorter than that of ordinary SCs. Araki et al. (2004) investigated the relationship between the amplitude and rise time of a DL using nighttime data at Guam. They showed that the amplitude and rise time of a DL have a negative correlation. They assumed that the rise time is determined by the passage time of the solar wind discontinuity across an effective length in the magnetosphere, which was estimated to be ~30 Re (Re denotes earth radius). If the effective length was also 30 Re in the anomalous SC, the solar wind speed would need to be 6000 km/s to explain the short rise time of 30 s. Such a high-speed solar wind velocity has never been observed and may not be realistic. An alternative explanation would be that the anomalous SC included the magnetic field variation of the PI in addition to the DL variation, even at low latitudes, since the rise time of the PI variation is typically shorter than that of the DL variation.

Fujita et al. (2003a, b) investigated the PI/MI magnetic field variation at auroral latitudes by performing a global MHD simulation based on magnetosphere-ionosphere coupling. They reproduced the PI/MI magnetic field variation via the FAC in the case of an ordinary SC. They clarified that the current system inducing the PI variation is generated by a dynamo, which converts flow kinetic energy into electromagnetic energy in the dayside magnetosphere when the solar wind discontinuity impinges on the magnetopause. The rate of energy conversion by the dynamo is expressed by

$$ \mathbf{J}\cdot \mathbf{E}=\frac{\rho }{2}\frac{\mathrm{d}{V}^2}{\mathrm{d}t}+\mathbf{V}\cdot \nabla p, $$

(1)

where the first term on the right side of the equation indicates the time rate of the change in flow kinetic energy, the second term describes the time rate of the change in thermal energy, and the left side of the equation describes the time rate of the change in the electromagnetic energy. When J ⋅ E has a negative sign, the dynamo converts flow kinetic energy or thermal energy into electromagnetic energy. Fujita et al. (2003a) suggested that the time rate of the change in flow kinetic energy is dominant for the dynamo when it generates the current system of the PI. Therefore, it is natural to expect that the PI variation depends on the solar wind speed and density.

To determine the mechanism controlling the amplitude of the PI and the rise time, we investigated the PI variation in detail by changing the solar wind speed and density by performing a high-resolution global MHD simulation based on magnetosphere-ionosphere coupling. Details of our simulation are described in “Simulation settings.” In “Simulation results,” we show the results for the solar wind parameter dependence of the PI amplitude and the rise time. Finally, in “Summary and discussion,” we summarize the results and discuss their application to extreme space weather events.

Simulation settings

The numerical global MHD model developed by Tanaka (2003) self-consistently solves the solar wind-magnetosphere-ionosphere coupling process (Moriguchi et al. 2008), which is necessary for the investigation of the SC. To achieve high resolution and accurately capture discontinuities, the MHD calculation employs a finite volume (FV) total variation diminishing (TVD) scheme with an unstructured triangular grid system. The number of triangular grid points is 30,722 in the horizontal direction and 240 grids in the radial direction. The outer and inner boundaries of the simulation are set at 200 Re and 3 Re, respectively. Details of the calculation of the inner boundary are given by Tanaka (2000). In this paper, the x-axis points toward the Sun, the z-axis points north, and the y-axis is chosen to satisfy the right-handed system. To examine the magnetic field variations on the ground, only the effect from the ionospheric Hall current is calculated using the Biot-Savart law. The distribution of the Hall current is calculated as the product of the ionospheric Hall conductivity and the electric field projected along the field line from the inner boundary to the ionosphere.

The outer boundary conditions correspond to the solar wind on the upstream side at x = 30 Re and a zero gradient on the downstream side at x = −200 Re. The simulation is started from a stationary state with solar wind parameters of N = 5/cc, V = 372 km/s, By = 2.5 nT, Bz = 4.3 nT, and T = 2 × 10⁵ K, which give the dynamic pressure P _dyn = 1.2 nPa. We applied four levels of dynamic pressure enhancement with a density jump only or a speed jump only as shown in Table 1.

Table 1 Solar wind parameters for the density jump and speed jump conditions used in this study

Simulation results

Global MHD simulations were executed using the eight sets of solar wind parameter settings shown in Table 1. Hereafter, the results for each set of solar wind parameters are denoted run02n, run05n, run10n, and run16n for the density jump conditions and run02v, run05v, run10v, and run16v for the speed jump conditions. Figure 1 shows the simulation results for run10n and run10v. The upper panels in each case show color contours of the FAC in the polar ionosphere. The lower panels in each case show color contours of the static pressure and velocity vectors in the equatorial plane of the magnetosphere. The sunward direction is to the left. Time progresses from left to right in each series. The contour interval is 2.22 μA/m² for the FAC and 41.6 × 10⁻¹¹ N/m² for the pressure. When an impulse with high dynamic pressure impinged on the magnetopause on the dayside, PI signatures, which consisted of a negative and positive pair of FACs, appeared in the dayside ionosphere in each simulation as shown in the leftmost panels. Subsequently, the MI signature, which consisted of an opposite-polarity pair of PI FACs, appeared in the dayside ionosphere in the second panel from the left. The MI signature moved from the dayside to nightside in the ionosphere following the vortices, which moved from the dayside to nightside in the equatorial plane of the magnetosphere as shown in the third panel from the left. These features are consistent with the results of the previous studies by Fujita et al. (2003a, b) and common to both run10n and run10v. On the other hand, there are some differences between these cases. The first is the intensity of the FAC. The intensity of the FAC in run10v is larger than that in run10n owing to the appearance of the PI and MI signatures. The peak intensity of the PI is 16.7 μA/m² in run10v and 5.57 μA/m² in run10n. The second difference is the temporal evolution of the development of the current systems, which is faster in run10v than in run10n. The period between the PI and MI signatures is 1 min in run10v and 2.6 min in run10n. Finally, an FAC with a complex structure appears after the PI and MI signatures in run10v at the time when some vortices are formed in the equatorial plane of the magnetosphere as shown in the rightmost panel in run10v. Figure 2 shows the horizontal component of the magnetic field variation derived from the ionospheric Hall current at three different longitudes (three panels) at four different latitudes (four lines in each panel) in the simulation results for run10n and run10v. One division in the vertical axis is 100 nT. The horizontal axis is time. These magnetic field variations are consistent with the PI and MI signatures of FACs in Fig. 1 and the results of Fujita et al. (2003a, b); Yu and Ridley (2009). The amplitude of magnetic field in run10v is larger than that in run10n owing to the appearance of the PI and MI signatures at any magnetic local time (MLT) and any latitude. Figure 3 shows the horizontal component of the magnetic field variation derived from the ionospheric Hall current in our simulation at 15 magnetic local time and 63° magnetic latitude. The upper panels show the results for the fivefold increase in solar wind dynamic pressure, and the lower panels show the results for the tenfold increase in solar wind dynamic pressure. The left panels show the results for the density jump condition, while the right panels show the results for the speed jump condition. There is a bipolar signature consisting of PI and MI magnetic field variations in each panel. First, we compare run10n with run10v. Although the solar wind dynamic pressure is the same in both cases, the PI magnetic field variation in run10v is more than three times that in run10n. The rise time of the PI in run10v is shorter than that in run10n. These features are consistent with the FAC signature of the PI shown in Fig. 1. To clarify how the solar wind dynamic pressure dependence of the PI variation differs between the density jump and velocity jump conditions, we compare the results for a tenfold increase in dynamic pressure with those for a fivefold increase. The PI magnetic field variation is almost the same in both run10n and run5n, while the variation is about twice as large in run10v as in run5v. The result indicates that the amplitude of the PI is larger for changes in solar wind velocity than for changes in solar wind density even though the dynamic pressure is the same, while the amplitude of the DL depends on the solar wind dynamic pressure only. It is therefore suggested that a high-speed solar wind is needed to create a large and sharp magnetic field variation of the PI.

Fig. 1

Overview of simulation results for density jump condition (run10n) and speed jump condition (run10v) with a tenfold increase in solar wind dynamic pressure. The upper panels of each case show color contours of the field-aligned current (FAC) in the polar ionosphere. The lower panels of each case show color contours of the static pressure and velocity vector unit at the equatorial plane in the magnetosphere. The sunward direction is to the left. The series of figures from left to right represents the time development of the SC. The contour interval is 2.22 μA/m² for the FAC and 41.6 × 10⁻¹¹ N/m² for the pressure

Fig. 2

The horizontal component of the magnetic field variation derived from the ionospheric Hall current at three different longitudes (three panels) at four different latitudes (four lines in each panel) in the simulation results for run10n and run10v. One division in the vertical axis is 100 nT. The horizontal axis is time. MLT magnetic local time

Fig. 3

Horizontal component of magnetic field variation derived from Hall current in the simulation at 15 magnetic local time and 63° magnetic latitude. The upper panels show the results for a fivefold increase, and the lower panels show the results for a tenfold increase. The left panels show the results for the density jump condition, and the right panels show the results for the speed jump condition. PI preliminary impulse, MI main impulse

To determine the solar wind density and velocity dependences of the PI amplitude, we investigate by simulation the PI amplitudes for 2–16-fold increases in the solar wind dynamic pressure for both speed jumps and density jumps as shown in the left panel of Fig. 4. The vertical axis indicates the magnetic field amplitude of the PI. The horizontal axis indicates the square root of the normalized solar wind dynamic pressure. Open triangles and open squares indicate the speed jump and density jump conditions, respectively. It is found that the PI amplitude for the speed jump condition has stronger dependence on the solar wind dynamic pressure than that for the density jump condition. The relative change in the PI amplitude for the speed jump condition is four times as large as that for the density jump condition. The maximum amplitude in our simulation is 500 nT in run16v.

Fig. 4

Dependences of PI amplitude (left) and rise time (right) on solar wind dynamic pressure. The horizontal axes are the square root of the normalized pressure (left) and the rise time (right). The vertical axis indicates the magnetic field amplitude of the PI. Open triangles and open squares indicate speed jump and density jump conditions, respectively

The right panel of Fig. 4 shows the relationship between the amplitude and the rise time of the PI. The vertical axis indicates the magnetic field amplitude of the PI. The horizontal axis indicates the rise time of the PI. The rise times in the speed jump simulations are shorter than those in the density jump simulations. The rise time is determined by the time in which the solar wind discontinuity passes through the effective length of the magnetosphere. The solar wind velocity in the speed jump simulations is higher than that in the density jump simulations. Therefore, the rise times in the speed jump simulations are shorter than those in the density jump simulations. In both the speed jump and density jump simulations, the PI amplitude increases as the PI rise time decreases. This tendency is the same as the relation between the rise time and the amplitude of the DL investigated by Araki et al. (2004). On the other hand, a difference is that the rise time of the PI is shorter than that of the DL. The rise time of the PI is at most 90 s, while the rise time of the DL is at least 1.5–9 min. We can estimate the effective length of the PI from the product of the solar wind velocity and the rise time in our simulation results, which is ~3.5 Re from the speed jump simulations. This indicates that the effective region of the PI is more localized than that of the DL, which has an effective length of 30 Re. The reason for the short PI rise time is that the dynamo generating the FAC is localized in the dayside magnetosphere. We next discuss the mechanism by which the dynamo generates the magnetic field disturbance on the ground via the magnetosphere-ionosphere current system.

The dynamo process in the dayside magnetosphere is investigated using our simulation results. The energy conversion rate by the dynamo is given by Eq. 1. The flow kinetic energy or thermal energy is converted into electromagnetic energy. Figure 5 shows the color contours of each term in the equatorial plane. The electric current vector is superposed on each color contour. The upper three panels show the results for run10v at 3.0 min. The top panel shows the color contour of J ⋅ E, the second panel shows the color contour of ρ/2 ⋅ dv ²/dt, and the third panel shows the color contour of V ⋅ ∇p. The bottom panel shows the result for run10n at 5.4 min, which is the color contour of J ⋅ E. Note that the color scale of J ⋅ E in run10n is one tenth of that in run10v. The contour interval is 8.9 × 10⁻¹¹ W/m³ for run10v and 8.9 × 10⁻¹² W/m³ for run10n. The white arrows indicate the positions of solar wind impulses. We find that J ⋅ E is strongly negative at the center of the top panel and that the electric current is intense in the same region. This region represents the dynamo region because J ⋅ E is negative. When we compare this region in the second panel with that in the third panel, it is found that the ρ/2 ⋅ dv ²/dt term is more negative than the V ⋅ ∇p term. This indicates that J ⋅ E is mainly determined not by the V ⋅ ∇p term but by the ρ/2 ⋅ dv ²/dt term. The dynamo converts the flow kinetic energy of the solar wind into electromagnetic energy through the deceleration of the flow. Next, it is useful to discuss the current in the equatorial plane to clarify the electric current system. The strong current in the dynamo region in the top panel decreases in the dusk direction. This means that the current converges in the equatorial plane. The converging current must flow out from the equatorial plane to the ionosphere through the magnetic field line because the divergence of the current should be zero. The yellow lines show the magnetic field lines drawn in the region with the enhanced PI in the ionosphere. The field lines are connected with the region where the current converges in the equatorial plane. Consequently, we have found that the dynamo current created by flow deceleration in the dayside magnetosphere is connected with the PI current in the ionosphere via the FAC. This is consistent with the results of Fujita et al. (2003a, b).

Fig. 5

Color contour of J ⋅ E in equatorial plane for run10v with run10n. The upper three panels show the results for run10v at 3.0 min in the equatorial plane. The electric current vector is superposed on each color contour. The top panel is the color contour of J ⋅ E, the second panel is the color contour of ρ/2 ⋅ dv ²/dt, and the third panel is the color contour of V ⋅ ∇p. The bottom panel shows the result for run10n at 5.4 min in the equatorial plane, which is the color contour of J ⋅ E. The color scale of J ⋅ E for run10n (bottom panel) is one tenth of that for run10v (top panel). The contour interval is 8.9 × 10⁻¹¹ W/m³ for run10v and 8.9 × 10⁻¹² W/m³ for run10n. The white arrows indicate the positions of the solar wind impulse. The yellow lines show the magnetic field lines drawn in the region with the enhanced PI in the ionosphere. The sunward direction is to the left

To clarify the difference in the dynamo between run10v and run10n, we compare the top panel with the bottom panel. The negative value of J ⋅ E and the current vector in run10v are clearly larger than those in run10n. The maximum value of J ⋅ E is 2.5 × 10⁻¹⁰ W/m³ for run10v and 1.5 × 10⁻¹¹ W/m³ for run10n. The large dynamo in run10v generates a stronger FAC than that in run10n even though the solar wind dynamic pressure is the same. This difference arises from the ρ/2 ⋅ dv ²/dt term. The time change of the velocity is effective for the dynamo under the speed jump condition while the density has a linear relation for the dynamo under the density jump condition. Therefore, the dynamo under the speed jump condition is larger than that under the density jump condition. Figure 6 shows the solar wind dynamic pressure dependence of the FAC for the speed jump and density jump conditions. The horizontal axis is the square root of the normalized pressure. The vertical axis is the peak strength of the FAC for the PI signature in the ionosphere. Open triangles and open squares indicate the speed jump and density jump conditions, respectively. The FAC under the speed jump condition has a stronger dependence on solar wind dynamic pressure than that under the density jump condition. This tendency is consistent with the amplitude of the PI shown in Fig. 4 (left). Moreover the peak of the FAC in run16n is almost saturated at 16 times the dynamic pressure. Although the peaks of the FAC in run10n and run16n have the same strength, the region of the FAC extends to a low latitude in run16n. Therefore, the PI magnetic field variation in run16n is larger than that in run10n. The FAC generated by the dynamo results in the variation of the PI magnetic field depending on the solar wind dynamic pressure.

Fig. 6

Dependence of FAC on solar wind dynamic pressure. The horizontal axis is the square root of the normalized pressure, and the vertical axis is the peak strength of the FAC for the PI signature in the ionosphere. Open triangles and open squares indicate the speed jump and density jump conditions, respectively

Summary and discussion

We have investigated the magnetic field variation of PI at auroral latitudes for several sets of solar wind speed and density jumps by performing global MHD ionosphere-magnetosphere coupling simulations. The new findings of this study are summarized as follows:

• 1) The amplitude of the PI has different values for a density jump and speed jump even though the solar wind dynamic pressure is the same.

• 2) The amplitude of the PI depends on the solar wind density and the speed. This is different from the DL response, which only depends on the solar wind dynamic pressure.

• 3) A high-speed solar wind is needed to create a large-amplitude PI because the dynamo that generates the FAC associated with the PI magnetic field variation is more effective at an enhanced speed than at an enhanced density.

• 4) The amplitude of the PI is larger when the rise time is shorter. This is consistent with the relationship between the rise time and amplitude of the DL, and the rise time of the PI is shorter than that of the DL because the effective length of the PI, i.e., the dynamo region, is shorter than that of the DL. Also, the rise time of the PI is shorter at an enhanced speed than at an enhanced density. It is therefore suggested that a high-speed solar wind is needed to create a rapid magnetic variation.

It is also found that a magnetic field variation similar to the so-called Psc appears after the PI/MI only under the speed jump condition as shown in the bottom right panel of Fig. 3. The oscillation has a period of a few minutes, corresponding to that of the vortices in the dayside equatorial plane in the magnetosphere, as shown in the bottom right panel of Fig. 1. When the high-speed solar wind impinges on the magnetosphere, vortices are repeatedly formed at the equatorial magnetopause, which is probably due to the K-H instability. It seems that the high pressure of the vortices plays an essential role as a current generator to drive the FAC and the magnetic field oscillation. The mechanisms of magnetic field oscillation driven by the dynamic pressure enhancement under the velocity jump condition will be discussed in detail in a separate paper.

We investigated the PI magnetic field variation at auroral latitudes based on the results of our simulation and that at low latitudes based on the results at auroral latitudes, because the validity of our simulation code is limited to the boundary of the ionosphere, corresponding to the inner boundary of the magnetosphere. To overcome this problem directly, we have to modify our simulation code to one with an inner boundary of less than 3 Re. However, considerable time is required to simulate the inner region because the Alfven speed in the model is greater than that in our model. In the future, we will modify our code and investigate the PI dependence on latitude.

Note that the parameters of the speed jump condition do not correspond to those of realistic solar wind because the shock conditions are not satisfied at the discontinuity between upstream and downstream. However, to investigate the effect of a speed jump, it was necessary to distinguish the effect of a density jump from that of a dynamic pressure jump to clarify which parameter should be varied to investigate PI magnetic variation. Here, we discuss a case in which the shock condition is satisfied. The solar wind discontinuity in the simulation must satisfy the shock condition to compare observations. Our simulation code is limited to high solar wind dynamic pressures of at least 16 times the solar wind dynamic pressure because the effect of a high solar wind pressure reaches the inner boundary in the simulation, which is at 3 Re. Therefore, we carried out the simulation while satisfying the shock condition by assuming 16 times the solar wind dynamic pressure, although the solar wind velocity is less than that predicted for the extreme event reported by Araki et al. (1997). The solar wind velocity in the simulation is 744 km/s and the density is 23.2/cc. According to the simulation, the amplitude of the PI is 350 nT. This amplitude is almost the same as that for the speed jump condition with a solar wind velocity of 744 km/s, which does not satisfy the shock condition. This suggests that it is not important for the amplitude and rise time of the PI to satisfy the shock condition. Therefore, the results in this paper may be reliable even though the shock condition is not satisfied.

Finally, we attempt to predict the amplitude of the super SC in the case of extreme solar wind. In this simulation, we investigated the magnetic field variation of the PI up to a solar wind velocity of 1500 km/s. The most well-known extreme event, called the Carrington event, may have had a solar wind velocity of above 2000 km/s (Tsurutani et al. 2003). Although it is not simple to extrapolate the solar wind velocity to obtain a value for an extreme event, the magnetic field variation of the PI was fitted by a linear function using our simulation results as shown in Fig. 7. The horizontal axis is the solar wind velocity, and the vertical axis is the horizontal component of the magnetic field. The solid line indicates the horizontal component of the magnetic field variation at 15 magnetic local time and 63° magnetic latitude, and the dashed line indicates that at 15 magnetic local time and 26° magnetic latitude. When the solar wind velocity is 2500 km/s, the magnetic field variation of the PI at 63° magnetic latitude is above 1000 nT, as shown by the solid line, and the rise time is 8.9 s, which is deduced from the effective length of 3.5 Re. The amplitude and rise time of the anomalously large SC reported by Araki et al. (1997) are 200 nT and 30 s at a latitude of 26.6°, respectively. As shown by the dashed line, the magnetic field variation of the PI at 26° magnetic latitude is 50 nT and the rise time is 8.9 s when the solar wind velocity is 2500 km/s. The results of Araki et al. (2004) suggest that the DL magnetic field variation is about 100 nT when the rise time is 90 s, which means that the solar wind velocity is about 2000 km/s. This is a smaller amplitude and a longer rise time than the values of the anomalously large SC. However, if the DL variation is superposed on the PI variation at a low latitude, we can explain the observed anomalously large SC because the amplitude of the DL is almost the same as that of the PI and the rise time of the DL is longer than that of the PI. Thus, it is inferred that a solar wind velocity of more than 2500 km/s in the downstream shock would be consistent with the super SC.

Fig. 7

Dependence of PI amplitude solar wind velocity. The horizontal axis is the solar wind velocity, and the vertical axis is the horizontal component of the magnetic field. The solid line indicates the horizontal component of the magnetic field variation at 15 magnetic local time and 63° magnetic latitude, and the dashed line indicates that at 15 magnetic local time and 26° magnetic latitude