Alfvénic Thermospheric Upwelling in a Global Geospace Model
HHogan et al. Alfvénic Thermospheric Upwelling 1
Key Points • A parameterized model for Alfv é n wave heating of the thermosphere was imple-mented in a global iono-sphere-thermosphere model • Alfv é n wave heating in-creases air density in the model by 20-30% in and near the cusp at CHAMP altitudes (
413 km) during an SIR event • Heating effects are not local-ized to the cusp during the event and improve model ac-curacy compared to orbit-av-eraged CHAMP-derived den-sity
Contents
1. Introduction 2. Methods 2.1 Overview of CMIT and the 26-27 March 2003 SIR event. 2.2 Parameterized Alfv é n wave model 2.3 Specification of back-ground parameters 2.4 Integrating Alfv é n wave heating into CMIT 3. Results and Discussion 4. Conclusions Acknowledgements References Supporting Information Correspondence
B. Hogan [email protected] W. Lotko [email protected] K. Pham [email protected]
Alfvénic Thermospheric Upwelling in a Global Geospace Model
Benjamin Hogan
William Lotko and Kevin Pham Thayer School of Engineering, Dartmouth College, Hanover, New Hampshire, USA, High Altitude Observatory, National Center for Atmospheric Research, Boulder, CO, USA, † Now at Laboratory for Atmospheric and Space Physics, University of Colorado Boul-der, Boulder, CO, USA and Department of Aerospace Sciences, University of Colorado Boulder, Boulder, CO, USA
Abstract
Motivated by low-altitude cusp observations of small-scale (~ 1 km) field-aligned currents (SSFACs) interpreted as ionospheric Alfvén resonator modes, we have investigated the effects of Alfvén wave energy depo-sition on thermospheric upwelling and the formation of air density enhancements in and near the cusp. Such den-sity enhancements were commonly observed near 400 km altitude by the CHAMP satellite. They are not predicted by empirical thermosphere models, and they are well-correlated with the observed SSFACs. A parameterized model for the altitude dependence of the Alfvén wave electric field, constrained by CHAMP data, has been devel-oped and embedded in the Joule heating module of the National Center for Atmospheric Research (NCAR) Cou-pled Magnetosphere-Ionosphere-Thermosphere (CMIT) model. The CMIT model was then used to simulate the geospace response to an interplanetary stream interaction region (SIR) that swept past Earth on 26-27 March 2003. CMIT diagnostics for the thermospheric mass density at 400 km altitude show: 1) CMIT without Alfvénic Joule heating usually underestimates CHAMP’s orbit-average density; inclusion of Alfvénic heating modestly im-proves CMIT’s orbit-average prediction of the density (by a few %), especially during the more active periods of the SIR event. 2) The improvement in CMIT’s instantaneous density prediction with Alfvénic heating included is more significant (up to 15%) in the vicinity of the cusp heating region, a feature that the MSIS empirical thermo-sphere model misses for this event. Thermospheric density changes of 20-30 % caused by the cusp-region Alf-vénic heating sporadically populate the polar region through the action of co-rotation and neutral winds.
1. Introduction
Thermospheric density anomalies with average values 20-30% larger than empirical model predictions are encountered by satellites orbiting near 400 km altitude in the high-latitude dayside thermosphere (Liu et al., 2005, Rentz & Lühr et al., 2008). About 40% of CHAMP satellite orbits exhibit such anomalies (Kervalishvili & Lühr, 2014). Soft electron precipitation, also referred to as broadband precipitation with energies of 100s eV (Živković et al, 2015), and intense, small -scale field-aligned currents (SSFACs) (Lühr, et al, 2004; Kervalishvili & Lühr, 2013) are well correlated with the anomalies in the cusp region. Enhanced Joule heating by quasistatic electric fields may contribute to the thermospheric upwelling that produces density anomalies, but observa-tions (Schlegel et al., 2005) and models (Zhang et al., 2012; Deng et al., 2013) indicate that this mechanism alone is insufficient. Simulation studies of the thermosphere (Brinkman et al., 2016), the coupled ionosphere-thermosphere (IT) (Deng et al, 2013) and the coupled magnetosphere-ionosphere-thermosphere (MIT) (Zhang et al., 2012, 2015a) show that soft electron precipitation in the cusp, combined with sufficiently intense, quasistatic electric fields, can produce density anomalies resembling those observed in the cusp region. In all of these studies, the soft electron precipitation augments the F -region Pedersen conductivity and the specific Joule heating rate at altitudes where the ambient air density is exponentially decreasing. The added heat and attendant increase in thermospheric temperature is effective in producing upwelling and a major density perturbation of the neutral gas (Jee et al., 2008; Clemmons et al., 2008; Brinkman et al., 2016). Lotko and Zhang (2018) showed that the Joule heating rate resulting from persistent and more or less contin-uously driven, small-scale Alfvénic perturbations representative of the SSFACs regularly observed by CHAMP in the cusp can be significant at F -region altitudes (200-400 km). Transient Alfvénic perturbations accompa-nying discontinuities in solar wind driving are also expected to augment F -region Joule heating (Tu et al., 2011) and contribute to enhancements in thermospheric density at these altitudes, but isolated transients in Alfvénic energy deposition cannot account for the relatively high probability of observing a thermospheric density anomaly on any given satellite pass through the cusp. In contrast with quasistatic electric fields, which are practically independent of altitude in the thermosphere, Alfvénic electric fields are more intense in the F -region than in the E -region owing to the influence of the ionospheric Alfvén resonator (IAR). The IAR traps Alfvén waves and intensifies their amplitudes near the F -region peak in plasma density where the wave speed lfvén Wave Occurrence Chaston et al. < 4200 km O cc u rr e n c e %
18 062412FAST “Broadband” Electron Precipitation
18 062412
Newell et al.
845 kmDMSP e l e c t r o n s / c m - s F || e ↓ “Quasistatic” E-field Variability Matsuo, Richmond 20 ( δ E ) ½ m V / m
18 062412DE-2
Liu et al. %ρ CHAMP / ρ MSIS90 − Air Density Anomalies ° ° ° Hogan et al. Alfvénic Thermospheric Upwelling 2 is slow (Trakhtengerts & Fel'Dshtein, 1984; Lysak, 1991). The results suggest that Alfvén wave energy depo-sition in the IAR may contribute to air density anomalies at CHAMP altitudes in the absence of Joule heating by intense quasistatic electric fields. The SSFACs recorded by CHAMP in association with density anomalies have been interpreted as IAR modes (Rother et al., 2007). However, distinguishing Alfvén waves with characteristic wave impedance and phase relations (e.g., Ishii et al., 1992) from quasistatic variability requires simultaneous electric and magnetic field measurements. CHAMP carried a sensitive magnetometer but no electric field sensor. Measurements from the Polar (Keiling et al., 2003) and FAST (Chaston et al., 2007; Hatch et al., 2017) satellites have been used to definitively identify Alfvénic fluctuations and their distributions in magnetic local time (MLT) and magnetic latitude (MLAT). The statistical distribution of the most probable occurrence of Alfvénic fluctuations from FAST shown in Figure 1 (lower right) correlate well with the MLT-MLAT distribution of CHAMP thermospheric density anomalies (Fig. 1, upper left), as does the distribution of soft electron precipitation derived from DMSP satellite measurements (Fig. 1, upper right). The MLT-MLAT distribution of quasistatic electric field variability is also shown in Fig. 1 (lower left) because it has been posited as a possible causal agent for the enhanced Joule heating needed to produce the CHAMP density anomalies (Deng et al., 2013; Brinkman et al., 2016). This particular distribution derived by Matsuo & Richmond (2008) from DE-2 satellite data for summer conditions and clock angles of 45 ° -135 ° of the inter-planetary magnetic field (IMF) also resembles the distribution of CHAMP anomalies. Distributions of electric field variability for winter and equinox conditions and other IMF clock angles are very different and, for some Figure 1.
Top left: One-year average difference between air density recorded by CHAMP in the northern hemisphere and estimated from MSIS. Illustrative CMIT convection streamlines superposed.Top right: Statistical number flux of broadband electron precipitation from DMSP. Feldstein statistical auroral oval superposed as dotted lines. Bottom left: Statistical root-mean-square electric field variability from DE-2, with same illustrative convection lines. Bottom right: Occurrence of Alfvén waves recorded by FAST, also with Feldstein oval superposed. All figures are in MLT-MLAT (magnetic local time, magnetic latitude) coordinates. Source of plot and nominal altitude of measurements are indicated at the bottom of each.
Supporting Information Text S1. CMIT treatment of IMF and Figure S1 compar-ing from the OMNI data set with the constrained used in the SIR simulation. Text S2. CMIT preconditioning. Text S3. Validity of neglecting neutral gas inertia in collisional Alfv é n wave dynamics and Fig-ure S2 showing the domain of validity for the model iono-sphere and thermosphere. Text S4. Implications of treating the cusp-region geomagnetic magnetic field as straight and uniform and adjustments of model parameters to accom-modate this simplification. Text S5. Illustrative figures of time series of Alfv é nic events constrained by CHAMP data (Figure S3) and amplitude-fre-quency spectra of modeled and observed field-aligned current events (Figure S4). Text S6 . Altitudes profiles for V, and (Figure S5) used to compute FDTD solutions. Text S7. Notes on Joule and frictional heating. ogan et al. Alfvénic Thermospheric Upwelling 3 conditions (e.g., summer and IMF clock angles between 135 ° -225) exhibit very low amplitudes where CHAMP anomalies occur. Matsuo & Richmond’s distributions of electric field variability assume the observed variation is spatial (no time variation) with contributions integrated over horizontal scales between 3 km and 500 km. The assumption of static electric field variability may limit the applicability of these results, particularly for contributions at the smaller scales (< few 100 km), for which Alfvénic rather than quasistatic variability is usually prominent (Ishii et al., 1992; Chaston et al., 2003; Lühr et al., 2015; Park et al., 2017). We focus in this paper on the effects of Alfvén wave energy deposition on the thermospheric mass density in and near the low-altitude cusp at the nominal CHAMP altitude of 400 km. The Alfvén wave model of Lotko & Zhang (2018) is adapted as a parameterized model in the NCAR Coupled Magnetosphere-Ionosphere model (CMIT), which is used to simulate and analyze the MIT response to the stream interaction region (SIR) that swept past Earth on 26-27 March 2003. The relative importance of Alfvénic Joule heating for this event is characterized as a difference in thermospheric mass densities at 400 km altitude with and without inclusion of Alfvénic heating. The CMIT model includes i) a cusp-finding algorithm (Zhang et al., 2013), which dynamically determines where the Alfvénic Joule heating is to be applied, and ii) a causally regulated empirical specifica-tion for soft electron precipitation. The soft precipitation is modeled as both direct-entry cusp precipitation and Alfvén wave-induced “broadband” precipitation (Zhang et al., 2015b). As mentioned above, soft precipitation and some type of electric field enhancement are important factors in producing thermospheric density en-hancement in the cusp region. We first provide an overview of the CMIT model and the interplanetary conditions for the SIR event to be simulated (Section 2.1). The parameterized Alfvén wave model to be embedded in CMIT is described (Section 2.2), along with specification of the ambient IT parameters that determine the characteristics of Alfvén propa-gation and absorption (Section 2.3) and integration of the Alfvén wave heating model into CMIT (Section 2.4). Discussion of the simulation results (Section 3) and principal conclusions (Section 4) follow.
2. Methods
NCAR’s CMIT model has been used previously to simulate the MIT response to solar wind and interplanetary driving (Wang et al., 2004; Wiltberger et al., 2004; Zhang et al., 2012, 2015a, Liu et al., 2018). The magneto-spheric component of CMIT is the Lyon-Fedder-Mobarry (LFM) global magnetosphere simulation model, which describes the coupling between the solar wind and magnetosphere by solving the three-dimensional equations of ideal magnetohydrodynamics (Lyon et al., 2004). The so-called “double-resolution” version (53x48x64 grid cells) of LFM is used here. Its physical domain is a distorted spherical grid extending 30 Earth radii (R E ) upstream from Earth into the solar wind, to 300 R E downstream on the nightside, and 100 R E on the sides. The LFM model is driven by time series of state variables for the solar wind and IMF obtained from interplanetary satellite measurements (the one-minute OMNI-combined data available at https://cdaweb.gsfc.nasa.gov/index.html/). The ionosphere-thermosphere (IT) component of CMIT is modeled by the Thermosphere-Ionosphere Electro-dynamics General Circulation Model (TIEGCM), which solves the three-dimensional equations of continuity, momentum and energy for neutral and ion gases of the IT system. The TIEGCM version used here determines state variables globally on a uniform geographic grid of 1.25 ° in both latitude and longitude with 57 pressure levels at different altitudes (Dang et al., 2018; Qian et al., 2013; Wang et al., 2008). In addition to its inputs from LFM – convection velocity and electron precipitation flux and energy – TIEGCM uses the F10.7 index for solar radio flux to parameterize its extreme ultraviolet ionization. The MIX solver couples TIEGCM and LFM on a uniform grid of 2 ° in MLAT and magnetic longitude and ex-tends from the poles down to 46 ° MLAT. MIX takes conductance information from TIEGCM and the field-aligned current (FAC) distribution from LFM to solve a Poisson equation for the electric potential at the iono-spheric height. The potential determines the ionospheric convection velocity in TIEGCM and a velocity bound-ary condition at LFM’s low-altitude boundary after mapping it along equipotential, magnetic dipole field lines to the boundary (Merkin & Lyon 2010). The time interval between LFM-MIX-TIEGCM data exchanges is 5 s. TIEGCM’s time step is set to be the same as the (5 s) exchange time. B x,GSE B y,GSM B z,GSM ρ V x2 n T n T n T n P a Start SI
Hogan et al. Alfvénic Thermospheric Upwelling 4
The low-altitude cusp finding algorithm developed by Zhang et al. (2013) from LFM state variables is used to define the cusp area. This simulated dynamic cusp is projected onto TIEGCM’s grid where low-energy, direct-entry cusp electron precipitation is specified, also using LFM state variables. Direct-entry cusp precipitation is distinct from broadband, monoenergetic and diffuse precipitation which MIX also determines from LFM state variables and passes to TIEGCM (Zhang et al., 2015b). Low-energy electron cusp precipitation has been shown to be important in thermospheric heating in the cusp because it enhances the bottomside F -region Pedersen conductivity and thermospheric heating there (Clemmons et al., 2008; Zhang et al., 2012, 2015a; Deng et al., 2013; Brinkman et al., 2016). We use LFM’s cusp finding algorithm to specify the TIEGCM cells in which Alfvénic heating is imposed (Section 2.4). We used CMIT to investigate the thermospheric heating and upwelling that occurred as a stream interaction region (SIR) swept past Earth on 26-27 March 2003 (Figure 2). We chose an SIR event because variability in the solar wind and IMF accompanying such events is expected to produce intense Alfvénic activity in the cusp. During this event the CHAMP mean altitude (413 km) is comparable to the altitude where the distribution of density anomalies in Fig. 1 occur (Liu et al., 2005). Hemispheric asymmetries in the ionospheric conductivities are minimized (though not entirely absent) because the event is near spring equinox. No SIRs occurred in the 11 days prior to this event (Jian et al., 2006), so a less disturbed MIT system is expected at model start time, making the geospace response to the SIR reasonably distinct. The sampled IT environment begins to respond to the SIR around 00:00 UT on the 27th (described in Section 3) which is the middle of the simulated period. The first day of the event (until 1800 UT) exhibits weak-to-moderate activity in the IMF; the second day exhibits more intense and persistent variability in the IMF. Using IMF data to drive LFM requires special treatment of a variable IMF B x in order to maintain ∇ · B=0 in the interplanetary medium (Lyon et al., 2004). As described in Text S1 of the Supporting Information, we used the method of Wiltberger et al. (2000) to specify the variable IMF B x for this event. To reduce memory of initial conditions in simulating the SIR event, TIEGCM and LFM were separately preconditioned as described in Text S2 of the Supporting Information. We consider energy deposition in the cusp ionosphere by low-frequency magnetohydrodynamic waves, pre-sumed to be generated at high altitudes near the magnetopause by interplanetary variability and/or magneto-pause dynamics associated with dayside reconnection and Kelvin-Helmholtz surface waves. In general, such disturbances stimulate slow and fast mode compressional waves and intermediate mode shear Alfvén waves. Slow mode waves encounter strong ion Landau damping in the magnetosphere where the ion to electron temperature ratio typically exceeds one. Consequently, very little slow mode power stimulated by high-altitude
Figure 2.
IMF and solar wind dynamic pressure from 00 UT on 26 March 2003 through 00 UT on 28 March ogan et al. Alfvénic Thermospheric Upwelling 5 disturbances reaches the ionosphere. The group velocity of fast mode waves is nearly isotropic, so fast mode power falls off as 𝑟𝑟 −2 with distance 𝑟𝑟 from the source. The power carried by nondispersive, shear Alfvén waves reaches the ionosphere most efficiently because (nondispersive) shear Alfvén waves are guided practically loss-free along magnetic field lines, and their energy is magnetically focused by the converging flux tube. The field-aligned Poynting flux 𝑆𝑆 ∥ of the Alfv é n wave varies approximately with magnetic field intensity 𝐵𝐵 as 𝑆𝑆 ∥ 𝐵𝐵⁄ . For these reasons, we consider energy deposition resulting only from shear mode Alfv é n waves and neglect signals and ionospheric energy deposition resulting from slow and fast mode waves. The linear model for collisional propagation and absorption of Alfvén waves in the ionosphere-thermosphere used in this study is described by Lotko & Zhang (2018). It neglects magnetic compressibility, ionospheric Hall currents and electron inertia, which, as discussed by Lotko and Zhang, is appropriate for angular wave fre-quencies ( 𝜔𝜔 = 2 𝜋𝜋𝜋𝜋 ) and perpendicular wavenumbers 𝑘𝑘 ⊥ satisfying 𝜇𝜇 𝑒𝑒 𝑛𝑛 𝑒𝑒 𝑚𝑚 𝑒𝑒 ⁄ ≫ 𝑘𝑘 ⊥2 ≫ 𝜇𝜇 𝜔𝜔𝜎𝜎 𝐻𝐻2 𝜎𝜎 𝑃𝑃 ⁄ for given Hall ( 𝜎𝜎 𝐻𝐻 ) and Pedersen ( 𝜎𝜎 𝑃𝑃 ) conductivities, the permeability of free space ( 𝜇𝜇 ) , and the electron charge ( 𝑒𝑒 ), mass ( 𝑚𝑚 𝑒𝑒 ) and density ( 𝑛𝑛 𝑒𝑒 ). We also neglect the effects of neutral-gas inertia, mediated by ion-neutral friction, on Alfvén wave propagation and dissipation. These effects are negligible when 𝜔𝜔 ≫ 𝜈𝜈 𝑖𝑖 𝜌𝜌 𝑖𝑖 𝜌𝜌 𝑛𝑛 ⁄ where 𝜌𝜌 𝑖𝑖 and 𝜌𝜌 𝑛𝑛 are the ion-species and neutral-gas mass densities and 𝜈𝜈 𝑖𝑖 is the ion-neutral collision frequency (see Support-ing Information Text S3 and Figure S2). The equations of Lotko & Zhang (2018) are: ( ) t P V ν ⊥ ⊥ ⊥ ∂ + = ∇ × E E B (1) ( ) t ∂ = − ∇ × B E (2) ( ) || || 0 ˆ E j σ η= = ⋅ ∇ × b B (3) 𝐄𝐄 ⊥ and 𝐁𝐁 ⊥ are the wave perpendicular electric and magnetic fields; 𝑗𝑗 ∥ is the wave field-aligned current derived from 𝐁𝐁 ⊥ ; 𝐛𝐛̂ is the unit vector in the local magnetic-field direction. In (1) and (3) the Pedersen dissipation rate ( 𝜈𝜈 𝑃𝑃 = 𝜎𝜎 𝑃𝑃 𝜀𝜀⁄ ) and magnetic diffusivity ( 𝜂𝜂 = 1 𝜎𝜎 𝜇𝜇 ⁄ ) are given in terms of the Pedersen and direct or parallel ( 𝜎𝜎 ) conductivities defined as νσ ν +Ω , s ss s n qm = ∑ (4) . ee e n em σ ν= . (5) The wave speed ( 𝑉𝑉 = 1 ( 𝜇𝜇 𝜀𝜀 ) / ⁄ = 𝑐𝑐 ( 𝜀𝜀 𝜀𝜀 ) ⁄ / ) in (2) depends on the speed of light in vacuum ( c ), the vacuum permittivity ( 𝜀𝜀 ) and the dielectric constant ( 𝜀𝜀 ) for Alfvén wave propagation including collisional effects
20 2 2s s ε ε 1 ν +Ω . pss ω = + ∑ (6) Plasma parameters in (4)-(6) include the number density ( 𝑛𝑛 𝑠𝑠 ) , charge ( 𝑞𝑞 𝑠𝑠 ) , mass ( 𝑚𝑚 𝑠𝑠 ) , plasma frequency ( 𝜔𝜔 𝑝𝑝𝑠𝑠 ) and collision frequency ( 𝜈𝜈 𝑠𝑠 ) of charged-particle species s ( e for electrons, i for ion species). The collision frequencies are taken from Schunk & Nagy (SN) (2009). For ions, 𝜈𝜈 𝑠𝑠 includes ion-neutral collisions (SN, equa-tion 4.146 and Table 4.4; Table 4.5); for electrons, it includes electron-ion (SN, equation 4.144) and electron-neutral collisions (SN, Table 4.6). Our focus is on low-altitude collisional Alfvénic energy deposition (below 1000 km and mostly below 500 km), so we follow the method of Lotko & Zhang (2018) and simplify the analysis by treating the geomagnetic field as straight and uniform. The implications and accommodations of this treatment are described in Text S4 of the Supporting Information. At this juncture we neglect spatial gradients in MLT in the background parameters and in Alfvén waveforms, which is also assumed in deducing SSFACs from low-altitude satellite data (see below). The 2D form of (1)-(3) in coordinates parallel ( z ) and perpendicular ( x ) to the background magnetic field ( x corresponds to the magnetic meridional or MLAT direction, positive toward north) become, after combining (2) and (3), ogan et al. Alfvénic Thermospheric Upwelling 6 , yP x BE Vt z ν∂ ∂+ = −∂ ∂ (7) η . y x y B E Bt z x x ∂ ∂ ∂ ∂= − +∂ ∂ ∂ ∂ (8) We neglect perpendicular (MLAT) spatial variations in background parameters 𝜈𝜈 𝑃𝑃 , V and η , which is valid when the length scales for such spatial variations are large compared to the perpendicular wavelength. Equa-tions (7) and (8) are solved by harmonic solutions of the form ( ) ( ) ( ) , , , ( ) cos cos , M Nx m n m m n nm n
E x z t E z k x t θ ω θ = = = + + ∑ ∑ , (9) ( ) ( ) ( ) , , , , ( ) cos cos . M Ny m n m m n n m nm n
B x z t B z k x t θ ω θ φ = = = + + + ∑ ∑ (10) 𝐸𝐸 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 ) and 𝐵𝐵 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 ) are altitude-dependent Fourier amplitudes; 𝜙𝜙 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 ) is the altitude-dependent phase dif-ference between the wave electric and magnetic field for each mode; 𝑘𝑘 𝑚𝑚 = 2 𝜋𝜋𝑚𝑚 𝜆𝜆 max ⁄ is the wavenumber; 𝜔𝜔 𝑛𝑛 = 2 𝜋𝜋𝑛𝑛 𝜏𝜏 max ⁄ is the angular frequency; and 𝜆𝜆 max and 𝜏𝜏 max are the maximum meridional perpendicular wave-length (100 km) and wave period (20 sec) considered. Lacking knowledge of any correlations or relationships between the different modes, we assume they are uncorrelated with random phases 𝜃𝜃 𝑚𝑚 and 𝜃𝜃 𝑛𝑛 on the interval [0, 2 π ]. The expansions include M × N = 400 ×
40 modes. The following considerations were made in choosing the upper and lower cutoffs on perpendicular wavenumbers and frequencies. − Modes with perpendicular wavelengths ≲ 𝜆𝜆 max 𝑀𝑀 = 250 m ⁄ , but thermospheric heating and upwelling are relatively insensitive to inclusion of these short wavelength modes. − The longest wavelength included in (9) and (10) ( 𝜆𝜆 max = 100 km) is an appropriate upper limit on the wavelength spectrum when modeling parameterized Alfvénic heating in CMIT when the TIEGCM cell size is 1.25 ° (140 km) and the MIX and LFM cell sizes are 220 km or larger when referenced to 100 km altitude. − Alfvén waves with frequencies comparable to or greater than the local proton gyrofrequency near the mag-netopause are locally absorbed by ion-cyclotron resonant interactions (Stawarz et al., 2016). For this rea-son, we do not include modes with frequencies > 𝑁𝑁 𝜏𝜏 max ⁄ = 2 Hz in (9) and (10) because they do not propagate from the near-magnetopause region to the low-altitude region of interest. − The lowest frequency mode ( 𝜋𝜋 min = 1 𝜏𝜏 max ⁄ ) in (9) and (10) is practically quasi-static. Including electric field variability with frequencies below 𝜋𝜋 min would augment the dc Joule heating rate in CMIT, which is treated in the Joule heating module of TIEGCM, as discussed in Section 2.4. In solving (7) and (8), the harmonic dependence in x is first introduced which converts 𝜕𝜕 𝜕𝜕𝑥𝑥 ⁄ in (8) to −𝑘𝑘 𝑚𝑚2 . The altitude dependent amplitudes 𝐸𝐸 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 ) and 𝐵𝐵 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 ) and phases 𝜙𝜙 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 ) , now parameterized by 𝑘𝑘 𝑚𝑚 , are then derived from the finite-difference time-domain (FDTD) numerical method described by Christ et al . (2002) with modes stimulated by a time-periodic driver with frequency 𝑛𝑛 𝜏𝜏 max ⁄ at altitude z = 4500 km and electric field amplitude 𝐸𝐸 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 𝑑𝑑 ) . Setting 𝐵𝐵 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 𝑑𝑑 ) = 𝐸𝐸 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 𝑑𝑑 ) 𝑉𝑉 ( 𝑧𝑧 𝑑𝑑 ) ⁄ at the driver stimulates a downward propagating Alfvén wave. The computational domain extends from the Earth’s surface ( z = 0 km altitude) up to z = 5000 km. The FDTD solver with Mur absorbing boundary conditions (Mur, 1981) is run until its solution converges to the harmonic forms (9) and (10) (Lotko & Zhang, 2018). We impose observational constraints on the FDTD solution by adjusting (through trial and error) the frequency and wavenumber dependence of 𝐸𝐸 𝑚𝑚 , 𝑛𝑛 ( 𝑧𝑧 𝑑𝑑 ) to produce an amplitude-frequency spectrum for the field-aligned current (FAC) similar to the average spectrum shown in Figure 3 of Rother et al. (2007), derived from CHAMP satellite measurements at an altitude of approximately 400 km. The procedure for imposing the constraints is ogan et al. Alfvénic Thermospheric Upwelling 7 as follows. − As is typical when determining the FAC from satellite measurements of the magnetic field (Lühr et al., 1996), Rother et al. assume (i) the FAC is effectively sheet-like (also consistent with the 2D approximation used in (7) and (8)) and (ii) the curl of the magnetic field is equivalent to the along-track Doppler derivative normal to the sheet, i.e., 𝜕𝜕 𝜕𝜕𝑥𝑥⁄ = 𝜕𝜕 𝑉𝑉 𝑠𝑠 𝜕𝜕𝜕𝜕⁄ . Here, x is the spatial coordinate normal to the sheet and the background magnetic field. 𝑉𝑉 𝑠𝑠 ≈ − In (10), we set 𝑥𝑥 = 𝑉𝑉 𝑠𝑠 𝜕𝜕 , differentiate the resulting expansion with respect to 𝑉𝑉 𝑠𝑠 𝜕𝜕 , and divide the result by 𝜇𝜇 to obtain the simulated, Doppler-derived FAC as a function of time. − Sample time series of the Doppler-derived FAC along a simulated CHAMP satellite track may then be constructed and its Fast Fourier Transform calculated to obtain a sample amplitude-frequency spectrum. Choosing the n and m dependence of the modal amplitudes at the driver to be (11) was found to approximate the spectral shape of the FACs reported by Rother et al. (2007). 𝐸𝐸 = 500 mV/m produced a peak FAC in the Doppler-derived time series comparable to the mean peak value of 350 μ A m ⁄ in Fig. 7 of Rother et al. (2007). The observed peak values range from 100 μ A m ⁄ (Rother’s threshold value) to 2000 μ A m ⁄ . 𝐸𝐸 was scaled up to 1500 mV/m to bracket a likely upper limit on thermospheric upwelling in the simulations with Alfvénic heating described in Sec. 3. Note that the actual electric field amplitude calculated from (9) is much less than 𝐸𝐸 owing to the mode number dependence in the denom-inator of (11). Model time series and spectra obtained from the above procedure, using profiles for the background param-eters specified in the next section, are described in the Supporting Information Text S5 and Figs. S3 and S4. FDTD solutions to (7) and (8) require specification of the altitude profiles ( z dependence) for parameters 𝜎𝜎 𝑃𝑃 , 𝜎𝜎 and ε . We determined these parameters from equations (4)-(6) using the most recent International Refer-ence Ionosphere (IRI-2016) (Bilitza 2018), the Navy Research Lab Mass Spectrometer and Incoherent Scatter radar (NRLMSISE-00, or MSIS) model (Picone et al., 2002) and an Earth-centered dipole model for the geo-magnetic field. The IRI model uses the date to look up historical data of daily F10.7, a measure of solar activity by measuring the solar radio flux at 10.7 cm as a proxy of solar activity. The MSIS model similarly uses F10.7, as well as user-input for the 3-hour Ap, a proxy of the planetary magnetic field for which historical data exists. Two different profiles for 𝜎𝜎 𝑃𝑃 , 𝜎𝜎 and ε were computed at the average cusp location (given by its MLT-MLAT center) at 12:00 UT on 26 March 2003 and on 27 March 2003 in a baseline CMIT run without Alfvénic heating. The first day exhibits weak-to-moderate solar wind/IMF driving. More intense SIR activity in the interplanetary medium is very prominent during the second day (Fig. 2). The electron density computed by IRI is scaled such that the peak electron density (NmF2) matches the peak in the average electron density profile in CMIT in the cusp region for the same 24-hour period. This rescaling of the plasma density from IRI values to the mean CMIT values for each day of the simulated interval is a step toward self-consistency of the parameterized model and CMIT. MSIS determines the thermosphere only up to 1000 km altitude. Above this altitude collisional effects are negligible and the collisionless Alfvén speed �𝐵𝐵 �𝜇𝜇 ∑ 𝑚𝑚 𝑠𝑠 𝑛𝑛 𝑠𝑠𝑠𝑠 ⁄ � is used up to 2000 km. Above 2000 km altitude, the Alfvén speed is constant and equal to its value at 2000 km. The Pedersen conductivity is interpolated on a log-linear scale between 1000 and 2000 km altitude using the conductivity calculated using MSIS-IRI at 1000 km and a specified small value (1.88x10 -12 S/m) at 2000 km altitude – the value in Lotko & Zhang (2018). The results are insensitive to this value as long as it is small. Above 2000 km we set 𝜎𝜎 𝑃𝑃 to its value at 2000 km. The parallel conductivity is interpolated on a log-linear scale from its computed value at 1000 km to a ( ) m n d E z E n n m = + A l t i t ude , k m quasistaticAlfvénic Hogan et al. Alfvénic Thermospheric Upwelling 8 large value at the top of the simulation space such that η becomes negligible; = 1x10 S/m was chosen. At altitudes below 90 km, below the range of CMIT, the wave speed is constant and equal to 12000 km/s. The wave speed is actually larger than 12000 km/s in the lower thermosphere, approaches the speed of light in the lower atmosphere and is artificially limited here for numerical efficiency. This choice does not influence the results in the IT region of interest. The conductivity profiles are interpolated on a log-linear scale between computed values at 90 km and 0 km altitude where both and are set to 6x10 -12 S/m (Lysak 1999, Lotko & Zhang, 2018). The wave solution at E - and F -region altitudes is relatively insensitive to the value specified for conductivities at Earth’s surface, provided it is small compared to the value above 90 km. The altitude profiles of , and wave speed used to compute FDTD solutions for the Alfvén wave fields are included in Figure S5 of the Supporting Information. The Alfvén wave solution, observationally constrained as described in the previous section, is used to augment CMIT’s Joule heating rate as a function of altitude. To this end, we use the fact that time scale for ion accel-eration by the Alfvénic fluctuations (wave periods ≤
20 sec) is short compared to the time scale for neutral-fluid acceleration ( − > 1 hour). This disparity in time scale renders the Alfvén wave-induced quiver velocity of the neutral gas—mediated by ion-neutral friction—negligible compared to the perturbed ion velocity for the considered range of wave frequencies (0.5-2 Hz) and IT conditions described in Sec. 2.3 (see SI Text S3 and Fig. S2). In this case, the effects of Alfvén wave acceleration on ion-neutral friction may be neglected in the neutral gas momentum equation, and the neutral gas may be considered stationary when evaluating the Alf-vénic contribution to Joule heating (SI Text S3). Also making use of the inequalities , ≪ Ω that underlie the derivation of (1)-(3), the Alfvénic Joule heating rate , may be calculated as (see SI Text S7). Rather than specifying a particular realization of the Alfvénic fluctuations encountered by the neutral fluid as it traverses the cusp region, we characterize the fluctuations in terms of an ensemble average over all possible random phases in (9) to obtain the Alfvénic Joule heating rate in terms of the RMS amplitude of the fluctua-tions: max max PJ A P x P x m n P A rmsx t m n
Q z E x z t dx dt E x z t E z E λ τ σσ σ σ δλ τ= = = ≡ ∑∫ ∫ (12) The average profile for the root-mean-square (RMS) Alfvénic electric field, , , is precomputed twice by the FDTD Alfvén wave solver, at 12:00 UT on 26 March 2003 and on 27 March 2003, at the average cusp locations where the profiles were determined from the TIEGCM-scaled, MSIS-IRI data as described in Sec. 2.3. These two representative , profiles are taken to be constant in x and t throughout the corresponding day and over the instantaneous cusp area determined by LFM’s cusp finding algorithm (Zhang et al., 2013). While the Alfvén wave power flowing into CMIT’s cusp area is actu-ally bursty in space and time, the thermosphere effectively inte-grates over the Alfvénic heating during its relatively slow transit through the cusp. Thus a space-time average for parameterized Alf-vénic heating of the thermosphere is appropriate. The profile ob-tained for , for day 2 of the event is shown in Figure 3 along with that of a static electric field with the same value in the E region for comparison. The Alfvénic electric field is ≈
50% larger at F -region altitudes than a static electric field of comparable intensity in the E -region ionosphere owing to the effect of the IAR. IAR modes at fre-quencies greater than ≈ F region (cf. Figs. 2 and 3 of Lotko & Zhang, 2018). The RMS electric field , from the Alfvén wave model is inte-grated into the Joule heating module of CMIT by augmenting the Figure 3.
Altitude profile of the rms Alfvén wave electric field , for day 2 of the event compared to that of a quasistatic electric field with the same amplitude in the E region. ogan et al. Alfvénic Thermospheric Upwelling 9 Joule heating per unit mass for the neutral �𝑄𝑄 𝐽𝐽 𝑇𝑇 𝑁𝑁 � and ion �𝑄𝑄 𝐽𝐽 𝑇𝑇 𝑖𝑖 � gasses according to the formulae ( ) N P A rmsTJ J E B n
B EQ C B σ δρ × ⊥ = − + v v (13) i n nT TJ Jn i mQ Qm m = + . (14) The Pedersen conductivity 𝜎𝜎 𝑃𝑃 in 𝑄𝑄 𝐽𝐽 𝑇𝑇 𝑁𝑁 is calculated by CMIT at each TIEGCM time step (5 seconds), 𝜌𝜌 is the thermospheric mass density from TIEGCM, B is the geomagnetic field (the International Reference Geomag-netic Field is used in TIEGCM), 𝐯𝐯 𝐸𝐸𝑥𝑥𝐸𝐸 and 𝐯𝐯 𝑛𝑛 ⊥ are, respectively, the ion E × B drift velocity and the neutral wind velocity perpendicular to B , both derived from CMIT, and 𝑚𝑚� 𝑛𝑛 and 𝑚𝑚� 𝑖𝑖 are the mean molecular mass of neutrals and ions, respectively. Equation (13) is the Joule heating rate in the reference frame of the neutral gas, which is equivalent to the frictional heating rate due to ion-neutral drag (see SI Text S7). For the reasons discussed in SI Text S3, the neutral gas may be considered stationary when evaluating the Alfvénic contribution to Joule heating, i.e., the second term on the right side of (13) The basic Joule heating module in TIEGCM includes the constant multiplicative factor ( 𝐶𝐶 𝐽𝐽 ) that augments the Joule heating specified by the first term in (13). This factor is intended to model the effects of quasistatic subgrid variability in the Joule heating (e.g., Figure 1) and effectively adds subgrid Joule heat where TIEGCM produces grid-resolved Joule heat. Although this simple scaling for the subgrid Joule heat is not likely to be accurate in detail, it improves TIEGCM’s baseline prediction for the thermospheric mass density. The addition of the parameterized Alfvénic heating model in CMIT also captures subgrid heating processes, so we reduced the nominal value 𝐶𝐶 𝐽𝐽 = 1.5 in the default TIEGCM code to 1.35. A time-dependent, assimilative adjustment of 𝐶𝐶 𝐽𝐽 , as in the assimilation method of Sutton (2018), may further improve the fidelity of the baseline prediction. Even though 𝛿𝛿𝐸𝐸 𝐴𝐴 , 𝑟𝑟𝑚𝑚𝑠𝑠 is constant during each day of the SIR event, the Alfvénic Joule heating rate evolves in CMIT as TIEGCM evolves because 𝑄𝑄 𝐽𝐽 , 𝐴𝐴 ( 𝑧𝑧 ) in (13) is calculated using TIEGCM’s 𝜎𝜎 𝑃𝑃 . This Alfvén wave coupling dynamically modulates Joule heating because the CMIT-determined 𝜎𝜎 𝑃𝑃 changes with dynamic variations in electron precipitation and ion-neutral chemistry, both of which are at play during the simulated SIR event. In theory and possibly in the future, both 𝛿𝛿𝐸𝐸 𝐴𝐴 , 𝑟𝑟𝑚𝑚𝑠𝑠 and 𝑄𝑄 𝐽𝐽 , 𝐴𝐴 ( 𝑧𝑧 ) could be regulated at the cadence of a CMIT time step by regulating the Alfvén wave driver amplitude in the FDTD solver using LFM’s Alfvén wave Poynting flux module (Zhang et al., 2015b) and TIEGCM’s IT dynamic state. The approach implemented here may be con-sidered a proof-of-concept to determine the value of continuing such future developments.
3. Results and Discussion
We first ran CMIT for 48 hours for the solar and interplanetary conditions during 26-27 March 2003 without the Alfvénic heating term in equation (13). This run provides a baseline for comparing changes when Alfvénic heating is included in CMIT. A one-orbit moving mean of the CHAMP data and CMIT results (Figure 4) was evaluated to gain an overview of cumulative effects and hour-scale trends in the thermospheric response to the SIR event. The IMF activity and solar wind driving are weak to moderate during the first 24 hours (26 March 2003) of the event, and the air density recorded along the CHAMP orbit increases modestly during the day, from about 1.6 × -12 kg/m to 2.3 × -12 kg/m . CMIT under-predicts the density by as much as 15% during the first day. Two curves are plotted for CMIT in Fig. 4 to bracket its prediction. The curve labeled CMIT is the baseline estimate. The curve labeled CMIT with Alfvénic heating is calculated for an RMS wave electric field ( 𝛿𝛿𝐸𝐸 𝐴𝐴 , 𝑟𝑟𝑚𝑚𝑠𝑠 ) that yields a peak Alfvénic FAC of 1 mA/m ( 𝐸𝐸 = 1500 mV/m in equation 11). The CMIT estimate for the orbit-average mass density of the neutral gas lies in the gray area between the two curves. The MSIS prediction is also plotted in Fig. 4 for comparison. Around 18:00 UT on 26 March, both the IMF activity and solar wind driving pick-up (cf. Fig. 2). The latency in the response of the thermosphere along this particular orbit is about 6 hours, after which the thermospheric density observed by CHAMP, simulated by CMIT and estimated by MSIS all begin a comparatively steep
20 MSIS246810-4-20246 CMIT CMIT Alfvén
26 Mar27 Mar
NHSH - ρ M o d e l / ρ C H A M P , % - ρ M o d e l / ρ C H A M P , %
00 12
26 Mar 2003 27 Mar 2003 ρ a i r , - k g / m CHAMP MSISCMIT withAlfvénic heatingCMITMoving Mean Density (averaging window = 1 orbit)
Hogan et al. Alfvénic Thermospheric Upwelling 10 increase shortly after 00:00 UT on March 27 (cf. Fig. 4). All three estimates of the air density exhibit a rel-ative maximum just after 18:00 UT on 27 March. CHAMP registers a density of 3.9 × -12 kg/m at this time; the corresponding MSIS estimate is near 3.0 × -12 kg/m while the CMIT prediction ranges from 3.6 to 3.8 × -12 kg/m . CMIT clearly tracks the CHAMP density more closely during the second day of the event. As expected the CMIT prediction is more dynamic than MSIS – a consequence of the minute-cadence interplanetary dynamics driving CMIT vs. the daily F10.7 and 3-hour Ap indices pa-rameterizing MSIS. The increase in air density due to Alfvénic heating relative to baseline is larger during the second day of the event than during the first day, with a maxi-mum difference of 0.08 × -12 kg/m during the first day, 0.15 × -12 kg/m during the second day. How-ever, the maximum difference over baseline is ≈
4% on both days because the baseline is greater on the second day. The larger increase in density from Alf-vénic heating on the second day is due in part to an
Figure One-orbit moving mean of CHAMP-derived thermospheric mass density for 26-27 March 2003, CMIT simulated thermospheric density along the CHAMP orbit, with and without Alfvénic heating, and the MSIS prediction of the density along the CHAMP orbit. The curve labeled CMIT is the baseline run. CMIT with Alfvénic heating is the run with average peak FAC of 1 mA/m . Gray area between them indicates the expected range for CMIT density enhancement with Alf-vénic heating effects. Figure 5.
Comparison of daily mean error (in %) inmodel predictions relative to CHAMP’s air density for MSIS, CMIT (baseline) and CMIT with Alfvénic heating for each hemisphere (NH, SH) and event day (26, 27 March). ogan et al. Alfvénic Thermospheric Upwelling 11 increase in TIEGCM’s Pedersen conductivity. To quantify the average daily cumulative error in model predictions, we calculated the instantaneous error ( − Model CHAMP ) ⁄ along the CHAMP orbit and then determined the daily average error for each hemisphere (Figure 5). By this measure, MSIS performs better than CMIT for weak-to-moderate activity (26 March), with Alfvénic heating decreasing CMIT’s relative error by ≈
2% in both hemispheres. For higher activity (27 March), CMIT performs better than MSIS. Including Alfvénic heating in CMIT decreases its error by ≈
3% in the NH on27 March, but it increases the error by 4% in the SH. We attribute the reduction in error when Alfvénic heating is included, in part, to the increased Joule heating resulting from the increase in CMIT’s Pedersen conductivity. The increase in error in the SH during 27 March may be due in part to inconsistencies in CMIT’s treatment of the geomagnetic field, which is based on the International Geomagnetic Reference Field in TIEGCM and on a Earth-centered dipole field in the LFM and MIX solvers. TIEGCM includes the intrinsic NH/SH asymmetries of the IGRF, but LFM and MIX drive TIEGCM with inputs from the magnetosphere that ignore this asymmetry. ρ a i r , - k g / m CHAMPLocation
CMITCMITAlfvén CHAMP
MSIS
N H
CHAMP ρ air , 10 -12 kg/m CHAMP ∆( CMIT) , %
75 60 45-75 -60 -45
DFE
27 Mar 2003 ρ a i r , - k g / m N H ρ air , 10 -12 kg/m ∆( CMIT) , %8 8.5 9UT, hours
CMIT
CMITAlfvén CHAMP
MSIS
CHAMPLocationCHAMP -75 -60 -4575 60 45CHAMP
ACB
Figure 6.
Air density along consecutive CHAMP orbits in the northern hemisphere (NH) on 27 March 2003 when CHAMP passes through the edges of high-latitude density enhancements. A: Simulated air density (color) at the CHAMP altitude ( ≈
413 km) vs MLT and MLAT at the UT corresponding to “CHAMP Location” in panel B. Black circle indicates location of CHAMP at that UT. B: Comparison of instantaneously observed, simulated and MSIS air density vs. UT. CMIT curve is the baseline run. Curve labeled “CMIT Alfvén” is a run with a peak Alfvénic FAC of 1 mA/m . As in Fig. 4 gray area indi-cates expected range of CMIT predictions with Alfvénic heating. C: % difference in CMIT air density with and without Alf-vénic heating. Δ (CMIT) = é ⁄ − . Panels D , E and F in same format. ogan et al. Alfvénic Thermospheric Upwelling 12 We also compared the instantaneous air density recorded on CHAMP with CMIT and MSIS predictions. The left panel of Figure 6 shows CMIT diagnostics for the air density on 27 March 2003 in the northern hemisphere (NH) at the CHAMP altitude (413 km). At 0815 UT CHAMP traverses the edge of a dayside density enhance-ment in the CMIT thermosphere including the effects of cusp-region Alfvénic heating (Fig. 6A). The point of closest approach at 0815 UT is indicated by the black circle located at 18.7 MLT and 81.0 ° MLAT. Fig. 6B compares the instantaneous air densities from CHAMP measurements, CMIT predictions and MSIS estimates along the CHAMP orbit vs UT. As in Fig. 4, the curve labeled CMIT is for the baseline run, while the curve labeled CMIT with Alfvénic heating is the CMIT run with a peak Alfvénic FAC of 1 mA/m . The gray area between these two curves indicates the expected range of CMIT predictions with Alfvénic heating. Fig. 6C is a NH map of the % difference in CMIT’s predictions with and without Alfvénic heating: Δ (CMIT) = 𝜌𝜌 𝐴𝐴𝐴𝐴𝐴𝐴𝐴𝐴 é 𝑛𝑛 𝜌𝜌 𝑏𝑏𝑏𝑏𝑠𝑠𝑒𝑒𝐴𝐴𝑖𝑖𝑛𝑛𝑒𝑒 ⁄ − . The maximum density difference in Fig. 6B is about 12%, but the main cusp density enhancement (not traversed by CHAMP) exhibits a difference closer to 20% in Fig. 6C. Changes in thermospheric density resulting from Alfvénic heating are not localized to the cusp or the dayside. The largest difference occurs between about 0200 and 0600 MLT near 80 ° MLAT. The effects of co-rotation and neutral winds evidently redistribute the additional cusp density enhancement produced by Alfvénic heating. The right panels of Fig. 6 continue the analysis in the same format for the subsequent CHAMP crossing of the NH. CHAMP traverses the edge of the density enhancement in CMIT at 9.67 MLT and -69.5 ° MLAT at 0955 UT (black circle in Fig. 6D). The comparison between CHAMP measurements, CMIT prediction and MSIS estimate along the orbit (Fig. 6E) shows that CMIT pre-dicts the location of CHAMP’s peak density in both hemispheres, and its accuracy improves with Alfvénic heating included. This feature is missed entirely by MSIS, which exhibits a relative minimum where both the CHAMP and MSIS densities peak. The maximum density difference for CMIT with and without Alfvénic heating is about 14% in Fig. 6E. As for the previous NH crossing, the main cusp density enhancement (which CHAMP does not traverse) exhibits a larger difference closer to 30% in Fig. 6F, and the redistribution and modification of the Alfvénic-enhanced density is prominent throughout the nightside polar region.
4. Conclusions • CMIT tracks orbit-averaged air density at CHAMP altitudes more closely than MSIS during intervals of solar wind and IMF activity characteristic of an SIR event (Fig. 4).
When comparing the IT response during the less active, first day of the simulated SIR event to the response during the more active second day, we find that CMIT replicates CHAMP-derived density measurements better than MSIS during the more active period. This finding is not surprising because the proxy-governed nature of MSIS relies on inputs of average parameters over longer periods (multi-hour to days) and is not as responsive to interplanetary variability as CMIT when driven by one-minute solar wind and IMF data. The MIT coupling intrinsic to CMIT allows effects such as soft electron precipitation, electric field variability and Alfvénic heating to modify the Pedersen conductivity, Joule heating rate and, therefore, air temperature and density in the cusp region on faster time scales than MSIS. • Alfvénic heating modestly improves CMIT’s daily mean prediction of orbit-averaged air density relative to CHAMP measurements (Fig. 5) and produces 20-30% regional enhancements in the high-latitude ther-mospheric density relative to a CMIT simulation without Alfvénic heating (Fig. 6 bottom panels).
Inclusion of Alfvénic heating effects in CMIT decreases its daily mean error by 2-3% in orbit-averaged air density relative to measurements along a particular CHAMP orbit, except in the southern hemisphere during the active second day of the SIR event when the error increased by 4%. This discrepancy may arise in part from the nonconforming geomagnetic field models used in CMIT’s three component models, which treat geo-magnetic hemispheric asymmetry differently. Enhancements of 20-30% in air density are seen throughout the high-altitude region (not just the cusp region) relative to a baseline CMIT run. This finding suggests that, during the course of the two-day long simulation, co-rotation and neutral winds redistribute the additional cusp density enhancements produced by Alfvénic heating. • Alfvénic heating produces significant thermospheric density enhancements at CHAMP altitudes in and near the cusp during times of intense Alfvénic Poynting flux (Fig. 6 middle panels). ogan et al. Alfvénic Thermospheric Upwelling 13
Alfvénic thermospheric heating is especially effective when magnetopause variability induced by solar wind and IMF activity stimulates intense Alfvénic Poynting fluxes flowing into the low-altitude cusp. The resulting Alfvén wave energy deposition in the ionospheric Alfvén resonator adds to the important effects of soft electron precipitation and quasistatic electric field variability on F -region Joule heating and thermospheric upwelling previously reported in the literature. These relatively fast dynamic effects operating in concert may explain the anomalous thermospheric densities CHAMP recorded relative to the MSIS empirical model (Liu et al., 2005). A detailed analysis showed that MSIS completely misses the air density enhancements registered by both CHAMP and CMIT along near cusp orbits; and CMIT with Alfvénic heating substantially improves CMIT’s baseline prediction of the enhancements (by 10-15%). During this SIR event, CHAMP did not actually traverse the central cusp region where CMIT’s density enhancements with Alfvénic heating are even larger (20-30%) than without it. Acknowledgments
CHAMP and MSIS data were provided by Eric Sutton. CHAMP data are archived at http://tinyurl.com/den-sitysets. IMF and Solar wind data were provided by J.H. King, N. Papatashvilli at AdnetSystems, NASA GSFC and CDAWeb (http://omniweb.gsfc.nasa.gov/). We thank Roger Varney for sharing his Fortran code for colli-sion frequencies, Jing Liu and Wenbin Wang for helpful insights into TIEGCM and thermospheric dynamics, and Wenbin for feedback on the manuscript. All code needed to produce the plots shown including CHAMP, MSIS, and CMIT data is available publicly: http://doi.org/10.5281/zenodo.3727080. We acknowledge support by the National Center for Atmospheric Research, a major facility sponsored by the National Science Foun-dation under Cooperative Agreement No. 1852977 and the Thayer School of Engineering at Dartmouth Col-lege. Computing resources were provided by NCAR’s Computational and Information Systems Laboratory (CISL). BH was supported by NASA NH Space Grant NNX15AH79; WL by NASA Grant 80NSSC19K0071, and WL and KP by NSF awards 1739188, 1522133.
References
Bilitza, D. (2018), IRI the International Standard for the Ionosphere.
Advanced Radio Science , , 1-11. https://doi.org/10.5194/ars-16-1-2018. Brinkman, D. G., R. L. Walterscheid, J. H. Clemmons, & J. H. Hecht (2016), High-resolution modeling of the cusp density anomaly: Response to particle and Joule heating under typical conditions. Journal of Geophys-ical Research: Space Physics , , 2645-2661. https://doi.org/10.1002/2015JA021658. Chaston, C. C., Bonnell, J. W., Carlson, C. W., McFadden, J. P., Ergun, R.E., & Strangeway, R. J. (2003). Properties of small-scale Alfvén waves and accelerated electrons from FAST. Journal of Geophysical Re-search, 108 (A4), 8003 https://doi.org/10.1029/2002JA009420. Chaston, C. C., Carlson, C. W., McFadden, J. P., Ergun, R. E., & Strangeway, R. J. (2007). How important are dispersive Alfvén waves for auroral particle acceleration?
Geophysical Research Letters, 34 , L07101. https://doi.org/10.1029/2006GL029144. Clemmons, J. H., Hecht, J. H., Salem, D. R., & Strickland D. J. (2008). Thermospheric Density in the Earth's Magnetic Cusp as Observed by the Streak Mission.
Geophysical Research Letters , , L24103. https://doi.org/10.1029/2008gl035972. Dang, T., Lei, J., Wang, W., Burns, A., Zhang, B., & Zhang, S.-R. (2018). Suppression of the polar tongue of ionization during the 21 August 2017 solar eclipse. Geophysical Research Letters, 45, 2918–2925. https://doi.org/10.1002/2018GL077328. Deng, Y., Fuller ‐ Rowell, T. J., Ridley, A. J., Knipp, D., & Lopez, R. E. (2013). Theoretical study: Influence of different energy sources on the cusp neutral density enhancement.
Journal of Geophysical Research: Space Physics , , 2340 ‐ Journal of Geophysical Research: Space Physics, 122 , 12,189–12,211. https://doi.org/10.1002/2017JA024175. Ish ii, M., Sugiura, M., Iyemori, T., & Slavin, J. A. (1992). Correlation between magnetic and electric field per-turbations in the field -aligned current regions deduced from DE 2 observations.
Journal of Geophysical Re-search, 97 (A9), 13,877–13,887. https://doi.org/10.1029/92JA00110. Jee, G., Burns, A. G., Wang, W., Solomon, S. C., Schunk, R. W., Scherliess, et al. (2008). Driving the TING model with GAIM electron densities: Ionospheric effects on the thermosphere.
J. Geophys. Res. , , ogan et al. Alfvénic Thermospheric Upwelling 14 A03305, https://doi.org/10.1029/2007JA012580. Jian, L., Russel, C. T., Luhmann, J. G., & Skoug, R. M. (2006). Properties of stream interactions at one AU during 1995–2004.
Solar Physics , 239(1-2), 337–392. https://doi.org/10.1007/s11207-006-0132-3. Keiling, A., Wygant, J. R., Cattell, C. A., Mozer, F. S., & Russell, C. T. (2003). The global morphology of wave
Poynting flux: P owering the aurora. Science, 299(5605), 383–386. https://doi.org/10.1126/science.1080073. Kervalishvili, G. N., & Lühr, H. (2013). The relationship of thermospheric density anomaly with electron tem-perature, small-scale FAC, and ion up- flow in the cusp region, observed by CHAMP and DMSP satellites.
Annales de Geophysique, 31 (3), 541–554. https://doi.org/10.5194/angeo-31-541-2013. Kervalishvili, G. N. & Lühr, H. (2014). Climatology of zonal wind and large-scale FAC with respect to the density anomaly in the cusp region: seasonal, solar cycle, and IMF B y dependence. Ann. Geophys., 32 , 249–261, https://doi.org/10.5194/angeo-32-249-2014. Liu, H., Lühr, H., Henize, V., & Köhler, W. (2005). Global distribution of the thermospheric total mass density derived from CHAMP.
J. Geophys. Res. , , A04301. https://doi.org/10.1029/2004JA010741. Liu, J., Wang, W., Zhang, B., Huang, C., & Lin, D. (2018). Temporal variation of solar wind in controlling solar wind ‐ magnetosphere ‐ ionosphere energy budget. Journal of Geophysical Research: Space Phys-ics , 123, 5862– 5869. https://doi-org.cuucar.idm.oclc.org/10.1029/2017JA025154. Lotko, W., & Zhang, B. (2018). Alfvénic heating in the cusp ionosphere-thermosphere.
Journal of Geophysical Research: Space Physics , 123. https://doi.org/10.1029/2018JA025990. Lyon, J.G., Fedder, J. A., & Mobarry, C. M., (2004). The Lyon–Fedder–Mobarry (LFM) Global MHD Magneto-spheric Simulation Code.
Journal of Atmospheric and Solar-Terrestrial Physics , (15-16), 1333–1350. https://doi.org/10.1016/j.jastp.2004.03.020. Lysak, R. L. (1991). Feedback instability of the ionospheric resonant cavity. Journal of Geophysical Research, 96 (A2), 1553–1568. https://doi.org/10.1029/90JA02154. Lysak, R. L. (1999). Propagation of Alfvén waves through the ionosphere: Dependence on ionospheric pa-rameters.
J. Geophys. Res. , A5), 10,017–10,030. https://doi.org/10.1029/1999JA900024. Lühr, H., Warnecke, J. F. & M. K. A. Rother, M. K. A. (1996). An algorithm for estimating field-aligned currents from single spacecraft magnetic field measurements: a diagnostic tool applied to Freja satellite data.
IEEE Transactions on Geoscience and Remote Sensing, 34(
Geophysical Research Letters, 31,
L06805. https://doi.org/10.1029/2003GL019314. Lühr, H., Park, J., Gjerloev, J. W., Rauberg, J., Michaelis, I., G. Merayo, J. M.,& Brauer, P. (2015). Field ‐ aligned currents' scale analysis performed with the Swarm constellation. Geophys. Res. Lett., 42 , 1– 8, https://doi.org/10.1002/2014GL062453. Park, J., Lühr, H., Knudsen, D. J., Burchill, J. K., & Kwak, Y. ‐ S. (2017). Alfvén waves in the auroral region, their Poynting flux, and reflection coefficient as estimated from Swarm observations.
Journal of Geophysical Research: Space Physics , , 2345– 2360, https://doi.org/10.1002/2016JA023527. Matsuo, T., & Richmond, A. D. (2008). Effects of high- latitude ionospheric electric field variability on global thermospheric Joule heating and mechanical energy transfer rate. Journal of Geophysical Research, 113,
A07309. https://doi.org/10.1029/2007JA012993. Merkin, V. G., & Lyon, J. G., (2010). Effects of the Low-Latitude Ionospheric Boundary Condition on the Global Magnetosphere.
J. Geophys. Res. , (A10). https://doi.org10.1029/2010ja015461. Mur, G. (1981). Absorbing boundary conditions for the finite -difference approximation of the time-domain elec-tromagnetic- field equations.
IEEE Transactions on Electromagnetic Compatibility , EMC-23 (4), 377–382. https://doi.org/10.1109/TEMC.1981.303970. Picone, J. M., A. E. Hedin, D. P. Drob, & A. C. Aikin, (2002). NRLMSISE-00 empirical model of the atmosphere: Statistical comparisons and scientific issues.
Journal of Geophysical Research: Space Physics , (A12), 1468, https://doi.org/10.1029/2002JA009430 Qian, L., et al. (2013). The NCAR TIE-GCM: A Community Model of the Coupled Thermosphere/Ionosphere System, in Modeling the Ionosphere-Thermosphere System, ed. by J.D Huba, R.W. Schunk & G.V. Khaza-nov, Geophysical Monograph Series 201, AGU, Washington D.C. https://doi.org/10.1029/2013GM001452. Rentz S., & Lühr H. (2008). Climatology of the cusp-related thermospheric mass density anomaly, as derived from CHAMP observations.
Annales Geophysicae. 26 (9), 2807-2823. https://doi.org/10.5194/angeo-26-2807. ogan et al. Alfvénic Thermospheric Upwelling 15
Rother, M., Schlegel, K., & Lühr, H. (2007). CHAMP observation of intense kilometer- scale field -aligned cur-rents, evidence for an ionospheric Alfvén resonator.
Annales de Geophysique, 25 (7), 1603–1615. https://doi.org/10.5194/angeo-25-1603-2007. Schunk, R., & Nagy, A. (2009).
Ionospheres, 2 nd Edition , Ch. 4. Cambridge, UK: Cambridge University Press. Stawarz, J. E., Eriksson, S., Wilder, F. D., Ergun, R. E., Schwartz, S. J., Pouquet, A., Burch, J. L., et al. (2016). Observations of turbulence in a Kelvin Helmholtz event on 8 September 2015 by the Magnetospheric Mul-tiscale mission.
J. Geophys. Res. , , 11,021–11,034. https://doi.org/10.1002/2016JA023458. Sutton, Eric K (2018). A New Method of Physics-Based Data Assimilation for the Quiet and Disturbed Ther-mosphere. Space Weather , (6), 736–753. https://doi.org/10.1002/2017sw001785. Trakhtengerts, V. Y. & Fel'Dshtein, A. Y. (1987). Turbulent regime of magnetospheric convection. Geomag-netism and Aeronomy , (2), 221-225. Tu, J., Song, P., & Vasyliūnas, V. (201
Journal of Geophysical Research, 116 , A09311. https://doi.org/10.1029/2011JA016620. Wang, W., Wiltberger, M., Burns, A. G., Solomon, S. C., Killeen, T. L., Maruyama, N., & Lyon, J. G. (2004). Initial results from the Coupled Magnetosphere–Ionosphere–Thermosphere Model: Thermosphere–iono-sphere responses.
Journal of Atmospheric and Solar-Terrestrial Physics , (15-16), 1425–1441. https://doi.org/10.1016/j.jastp.2004.04.008. Wang, W., Lei, J., Burns, A. G., Wiltberger, M., Richmond, A. D., Solomon, S. C., Killeen, T. L., Talaat, E. R., & Anderson, D. N. (2008). Ionospheric electric field variations during a geomagnetic storm simulated by a coupled magnetosphere ionosphere thermosphere (CMIT) model. Geophysical Research Letters , , L18105. https://doi.org/10.1029/2008GL035155. Wiltberger, M., Pulkkinen, T. I., Lyon, J. G., & Goodrich, C. C. (2000). MHD simulation of the magnetotail during the December 10, 1996, substorm. Journal of Geophysical Research , (A12), 27649– 27663, https://doi.org/10.1029/1999JA000251. Wiltberger, M., Wang, W., Burns, A.G., Solomon, S.C., Lyon, J. G., & Goodrich, C. C. (2004). Initial results from the coupled magnetosphere-ionosphere-thermosphere model: Magnetospheric and ionospheric re-sponses. J. Atmos. Sol.-Terr. Phys. , (15-16), 1411-1423, https://doi.org/10.1016/j.jastp.2004.03.026. Zhang, B., Lotko, W., Brambles, O., Wiltberger, M., Wang, W., Schmitt, P., & Lyon, J. (2012). Enhancement of thermospheric mass density by soft electron precipitation. Geophysical Research Letters, 39,
L20102. https://doi.org/10.1029/2012GL053519. Zhang, B., Brambles, O., Lotko, W., Dunlap ‐ Shohl, W., Smith, R., Wiltberger, M., & Lyon, J. (2013). Predicting the location of polar cusp in the Lyon ‐ Fedder ‐ Mobarry global magnetosphere simulation.
Journal of Geo-physical Research: Space Physics , , 6327–6337. https://doi.org/10.1002/jgra.50565. Zhang, B., Varney, R. H., Lotko, W., Brambles, O. J., Wang, W., Lei, J., Wiltberger, M., et al. (2015a). Path-ways of F region thermospheric mass density enhancement via soft electron precipitation. Journal of Geo-physical Research: Space Physics, 120,
Journal of Geophysical Research: Space Physics , Živković T, Buchert S, Ritter P, Palin L, & Opgenoorth, H. (2015). Investigation of energy transport and ther- mospheric upwelling during quiet magnetospheric and ionospheric conditions from the studies of low- and middle-altitude cusp.
Annales Geophysicae , (6), 623-635. https://doi.org/10.5194/angeo-33-623-2015. ogan et al. Alfvénic Thermospheric Upwelling 16 SUPPORTING INFORMATION
Alfvénic Thermospheric Upwelling in a Global Geospace Model
Benjamin Hogan , William Lotko and Kevin Pham Thayer School of Engineering, Dartmouth College, Hanover, NH, USA High Altitude Observatory, National Center for Atmospheric Research, Boulder, CO, USA † Now at Laboratory for Atmospheric and Space Physics, University of Colorado Boulder, Boulder, CO, USA and Department of Aerospace Sciences, University of Colorado Boulder, Boulder, CO, USA
Contents of Supporting Information
Text S1 to S7 Figures S1 to S5
Introduction
This Supporting Information includes the following sections and figures: Text S1. CMIT treatment of IMF 𝐵𝐵 𝑥𝑥 and Figure S1 comparing 𝐵𝐵 𝑥𝑥 from the OMNI data set with the con-strained 𝐵𝐵 𝑥𝑥 used in the SIR simulation. Text S2. CMIT preconditioning. Text S3. Validity of neglecting neutral gas inertia in collisional Alfvén wave dynamics and Figure S2 showing the domain of validity for the model ionosphere and thermosphere. Text S4. Implications of treating the cusp-region geomagnetic magnetic field as straight and uniform and adjustments of model parameters to accommodate this simplification. Text S5.
Illustrative figures of time series of Alfvénic events constrained by CHAMP data (Figure S3) and amplitude-frequency spectra of modeled and observed field-aligned current events (Figure S4). Text S6. Altitudes profiles for V, 𝜎𝜎 𝑃𝑃 and 𝜎𝜎 (Figure S5) used to compute FDTD solutions. Text S7. Notes on Joule and frictional heating. ogan et al. Alfvénic Thermospheric Upwelling 17 Text S1. CMIT Treatment of IMF B x Single point measurements acquired by a solar wind monitor are typically linearly advected in the direction of the solar wind flow, which is very nearly the GSE x direction. Consequently, the time varia-bility in IMF implies ≠ ⁄ . Since information on the y and z dependence of IMF and is not available from single-point upstream measurements, the solenoidal condition, ∇ ∙ = 0 , cannot be maintained without introducing additional constraints on the IMF. We follow the procedure of Wilt-berger et al. (2000) and assume the variability in IMF can be represented as ( ) ( ) ( ) x y z B t a bB t cB t = + + . (S1.1) The coefficients a , b , and c are determined by a multiple linear regression fit to all one-minute IMF data samples for the two-day SIR event (60 samples per hour ×
48 hours). The fit establishes a plane oriented at a fixed angle relative to the GSE x direction with normal direction b cx y za a a = − − n (S1.2) along which the normal magnetic field is constant. This constraint allows a time-dependent B x to be introduced into the simulation that maintains ∇ ∙ = 0 . For the 26-27 March 2003 event, the fitting procedure gives: = (0) = 0.99528304; = = − − Figure S1 compares B x from the OMNI dataset with the value used in the CMIT simula-tion and that calculated from S1.1. Wiltberger, M., Pulkkinen, T. I., Lyon, J. G., & Goodrich, C. C. (2000). MHD simulation of the magnetotail during the December 10, 1996, substorm.
Journal of Geophysical Research , (A12), 27649– 27663, https://doi.org/10.1029/1999JA000251. Figure S1.
Comparison of measured IMF B x from the OMNI dataset with the values resulting from the multiple linear regression fit for the 26-27 March 2003. ogan et al. Alfvénic Thermospheric Upwelling 18 Text S2. CMIT preconditioning
CMIT is preconditioned for runs with real-time solar wind inputs as follows. The magnetosphere portion of CMIT, LFM, is typically preconditioned for four hours and fifty minutes before the start of the simulation time (in this study 00 UT March 26 2003). 1.
LFM is first initialized with an Earth-centered dipole magnetic field permeating the entire simu-lation domain, populated by a low density fluid ( ≈ -3 ) composed of protons and alpha particles with number density ratio 𝑛𝑛 𝐻𝐻𝑒𝑒 𝑛𝑛 𝑝𝑝 ⁄ = 0.04 . 2. For the first fifty minutes of the LFM preconditioning, the upstream boundary conditions on the IMF and solar wind are specified with a 0 nT IMF, solar wind velocity ( V x , V y , V z ) = ( − 𝑛𝑛 𝑆𝑆𝑆𝑆 = 5 cm -3 and sound speed 𝑐𝑐 𝑆𝑆𝑆𝑆 = 40km/s. As the resulting solar wind dy-namic pressure interacts with the dipole magnetic field, it compresses the magnetic cavity on the dayside and stretches the nightside field to form the expected magnetospheric shape. This first phase of initialization prevents the zero-IMF solar wind front from advecting onto a grid cells containing dipole magnetic field lines so the resulting magnetosphere is in the initial low-density, near-vacuum state. 3.
At the end of this 50-minute period, the IMF at the upstream boundary is set to ( B x ,B y ,B z ) = (0,0, −
5) nT for the next two hours with the same solar wind parameters. When the southward IMF contacts the magnetopause, the resulting magnetic reconnection opens the magneto-sphere and allows solar wind fluid to enter the magnetosphere. 4.
At 2 hours, 50 minutes into the preconditioning, IMF B z is changed to +5 nT, which produces a quiet state and serves to expand and fully populate the magnetosphere with fluid from the so-lar wind. 5. The conductances required by MIX during the standalone LFM-MIX preconditioning are calcu-lated from empirical models for precipitation and EUV-induced ionization. These conduct-ances are obtained from TIEGCM after the models are coupled starting at 00 UT 26 March 2003. CMIT’s ionosphere-thermosphere model (TIEGCM) is preconditioned for five days (starting at 00 UT 21 March 2003) before starting the CMIT simulation interval at 00 UT 26 March 2003. For the pre-conditioning, standalone TIEGCM uses three empirical models as inputs, two of which (high-latitude electric field and electron precipitation models) are changed at 00 UT 26 March 2003 to LMF-MIX in-puts. The high-latitude electric field and resulting ionospheric drift velocity during preconditioning are derived from the Kp-parameterized, Heelis empirical model . The extreme ultraviolet (EUV) solar irradi-ance required by TIEGCM is obtained from the F10.7-parameterized, empirical EUVAC model . TIEGCM uses its default Kp-parameterized, empirical auroral precipitation model during the five-day initializa-tion. TIEGCM retains a history file of Kp and F10.7 values by date, which regulate the empirical models during the five-day preconditioning. EUVAC continues to be regulated by the historical F10.7 during 00 UT 26 March 2003 to 00 UT 28 March 2003. LFM-MIX and TIEGCM initializations end at 00 UT 26 March 2003, when they are coupled and run as the CMIT model. In CMIT mode, LFM-MIX provides high-latitude convection and electron precipitation inputs to TIEGCM, while TIEGCM provides the Pedersen and Hall conductances required by the MIX Zhang, B., Lotko, W., Brambles, O., Wiltberger, M. & Lyon, J. (2015b). Electron precipitation models in global magnetosphere simulations. Journal of Geophysical Research: Space Physics, 120, 1035–1056. https://doi.org/10.1002/2014JA020615. Wiltberger, M., Wang, W., Burns, A. G., Solomon, S. C., Lyon, J. G., & Goodrich, C. C. (2004). Initial results from the coupled magnetosphere ionosphere thermosphere model: magnetospheric and ionospheric responses,
Journal of Atmospheric and Solar-Terrestrial Physics, 66 (15-16), 1411-1423. https://doi.org/10.1016/j.jastp.2004.03.026. Heelis, R. A., Lowell, J. K., & Spiro, R. W. (1982). A model of the high‐latitude ionospheric convection pattern.
J. Geophys. Res. , 87(A8), 6339– 6345. https://doi.org/10.1029/JA087iA08p06339. Richards, P. G., Fennelly, J. A., & Torr, D. G. (1994). EUVAC: A solar EUV Flux Model for aeronomic calculations.
J. Geophys. Res. , 99(A5), 8981– 8992. https://doi.org/10.1029/94JA00518. ogan et al. Alfvénic Thermospheric Upwelling 19 potential solver. As described in Sec. 2.4, the pre-computed electric field in the Alfvén wave subgrid model is also included in TIEGCM at this time, first with the precomputed value for 26 March 2003, then with the precomputed value for 27 March 2003. When LFM-MIX and TIEGCM are first coupled, immediately after the preconditioning as standalone models, an impulsive change of state occurs in the simulated geospace environment due to discontinuous changes in the conductances used by MIX and the high-latitude convection and precipitation in TIEGCM as well as the addition of Alfvénic Joule heating to TIEGCM. To facilitate relaxation of the impulsive change before the onset of magnetic activ-ity in the SIR event at approximately 1800 UT on 26 March 2003 (cf. Fig. 2) we chose to start running CMIT 18 hours earlier at 00 UT 26 March 2003. This 18-hour buffer period serves to facilitate relaxation of the impulsive change of state in CMIT. Preconditioning intervals longer than five days for TIEGCM yield modestly different initial states at the end of the preconditioning interval, e.g., differences in neutral density at CHAMP altitudes for twenty-day preconditioning differ from values obtained for five-day preconditioning by less than 20%. Our controlled simulation experiments are concerned primarily with differences in CMIT predictions with and without parameterized Alfvénic heating, so the results of interest are relatively insensitive to moderate differences in IT states obtained from different preconditioning intervals. The buffer interval starting at 00 UT on 26 March 2003 before the onset of intense IMF activity of the SIR event 18 hours later further moderates any sensitivities to preconditioning after 18 UT on 26 March 2003.
Text S3. Validity of neglecting neutral gas inertia in collisional Alfvén wave dynamics
Neutral-gas inertia, mediated by ion-neutral collisions, has been shown to modify Alfvén wave characteristics at sufficiently low wave frequencies and high ion-neutral collision frequencies.
We derive the conditions of validity that allow this effect to be neglected for the model ionosphere and Alfvén wave spectrum considered in the paper. The momentum equation for a cold neutral gas with velocity 𝐯𝐯 𝑛𝑛 , mass 𝑚𝑚 𝑛𝑛 and number density 𝑛𝑛 𝑛𝑛 interacting with ion species i with velocity 𝐯𝐯 𝑖𝑖 , mass 𝑚𝑚 𝑖𝑖 and number density 𝑛𝑛 𝑖𝑖 via friction with collision frequency 𝜈𝜈 𝑛𝑛𝑖𝑖 is ( ) ( ) nn m n n ni n i i i in i n n m n m n mt ν ν∂ = − − = −∂ v v v v v . (S3.1) The last equality follows from momentum conservation during collisions: 𝑚𝑚 𝑛𝑛 𝑛𝑛 𝑛𝑛 𝜈𝜈 𝑛𝑛𝑖𝑖 = 𝑚𝑚 𝑖𝑖 𝑛𝑛 𝑖𝑖 𝜈𝜈 𝑖𝑖𝑛𝑛 . Mo-mentum transfer from electrons is less by a factor of order 𝑚𝑚 𝑒𝑒 𝜈𝜈 𝑒𝑒𝑛𝑛 𝑚𝑚 𝑖𝑖 𝜈𝜈 𝑖𝑖𝑛𝑛 ⁄ and has been neglected. Summing (S3.1) over n , assuming all neutral species have the same wind velocity u n , gives ( ) ( ) nn n n ni n i i i i i nn n m n mt ρ ν ν∂ = − − = −∂ ∑ u u v v u (S3.2) where n n nn n m ρ ≡ ∑ and i inn ν ν≡ ∑ . Define 𝛿𝛿 𝑖𝑖 ≡ 𝑛𝑛 𝑖𝑖 𝑚𝑚 𝑖𝑖 𝜌𝜌 𝑛𝑛 ⁄ , assume time-harmonic variation, and Fourier Transform (S3.2) with respect to time to obtain the following relationship between Fourier amplitudes (designated as 𝐮𝐮� 𝑛𝑛 and 𝐯𝐯� 𝑖𝑖 ) : i in ii i i δνω δν= − + u v . (S3.3) Hasegawa, Akira & Uberoi, Chanchal (1982).
The Alfven Wave Song, P., Vasyliūnas, V. M., & Ma, L. (2005). Solar wind‐magnetosphere‐ionosphere coupling: Neutral atmosphere effects on signal propagation.
J. Geophys. Res., 110 , A09309. https://doi.org/10.1029/2005JA011139. ogan et al. Alfvénic Thermospheric Upwelling 20
The modulus of (S3.3) is . i i in i i u δνω δν= + v (S3.4) Lysak’s (1999) wave equations, from which equations (1)-(3) of Sec. 2.2 are derived, assume the neutrals are effectively stationary ( ≪ v ) . Equation (S3.4) shows that this assumption is valid when ≫ . In this case, the inertia of the neutral gas impedes its response to rapid oscillations in the ion velocity. Then ( ) ( )( ) . ii i i i i i n i ii i i i i i m n m n ent m n en νν∂ = − − + + ×∂ ≈ − + + × v v u E v Bv E v B (S3.5) This approximation is well-satisfied (Fig. S2) for the model ionosphere-thermosphere (described in SI Text S6 and Figure S5) and range of wave frequencies treated in the paper, especially in the lower ionosphere. Collisional momentum exchange is conserved between ion species and the neutral gas, which is enslaved to the (Alfvénic) oscillating ion gas when ≫ : ≅⁄ . From this rela-tion and (S3.4) with ≫ , we estimate the magnitude of the neutral-gas “quiver” velocity induced by Alfvén-wave ion oscillations as ~( ⁄ ) v ,. The maximum Alfvén wave electric field in the cusp is ~ 100 mV/m (cf. Fig. S3) with an implied maximum E × B velocity of ~ 2 km/s. Using the mini-mum frequency of 0.05 Hz, the factor ⁄ from Fig. S2, is 0.04, 10 -4 , and 10 -6 at altitudes of 600 km, 300 km and 150 km, respectively. The maximum quiver velocity in the neutral wind is then ~ 80 m/s near the exobase, ~ 4 m/s in the F region, and 0.04 m/s in the E region. The average ion velocity en-countered by the neutral gas in its transit through the cusp is about 50x less for the example in Fig. S3, and it is less by an additional 10x for Alfvén waves with frequencies (> 0.5 Hz), which are most effective in augmenting the Joule heating rate and thermospheric upwelling at F -region altitudes. From these Lysak, R. L. (1999). Propagation of Alfvén waves through the ionosphere: Dependence on ionospheric parame-ters.
Journal of Geophysical Research, 104(
A5), 10,017–10,030. https://doi.org/10.1029/1999JA900024. Lotko, W., & Zhang, B. (2018). Alfvénic heating in the cusp ionosphere-thermosphere.
Journal of Geophysical Re-search: Space Physics, 123 . https://doi.org/10.1029/2018JA025990. A l t i t u d e , k m m i n i ν i /2 πρ n , s -1 O + O NO + Figure S2.
Altitude pro-file of /2 for O + , NO + , and O (domi-nant ion composition) vs altitude and comparison with the range of Alfvén wave frequencies (gray region, 0.05 Hz to 2 Hz) considered in the paper. Ion and neutral densities and ion-neutral collision frequencies are derived from IRI and MSIS mod-els evaluated in the cusp region for 27 March 2003 (see SI Text S6). ogan et al. Alfvénic Thermospheric Upwelling 21 estimates, we conclude that Alfvénic ion forcing of neutral-gas momentum in the cusp is negligible compared to other thermospheric forcings (e.g., frictional quasi-static convection, pressure gradients, advection, etc.). Text S4. Implications of treating the geomagnetic field as straight and uniform
The actual spatial variation in the background magnetic field introduces the following effects for an Alfvén wave propagating from 1000 km to 90 km altitude at the nominal magnetic latitude (75°) of the cusp: (i) The magnetic field increases by a factor of 1.5, which increases the speed of wave propa-gation; (ii) the flux tube convergence reduces the perpendicular length scale of the wave structure by a factor of 0.82; and (iii) the field line curvature shifts the geographic location of the magnetic foot point to higher latitude by about 1°. Effect (i) is included in the calculation of the dielectric constant ε in equation (6) and the wave speed V . Effect (ii) tends to reduce the Joule heating rate and increase the Ohmic heating rate, in going from high to low altitude and is modeled as described below. The small effect (iii) influences the interpretation of the geomagnetic latitude of wave energy deposition. To include effect (ii) in the Alfvén wave solution, we follow Lotko & Zhang (2018) and note that the radial variation of the square of the effective perpendicular wavenumber 𝑘𝑘 ( 𝑟𝑟 ) due to magnetic focusing is 𝑘𝑘 ( 𝑟𝑟 ) = 𝑘𝑘 ( 𝑟𝑟 𝑟𝑟 ) ⁄ = (2 𝜋𝜋 𝜆𝜆 ⁄ ) ( 𝑟𝑟 𝑟𝑟 ) ⁄ where 𝑘𝑘 is the wavenumber specified at 400 km altitude, 𝑟𝑟 ≡ R E + 400 km = 6800 km and r = 1 R E + z is the radial distance at altitude z . The altitude dependence of the wave solution due to the altitude dependence of the wavenumber may be absorbed into an effective magnetic diffusivity defined as 𝜂𝜂 𝑒𝑒𝐴𝐴𝐴𝐴 ( 𝑧𝑧 ) ≡ [ 𝑟𝑟 (R 𝐸𝐸 + 𝑧𝑧 ) ⁄ ] 𝜇𝜇 ⁄ 𝜎𝜎 ( 𝑧𝑧 ) . With this definition, the solutions for the Alfvén wave electric field are pa-rameterized by the perpendicular wavelength at 400 km altitude. Text S5. Alfvénic time series and observed and model spectra for FACs
Choosing modal amplitudes of the Alfvén wave electric field at the driver to be , 2 2
500 mV/m( ) 2.4 m n d
E z n n m = + (S4.1) yields a peak amplitude for the FAC of approximately 350 µ A/m at 400 km altitude, the mean peak value in Fig.7 of Rother et al. (2007) . A sample time series for the Doppler-derived FAC is shown in Fig. S3, calculated from the Doppler derivative 𝜕𝜕 𝑉𝑉 𝑠𝑠 𝜕𝜕𝜕𝜕⁄ of equation (10). Time series for the Alfvén wave electric and magnetic fields, curl-derived FAC (spatial rather than satellite Doppler derivative) and Joule heating rates are also shown in Fig. S3. The amplitude-frequency spectrum from the octave averaged FFT of the FAC in this time series is shown in Fig. S4, superposed on the average spectrum from Fig. 3 of Rother et al. (2007) . The shape of the model spectrum approximates that of observed average spectrum reasonably well at low-frequencies ( f < 5 Hz), but it rolls-off more gradually above 5 Hz. The model spectrum at frequencies > 5 Hz is due to Doppler-shift of high wavenumber ( k ) modes for which the frequency in the spacecraft frame is f = kV s. with V s ≈ ). Our results for ther-mospheric Joule heating are relatively insensitive to this high- k dissipation. The following analysis is for an Earth-centered dipole field with an equatorial surface value of 3.1 × T. Rother, M., Schlegel, K., & Lühr, H. (2007). CHAMP observation of intense kilometer-scale field-aligned currents, evidence for an ionospheric Alfvén resonator.
Annales de Geophysique, 25 (7), 1603–1615. https://doi.org/10.5194/angeo-25-1603-2007. n T Time, s0 0 5 10 15 20 255 n W / m µ A / m -2000200 -50050 m V / m E x B y σ P E x2 j z-curl j z-Doppler rms19.9-2.2avgrms86.818.7avg -2000200 µ A / m rms54.40.9avgrms55.20.4avg0.4avg -334-338 359-939 Hogan et al. Alfvénic Thermospheric Upwelling 22
Figure S3.
Illustrative synthetic time series of Alfvénic fields along a simulated CHAMP orbital segment through a 200-km wide (MLAT) cusp. The fields are calculated from (9) and (10) with = . From top to bottom, MLAT component of the electric field , MLT component of the magnetic field , FAC calculated as the curl of , FAC calculated as the Doppler derivative of , and volumetric Joule heating rate . Maximum value of each field is marked near the peak; rms and average values to the right. Figure S4.
Octave-average spectra for the field-aligned current determined from the curl ( ) and the Doppler-derivative ( , ) of ( , , ) for the model ensemble of Alf-vén waves and from CHAMP data averaged over 2883 events (Rother et al., 2007) !1 . The model spectra are calculated for average background parameters on 27 March 2003. A m p l i t u d e , µ A / m j z j z,D CHAMP Average Spectrum, 2883 Events
ModelSpectra ogan et al. Alfvénic Thermospheric Upwelling 23
Text S6. Altitudes profiles of V, and used to compute FDTD solutions
Figure S5 shows profiles for the wave speed V , Pedersen conductivity , and parallel conductivity up to 2000 km altitude on 27 March 2003, calculated using IRI and MSIS to determine the iono-spheric and thermospheric variables in equations (4)-(6). These parameters are documented in the re-pository as “iriday2parameters.txt” and “msisday2parameters.txt.” The electron density profile ob-tained from IRI is subsequently multiplied by a constant so that the peak electron density, or NmF2, matches that of the average electron density profile in the central cusp from the baseline CMIT simula-tion for each day of the study. These profiles are essentially the same for 26 March 2003. Above 2000 km, V and are constant and equal to their respective values at 2000 km up to 5000 km. To facilitate convergence of the wave solver, is interpolated to a large value specified at the 5000 km top of sim-ulation space, here 10 S/m. The results are not sensitive to this value as long as it is sufficiently large. The resulting profiles are used to calculate FDTD wave solutions as described in section 2.3 of the main text.
Text S7. Notes on electromagnetic energy deposition and frictional heating
The term Joule heating as conventionally used in the IT literature (and in most of the papers refer-enced in the Introduction to this paper) refers to the rate of IT electromagnetic energy dissipation cal-culated as P E E σ σ ⊥ ′ ′⋅ = + J E . (S7.1) ′ = + × is the electric field in the reference frame of the neutral gas moving with velocity , is the background (geomagnetic) field, and and are the Pedersen and direct conductivities given in Sec. 2.2. The contribution ∥ to Joule heating is sometimes distinguished as Ohmic heating owing to the functional similarity between and the electrical conductivity of a simple solid conduc-tor. In most applications the dominant term in (S7.1) is ⊥ , which alone is often taken to be the Joule heating rate in IT studies. Figure S5.
Altitude profiles up to 2000 km for the Alfvén speed ( V) , Pedersen ( ) and parallel ( ) conductivi-ties on 27 March 2003 used by the FDTD solver. H e i g h t , k m V km/s S/m10 -10 -5 σ P S/m10 -10 σ ogan et al. Alfvénic Thermospheric Upwelling 24 Equation (S7.1) assumes the IT variability is quasistatic, i.e., it occurs on time scales much larger than both the gyroperiod and the mean time between collisions. This assumption allows an iono-spheric Ohm’s law to be expressed as P H E σ σ σ ⊥ ′ ′= + × + J E b E b . (S7.2) From the perspective of kinetic theory, IT Joule and Ohmic heating arise from friction between the electron, ion and neutral gases. When the energy transfer integral for collisional interactions can be given in terms of Maxwell molecule collisions, the heating rate for species 𝑟𝑟 due to collisions with spe-cies 𝑠𝑠 with frequency 𝜈𝜈 𝑟𝑟𝑠𝑠 is ( ) rs r r rs b s r s r sr s Q m n k T T mt m m δ νδ = − + −+ v v . (S7.3) 𝑇𝑇 𝑠𝑠 , 𝐯𝐯 𝑠𝑠 , 𝑚𝑚 𝑠𝑠 and 𝑛𝑛 𝑠𝑠 are the species temperature, velocity, atomic mass and number density; 𝑘𝑘 𝐸𝐸 is Boltz-mann’s constant. At altitudes below 400-500 km, the length scale 𝐻𝐻 for spatial variations and the time scale 𝜏𝜏 due to ion gyro motion Ω 𝑖𝑖−1 and ion-neutral collisions ν 𝑖𝑖𝑛𝑛−1 satisfy | 𝐯𝐯 𝑖𝑖 − 𝐮𝐮 𝑛𝑛 | 𝜏𝜏 ≪ 𝐻𝐻 . In this case, collisional ion-neutral energy exchange achieves an approximate balance between temperature and friction cou-pling with 𝑘𝑘 𝐸𝐸 ( 𝑇𝑇 𝑖𝑖 − 𝑇𝑇 𝑛𝑛 ) ≅ 𝑚𝑚 𝑛𝑛 | 𝐯𝐯 𝑖𝑖 − 𝐮𝐮 𝑛𝑛 | . The heating rate of the neutral gas due to ion friction, 𝑄𝑄 𝐽𝐽𝑇𝑇 𝑛𝑛 in the notation of Sec. 2.4, is then n niTJ i i in i ni i QQ m nt δ νδ= ≅ − ∑ ∑ v u . (S7.4)
Substituting the ion-neutral velocity difference in the quasi-static limit, i n i i in e Em ν′− = Ω + v u , (S7.5) in (S7.4), the frictional heating rate of the neutral gas can then be equivalently calculated as the Joule heating rate in the neutral-wind frame:
22 22 22 2 . n i inTJ i i in i ni i i i inPi E B n Pi e nQ m n EmE νν νσ σ × ′= − = Ω +′= − = ∑ ∑∑ v uv u (S7.6) The partial Pedersen conductivities 𝜎𝜎 𝑃𝑃𝑖𝑖 are defined implicitly by the third equality in (S7.6). This particular form is used to calculate the Joule heating rate (né ion-neutral frictional heating) in the TIEGCM . This result can also be derived in the MHD approximation in the center-of-mass frame for plasma Schunk, R., & Nagy, A. (2009). Ionospheres, 2nd Edition, Ch. 4. Cambridge, UK: Cambridge University Press. Schunk, R. W. (1975). Transport equations for aeronomy.
Planetary and. Space Science, 23 , 437-485. https://doi.org/10.1016/0032-0633(75)90118-X. St.‐Maurice, J.‐P., & Schunk, R. W. (1981). Ion‐neutral momentum coupling near discrete high‐latitude iono-spheric features.
J. Geophys. Res., 86 (A13), 11299– 11321. https://doi.org10.1029/JA086iA13p11299. Thayer, Jeffrey P. & Semeter, Joshua (2004). The convergence of magnetospheric energy flux in the polar at-mosphere.
Journal of Atmospheric and Solar-Terrestrial Physics, 66 , 807–824. https://doi.org/10.1016/j.ja-stp.2004.01.035. ogan et al. Alfvénic Thermospheric Upwelling 25 and neutral constituents (essentially the neutral-gas frame), when small corrections for Ohmic heat-ing from the generalized Ohm’s law and the work done on the neutrals in this frame are neglected. Alternatively, the Joule heating rate in the ion reference frame is essentially the Ohmic heating rate in MHD except at altitudes from about 120-300 km, where it differs by the additional friction term 𝑛𝑛 𝑒𝑒 𝑚𝑚 𝑒𝑒 𝜈𝜈 𝑒𝑒𝑛𝑛 | 𝐮𝐮 𝑀𝑀𝐻𝐻𝐷𝐷 − 𝐮𝐮 𝑛𝑛 | . Ion-neutral friction appears as separate mechanical energy dissipation in the ion reference-frame formulation and is not associated with electromagnetic energy dissipation in that frame. These formal distinctions regarding neutral gas heating, whether considered as electromag-netic or mechanical energy dissipation, do not change the calculated rate of neutral gas heating or its physical origin arising from the collisional interaction between charged particles and the neutral gas. In our treatment of Joule dissipation (Sec. 2.4), Ohmic dissipation is unimportant for the perpen-dicular wavelengths of interest and only conventional Joule dissipation is evaluated for Alfvénic en-ergy deposition 𝑄𝑄 𝐽𝐽 , 𝐴𝐴 . We now show that (S7.6) can be used to evaluate collisional Alfvénic energy deposition with the result 𝑄𝑄 𝐽𝐽 , 𝐴𝐴 = 𝜎𝜎 𝑃𝑃 𝐸𝐸 𝐴𝐴⊥2 , even though Alfvén wave dynamics do not fully satisfy the quasistatic approximation assumed in deriving (S7.6). This simple Joule heating formula is applicable to Alfvén dynamics because, as discussed in SI Text S3, at sufficiently high wave frequency and low ion-neutral collision frequency, the ratio of the perturbed neutral gas velocity to Alfvén wave ion ve-locity is negligibly small. Since Lysak’s formulation is the basis for equations (1)-(3) in Sec. 2.2, we start with his equation (A3) for the linear velocity perturbation for ionized species s : s s ss s em t νϖ ⊥ ⊥ ⊥ ∂= + + ×∂ v E E Ω . (S7.7) We defined 𝜛𝜛 𝑠𝑠2 ≡ 𝜈𝜈 𝑠𝑠2 + Ω 𝑠𝑠2 . As discussed by Lysak, this expression with ∂ / ∂ t → ω assumes 𝜔𝜔 , 𝜔𝜔𝜈𝜈 𝑠𝑠 ≪Ω 𝑠𝑠2 . Retaining corrections with 𝜔𝜔𝜈𝜈 𝑠𝑠 ≥ Ω 𝑠𝑠2 improves the accuracy of the polarization drift in the E layer (term proportional to ∂ / ∂ t ), but the resulting differences in Alfvén wave propagation are minor be-cause the polarization current in the E layer is much less than the Pedersen current. Equation (S7.7) also neglects the small, wave-induced perturbation in neutral velocity stimulated by ion-neutral fric-tion (see SI Text S3). From (S7.7), we obtain s s i s ss s s sss e e Em t m t tEB ν ν ϖϖ ϖϖ ⊥ ⊥⊥ ⊥ ⊥ ⊥ ⊥⊥ ∂ ∂ ∂= + + × = + ⋅ +∂ ∂ ∂Ω E EE E
Ω E v (S7.8) where 𝜔𝜔 , 𝜔𝜔𝜈𝜈 𝑠𝑠 ≪ Ω 𝑠𝑠2 has been used again in the last equality. In light of the analysis in SI Text S3 we can neglect 𝐮𝐮 𝒏𝒏 compared with 𝐯𝐯 𝒊𝒊 in 𝐄𝐄′ ⊥ in (S7.6) to obtain for the Joule heating rate due to Alfvén-wave energy deposition A A A i i iJ A i i i iA i i i Pi i ii i i
E n eQ m n m n E EB m νν ν σϖ ϖ ⊥⊥ ⊥ ⊥
Ω= = = ∑ ∑ ∑ v . (S7.9) Vasyliūnas, V. M., & Song, P. (2005). Meaning of ionospheric Joule heating.
J. Geophys. Res., 110,
A02301, https://doi.org/10.1029/2004JA010615. Strangeway, R. J. (2012). The equivalence of Joule dissipation and frictional heating in the collisional iono-sphere.
J. Geophys. Res., 11
7, A02310, https://doi.org/10.1029/2011JA017302. Gombosi, T. I. (1998). Physics of the Space Environment, Ch, 4. Cambridge, UK: Cambridge University Press. Sydorenko, D. & Rankin, R. (2017). The stabilizing effect of collision-induced velocity shear on the ionospheric feedback instability in Earth’s magnetosphere.