Interaction of Small-Scale Gravity Waves with the Terdiurnal Solar Tide in the Mesosphere and Lower Thermosphere

Abstract. Implementing a nonlinear whole atmosphere gravity wave (GW) parameterization into the Middle and Upper Atmosphere Model extending to the lower thermosphere (160 km), we study the response of the atmosphere in terms of the circulation patterns, temperature distribution, and migrating terdiurnal solar tide activity to the upward propagating small-scale internal GWs originating in the lower atmosphere. We perform three test simulations for the Northern Hemisphere winter conditions in order to assess the effects of variations in the initial GW spectrum on the dynamics of the mesosphere and lower 5 thermosphere. We find that the overall strength of the source level momentum flux has a comparatively small impact on zonal mean dynamics. The tails of the GW source level spectrum, however, are crucial for the lower thermosphere dynamics. With respect to the terdiurnal tide, we find a strong dependence of tidal amplitude on the induced GW drag, generally being larger when GW drag is increased.

respective stationary planetary waves (SPWs) with wavenumbers 1-3 (see also Lilienthal et al., 2017). The horizontal resolution is 5 • × 5.625 • (latitudes × longitudes). There are 56 vertical levels, which are evenly spaced in logarithmic pressure height with p s = 1000 hPa as the reference pressure level and H = 7 km as the scale height. Vertical spacing is about 2.8 km, and consequently, the upper boundary is located at z 56 ≈ 160 km. The lower levels up to 30 km height are nudged with 2000-2010 60 mean monthly mean ERA-Interim zonal mean temperature fields to correctly represent the dynamics in the lower atmosphere (Jacobi et al., 2015;Lilienthal et al., 2018). Solar and infrared radiative processes are parameterized according to the works by Strobel (1978) and Fomichev and Shved (1985), respectively. These parameterizations focus on (i) the absorption and emission processes of the most important atmospheric constituents like H 2 O (troposphere), CO 2 and O 3 (stratosphere) as well as on (ii) absorption bands like the extrem ultra violet (EUV) band in the thermosphere. The H 2 O, CO 2 , and O 3 distributions are prescribed. Further parameterizations deal with thermospheric processes such as Rayleigh friction, ion drag, and Newtonian cooling.
Solar tides are generated self-consistently in the model by the absorption of solar radiation, mainly due to water vapor and ozone. Unlike other mechanistic models, there is no explicit tidal forcing at the lower boundary. The sources of TDTs within MUAM were demonstrated in the work by Lilienthal et al. (2018). These sources are predominantly solar heating in the 70 troposphere and stratosphere, nonlinear interactions between the DT and SDT in the mesosphere, and GW-tide interactions in the thermosphere.
In contrast to earlier MUAM versions (e.g., Jacobi et al., 2006), where a linear, Lindzen-type GW scheme was applied, we now use a nonlinear whole atmosphere GW scheme according to the work by Yiǧit et al. (2008). An increasing number of whole atmosphere models are being developed around the world, which further indicate the necessity of the use of GW schemes that 75 are suitable for the whole atmosphere region. The whole atmosphere scheme describes the upward propagation of small-scale GWs and their dissipation due to various realistic atmospheric dissipation processes, such as nonlinear interactions, molecular diffusion and thermal conduction, turbulent viscosity, ion drag, and radiative damping. The most dominant dissipative processes in the middle and upper atmosphere are the nonlinear interactions between GWs and the dissipation due molecular processes.
This physical rationale and the detailed description of the scheme are given in a number of publications (e.g., Yiǧit et al., 2008;80 Yigit and Medvedev, 2013;Miyoshi and Yiǧit, 2019). It has been also used in studies of GW effects with Martian GCMs (e.g., Yigit et al., 2018).
Vertical profiles of the Newtonian cooling coefficient, the eddy diffusion coefficient and electron density are currently not used in our implementation of the parameterization, i.e. these parameters are set to zero. The source level of GWs is defined at about 15 km similar to previous implementations of the whole atmosphere scheme (e.g., Yigit et al., 2014;Yiǧit and Medvedev, 85 2017). The GW spectrum specifies the GW momentum fluxes as a function of horizontal phase speeds.
We first conduct a reference (benchmark) simulation, which will later facilitate an assessment of the changes induced by variations in the GW source spectrum. The benchmark case (called EXP1 hereafter) is generated by spinning up the model for a period of 390 days, in which the mean circulation is built up and different waves such as GWs (after day 60), planetary waves and tides (after day 180) are included. After the spin-up, we run the model for 30 days, which are used for our analysis of the 90 background dynamics (zonal/meridional wind and temperature) and wave parameters shown in sections 3-5. Because MUAM is driven only by monthly mean boundary conditions and reaches almost a steady state with small day-to-day variations after the spin-up period, the average of these last 30 days represents the monthly mean state of the atmosphere.
For the GW spectrum of EXP1, we adapt the original spectrum by Yiǧit et al. (2008). It includes a total of nh = 30 harmonics with the horizontal phase speeds c i ranging between ±2 and ±80 m s −1 . The peak momentum flux at the source level is 95 u w (z 0 ) = 0.00025 m 2 s −2 and the full width at half maximum (FWHM) of the spectrum function is located at c w = 35 m s −1 .
The horizontal wavelength of GWs is assumed to be λ H = 300 km. This spectrum has also been used in a number of recent publications (e.g., Yiǧit et al., 2009Yiǧit et al., , 2012Miyoshi and Yiǧit, 2019).
In an additional experiment (EXP2) we retained the properties of the GW spectrum, except for increasing the peak momentum flux at the source level to u w (z 0 ) = 0.00035 m 2 s −2 . For a third experiment (EXP3), the peak momentum flux was 100 the same like in EXP1, but the spectrum was modified now including nh = 34 harmonics and a FWHM of c w = 26 m s −1 .
Thereby, the total momentum flux at the source level, i.e. i u w i (z 0 ) is the same for EXP3 like for EXP1. The GW spectra for EXP1, EXP2, and EXP3 are presented in Fig. 1.

Background Circulation
We first study the results of the benchmark simulation (EXP1) based on the standard GW spectrum. Figure 2 shows altitude-105 latitude cross sections of the monthly mean zonal mean (a) zonal (u) and (b) meridional wind (v) as well as the (c) neutral temperature (T ) for Northern Hemisphere (NH) winter conditions (January). Figure 3 shows the GW effects in the same manner for (a) zonal GW drag, (b) meridional GW drag, and (c) GW heating/cooling. A strong westerly wind system exceeding 80 m s −1 is prevalent in the middle atmosphere in the NH, while summer easterlies dominate the Southern Hemisphere (SH).
These middle atmosphere jets extend up to an altitude of about 90 km (Fig. 2a), and significantly influence the upward propaga-110 tion condition of small-scale GWs via wave filtering and critical level interactions. Thus, in the NH (winter) mainly westward directed GWs can propagate into the upper atmosphere, while the eastward directed GWs are substantially damped or largely filtered out. The opposite phenomenon prevails in the summer Southern Hemisphere (SH). This distribution of GW drag is clearly seen by the eastward GW drag in the SH and westward GW drag in the NH between 80-100 km (Fig. 3a), which is primarily responsible for the reversal of the zonal mean flow above 90 km, from eastward to westward direction (∼ −40 m s −1 ) cells leading to the reversals in the zonal and meridional winds in the MLT are driven by GW dynamics. The associated changes in the residual mean circulation lead to the adiabatic cooling and warming of the summer and winter hemispheres, respectively, as can be seen by the temperature distribution in the MLT (Fig. 2c). Owing to the strong northward winds at around 80 km, air descends and warms up the winter mesopause adiabatically (∼ 180 K), while it ascends and cools down adiabatically the summer mesopause (< 150 K). During polar night, when the polar vortex establishes, there is a lack of incoming solar radiation 130 so that the temperature is decreasing in the NH in the polar stratosphere. However, in the SH, the temperature rises up to 270 K in the polar stratosphere due to the absorption of solar radiation by ozone. Above 120 km in the thermosphere the temperature gradually increases, exceeding, e.g., 900 K in summer due primarily to the enhanced absorption of solar UV and EUV by molecular and atomic oxygen. GW heating and cooling [in K d −1 ] presented in Fig. 3c suggests that the thermal effects of GWs increase with increasing altitude in the thermosphere. The primary thermal effect of GWs is to cool the thermosphere 135 above 120 km altitude, which is in agreement with previous studies (Yigit and Medvedev, 2009;Yiǧit and Medvedev, 2017).
The thermal effects of GWs are produced by the combination of the frictional heating and the differential (dynamic) cooling by dissipating GWs (Medvedev and Klaassen, 2003;Yiǧit et al., 2008).
We compare our results with reference climatologies such as the Upper Atmosphere Research Satellite (UARS) Reference Atmosphere Project (URAP; Swinbank and Ortland, 2003) or the Committee on Space Research (COSPAR) International being comparable to CIRA-86 and HWM-14, but about 20 m s −1 larger than the jet maximum proposed by GEWM and URAP. Note, however, that GEWM is only available for altitudes between 70 and 100 km (slightly above our jet maximum) and URAP is interpolated for large areas near the jets. Due to the GW scheme after Yiǧit et al. (2008), included in the present simulations, the mesospheric jet has strengthened by additional 20 m s −1 and the jet maximum is slightly shifted towards the 150 North compared to the MUAM simulations presented in the work by Samtleben et al. (2019). The zonal wind jet is of similar magnitude as the one simulated by Miyoshi and Yiǧit (2019), but again stronger than the one predicted by the WACCM-X model (Qian et al., 2019). As in WACCM-X, the magnitude and altitude-latitude structure of the meridional MLT wind jet in MUAM is comparable to the ones predicted by the radar-based GEWM (Portnyagin et al., 2004;Jacobi et al., 2009). Our results can also be compared with meteor radar retrievals of mesospheric winds in the ∼ 85-90 km range. For example, in the 155 zonal wind climatology determined by Pramitha et al. (2019) it is seen overall that weak eastward winds reverse direction to westward at higher altitudes during January at equatorial and low-latitudes in the NH, which is also seen in our simulations. GW generation processes in the lower atmosphere are complex and a number of processes contribute to the formation of the GW spectrum. We next would like to test the response of the GCM to variations in the initial GW spectrum. Two more 160 experiments have been performed as has been described in section 2. The same mean fields are analyzed and presented for these simulations as has been done for EXP1. The field of a respective experiment is shown as contour lines, while the differences with respect the benchmark run (i.e., EXP2-EXP1 and EXP3-EXP1) are given in color shading.

Modified GW Spectrum: Increased Flux at the Source Level
The impact of an increased GW flux at the source level on the mean circulation is shown in the difference plots in Fig. 4. 165 Fig. 4a indicates that right below the region of eastward wind reversal in the SH and around the westward wind reversal in the NH, the mean zonal wind has become relatively westerly and easterly, respectively. Studying the associated GW drag results in Figure 5a can provide some insight into this result. Increasing the source maximum momentum strength affects primarily the dissipation altitude, thus the saturation level of individual GW harmonics as well as the associated drag produced by them.
Larger momentum flux means that GWs dissipate at lower altitude due to increased nonlinear interactions in the MLT, which 170 then enhances the mesospheric reversals in both hemispheres, as a consequence of the increased eastward and westward GW drag in the SH and NH around 80-100 km. The secondary enhancement of the westward GW drag takes place in the lower thermosphere, for example, in the SH around 120 km, owing to the enhanced dissipation of the surviving faster GW harmonics due to increasing molecular viscosity with height, which has been previously discussed.
The meridional circulation presented in Fig. 4b shows a strengthening of the mesospheric circulation by 1-2 m s −1 , i.e. a stronger southward wind as well as an enhancement of the lower thermospheric mean northward circulation, which are primarily driven by the intensification of the mesospheric and lower thermospheric GW momentum deposition. Owing to the enhanced mesospheric meridional circulation the upward (downward) movement in the polar region in the SH (NH) intensifies, which leads to a stronger adiabatic cooling (warming). This effect can be seen in the mesospheric temperature differences ( Fig. 4c), which are negative (positive) in the SH (NH) around 80 km in the polar region with up to −3 K (+5 K). A similar 180 effect but in the opposite sense happens in the lower thermosphere (around 110 km), where there is a relative adiabatic warming in the SH polar latitudes but a relative cooling in the NH middle-to high-latitudes. Higher up in the thermosphere (> 120 km) the combined effect of the changes in the GW-induced mean meridional circulation and GW heating/cooling lead to a slight relative cooling of a few K d −1 with respect to the control simulation, while the changes in the GW induced heating/cooling are overall relatively small (Figure 5c), despite the significant increase in the initial source strength. The only exception is the NH 185 polar thermosphere below 100 km, where the GW heating/cooling characteristic dipole structure with respect to the altitude shows some intensification.

Modified GW Spectrum: Same Total Flux
We next present the simulation results of EXP3, in which we have increased the number of wave harmonics from 30 to 34, keeping the maximum phase speed, the total momentum flux, and the peak momentum flux the same as in the benchmark case, 190 EXP1. For this purpose, the FWHM was reduced to c w = 26 m s −1 . This adjustment has shifted the phase speeds to slightly larger values in the tail of the spectrum, while significantly decreasing the momentum flux they carry (Fig. 1). The impact of the spectral changes, primarily in the tail of the spectrum, which is now populated with slightly faster waves but with smaller wave fluxes can be seen in the zonal wind (Fig. 6a). It shows that the height of the wind reversal is slightly shifted upward and its magnitude is reduced, indicated by the relative negative/positive zonal wind differences, EXP3-EXP1, around 100 km 195 in the SH and NH, respectively. Overall, the strengths of the GW momentum flux controls the magnitude of the zonal wind reversal in the MLT. Higher up in the NH lower thermosphere, the easterlies become weaker by more than 20 m s −1 compared to EXP1. These changes can be explained by the reduction of the momentum flux in the tail of the spectrum in EXP3 relative to EXP1.
The meridional circulation (Fig. 6b) in the mesosphere as well as in the thermosphere also shows a weakening. The south-200 ward wind (northward wind) between 80 and 100 km (100 and 120 km) is reduced by more than 2 m s −1 (4 m s −1 ). The weakening of both circulation patterns affects the intensity of the downward and upward movements in the polar region, which is therefore also less pronounced. This leads to a colder (warmer) winter (summer) mesopause, which can be seen in the negative (positive) temperature anomalies of −3 K (+3 K) shown in Fig. 6c. In the polar thermosphere, the effect is even stronger with more than +5 K temperature difference. These changes are primarily controlled by the associated changes in the zonal 205 GW drag, rather than by the changes in the meridional GW drag. The meridional GW drag anomalies are confined to the polar latitude in the NH MLT (Fig. 7b). Up to 8-10 K d −1 relative reduction in the resultant GW cooling controls the thermal budget of the lower thermosphere to a much lesser extent than the dynamical changes.
Between 100 and 120 km the zonal GW drag difference (Fig. 7a) shows a strong positive anomaly of more than +50 m s −1 d −1 , which may be an effect of the vertical shift of the eastward directed zonal GW drag on the SH. While the zonal GW drag anoma-210 lies are stronger on the SH, the meridional GW drag anomalies (Fig. 7b) are more pronounced in the NH, especially in the polar thermosphere.We observe a weakening of the meridional GW drag of more than 20 m s −1 d −1 . Thereby, the region of southward directed GW drag near 120 km in EXP1 (Fig. 3b) disappears in EXP3 (Fig. 7b). The northward directed GW drag near 110 km, however, persists. With respect to the GW heating anomalies (Fig. 7c), the induced cooling in the thermosphere is strongly reduced by more than 10 K d −1 . Thus, the EXP3 cooling is half the EXP1 cooling.

215
Note that in EXP2, we had increased the source peak momentum flux by more than a factor of 1.5. This corresponds to increasing the momentum flux of each GW harmonic by 50% as well as the total momentum flux by 50%. In EXP3, however, we have kept the total momentum flux constant, as well as the peak momentum flux and the phase speed range, varying only the FWHM of the spectrum and the number of harmonics. Interestingly, the response of the mesosphere and lower thermosphere in terms of changes in circulation patterns and temperature distributions are stronger in EXP3 than in EXP2, which emphasizes 220 the dynamical significance of the faster (nonorographic) GWs for the structure of the mesosphere and lower thermosphere. This means that uncertainty in the peak momentum flux of the GW spectrum has less impact than the uncertainty in the spectral shape of the spectrum on the structure of the atmosphere up to the lower thermosphere.

Relation between Gravity Waves and Terdiurnal Tides
The migrating TDT amplitudes and phases in the zonal wind, meridional wind, and temperature are shown in Fig. 8 for the 225 benchmark simulation EXP1. The amplitudes of all components have a maximum in the lower thermosphere between 120 and 140 km, being larger in the summer SH than winter NH. Ground based radar observations and satellite measurements of the TDT activity have overall focused on the MLT region between 80 to 110 km and have reported an autumn and winter maximum of the TDT amplitudes (e.g., Beldon et al., 2006;Jacobi, 2012;Liu et al., 2019;Pancheva et al., 2013;Yue et al., 2013). The TDT amplitudes simulated by the GCM are also larger in the winter hemisphere than the summer hemisphere in the MLT,230 in relatively good agreement with these measurements. Vertical wavelengths in the MLT and above, taken from the vertical TDT phase gradients in Fig. 8d-f, are on the order of 30 km in the summer hemisphere, but larger in winter. This agrees with observations (Thayaparan, 1997;Jacobi, 2012).
The most important dynamical feature in the lower thermosphere in our modeling study are GWs of lower atmospheric origin parameterized by the whole atmosphere GW parameterization, and therefore they are the most obvious candidate for 235 tidal modulation at those heights. We next want to examine their influence on the TDT. The terdiurnal components of the GW parameters presented in Fig. 9 can be used as a proxy for GW-TDT interactions. The amplitudes of terdiurnal GW drag and heating maximize near 120-130 km, where TDT amplitudes in temperature and wind maximize as well. This indicates that the TDT in the thermosphere is strongly influenced by the momentum deposition of dissipating GWs. We also present phases of terdiurnal GW drag and heating in Fig. 9d-f. Similar to the wind and temperature TDT (Fig. 8d-f), the phases of the terdiurnal 240 GW effects are rather irregular in the mesosphere but become more organized in the lower thermosphere with longer vertical wavelengths.
Figures 10a-c show the TDT amplitudes for the simulation EXP2 (contour lines), i.e. the one with increased GW source flux, and their differences with respect to EXP1 (color shading). Their general structure in wind and temperature is similar to that of EXP1 (Fig. 8a-c). Differences EXP2-EXP1 amount to ±1 to ±2 m s −1 and K. The terdiurnal signature in GW 245 parameters ( Fig. 10d-f) are mainly increased near their maxima, which can be interpreted as a direct result of the increased source momentum flux in simulation EXP2. There are also some negative changes between EXP2 and EXP1 but these are rather irregular and most likely a result of the slightly altered background dynamics described in section 4.1, influencing wave propagation conditions in general.
Furthermore, there is good agreement between the TDT amplitude changes in the zonal wind (Fig. 10a) and zonal GW drag 250 (Fig. 10d). For example, the latter one is increased at northern low latitudes at about 120-140 km and decreased at southern low latitudes at a similar altitude. This pattern is visible in TDT zonal wind amplitude changes, as well. In the meridional components (Fig. 10b,e), the positive/negative change patterns also largely agree. In the thermal component (Fig. 10c,f), however, the TDT temperature amplitude seems to be damped where the difference in terdiurnal GW heating is positive. Figure 11 is similar to Fig. 10, but refers to the differences between EXP3 and EXP1. This means that differences are based 255 on a GW spectrum with more harmonics in EXP3 than in EXP1, but a smaller momentum flux for each individual harmonic, in particular in the mid-range phase speeds (see Fig. 1). As a result, the sign in the EXP3-EXP1 differences of terdiurnal GW parameters is opposite to that of EXP2-EXP1 differences, i.e. it is negative almost everywhere (see Fig. 11d-f). Accordingly, the TDT wind amplitudes (Fig. 11a,b) are also mostly smaller in the EXP3 simulation compared to EXP1, decreasing by about 1-2 m s −1 in the area of maximum amplitudes. They partly increase by a similar magnitude, but the decrease dominates, 260 especially in the zonal wind component. Similar to the differences of EXP2-EXP1, the relation between terdiurnal GW heating ( Fig. 11f) and TDT temperature amplitude (Fig. 11c) is less clear and we also observe large patches of positive amplitude changes up to 1.5 K in the temperature component that seem to be corresponding to negative GW heating.

Summary and Conclusions
We have successfully implemented a whole-atmosphere GW paramterization according to Yiǧit et al. (2008) into the mech-265 anistic MUAM GCM. The zonal mean horizontal wind patterns in the middle atmosphere as well as the global temperature distribution is realistically reproduced by the model, and reasonably agree with respect to established climatologies and other GCM predictions. We have successfully simulated the self-consistent direct GW penetration into the thermosphere as has been previously done by other GCMs using the whole atmosphere scheme (Yiǧit et al., 2009;Miyoshi and Yiǧit, 2019). Compared to earlier MUAM versions that use a highly tuned linear Lindzen-type GW parameterization for the middle atmosphere (e.g. In the lower thermosphere, the main impact of the whole atmosphere GW parameterization is to drive a vertical/meridional circulation, mainly by the GWs with high phase speeds, which are able to propagate through the MLT wind jets. We performed two further experiments to investigate the response of the middle and upper atmosphere dynamics to changes 275 of the initial GW spectrum and horizontal momentum flux. In the first experiment (EXP1), we increased the momentum flux at the source level by roughly 50% of the original value, while keeping the spectral shape unchanged. As expected, these modifications result in a stronger GW dissipation at lower levels (near 80 km), connected with a slightly intensified wind reversal in the MLT. This intensifies both the upper mesosphere and lower thermosphere meridional winds, related to stronger adiabatic warming/cooling of the winter mesosphere/lower thermosphere and reverse effects in the summer hemisphere. However, it 280 turns out that the relatively drastic change of the peak source momentum flux does not influence the global dynamical patterns by such a large degree. One may conclude that an adjustment of this parameter is not crucial for the implementation of this GW parameterization into a GCM like MUAM.
In the final experiment (EXP3), the total momentum flux was kept constant with respect to the benchmark case EXP1, but the total number of harmonics was increased and the width of the spectrum was decreased, yielding smaller fluxes for 285 the high phase speed tail of the spectrum. In this experiment, the mesospheric wind system was less affected, but the lower thermospheric jets were weakened, connected with cooling/warming of the winter/summer mesosphere, and reversed thermal effects above. Considering EXP2 and EXP3 with respect to the benchmark case, in general, the lower thermospheric circulation and temperature distributions responded more strongly to the changes in the spectral shape of the GW spectrum than to the increase in the peak momentum flux. This emphasizes the dynamical significance of the faster GWs for the structure of the 290 lower thermosphere, as these waves are less affected by dissipation and filtering processes in the stratosphere and mesosphere, and thus penetrate into the thermosphere (e.g. Yiǧit et al., 2008Yiǧit et al., , 2009Yiǧit and Medvedev, 2017).
We also investigated the effect of the GW source distribution on the amplitudes of the TDT. The simulated latitudinal-vertical distribution of TDT amplitudes, and their vertical wave structure was found to be realistic when compared with observations.
Modifications of the GW source spectrum or total momentum flux mainly influences the TDT in the lower thermosphere, 295 and to a much lesser degree in the upper mesosphere. We found that increasing the GW momentum flux essentially leads to an increased terdiurnal variation in the GW drag and increased amplitudes in the lower thermosphere, which we interpret as a direct result of increased momentum flux of fast GWs. In turn, narrowing of the width of the GW spectrum mostly leads to a reduced terdiurnal variation in the GW drag and lower TDT amplitudes, as a consequence of the reduced momentum flux of fast GWs in the lower thermosphere. Miyahara and Forbes (1991) already demonstrated in their simulations that an 300 interaction between the DT and GWs can generate a secondary TDT. Recently, GW-tide interactions were again underlined as an important excitation mechanism of the TDT (Lilienthal et al., 2018;Lilienthal and Jacobi, 2019). With respect to DT and SDT amplitudes, Ribstein and Achatz (2016) pointed out that tidal amplitudes are very sensitive to the respective implemented GW parameterization. Thus, it is reasonable that a modified GW drag also strongly influences the TDT amplitudes as has been shown in our simulations.

305
Our modeling experiments highlight the importance of taking into account self-consistent propagation of GWs and their dissipation in the mesosphere and lower thermosphere. GWs play a crucial role for the mesospheric and thermospheric circulation patterns and temperature variations as well as for the TDTs. It is noteworthy that small modifications of momentum fluxes of fast GWs have substantial impact on the lower thermosphere mean circulation, and also on tidal amplitudes in particular with respect to the TDT.   10 −4 K d −1 Figure 11. Same as Fig. 10 but for EXP3 (contour lines) and differences EXP3-EXP1 (color shading).