Zhaoyang Chai a , b , Minghua Zhang c , * , Qingcun Zeng a , b , Jinbo Xie d , Ting You b , e , He Zhang a , b

a International Center for Climate and Environment Sciences, Institute of Atmospheric Physics, Chinese Academy of Sciences, Beijing, China

b College of Earth and Planetary Science, University of Chinese Academy of Sciences, Beijing, China

c School of Marine and Atmospheric Sciences, State University of New York at Stony Brook, NY, USA

d State Key Laboratory of Numerical Modeling for Atmospheric Sciences and Geophysical Fluid Dynamics, Institute of Atmospheric Physics, Chinese Academy of Sciences,Beijing, China

e Center for Monsoon System Research, Institute of Atmospheric Physics, Chinese Academy of Sciences, Beijing, China

Keywords:QBO Gravity waves Momentum budget Zonal wind

ABSTRACT The quasi-biennial oscillation (QBO), a dominant mode of the equatorial stratospheric (～100—1 hPa) variability,is known to impact tropospheric circulation in the middle and high latitudes. Yet, its realistic simulation in general circulation models remains a challenge. The authors examine the simulated QBO in the 69-layer version of the Institute of Atmospheric Physics Atmospheric General Circulation Model (IAP-AGCML69) and analyze its momentum budget. The authors find that the QBO is primarily caused by parameterized gravity-wave forcing due to tropospheric convection, but the downward propagation of the momentum source is significantly offset by the upward advection of zonal wind by the equatorial upwelling in the stratosphere. Resolved-scale waves act as a positive contribution to the total zonal wind tendency of the QBO over the equator with comparable magnitude to the gravity-wave forcing in the upper stratosphere. Results provide insights into the mechanism of the QBO and possible causes of differences in models.

1. Introduction

The quasi -biennial oscillation (QBO) of the zonal-mean zonal wind is the primary mode of variability of the tropical stratosphere. The oscillation consists of easterly and westerly wind regimes alternating in time and propagating downward from about 1 to 100 hPa, with an average period of 28 months ( Naujokat, 1986 ; Baldwin et al., 2001 ) and a downward propagating speed of ～1 km/month ( Yu et al., 2017 ). It is well known that the QBO is produced by the interactions between waves and mean flow (e.g., Lindzen and Holton, 1968 ) and mainly driven by the deposition of momentum in the shear zones, which is carried by vertically propagating tropical eastward and westward propagating waves (e.g., Holton and Lindzen, 1972 ). The QBO has been shown to impact the strength of the stratospheric polar vortex (the polar stratospheric variability) by modulating the propagation of larger-scale waves into the polar region, which is well-known as the Holton—Tan mechanism ( Holton and Tan, 1980 ). The QBO also impacts the distribution of stratospheric chemical constituents such as ozone and water vapor( Randel and Wu, 1996 ). Observational studies have suggested that the QBO can impact the extratropical storm tracks ( Wang et al., 2018 ) and modulate the Madden—Julian Oscillation ( Kim et al., 2020 ). Also, the QBO affects mesospheric variability by selectively filtering waves that propagate upward through the tropical stratosphere ( Baldwin et al.,2001 ).

Although the QBO mechanism has been understood for many years,it remains a challenge to simulate the QBO in general circulation models(GCMs), with only a few GCMs being able to reproduce it. In phase 5 of the Coupled Model Intercomparison Project, only five models could generate the QBO internally ( Schenzinger et al., 2017 ; Butchart et al.,2018 ). As the QBO results from the wave—mean flow interactions, the waves —including tropical Kelvin waves, mixed Rossby gravity waves,inertial-gravity waves, and small-scale gravity waves (GWs) —need to be correctly represented to simulate a realistic QBO. However, many GCMs still cannot simulate a realistic spectrum of tropical waves because of their low resolution and their deficiencies in the parameterization of small-scale GW forcing ( Ricciardulli and Garcia, 2000 ; Lott et al.,2014 ). Studies have suggested that an adequately fine vertical resolution(vertical grid spacing of ～500—700 m) of the troposphere and lower stratosphere is also necessary to simulate the QBO due to the forcing of some resolved waves with small vertical wavelength and the need to capture the wind shear ( Giorgetta et al., 2006 ; Richter et al., 2014 ;Geller et al., 2016 ). Garcia and Richter (2019) showed that both the 70-layer and 110-layer versions of the Whole Atmosphere Community Climate Model (WACCM) can simulate the QBO, but the simulated QBO in the 70-layer model is relatively weak and fails to propagate to below 50 hPa. It has been shown that the period of the simulated QBO relates to both the characteristic vertical wavenumber ( Scaife et al., 2000 ) and the GW source strength ( Scaife et al., 2000 ; Giorgetta et al., 2006 ).

The purpose of this study is to understand the momentum budget in generating the QBO by using simulations from the Institute of Atmospheric Physics Atmospheric General Circulation Model, version 4.1(IAP-AGCM 4.1). The dominant zonal momentum sources and sinks are diagnosed to highlight the roles of the Brewer—Dobson circulation, nonlinear transport, and resolved-scale waves in addition to parameterized GWs. A companion paper will report sensitivities of the model’s QBO to model layers, vertical resolution, model top height, and forcing configuration based on the understanding from the momentum budget analysis.The paper is organized as follows. Section 2 describes the model and its configuration. Section 3 presents the results, including an evaluation of the simulated QBO and analysis of the momentum budget. Section 4 contains a summary and discussion.

2. Model and configuration

2.1. Model description

Full details about the model used here, IAP-AGCM 4.1, are given in Zhang (2009) , Sun et al. (2012) , and Zhang et al. (2013) . Briefly,IAP-AGCM 4.1 uses a horizontal resolution of ～1.4°×1.4° and 30 vertical layers, with the model top at ～2.26 hPa. While the newer version, IAP-AGCM 5.0, uses new physical parameterizations of convection, clouds, and other updates ( Zhang et al., 2020 ), IAP-AGCM 4.1 uses the full physical package from the Community Atmosphere Model, version 5 (CAM5; Neale et al., 2012 ), with only minor adjustments to certain parameters. The parameterization of GWs specifies three dominant wave sources separately —orography, fronts, and convection —following Richter et al. (2010) . The regimes of wave propagation, saturation,breaking, and momentum deposition to the mean flow are based on Lindzen (1987) . The orographic GW drag source parameterization follows McFarlane (1987) and non-orographic GW source parameterizations follow Richter et al. (2010) and Beres et al. (2004 , 2005 ) as implemented in the WACCM version of CAM5.2. Convectively generated GWs are parameterized to be launched at the top of the convection whenever the deep convection parameterization ( Zhang and McFarlane, 1995 ) is activated. The amplitude of convective GWs is parameterized to be proportional to the square of the maximum tropospheric heating rate due to convection. In particular, the tunable parameter called the “efficiency factor ”in the convective GW parameterization can be interpreted as the fraction of time during which the model grid-box is covered by GWs.Richter et al. (2010 , 2014 ) used a varying efficiency factor from 0.1 to 0.55 for different grid-box sizes of the model, and it is set to 0.5 in this study.

Fig. 1. The (a) vertical grids and (b) vertical grid spacings used in the 69-layer version of IAP-AGCM (IAPL69). The vertical spacings are about 0.5 km between 800 hPa and 40 hPa, and gradually increase to about 4.6 km near the model top(～0.01 hPa).

2.2. Model configuration

In this study, the default 30-layer IAP-AGCM is reconfigured to a 69-layer model —namely, IAP-AGCML69 (referred to simply as IAPL69 hereafter). Fig. 1 shows the vertical grid distribution and its grid spacings for IAPL69. The vertical resolution is about 0.5 km in the free troposphere and lower stratosphere. The grid spacing then gradually increases to about 4.6 km near the model top at ～0.01 hPa. The choices of 69 levels and the model top are somewhat arbitrary but are guided by the model’s ability to simulate the QBO. Sensitivities to these choices,including the performance of a 91-layer model ( Chai et al., 2020 ), will be reported in a separate study. The experimental simulation is carried out for the AMIP run type from 1979 to 2005 with prescribed observed global Hadley Centre Sea Surface Temperature and Sea Ice ( Taylor et al.,2012 ).

2.3. Data

All analyses are based on the monthly mean output from IAPL69 and the monthly mean reanalysis from the ECMWF’s third-generation reanalysis product, ERA-Interim (ERAI; Dee et al., 2011 ). The ERAI data have 37 pressure levels up to 1 hPa and have a horizontal resolution of 1.5°×1.5°. Zonal wind, meridional wind, temperature, and surface pressure of ERAI are used in our analyses.

3. Results

3.1. The simulated QBO in the IAPL69

The observed and simulated zonal winds are compared in Fig. 2 by using a 27-year time—pressure cross section of the monthly mean averaged from 2°S to 2°N. Fig. 2 (a) shows that in ERAI the typical QBO has an average period of about 28 months. The average of the maximum easterly and westerly is around ? 40 m sand 20 m sat ～20 hPa, respectively, implying an asymmetric oscillation of the zonal wind. Fig. 2 (b)shows that in IAPL69 the period of the zonal wind oscillation is about 27 months, which is very close to that of ERAI (～28 months), but the westerly of the QBO is much stronger (～35 msat 7 hPa) than that of ERAI (～20 msat 20 hPa), so the easterly is squeezed to a thin band (persisting for ～9 months at 7 hPa vs. ～16 months at 20 hPa in ERAI) with insufficient downward penetration. This stronger westerly phase also exists in the 60-layer version of CAM5 ( Richter et al., 2014 ),which may be associated with the stronger forcing from the parameterized convectively generated eastward-propagating GWs. The peak value of the amplitude of the simulated oscillation is above 10 hPa, which is mixed with semi-annual oscillation near 1 hPa, while the peak value for ERAI is near ～20 hPa.

The similarities and differences of the simulated QBO amplitude with that in ERAI are more clearly shown in Fig. 3 (a, b). The model captures the overall structure of the zonal wind amplitude ( Fig. 3 (a)), including the poleward extension to ～15°S/N, but this extension is a little narrower at lower levels relative to observation, as seen in many GCMs( Schenzinger et al., 2017 ). The model has a larger amplitude, with a maximum of ～33 msvs. ～28 msin ERAI ( Fig. 3 (b)). In addition,the zonal wind amplitude peaks between ～3 and 10 hPa in IAPL69 vs. between ～10 and 25 hPa in ERAI. The amplitude distribution of the accompanying temperature in the model is compared with ERAI in Fig. 3 (c, d). Consistent with the zonal wind difference, IAPL69 simulates a stronger amplitude in the upper stratosphere and weaker amplitude in the lower stratosphere than ERAI. The maximum amplitude in the model is ～4.6 K vs. ～3.2 K in ERAI. Additional experiments indicate drastic sensitivity of these amplitudes to the configuration of model vertical resolution, top height, and the parameterized strength of the wave sources.The higher position of the simulated QBO may be linked with the coarse vertical resolution above 40 hPa in the model, which is too coarse to accurately represent the magnitude of the vertical shear of the zonal wind.Garcia and Richter (2019) showed the difference in the simulated QBO between the 70-layer and 110-layer versions of WACCM. The higherresolution model better captured the position of the QBO center. This is also consistent in other GCMs ( Schenzinger et al., 2017 ).

Fig. 2. Time—pressure cross section of the monthly mean zonal wind (filled contours; ms?1 ) averaged from 2°S to 2°N for the years 1979—2005 for (a) ERAI and(b) IAPL69. Red represents positive (westerly) winds and green negative (easterly). Bold solid, normal solid and dashed black contour lines are 0, 20, and? 20 ms?1 , respectively. Maximum and minimum values (rounded) of zonal wind are given in the square brackets on the top-left side of each plot.

Fig. 3. Amplitudes of de-seasonal (a, b) zonal wind (ms?1 ) and (c, d) temperature (K) signals of the QBO as a function of latitude and pressure in (a, c) ERAI and (b, d) IAPL69. Maximum and minimum values (rounded) are indicated in the square brackets on the top side of each plot.

3.2. Momentum forcing of the QBO

The zonal momentum equation with a pressure coordinate for IAPAGCM can be written as

where Fand Fare u tendencies due to model dynamics and physics,respectively, and t is time. The dynamic term Fcan be written as

Fig. 4. Time—pressure cross section of simulated QBO-filtered (27-month period) zonal wind tendencies (m s?1 /month; contours shaded red and blue for positive and negative respectively), averaged from 2°S to 2°N in (a—d) and from 10°S to 10°N in (f—i), due to (a, f) GWD (parameterized convective GWs), (b, g) WAV (resolved waves, including vertical and horizontal wave momentum transports), (c, h) ADV (advection, including vertical and horizontal components), and (d, i) SUM (the sum of all the tendencies). Panels (e, j) are the residual terms. The QBO-filtered zonal wind is also plotted (unfilled contour lines; ms?1 ; contour interval of 5 ms?1 ;solid and dashed lines are positive and negative respectively; zero-lines in heavy bold type). Maximum and minimum values (rounded) of each tendency are given in the square brackets on the top side of each plot.

The zonal wind tendency due to resolved-scale waves ( F) is written as

where F, F, F, F, and Fare tendencies due to deep convection(dc), shallow convection (sc), vertical diffusion (vd), Rayleigh Friction(rf), and gravity waves (gw). Here, Fcan be further split into three parts as

where F, F, and Fare u tendencies due to convective, frontal,and orographic GWs correspondingly.

In the following, we show the contribution of these processes to the total zonal wind tendency in generating the QBO.

The forcing from GWD (parameterized convective GWs), ADV (dynamical advection), and WAV (resolved waves) play different roles in the evolution of the QBO. Fig. 4 shows these terms averaged from 2°S to 2°N in the upper row and from 10°S to 10°N in the bottom row. The total tendency (SUM: Fig. 4 (d, i)) is plotted over the zonal wind. Since only the peak Fourier component and phase (period of 27 months) are plotted, an identically repeating pattern is formed, which is shown for five years. The SUM tendency peaks in the shear regions as expected,which is relatively homogeneous from 50 to 1 hPa ( Fig. 4 (d, i)), indicating a phase shift of π∕2 between the SUM tendency (filled contours) and the zonal wind (unfilled contours). It drives the evolution of the QBO wind, forcing downward propagation of the westerly wind and easterly wind alternately. The SUM tendency averaged between 2°S and 2°N is a little larger, with a maximum of 5.4 ms/month vs. 4.7 ms/month averaged between 10°S and 10°N.

The GWD tendency ( Fig. 4 (a, f)) peaks at a phase lag of less than π∕4 behind the SUM tendency and dominates above 10 hPa. The magnitude averaged between 10°S and 10°N ( Fig. 4 (f)) is much larger, with a maximum of 24 ms/month, than that averaged between 2°S and 2°N, with a maximum of 15 m s/month ( Fig. 4 (a)), which may be associated with the meridional double-peak structure of convection in the tropics.The GWD forcing contributes a lot to the SUM tendency and is partly canceled out by other processes.

For the WAV tendency, there are some differences between the 2°S—2°N average ( Fig. 4 (b)) and the 10°S—10°N average ( Fig. 4 (g)). The former ( Fig. 4 (b)) peaks at a phase lead of less than π∕4 in advance of the SUM tendency, but the latter ( Fig. 4 (g)) peaks a phase lead of almost π∕2 before the SUM tendency ( Fig. 4 (i)). This phase shift between Fig. 4 (b)and (g) results from the meridional structure of the WAV (not shown).For the average between 2°S and 2°N, the WAV forcing makes a positive contribution to the SUM tendency, but for the average between 10°S and 10°N it makes no contribution to the SUM but decelerates the zonal wind. As a result, the WAV forcing drives the descent of the QBO phase together with the GWD forcing over the equator and spreads out the easterly or westerly momentum poleward to the extratropics.

Table 1 Pearson linear cross-correlations among GWD, WAV, ADV, SUM, and U, averaged for 2°S—2°N and 10°S—10°N, respectively. For example, the figure 0.760 represents the correlation between SUM and WAV of 2°S—2°N; the figure ? 0.233 represents the correlation between SUM and WAV of 10°S—10°N. The figure in parentheses is the original amplitude of GWD, WAV, ADV, SUM, and U. For GWD, WAV, ADV, and SUM, the unit is m s ?1 ∕ month , and for U it is m s ?1 .The figure in square brackets is the corresponding phase difference, calculated as the inverse cosine of the correlation.

A more quantitative description of the forcing contribution to the simulated QBO at 5 hPa and phase relations is given by using the projection of the forcing components onto the SUM in Table 1 . Since the Pearson linear correlation coefficient is the cosine of the phase difference, the projection is simply the amplitude of the forcing multiplied by the correlation coefficient. At the equator, GWD ( 10 . 4 ×0 . 909 =9 . 5 m s∕ month ) and WAV ( 8 . 6 ×0 . 760 = 6 . 5 m s∕ month ) are the leading contributions ( 6 . 5 + 9 . 5 = 16 . 0 m s∕ month ) to drive the QBO,forcing the descent of the zonal wind shear, while ADV offsets these contributions ( 11 . 2 ×( ?0 . 999 ) = ?11 . 2 m s1 ∕ month ) with a remainder of 4.8 m s1 ∕ month , which is close to 5.3 m s1 ∕ month shown in Table 1 . However, in the tropics between 10°S and 10°N, GWD becomes the sole leading contribution ( 15 . 8 ×0 . 957 = 15 . 1 m s1 ∕ month )to drive the QBO, while ADV offsets this contribution ( 10 . 4 ×( ?0 . 992 ) =?10 . 3 m s1 ∕ month ) with a remainder of 4.8 m s1 ∕ month , which is close to 4.6 m s1 ∕ month shown in Table 1 . As a result, the phase shift between ADV and SUM is mainly due to the upwelling in the Brewer—Dobson circulation; the phase shift between GWD and SUM is relatively small, which is caused by the leading contribution of convectively generated GWs; the phase shift between WAV and SUM is between～0 and π∕2 , which is mainly linked with the larger convergence of meridional zonal momentum flux by the resolved-scale transient waves in the late half of the easterly and westerly phase than in their early half.

As a qualitative reference, we have made comparisons with the reanalysis results (including ERAI, NCEP-2, and MERRA-2; not shown).Although there are relatively large uncertainties in the forcing terms, in particular the forcing term of the resolved-scale wave (WAV), the basic features of the QBO forcing terms are qualitatively consistent with the model’s simulation as follows: (1) GWD forcing is about 1.1—2.5 times that of SUM in an almost synchronized manner; (2) ADV forcing offsets about 30%—60% of GWD forcing; (3) WAV forcing may have a positive contribution to SUM with a phase shift of 0—π∕2 , with relatively large uncertainties in both phase difference and magnitude.

4. Summary and discussion

The QBO is reasonably simulated in IAPL69, showing downward propagation of easterlies and westerlies with comparable period and meridional structure relative to ERAI reanalysis. The amplitude is larger than observation (～33 msvs. ～28 msin zonal wind) owing largely to the stronger westerly in the model. Also, the maximum of the simulated QBO zonal wind is higher than in ERAI, as seen in many GCMs( Schenzinger et al., 2017 ).

The zonal wind time tendency of the QBO is diagnosed by showing four leading terms: convectively generated gravity wave drag (GWD),advection tendency (ADV) by the zonally averaged meridional circulation, resolved-scale waves (WAV), and residual tendency (RES). We have found that, in the formation of the simulated QBO, GWD, and ADV dominate the zonal wind tendency. GWD contributes directly to the zonal wind tendency and forces the downward propagation of the zonal wind shear. ADV cancels out most of the contribution of GWD forcing via the equatorial upwelling of the Brewer—Dobson circulation in the stratosphere. WAV presents as a positive contribution to the zonal wind tendency of the QBO over the equator together with GWD.

In terms of the phase relationships between the processes, relative to the total tendency SUM, the waveform of GWD lags by a phase of less than π∕4 . The waveform of WAV leads by a phase of less than π∕4 to SUM. The combination of these two processes makes SUM lead by a phase of π∕2 to the zonal wind, which steers the evolution of zonal wind with time. This combination is offset more than half by ADV.

The above analyses suggest that multiple processes are at play in the evolution of the QBO. Differences in any of the leading processes discussed above can result in different features of the simulated QBO or the disappearance of the QBO in the case of sufficiently strong equatorial upwelling, consistent with the previous work of Dunkerton (1991) .Simulation of a more realistic QBO in models likely requires good simulation in all these processes. It would be interesting to examine how the momentum budget in different models compares. A follow-up study will discuss the sensitivities of the forcing processes to different model resolutions and other configurations along with the accompanying QBO.

Funding

This research was supported by the National Major Research High Performance Computing Program of China [grant number 2016YFB0200800 ], the National Natural Science Foundation of China[grant numbers 41630530 and 41706036 ], and the National Key Scientific and Technological Infrastructure project “Earth System Science Numerical Simulator Facility ”(EarthLab).

Acknowledgments

We gratefully acknowledge the assistance of Dongling Zhang in supporting computing resources to facilitate the development of IAP-AGCM.We also thank Jiangbo Jin for providing useful experimental design advice for evaluating the simulations in IAP-AGCM.

Supplementary materials

Supplementary material associated with this article can be found, in the online version, at doi:10.1016/j.aosl.2020.100021 .