Climatological impact of the Brewer–Dobson circulation on the N₂O budget in WACCM, a chemical reanalysis and a CTM driven by four dynamical reanalyses

The Brewer–Dobson circulation (BDC) is a stratospheric circulation characterized by upwelling of tropospheric air in the tropics, poleward flow in the stratosphere, and downwelling at mid and high latitudes, with important implications for chemical tracer distributions, stratospheric heat and momentum budgets, and mass exchange with the troposphere. As the photochemical losses of nitrous oxide (N₂O) are well known, model differences in its rate of change are due to transport processes that can be separated into the mean residual advection and the isentropic mixing terms in the transformed Eulerian mean (TEM) framework. Here, the climatological impact of the stratospheric BDC on the long-lived tracer N₂O is evaluated through a comparison of its TEM budget in the Whole Atmosphere Community Climate Model (WACCM), in a chemical reanalysis of the Aura Microwave Limb Sounder version 2 (BRAM2) and in a chemistry transport model (CTM) driven by four modern reanalyses: the European Centre for Medium-Range Weather Forecasts Interim reanalysis (ERA-Interim; Dee et al., 2011), the Japanese 55-year Reanalysis (JRA-55; Kobayashi et al., 2015), and the Modern-Era Retrospective analysis for Research and Applications version 1 (MERRA; Rienecker et al., 2011) and version 2 (MERRA-2; Gelaro et al., 2017). The effects of stratospheric transport on the N₂O rate of change, as depicted in this study, have not been compared before across this variety of datasets and have never been investigated in a modern chemical reanalysis. We focus on the seasonal means and climatological annual cycles of the two main contributions to the N₂O TEM budget: the vertical residual advection and the horizontal mixing terms.

The N₂O mixing ratio in the CTM experiments has a spread of approximately ∼20 % in the middle stratosphere, reflecting the large diversity in the mean age of air obtained with the same CTM experiments in a previous study. In all datasets, the TEM budget is closed well; the agreement between the vertical advection terms is qualitatively very good in the Northern Hemisphere, and it is good in the Southern Hemisphere except above the Antarctic region. The datasets do not agree as well with respect to the horizontal mixing term, especially in the Northern Hemisphere where horizontal mixing has a smaller contribution in WACCM than in the reanalyses. WACCM is investigated through three model realizations and a sensitivity test using the previous version of the gravity wave parameterization. The internal variability of the horizontal mixing in WACCM is large in the polar regions and is comparable to the differences between the dynamical reanalyses. The sensitivity test has a relatively small impact on the horizontal mixing term, but it significantly changes the vertical advection term and produces a less realistic N₂O annual cycle above the Antarctic. In this region, all reanalyses show a large wintertime N₂O decrease, which is mainly due to horizontal mixing. This is not seen with WACCM, where the horizontal mixing term barely contributes to the TEM budget. While we must use caution in the interpretation of the differences in this region (where the reanalyses show large residuals of the TEM budget), they could be due to the fact that the polar jet is stronger and is not tilted equatorward in WACCM compared with the reanalyses.

We also compare the interannual variability in the horizontal mixing and the vertical advection terms between the different datasets. As expected, the horizontal mixing term presents a large variability during austral fall and boreal winter in the polar regions. In the tropics, the interannual variability of the vertical advection term is much smaller in WACCM and JRA-55 than in the other experiments. The large residual in the reanalyses and the disagreement between WACCM and the reanalyses in the Antarctic region highlight the need for further investigations on the modeling of transport in this region of the stratosphere.

1 Introduction

The Brewer–Dobson circulation (BDC; Dobson et al., 1929; Brewer, 1949; Dobson, 1956) in the stratosphere is characterized by upwelling of tropospheric air to the stratosphere in the tropics, followed by poleward transport in the stratosphere and extra-tropical downwelling. For tracer-transport purposes the BDC is often divided into an advective component, the residual mean meridional circulation (hereafter residual circulation) and a quasi-horizontal two-way mixing, which causes net transport of tracers but not of mass (Butchart, 2014).

The BDC is driven by tropospheric waves breaking into the stratosphere (Charney and Drazin, 1961), which transfer angular momentum and force the stratosphere away from its radiative equilibrium. This is balanced by a poleward displacement of air masses, which implies tropical upwelling and extra-tropical downwelling (Holton, 2004). The residual circulation can be further separated into three branches: the transition branch, the shallow branch and the deep branch (Lin and Fu, 2013). The transition branch encompasses the upper part of the transition layer between the troposphere and the stratosphere (the tropical tropopause layer; Fueglistaler et al., 2009); the shallow branch is a year-round, lower stratospheric, two-cell system driven by breaking of synoptic-scale waves; and the deep branch is driven by Rossby and gravity waves breaking in the middle and high parts of the stratosphere during winter (Plumb, 2002; Birner and Bönisch, 2011). The contributions of different wave types to driving the BDC branches has been quantified using the downward control principle, which states that the poleward mass flux across an isentropic surface is controlled by the Rossby or gravity waves breaking above that level (Haynes et al., 1991; Rosenlof and Holton, 1993), and using eddy heat flux calculations as an estimate of the wave activity from the troposphere (e.g., Newman and Nash, 2000).

The quasi-horizontal two-way mixing is generated by two-way transport due to the adiabatic motion of Rossby waves. In the stratosphere this motion is ultimately combined with molecular diffusion which makes the total process irreversible (Shepherd, 2007). The two-way mixing is stronger in a specific latitudinal region of the winter stratosphere, the “surf zone” (McIntyre and Palmer, 1983), and in the subtropical lower stratosphere all year round (e.g., Fig.1 of Bönisch et al., 2011). The mixing process homogenizes the tracer concentration in the surf zone and creates sharp tracer and potential vorticity (PV) gradients on its edges (in the subtropics and at the polar vortex edge), indicating an inhibition of mixing. For this reason, the subtropics and the polar vortex edge are often called transport barriers (Shepherd, 2007).

The BDC plays a major role in controlling the spatial and temporal distributions of chemical tracers, such as ozone, water vapor, aerosols and greenhouse gases, as well as in coupling stratospheric processes with the climate system (Riese et al., 2012; Butchart, 2014; Tweedy et al., 2017). The natural variability of the atmosphere largely influences the BDC (Hardiman et al., 2017). All three branches of the BDC are affected by changes in sea surface temperatures and the El Niño–Southern Oscillation (Yang et al., 2014; Diallo et al., 2019) as well as the phase of the quasi-biennial oscillation (QBO, Diallo et al., 2018) and the Arctic oscillation (Salby and Callaghan, 2005).

Modeling studies have predicted an acceleration of the BDC over the last few decades and the 21st century due to the increase in well-mixed greenhouse gases (Butchart et al., 2010; Hardiman et al., 2014; Palmeiro et al., 2014) and ozone-depleting substances (Polvani et al., 2018), but these results cannot be evaluated easily because the BDC cannot be observed directly (Butchart, 2014). Observational studies over short periods (typically 2003–2012) show significant evidence of a changing BDC in the boreal lower stratosphere (Schoeberl et al., 2008; Stiller et al., 2012; Hegglin et al., 2014; Mahieu et al., 2014; Haenel et al., 2015), but balloon-borne observations of SF₆ and CO₂ in the northern midlatitudes show a nonsignificant trend of the deep branch of the BDC over the past few decades (Engel et al., 2009, 2017). This difficulty in deriving observational trends in the BDC can be partly attributed to the spatial and temporal sparseness of the observations as well as its large dynamical variability and the uncertainty of trends derived from nonlinearly increasing tracers (Garcia et al., 2011; Hardiman et al., 2017; Fritsch et al., 2020). Before investigating multi-decadal changes of the BDC, it is important to perform an accurate evaluation of its climatological state and interannual variability, which is the aim of this paper.

In this study, we use N₂O as a tracer to study the BDC. N₂O is continuously emitted in the troposphere (with larger abundances in the Northern Hemisphere, NH) and transported into the stratosphere where it is destroyed by photodissociation and, to a lesser extent, by reaction with O(¹D). The estimated lifetime of N₂O is approximately 120 years, which makes it an excellent long-lived tracer for transport studies in the middle atmosphere (Brasseur and Solomon, 2006; Seinfeld and Pandis, 2016).

We use the transformed Eulerian mean (TEM, Andrews et al., 1987) analysis to separate the local rates of change of N₂O due to transport and chemistry (Randel et al., 1994). The transport term can be further separated into the contribution of isentropic mixing and residual advection, as done previously for O₃ and CO (Abalos et al., 2013). The isentropic mixing and the residual advection can be additionally separated in their horizontal and vertical contributions. In the tropical lower stratosphere, the distinction between vertical and horizontal transport is important, as they impact the seasonality of N₂O differently in the northern and southern tropics (Tweedy et al., 2017). We choose to focus our study on the horizontal mixing and vertical advection, because their magnitudes are larger than the vertical mixing and the meridional residual advection in most of the stratosphere.

Chemistry climate models (CCMs) include the full representation of dynamical, radiative and chemical processes in the atmosphere as well as their interactions. In particular, they combine the feedbacks of the chemical tracers on the heat budget and dynamics, which ultimately affects tracer transport. We use version 4 of the Whole Atmosphere Community Climate Model (WACCM; Garcia et al., 2017) to simulate the N₂O TEM budget in the stratosphere for the 2005–2014 period. WACCM has been widely used for studies of tracers' transport in the stratosphere and upper troposphere based on the TEM analysis (e.g., Abalos et al., 2017). WACCM simulations of the climatological N₂O over the 2005–2014 period have also been evaluated favorably with satellite observations in the stratosphere (Froidevaux et al., 2019).

In order to assess their representation of the atmospheric processes, CCMs are often compared to reanalyses (e.g., Gerber et al., 2010). Reanalysis products merge dynamical atmospheric observations (e.g., surface pressure, wind and temperature) with a global forecast model using an assimilation scheme to offer the best reproduction of the past climate. They provide a multivariate, consistent record of the global atmospheric state. Reanalyses are made using different assimilation methods and forecast models (Cameron et al., 2019), and they are often compared among each other and with CCMs (Rao et al., 2015). The SPARC (Stratosphere-troposphere Processes And their Role in Climate) Reanalysis Intercomparison Project (S-RIP) coordinates the intercomparison of all major global atmospheric reanalyses and provides reports to document these results (Fujiwara et al., 2017; Long et al., 2017).

Meteorological fields from reanalyses are often used to drive chemistry transport models (CTMs) in order to study the BDC through a common diagnostic, namely the age of air (AoA; Waugh and Hall, 2002), and simulate realistic distributions of chemical tracers (Monge-Sanz et al., 2012; Ménard et al., 2020). Thanks to their simplicity, CTMs are useful to compare different reanalyses within the same transport framework, thereby contributing to the study of the BDC in S-RIP (chapter 5 of the S-RIP report; see Fig. 1 in Fujiwara et al., 2017). CTMs may use either sigma-pressure levels with a kinematic transport scheme and vertical velocities simply derived from mass conservation or isentropic levels with a diabatic transport scheme (Chipperfield, 2006). Recent intercomparisons showed that the AoA depends to a large extent on the input reanalysis, both using the kinematic approach (Chabrillat et al., 2018) and the diabatic approach (Ploeger et al., 2019).

Here, we use the same CTM as for the kinematic AoA study, i.e., the Belgian Assimilation System of Chemical ObsErvation (BASCOE) CTM. Observations of another long-lived stratospheric tracer, HCFC-22, were recently interpreted with WACCM and BASCOE CTM simulations, showing the interest of this model intercomparison (Prignon et al., 2019). In order to contribute further to the S-RIP BDC activity, four different dynamical reanalyses are used here to drive the BASCOE CTM simulations, compute the N₂O TEM budget and compare its components with the results derived from WACCM, namely the European Centre for Medium-Range Weather Forecasts Interim reanalysis (ERA-Interim; Dee et al., 2011), the Japanese 55-year Reanalysis (JRA-55; Kobayashi et al., 2015), and the Modern-Era Retrospective analysis for Research and Applications version 1 (MERRA; Rienecker et al., 2011) and version 2 (MERRA-2; Gelaro et al., 2017).

While dynamical reanalyses do not assimilate observations of chemical compounds, chemical reanalyses achieve this step and can be used to evaluate CCMs or study differences between instruments using the reanalysis as a transfer tool (Errera et al., 2008; Lahoz and Errera, 2010; Davis et al., 2016). To our knowledge, chemical reanalyses driven by meteorological fields from modern dynamical reanalyses have not been used to study tracer transport in the stratosphere using the TEM framework. WACCM and the CTM experiments are compared with a chemical reanalysis of the Aura Microwave Limb Sounder (MLS) using the BASCOE data assimilation system (DAS) driven by the ERA-Interim reanalysis (BRAM2; Errera et al., 2019).

To summarize, in this study we analyze the representation of the BDC in WACCM through an analysis of the TEM budget of N₂O, and we evaluate the simulation of this budget through comparisons with the BASCOE CTM (driven by four dynamical reanalyses) and the BRAM2 chemical reanalysis. In Sect. 2, we describe the datasets used in the study and the TEM analysis of N₂O. In Sect. 3, we analyze the seasonal mean patterns of the TEM N₂O budget in each dataset and their differences. Sections 4 and 5 investigate the mean annual cycle and the variability of the N₂O TEM budget terms, respectively, with a focus on the differences between the datasets. Section 6 concludes the study with a summary of our findings and directions for possible future research.

Marsh et al. (2013)Garcia et al. (2017)Chabrillat et al. (2018)(Dee et al., 2011)Chabrillat et al. (2018)(Kobayashi et al., 2015)Chabrillat et al. (2018)(Rienecker et al., 2011)Chabrillat et al. (2018)(Gelaro et al., 2017)Errera et al. (2019)(Dee et al., 2011)(Livesey et al., 2015)

Table 1 Overview of the datasets used in this study.

2 Data and method

This work uses seven datasets that were generated by WACCM, the BASCOE CTM and the BASCOE DAS. Table 1 provides an overview of these datasets and their main differences, and the next three subsections provide details about the models and systems that generated them.

2.1 WACCM

WACCM (Garcia et al., 2017) is the atmospheric component of the Community Earth System Model version 1.2.2 (Hurrell et al., 2013), which has been developed by the U.S. National Center for Atmospheric Research. It is the extended (whole atmosphere) version of the Community Atmosphere Model version 4 (CAM4; Neale et al., 2013).

We ran one realization of the public version of WACCM (hereafter WACCM4; Marsh et al., 2013) with a similar setup (e.g., lower boundary conditions) to the CTM experiments; the source code of WACCM4 is available for download at https://svn-ccsm-models.cgd.ucar.edu/cesm1/release_tags/cesm1_2_2/, last access: 28 October 2020. In this study we also use three realizations of the REF-C1 simulation used in the SPARC Chemistry-Climate Model Initiative (CCMI; Morgenstern et al., 2017). The CCMI experiments, hereafter WACCM-CCMI, differ from WACCM4 with respect to the modified gravity wave parameterization and the updated heterogenous chemistry (Garcia et al., 2017). The inclusion of WACCM4 allows us to carry out a sensitivity test regarding the impact of the modified gravity wave parameterization on the simulation of the N₂O transport (see Sect. 4 for detailed analysis). We use 3-D daily mean output over the 2005–2014 period to allow a fair comparison with the BRAM2 dataset (see Sect. 2.3 for detailed analysis). WACCM has a longitude–latitude grid of ${\mathrm{2.5}}^{\circ }\times {\mathrm{1.9}}^{\circ }$ and 66 vertical levels ranging from the surface to an altitude of about 140 km. The vertical coordinate is hybrid-pressure, i.e., terrain-following below 100 hPa and purely isobaric above this level. The vertical resolution depends on the height: it is approximately 3.5 km above 65 km, 1.75 km around the stratopause (50 km), 1.1–1.4 km in the lower stratosphere (below 30 km) and 1.1 km in the troposphere. The time step for the physics in the model is 30 min.

The physics of WACCM is the same as CAM4, and the dynamical core is a finite volume with a horizontal discretization based on a conservative flux-form semi-Lagrangian (FFSL) scheme (Lin, 2004). The gravity wave parameterization accounts for momentum and heat deposition separating orographic and non-orographic sources. The orographic waves are modified according to Garcia et al. (2017), whereas non-orographic waves are parameterized depending on the convection and the frontogenesis occurrence in the model (Richter et al., 2010).

In this study, the WACCM versions considered are not able to internally generate the quasi-biennial oscillation (QBO; see, e.g., Baldwin et al., 2001). Thus, the QBO is nudged by a relaxation of stratospheric winds to observations in the tropics (Matthes et al., 2010). The solar forcing uses the Lean et al. (2005) approach.

WACCM includes a detailed coupled chemistry module for the middle atmosphere based on the Model for Ozone and Related Chemical Tracers, version 3 (MOZART-3; Kinnison et al., 2007; Marsh et al., 2013). The species included within this mechanism are contained within the O_x, NO_x, HO_x, ClO_x and BrO_x chemical families, along with CH₄ and its degradation products. In addition, 20 primary non-methane hydrocarbons and related oxygenated organic compounds are represented along with their surface emission. There are a total of 183 species and 472 chemical reactions; this includes 17 heterogeneous reactions on multiple aerosol types, i.e., sulfate, nitric acid trihydrate and ice. In WACCM-CCMI the heterogeneous chemistry is updated by Solomon et al. (2015).

2.2 BASCOE CTM

The BASCOE data assimilation system (Errera et al., 2019) is built on a chemistry transport model, which consists in a kinematic transport module with the FFSL advection scheme (Lin and Rood, 1996) and an explicit solver for stratospheric chemistry, comprising 65 species and 243 reactions (Prignon et al., 2019). Chabrillat et al. (2018) explain in detail the preprocessing procedure that allows the BASCOE CTM to be driven by arbitrary reanalysis datasets as well as the setup of model transport. As usual for kinematic transport modules, the FFSL scheme only needs the surface pressure and horizontal wind fields from reanalyses as input, because it is set on a coarser grid than the input reanalyses, and relies on mass continuity to derive vertical mass fluxes corresponding to its own grid. Similar to Chabrillat et al. (2018), the model is driven by four different reanalysis datasets on a common, low-resolution, latitude–longitude grid (${\mathrm{2.5}}^{\circ }\times {\mathrm{2}}^{\circ }$), but it retains their native vertical grids. In this way, we avoid any vertical regridding, and the intercomparison explicitly accounts for the different vertical resolutions.

The four input reanalyses, the European Centre for Medium-Range Weather Forecasts Interim reanalysis (ERA-Interim, hereafter ERAI;  Dee et al., 2011), the Japanese 55-year Reanalysis (JRA-55;  Kobayashi et al., 2015), and the Modern-Era Retrospective analysis for Research and Applications (MERRA;  Rienecker et al., 2011) and its version 2 (MERRA-2;  Gelaro et al., 2017), are part of the SPARC Reanalysis Intercomparison Project (S-RIP), which is a coordinated intercomparison of all major global atmospheric reanalyses, and they are described in Fujiwara et al. (2017). ERAI and JRA-55 have 60 levels up to 0.1 hPa while MERRA and MERRA-2 have 72 levels up to 0.01 hPa. The CTM time step is set to 30 min. As for the WACCM experiment, we used the daily mean outputs from the BASCOE CTM over the 2005–2014 period.

2.3 BASCOE reanalysis

BRAM2 is the BASCOE reanalysis of Aura MLS, version 2, which covers the period from August 2004 to August 2019 (Errera et al., 2019). For BRAM2, BASCOE is driven by dynamical fields from ERA-Interim, with a ${\mathrm{3.75}}^{\circ }\times {\mathrm{2.5}}^{\circ }$ (longitude × latitude) horizontal resolution. The vertical grid is represented by 37 hybrid-pressure levels which are a subset of the ERA-Interim 60 levels.

In BRAM2, N₂O profiles from the MLS version 4 standard product have been assimilated within the 0.46–68 hPa pressure ranges (Livesey et al., 2015). This dataset is retrieved from the MLS 190 GHz radiometer instead of the 640 GHz radiometer in an earlier MLS version. The 640 GHz radiometer, which provided a slightly better quality retrieval down to 100 hPa, ceased to be delivered after August 2013 because of instrumental degradation in the band used for that retrieval. To avoid any artificial discontinuity due to switching from one product to the other in August 2013, BRAM2 has assimilated the 190 GHz N₂O during the whole reanalysis period.

BRAM2 N₂O has been validated between 3 and 68 hPa against several instruments with a general agreement within 15 % depending on the instrument and the atmospheric region (the middle stratosphere or the polar vortex; see Errera et al., 2019). It is not recommended to use the BRAM2 N₂O reanalysis outside of this pressure range. BRAM2 N₂O is also affected by a small drift of around −4 % between 2005 and 2015 (see also Froidevaux et al., 2019).

2.4 TEM diagnostics

For atmospheric tracers, the TEM analysis (Andrews et al., 1987) allows for the separation of the local change in a tracer with the volume mixing ratio χ in terms of transport and chemistry (Eq. 1).

$$\begin{array}{}\text{(1)}& {\stackrel{\mathrm{‾}}{\mathit{\chi }}}_t=-{\stackrel{\mathrm{‾}}{v}}^{\ast }{\stackrel{\mathrm{‾}}{\mathit{\chi }}}_y-{\stackrel{\mathrm{‾}}{w}}^{\ast }{\stackrel{\mathrm{‾}}{\mathit{\chi }}}_z+e^{z/H}\mathrm{\nabla }\cdot M+\stackrel{\mathrm{‾}}{S},\end{array}$$

where χ is the volume mixing ratio of N₂O, and M is the eddy flux vector, which is defined as

$$\begin{array}{}\text{(2a)}& {\displaystyle }& {\displaystyle }{M}^{\left(y\right)}\equiv -e^{-z/H}\left(\stackrel{\mathrm{‾}}{v^{\prime }{\mathit{\chi }}^{\prime }}-\stackrel{\mathrm{‾}}{v^{\prime }{\mathit{\theta }}^{\prime }}{\stackrel{\mathrm{‾}}{\mathit{\chi }}}_z/{\stackrel{\mathrm{‾}}{\mathit{\theta }}}_z\right),\text{(2b)}& {\displaystyle }& {\displaystyle }{M}^{\left(z\right)}\equiv -e^{-z/H}\left(\stackrel{\mathrm{‾}}{w^{\prime }{\mathit{\chi }}^{\prime }}+\stackrel{\mathrm{‾}}{v^{\prime }{\mathit{\theta }}^{\prime }}{\stackrel{\mathrm{‾}}{\mathit{\chi }}}_y/{\stackrel{\mathrm{‾}}{\mathit{\theta }}}_z\right).\end{array}$$

${\stackrel{\mathrm{‾}}{v}}^{\ast }$ and ${\stackrel{\mathrm{‾}}{w}}^{\ast }$ are the meridional and vertical components of the residual mean meridional circulation, respectively, and they are defined as

$$\begin{array}{}\text{(3a)}& {\displaystyle }& {\displaystyle }{\stackrel{\mathrm{‾}}{v}}^{\ast }\equiv \stackrel{\mathrm{‾}}{v}-e^{z/H}{\left(e^{-z/H}\stackrel{\mathrm{‾}}{v^{\prime }{\mathit{\theta }}^{\prime }}/\stackrel{\mathrm{‾}}{{\mathit{\theta }}_z}\right)}_z,\text{(3b)}& {\displaystyle }& {\displaystyle }{\stackrel{\mathrm{‾}}{w}}^{\ast }\equiv \stackrel{\mathrm{‾}}{w}+{\left(a\cos \mathit{\varphi }\right)}^{-\mathrm{1}}{\left(\cos \mathit{\varphi }\stackrel{\mathrm{‾}}{v^{\prime }{\mathit{\theta }}^{\prime }}/\stackrel{\mathrm{‾}}{{\mathit{\theta }}_z}\right)}_{\mathit{\varphi }}.\end{array}$$

Here, $\stackrel{\mathrm{‾}}{v}$, $\stackrel{\mathrm{‾}}{w}$ and $\stackrel{\mathrm{‾}}{\mathit{\theta }}$ are the Eulerian zonal-mean meridional and vertical velocities and the potential temperature, respectively; ϕ is the latitude; and S is the net rate of change due to chemistry, i.e., $\stackrel{\mathrm{‾}}{S}=\stackrel{\mathrm{‾}}{P}-\stackrel{\mathrm{‾}}{L}$, where $\stackrel{\mathrm{‾}}{P}$ and $\stackrel{\mathrm{‾}}{L}$ are the respective zonal-mean chemical production and loss rates. Overbar quantities represent zonal-mean fields, primed quantities represent the departures from the zonal mean and subscripts denote derivatives. Meridional derivatives are evaluated in spherical coordinates and vertical derivatives with respect to log-pressure altitude $z\equiv -H{\log }_e(p/p_{\mathrm{s}})$, with p_s=10⁵ Pa and H=7 km.

Hence, transport is separated into advection due to the residual circulation (first two terms on the right-hand side, RHS, of Eq. 1) and irreversible quasi-horizontal isentropic eddy mixing, $e^{z/H}\mathrm{\nabla }\cdot \mathit{M}$.

In order to better understand the role of each term in the tracer balance, it is useful to separate the components of the vector M and rearrange the terms of Eq. (1) as follows:

$$\begin{array}{}\text{(4)}& {\stackrel{\mathrm{‾}}{\mathit{\chi }}}_t=A_y+M_y+A_z+M_z+\left(\stackrel{\mathrm{‾}}{P}-\stackrel{\mathrm{‾}}{L}\right)+\stackrel{\mathrm{‾}}{\mathit{ϵ}},\end{array}$$

where

$$\begin{array}{}\text{(5a)}& {\displaystyle }& {\displaystyle }{A}_y=-{\stackrel{\mathrm{‾}}{v}}^{\ast }{\stackrel{\mathrm{‾}}{\mathit{\chi }}}_y,\text{(5b)}& {\displaystyle }& {\displaystyle }{M}_y=e^{z/H}\cos {\mathit{\varphi }}^{-\mathrm{1}}{\left(M^{\left(y\right)}\cos \mathit{\varphi }\right)}_y,\text{(5c)}& {\displaystyle }& {\displaystyle }{A}_z=-{\stackrel{\mathrm{‾}}{w}}^{\ast }{\stackrel{\mathrm{‾}}{\mathit{\chi }}}_z,\text{(5d)}& {\displaystyle }& {\displaystyle }{M}_z=e^{z/H}{\left(M^{\left(z\right)}\right)}_z.\end{array}$$

Here, A_y represents the meridional residual advection, M_y represents the horizontal transport due to eddy mixing, A_z represents the vertical residual advection and M_z represents the vertical eddy mixing (all expressed in ppbv d⁻¹). It is important to note that the total mixing term (M_y+M_z) includes not only the effects of irreversible mixing but also some effects of the advective transport that are not resolved by the residual advection (Andrews et al., 1987; Holton, 2004).

Before any TEM calculation, all of the input fields are interpolated to constant pressure levels from the hybrid-sigma coefficients, which retain the same vertical resolution as the original vertical grid of each dataset (Table 1). Each derivative is computed using a centered differences method.

In addition to the physical TEM terms (Eq. 1), it is necessary to include an additional term on the RHS of Eq. (4): the residual term ϵ. It is the difference between the actual rate of change of χ (LHS of Eq. 4) and the sum of all of the transport and chemical terms of the TEM budget. This nonzero residual has several causes (Abalos et al., 2017). The TEM calculations for WACCM rely on the diagnostic variable w, which is not used to advect the tracers, because the model is based on a finite volume dynamical core (Lin, 2004). Furthermore, an implicit numerical diffusion is added to the transport scheme in WACCM in order to balance the small-scale noise without altering the large-scale features. This numerical diffusion is not included in the TEM budget and is larger in regions with large small-scale features, i.e., regions where gradients are stronger (Conley et al., 2012). All TEM calculations are done using daily mean data, even though WACCM and BASCOE both run with a much smaller time step of 30 min. The daily mean fields are interpolated from their native hybrid-sigma levels to constant pressure levels prior to the TEM analysis, leading to numerical errors in the lower stratosphere. The BASCOE datasets have a coarser horizontal resolution than their input reanalyses (especially BRAM2; see Table 1). This affects the accuracy of the vertical and horizontal derivatives, with possible implications for the residual. The possible causes of the residual in the five reanalysis datasets are discussed in more detail in Sect. 3. For WACCM-CCMI, the TEM budget is computed for each realization, allowing for the examination of both the ensemble mean (e.g., for seasonal means) or the model envelope (e.g., for line plots). In order to validate our N₂O TEM budget, we reproduced the findings reported in Tweedy et al. (2017, Fig. 7) with WACCM-CCMI in the tropical lower stratosphere, and we noticed similar results (not shown).

In order to interpret the TEM analysis of the N₂O budget, we also compute the Eliassen–Palm flux divergence (EPFD). The Eliassen–Palm flux is a 2-D vector defined as $\mathit{F}\equiv (F^{\left(\mathit{\varphi }\right)},F^{\left(z\right)})$ (Andrews et al., 1987), with its meridional and vertical components, respectively, given by

$$\begin{array}{ll}\text{(6a)}& {\displaystyle }& {\displaystyle }{F}^{\left(\mathit{\varphi }\right)}\equiv e^{-z/H}a\cos \mathit{\varphi }\left({\stackrel{\mathrm{‾}}{u}}_z\stackrel{\mathrm{‾}}{v^{\prime }{\mathit{\theta }}^{\prime }}/\stackrel{\mathrm{‾}}{{\mathit{\theta }}_z}-\stackrel{\mathrm{‾}}{v^{\prime }{u}^{\prime }}\right){\displaystyle }& {\displaystyle }{F}^{\left(z\right)}\equiv e^{-z/H}a\cos \mathit{\varphi }\mathit{\{}\left[f-\left(acos\mathit{\varphi }{)}^{-\mathrm{1}}\right(\stackrel{\mathrm{‾}}{u}\cos \mathit{\varphi }{)}_{\mathit{\varphi }}\right]\\ \text{(6b)}& {\displaystyle }& {\displaystyle }\stackrel{\mathrm{‾}}{v^{\prime }{\mathit{\theta }}^{\prime }}/\stackrel{\mathrm{‾}}{{\mathit{\theta }}_z}-\stackrel{\mathrm{‾}}{w^{\prime }{u}^{\prime }}\mathit{\}}.\end{array}$$

The EPFD reflects the magnitude of the eddy processes and provides a direct measure of the dynamical forcing of the zonal-mean state by the resolved eddies (Edmon et al., 1980).

The four dynamical reanalyses used in this study provide overall consistent temperature and winds in the stratosphere, but they can lead to a different representation of large-scale transport (e.g., Chabrillat et al., 2018) due to the biases in the temperature and wind fields (Kawatani et al., 2016; Tao et al., 2019). Note that the TEM quantities are not directly constrained by observations, especially the upwelling velocity ${\stackrel{\mathrm{‾}}{w}}^{\ast }$, that can vary considerably in the dynamical reanalyses, as it is a small residual quantity (Abalos et al., 2015).

In the rest of the paper, we will assume that the BRAM2 product provides the best available approximation of the TEM budget for N₂O – at least where the residual is smaller than the vertical advection and horizontal mixing terms. This assumption relies on the combination in BRAM2 of dynamical constraints from ERA-Interim and chemical constraints from MLS (Errera et al., 2019).

Figure 1 Latitudinal profiles of the N₂O TEM budget terms at 15 hPa averaged in DJF (2005–2014) with WACCM-CCMI (a), JRA-55 (b), MERRA-2 (c), MERRA (d), ERAI (e) and BRAM2 (f). The color code is shown in the legend. Units are in parts per billion by volume per day (ppbv d⁻¹).

Figure 2 Same as for Fig. 1 but for JJA.

Figures 1 and 2 show the N₂O TEM budget terms at 15 hPa for all of the datasets for the boreal winter (December–January–February, DJF, mean) and summer (June–July–August, JJA, mean), respectively. The 15 hPa level (around 30 km altitude) was chosen because large differences can be found between WACCM-CCMI, BRAM2 and the CTM runs at this level and because the dynamical reanalyses are not constrained as well by meteorological observations at higher levels (Manney et al., 2003). Figures 1 and 2 aim to show how the dynamical and chemical terms of the budget balance each other to recover the tendency ${\stackrel{\mathrm{‾}}{\mathit{\chi }}}_t$ at different latitudes. A discussion on the differences between the datasets and their possible physical causes is given in the next sections.

The vertical advection term A_z shows how the upwelling contributes to increasing the N₂O abundances in the tropics and summertime midlatitudes and how polar downwelling contributes to decreasing the N₂O abundances in the winter hemisphere. The horizontal transport out of the tropics due to eddies, as represented by M_y, reduces the N₂O abundance in the tropical latitudes of the wintertime hemisphere and increases the N₂O mixing ratio at high latitudes in the winter hemisphere. The other terms of the TEM budget are weaker than A_z and M_y: the meridional advection term A_y tends to increase the N₂O abundance in the winter subtropics and extra-tropics; the vertical transport term due to eddy mixing, M_z, decreases the N₂O mixing ratio over the northern polar latitudes; and the chemistry term P−L shows that N₂O destruction by photodissociation and O(¹D) oxidation contributes to the budget in the tropics and also in the summertime hemisphere. All budget terms are weaker in the summer hemisphere than in the winter hemisphere. Over the southern polar winter latitudes, the reanalyses deliver negative M_y that are balanced by large positive residuals, which implies a less robust TEM balance (Fig. 2). This is not the case with WACCM, where M_y tends to increase the N₂O abundance in the polar vortex. Such differences between the datasets are highlighted and discussed in the next sections.

Figure 3 Climatological (2005–2014) latitude–pressure cross sections of three N₂O TEM budget terms averaged in DJF (ppbv d⁻¹) for the horizontal mixing term (a, d, g, j), the vertical residual advection term (b, e, h, k) and the residual term (c, f, i, l). The datasets are, from top to bottom, WACCM-CCMI, JRA-55, MERRA-2 and BRAM2. The residual term for WACCM-CCMI is from a single realization of the model. The thin black lines show the zonal-mean zonal wind (from 0 to 40 m s⁻¹ every 10 m s⁻¹), the black thick line represents the dynamical tropopause for the considered season and the green thin lines show the climatological mixing ratio of N₂O (from 20 to 300 ppbv with 40 ppbv spacing).

Figure 4 Same as Fig. 3 but for JJA and with a different color scale. The thin black contours show the zonal-mean zonal wind (from 0 to 100 m s⁻¹ every 20 m s⁻¹).

3 Latitude–pressure cross sections

Figures 3 and 4 show the DJF and JJA means of three contributions to the N₂O TEM budget, respectively, namely the horizontal mixing M_y, the vertical advection A_z and the residual terms ϵ, for WACCM-CCMI, JRA-55, MERRA-2 and BRAM2. For those datasets, the remaining terms of the TEM budget (A_y, M_z and P−L) for DJF and JJA are shown in Supplement Figs. S1 and S2, respectively. The full N₂O TEM budgets obtained with MERRA and ERAI for DJF and JJA are shown in Figs. S3 and S4, respectively. In the case of WACCM-CCMI, the seasonal means were computed separately for each realization, and we verified that the ensemble means show the same features as the individual realizations. Large differences arise in the dynamical terms of the TEM budget between summer and winter for both hemispheres in the extra-tropics. The strong seasonality of the deep branch of the BDC and of the transport barriers are the causes of these differences, as for the seasonal variations of the age of air spectrum (Li et al., 2012).

We also reproduced the results of Randel et al. (1994, Fig. 8) for the WACCM-CCMI multi-model mean and the reanalysis mean in DJF (Figs. S5 and S5, respectively). The WACCM-CCMI and the reanalysis means agree with the Community Climate Model version 2 of the early 1990s with regard to the general pattern of the TEM terms, although both deliver stronger contributions, especially the reanalyses mean.

We first compare the contribution of A_z across the datasets in Figs. 3 and 4. The tropical upwelling increases the abundance of N₂O mostly in the mid–high stratosphere (between 1 and 15 hPa) with the maximum contribution in the summer tropics, whereas the downwelling decreases it mostly in the wintertime extra-tropics in the middle and low stratosphere (between 5 and 100 hPa). This reflects the path followed by the deep branch of the BDC (Birner and Bönisch, 2011). During boreal winter, these features are very similar across all datasets (Fig. 3), but noticeable differences appear during the austral winter (Fig. 4): the tropical upwelling has a larger secondary maximum in the southern tropics with JRA-55 and MERRA-2 than with the other datasets, and the extra-tropical downwelling extends to the South Pole in WACCM-CCMI and JRA-55 whereas it is mostly confined to the midlatitudinal surf zone in the other reanalyses. In the lower stratosphere, A_z shows the contribution of the residual advection by the shallow branch of the BDC to the N₂O abundances in the winter and summer hemispheres. The two-cell structure, consisting of upwelling of N₂O in the subtropics and downwelling in the extra-tropics, consistently agrees across all datasets. The meridional residual advection term A_y contributes to the poleward transport of air masses in the middle stratosphere, mostly during the winter, and its contribution to the N₂O TEM budget is weaker than A_z. A_y agrees well among the datasets in boreal winter (Figs. S1, S3), while during austral winter WACCM-CCMI overestimates it by around 30^∘ S compared with the reanalyses (Figs. S2, S4).

We move now to the mixing contributions to the N₂O budget. The horizontal mixing is the predominant contribution to the poleward tracer transport in the middle and lower stratosphere (Abalos et al., 2013), as it flattens the tracer gradients generated by the residual advection. In the N₂O TEM budget during boreal winter, M_y mostly balances the extra-tropical downwelling and part of the tropical upwelling (Figs. S3, S4). The surf zone is characterized by strong horizontal mixing, depicted here as large positive M_y contributions, and delimited by transport barriers which appear as intense gradients of M_y in the winter hemispheres (panels b, e, h and k in Figs. 3 and 4). In the wintertime NH, the patterns of M_y are similar in all datasets (Fig. 3), but the effect of horizontal eddy mixing on N₂O is stronger in the reanalyses than in WACCM-CCMI. In Sect. 4, we quantitatively analyze the differences in the mid-stratospheric M_y between datasets. The residual terms in the reanalyses (Fig. 3c, f, i, l) are largest in the middle stratosphere at the latitudes of the transport barriers, and their signs are opposite to M_y.

In the austral winter, over the Antarctic polar cap and below 30 hPa, M_y agrees remarkably well in all datasets (Fig. 4). Closer to the vortex edge and above 30 hPa, the wintertime decrease in N₂O is mainly due to downwelling in WACCM-CCMI, whereas the reanalyses, especially BRAM2, show that the horizontal mixing plays a major role (Fig. 4). The impact of horizontal mixing on N₂O inside the wintertime polar vortex is not negligible (e.g., de la Cámara et al., 2013; Abalos et al., 2016a), as Rossby wave breaking occurs there as well as in the surf zone. In contrast to the reanalyses, in WACCM-CCMI the M_y contribution is close to zero in the Antarctic vortex and maximum along the vortex edge (Fig. 4). This disagreement can be related to differences in the zonal wind: it is overestimated in WACCM above 30 km at subpolar latitudes compared with MERRA (Garcia et al., 2017), and the polar jet is not tilted equatorward as in the reanalyses (see black thin lines in Figs. 4 and 3 of Roscoe et al., 2012). However, the differences in M_y and A_z above the Antarctic in winter should be put into perspective with the relatively large residual terms that point to incomplete TEM budgets in the reanalyses (Figs. 4 and S4 right columns). Near the Antarctic polar vortex, the assumptions of the TEM analysis (such as small amplitude waves) are less valid, leading to larger errors in the evaluation of the mean transport and eddy fluxes (Miyazaki and Iwasaki, 2005).

As the relative importance of the residual is considerable above the Antarctic in the reanalyses (Fig. 4), it is necessary to better understand its physical meaning. Dietmüller et al. (2017) applied the TEM continuity equation to the age of air in CCM simulations. Computing the “resolved aging by mixing” (i.e., the AoA counterpart of M_y+M_z) as the time integral of the local mixing tendency along the residual circulation trajectories and the “total aging by mixing” as the difference between the mean AoA (mAoA) and the residual circulation transit time, they defined the “aging by mixing on unresolved scales” (i.e., by diffusion) as the difference between the latter and the former. This “aging by diffusion”, which can be related by construction to our residual term, arises around 60^∘ S from the gradients due to the polar vortex edge. Even though we use a real tracer (N₂O), we find a qualitative agreement with this analysis based on AoA: our residual term is larger in regions characterized by strong gradients such as the Antarctic vortex edge and larger with dynamics constrained to a reanalysis than with a free-running CCM (see EMAC results in Fig. 1d by Dietmüller et al., 2017). Thus, we interpret the residual as the sum of mixing at unresolved scales and numerical errors (Abalos et al., 2017).

In the summertime lower stratosphere, we note a stronger contribution of M_y to the N₂O abundances above the subtropical jets in both hemispheres and for all datasets compared with higher levels in summer (panels b, e, h and k in Figs. 3 and 4). This behavior is consistent with calculations of the effective diffusivity and age spectra (Haynes and Shuckburgh, 2000; Ploeger and Birner, 2016). It is due to transient Rossby waves that cannot travel further up into the stratosphere owing to the presence of critical lines, i.e., where the phase velocity of the wave matches the background wind velocity, generally leading to wave breaking (Abalos et al., 2016b). In particular, above the northern tropics during the boreal summer (Figs. 4, S2, S4), horizontal mixing is primarily associated with the Asian monsoon anticyclone, and it causes a decrease in N₂O (Konopka et al., 2010; Tweedy et al., 2017). In the lower stratosphere, the contributions from M_y combine with that from A_z in the total impact of the shallow branch of the BDC on N₂O all year round (Diallo et al., 2012).

The vertical mixing contribution M_z is very small during boreal winter, except in the middle and lower stratosphere poleward of 60^∘ N, where it tends to balance the M_y contribution (Figs. S1, S3). In austral winter, there is a strong disagreement between WACCM-CCMI and the reanalyses around 60^∘ S between 5 and 15 hPa (Fig. S2). WACCM-CCMI simulates a strong M_z contribution at the polar jet core that decreases the N₂O abundances and tends to balance M_y, whereas M_z is weaker and increases N₂O in the higher stratosphere in the reanalyses.

In the next section, we focus on a single level in the middle stratosphere to quantitatively study the disagreement between WACCM-CCMI and the reanalyses.

4 Climatological seasonal cycles

After investigating the seasonal means of A_z and M_y, it is interesting to examine their climatological mean annual cycles in order to study the month-to-month variations over the year and their dependance on the latitude in the middle stratosphere. The cycles are shown for three latitude bands in each hemisphere corresponding to the tropics (0–20^∘), the surf zones (40–60^∘) and the polar regions (60–80^∘). For WACCM-CCMI, we examine the envelope (the yellow shaded areas in Figs. 5–12) of the three model realizations in order to evaluate the role of the internal variability and its relative importance for each month and latitude band. In the following, we will consider BRAM2 as the reference when comparing N₂O mixing ratios between datasets, because its dynamics and chemistry are both constrained by observational datasets.

Figure 5 Monthly mean annual cycle of (a, b) the EPFD ($\mathrm{m}{\mathrm{s}}^{-\mathrm{1}}{\mathrm{d}}^{-\mathrm{1}}$) averaged between 3 and 50 hPa and (c, d) ${\stackrel{\mathrm{‾}}{w}}^{\ast }$ (mm s⁻¹) at 15 hPa for (a, c) the Antarctic region (60–80^∘ S) and (b, d) the Arctic region (60–80^∘ N). The color code is shown in the legend. The yellow envelope shows the three realizations of the WACCM-CCMI simulation.

4.1 Polar regions

The EPFD is often used to quantify the forcing of the wave drag due to resolved (planetary) waves (e.g., Gerber, 2012; Konopka et al., 2015). We first show the monthly mean climatological annual cycles of EPFD averaged between 3 and 50 hPa and the residual vertical velocity ${\stackrel{\mathrm{‾}}{w}}^{\ast }$ at 15 hPa for the polar regions (60–80^∘ S and N, Fig. 5). We arbitrarily average the EPFD between 3 and 50 hPa in order to identify the wave forcing for the deep branch of the BDC (Plumb, 2002; Konopka et al., 2015). However, the qualitative results do not depend on the choice of the lower boundary level. We also show one realization of the earlier version WACCM4 that suffered from a larger cold bias above the Antarctic (see Sect. 2.1). In WACCM-CCMI, the parameterization of gravity waves was adjusted in order to reduce this issue while not significantly changing the dynamics in the NH, which results in an enhanced polar downwelling above the southern polar region (Garcia et al., 2017). Above the Antarctic, the forcing from resolved waves peaks in October in the reanalyses as a result of the vortex breakup that allows enhanced wave activity compared with austral winter (Randel and Newman, 1998). The WACCM simulations miss this strong springtime peak, and they are in good agreement with the reanalyses during the rest of the year (Fig. 5a). The residual vertical velocity w^∗ above the Antarctic is shown in Fig. 5c. This comparison between the WACCM versions was already shown in Garcia et al. (2017, Fig. 10), but we repeat it here adding the dynamical reanalyses. In November–December the weaker downwelling in WACCM-CCMI agrees well with the reanalyses. Throughout the rest of the year WACCM-CCMI simulates a stronger downwelling than all reanalyses (also at lower levels, not shown). This difference raises the question of whether the residual vertical velocity is correctly represented in WACCM-CCMI or in the dynamical reanalyses. Above the Arctic, the WACCM simulations underestimate the EPFD contribution during boreal winter compared with the reanalyses (Fig. 5b), and the downwelling velocities simulated by WACCM are weaker than the reanalyses during that period, with no significant differences between the WACCM versions (Fig. 5d). The differences between WACCM and the reanalyses in EPFD and w^∗ in the polar regions will help with the interpretation of the differences in A_z and M_y.

Figure 6 Annual cycles over the 2005–2014 period at 15 hPa for (a, b) the N₂O volume mixing ratio (ppbv), (c, d) M_y (ppbv d⁻¹) and (e, f) A_z (ppbv d⁻¹) in the (a, c, e) Antarctic region (60–80^∘ S) and the (b, d, f) Arctic region (60–80^∘ N). Note that the vertical scale differs for M_y and A_z. The olive envelope shows the three realizations of the WACCM-CCMI simulation. The color codes for the four CTM simulations follow the conventions of S-RIP (Fujiwara et al., 2017). BRAM2 is depicted using a black line and symbols, as usually done for observations, because it is constrained by both dynamical and chemical observations. As the N₂O mixing ratio in BRAM2 has been evaluated with a 15 % uncertainty (1σ standard deviation) at 15 hPa (Errera et al., 2019), this is highlighted by a dark gray region in panels (a) and (b). The light gray shading around the BRAM2 cycles represents the uncertainty arising from the residual term in the TEM budget, i.e., it is entirely interpreted first as an uncertainty on A_z and then as an uncertainty on M_y in order to remain cautious.

Figure 6 shows the monthly mean climatological annual cycle of the N₂O mixing ratio, A_z and M_y for the polar regions (60–80^∘ S and N) at 15 hPa for all of the datasets. First, we investigate the N₂O mixing ratio in the Antarctic region (Fig. 6a). During winter, the N₂O abundances are smaller than the rest of the year due to the suppressed transport from the lower latitudes caused by the onset of the polar barrier. After the vortex breakup, the N₂O increase during spring and early summer is smaller in all of the simulations than in BRAM2. In WACCM-CCMI, the modification of the parameterization of gravity waves also results in a shift towards earlier vortex breakup dates in the austral spring compared with WACCM4 (Garcia et al., 2017). The earlier vortex breakup in WACCM-CCMI allows for the transport of N₂O-rich air from lower latitudes for a longer period compared with WACCM4, resulting in larger and more realistic simulations of the N₂O mixing ratios during austral spring and early summer (Fig. 6a).

In the Antarctic region, the downwelling decreases N₂O during most of the year (A_z term in Fig. 6c). Here, JRA-55 and WACCM-CCMI are outliers: both present stronger A_z contributions in fall and winter, especially WACCM-CCMI, which reaches values that are 3 times stronger than BRAM2 in early winter, as a result of the stronger downwelling velocity simulated by WACCM-CCMI in that region. While this strong disagreement is questioned by the large residuals, we note that all of the reanalyses confirm it except JRA-55. During fall and summer, A_z is stronger in WACCM-CCMI than in WACCM4 as a consequence of the stronger downwelling in WACCM-CCMI that results from the modification of the gravity wave parameterization.

We now turn to the contribution from M_y. In the Antarctic region, M_y is very different among the datasets during winter: in BRAM2 it contributes to the N₂O decrease during fall and winter, with the strongest contribution in July; however, this contribution is 2 times weaker with the CTM simulations, and the horizontal mixing has almost no effect on N₂O in WACCM-CCMI (Fig. 6e). As already mentioned, the TEM analysis suffers from large residuals in the wintertime Antarctic region. Nevertheless, we note that the disagreement between WACCM-CCMI and BRAM2 is significant, because the envelope of WACCM-CCMI realizations falls completely outside of the possible BRAM2 values in fall and winter when accounting for the residual. During the austral spring, the vortex breakup leads to an increased wave activity reaching the Antarctic, and M_y is in better agreement among all datasets compared with austral winter. Note that WACCM-CCMI exhibits large internal variability in this season (Fig. 6e).

It is interesting to highlight the differences between the wintertime Arctic and Antarctic regions, because the hemispheric differences in wave activity (Scaife and James, 2000; Kidston et al., 2015) play a crucial role in the N₂O abundances and TEM budget. Above the Arctic, the N₂O abundances simulated by WACCM agree with the BRAM2 reanalysis, except in December and January, and the CTM experiments driven by MERRA and MERRA-2 deliver smaller N₂O mixing ratios compared with BRAM2 (Fig. 6b). A_z is also in good agreement between the datasets above the Arctic, with the exception of ERAI and JRA-55 which provide stronger contributions (Fig. 6d). The Arctic is characterized by a very variable polar vortex with a shorter life span than the Antarctic vortex (Randel and Newman, 1998; Waugh and Randel, 1999), resulting in an enhanced contribution of the horizontal mixing to the N₂O budget during winter compared with the Antarctic (Fig. 6f). In contrast to the dynamical reanalyses and BRAM2, WACCM shows a 2-fold underestimation of the N₂O changes in the Arctic due to horizontal mixing during winter. Note that the Arctic extended winter presents the largest internal variability compared with the other regions, as shown by the spread in the WACCM realizations. The weaker contribution from M_y in WACCM is meaningful because the relative importance of the residual term is small in the NH. The horizontal mixing is predominately influenced by the forcing from the breaking of resolved (planetary) waves (Plumb, 2002; Dietmüller et al., 2018). In the Arctic region, WACCM underestimates the forcing from resolved waves compared with the dynamical reanalyses in the middle stratosphere (see Fig. 5d). This discrepancy in the resolved wave driving could contribute to the large differences in the wintertime M_y between the CCM simulations and the CTM experiments above the Arctic. On the other hand, the role of different waves in driving mixing processes is an ongoing research topic, and additional data and sensitivity tests are needed in order to establish a clear separation of the waves' contributions (e.g., gravity wave parameterization, spatial resolution, etc.; Dietmüller et al., 2018). It should also be emphasized that WACCM is among the CCMI models with the lowest contribution of aging by mixing to the age of air (Fig. 2 in Dietmüller et al., 2018).

Figure 7 As for Fig. 5 but for the middle latitudes: (a, c) southern midlatitudes (40–60^∘ S); (b, d) northern midlatitudes (40–60^∘ N).

4.2 Middle latitudes

Figure 7 shows the monthly mean climatological annual cycle of w^∗ at 15 hPa and the EPFD averaged between 3 and 50 hPa over the surf zones (40–60^∘ S and N), and Fig. 8 shows the monthly mean climatological annual cycle of the N₂O mixing ratio, A_z and M_y at 15 hPa averaged over the same latitudes. The subtropical barriers are not shown because M_y and A_z change sign in these regions and averaging across the regions would hinder the interpretation of their means.

Figure 8 As for Fig. 6 but for the middle latitudes: (a, c, e) southern midlatitudes (40–60^∘ S); (b, d, f) northern midlatitudes (40–60^∘ N).

In the southern midlatitudes, the EPFD peaks in austral spring in the reanalyses, due to the enhanced wave activity in the Southern Hemisphere (SH) during austral spring compared with winter (Konopka et al., 2015), whereas the WACCM simulations deliver an earlier and weaker peak during austral winter (Fig. 7a). The downwelling velocity w^∗ shows a similar pattern to the EPFD (Fig. 7c), as it is also driven by the breaking of resolved waves (Abalos et al., 2015). In the northern midlatitudes, the EPFD peaks in winter in all of the datasets, reflecting the stronger wave forcing in the surf zone in this season, and WACCM simulates lower EPFD values compared with the reanalyses (Fig. 7b) that leads to a weaker downwelling velocity in the WACCM simulations (Fig. 7d). As for the polar regions, the differences in EPFD and w^∗ between the WACCM simulations and the reanalyses will help with the interpretation of the differences in A_z and M_y.

With regard to the N₂O mixing ratio in both hemispheres, the CTMs driven by JRA-55 and ERAI are in good agreement with BRAM2, whereas MERRA-2 and MERRA underestimate it (Fig. 8a, b). The WACCM-CCMI simulations agree well with the BRAM2 chemical reanalysis, confirming the results obtained through the direct comparison with MLS observations (Froidevaux et al., 2019).

We now investigate the contribution from A_z and M_y. In the southern midlatitudes, A_z is negative in all seasons except during summer and there is again a good agreement among the datasets except for WACCM-CCMI and JRA-55 (Fig. 8c). These two datasets appear to have a purely annual cycle in this region, whereas the other four show a semiannual component. The peak in the A_z contribution in the reanalyses in September results from the increased forcing from the resolved waves (see Fig. 7a) and from the stronger contribution from gravity waves to the mass flux during spring (Sato and Hirano, 2019, Fig. 11). In the same region, M_y increases throughout the winter, reflecting the mixing associated with the surf zone, and also peaks in early spring in the reanalyses (Fig. 8e). During summer and early fall, M_y does not contribute significantly to the TEM budget, and M_y reaches negative values in November which are comparable to the residual term. Both A_z and M_y peak in mid-winter in the WACCM-CCMI simulations, whereas these maxima are reached 3 months later in the reanalyses. This difference is related to the earlier minimum in the downwelling velocity ${\stackrel{\mathrm{‾}}{w}}^{\ast }$ simulated by WACCM-CCMI (see Fig. 7c), which directly affects A_z (Fig. 8c) and, due to compensation, M_y (Fig. 8e). Among the reanalyses, the compensating contributions of A_z and M_y are stronger for JRA-55 than for the other reanalyses (up to 2 times larger in September, see Fig. 8c, e). This reflects the strong BDC in JRA-55 that resulted in the youngest mean AoA in the whole stratosphere (Chabrillat et al., 2018).

In the northern middle latitudes, A_z shows the effect of the wintertime downwelling to lower levels on N₂O, with the WACCM experiments simulating a slightly weaker contribution than the reanalyses (Fig. 8d). Such disagreement mostly originates from the weaker downwelling velocity in the CCM compared with the reanalyses shown in Fig. 7b. In the northern midlatitudes, the strong M_y contribution tends to increase the N₂O abundances in the surf zone during winter (Fig. 8f). The reanalyses show a large spread, with values reaching ∼1.5 in BRAM2 and ∼0.9 ppbv d⁻¹ in the MERRA runs, and WACCM-CCMI presents a large underestimation with respect to the reanalyses. While the spread across the reanalyses cannot be explained by the forcing from the resolved waves, the weaker M_y contribution simulated by WACCM could be partly attributed to the weaker EPFD in the CCM compared with the reanalyses (see Fig. 7d).

Figure 9 As for Fig. 6 but for the tropics: (a, c, e) southern tropics (0–20^∘ S); (b, d, f) northern tropics (0–20^∘ N).

4.3 Tropics

Figure 9 shows the climatological annual cycle for the N₂O mixing ratio, A_z and M_y for the southern and northern tropics (0–20^∘ S and N) at 15 hPa across all of the datasets. The same latitude bands for the cycles of w^∗ and EPFD are shown in the Supplement (Fig. S7). In the tropical regions, the N₂O mixing ratio in WACCM-CCMI agrees well with the reanalysis of the Aura MLS, whereas the CTM results show large differences in the N₂O abundances depending on the input reanalysis (Fig. 9a, b). In regions where the AoA is less than 4.5 years and N₂O is greater than 150 ppb, i.e., in the tropical regions and lower stratospheric middle latitudes (Strahan et al., 2011), the N₂O mixing ratio is inversely proportional to the mAoA, because faster upwelling (younger air) implies more N₂O transported from lower levels, decreasing its residence time and resulting in a limited chemical destruction (Hall et al., 1999; Galytska et al., 2019). The dynamical reanalyses also produce large differences in mAoA at 15 hPa: MERRA delivers the oldest mAoA, and MERRA-2, ERAI and JRA-55 progressively show a younger mAoA (Fig. 4b in  Chabrillat et al., 2018). Hence, the large discrepancies in the N₂O mixing ratio can be explained by the large differences in the mAoA, while M_y and A_z contribute to rates of change of N₂O.

We continue by investigating the contribution from A_z. In both of the tropical regions, the upwelling term A_z is positive all year round, showing the effect of tropical upwelling, and agrees very well in the reanalyses (Fig. 9c, d), as a result of the good agreement in the tropical upwelling velocity at 15 hPa (Fig. S7 bottom row), and also as depicted by mAoA diagnostics (Fig. 4d in Chabrillat et al., 2018). Large interhemispheric differences arise in the seasonality of A_z between the tropical regions. The largest values of A_z in the southern tropics are in the boreal late-fall and winter (Fig. 9c), whereas no large seasonal variations can be detected in the annual cycle of the A_z in the northern tropics (Fig. 9d). This is the result of the more pronounced seasonality of the upwelling velocity in the southern tropics compared with the northern tropics (Fig. S7 bottom row).

We now turn to the contribution from M_y. In the southern tropics, M_y causes a decrease in the N₂O abundances from May to October (when N₂O is transported to the middle latitudes), and it has a near-zero contribution over the rest of the year, which is generally common in all of the datasets considered (Fig. 9e). The BRAM2 uncertainty is smaller than for the polar region and middle latitudes, indicating better performance of the TEM analysis outside of the high latitudes. In the northern tropics, M_y is negative from November to April and presents a marked seasonality in the reanalyses that is much weaker in WACCM-CCMI (Fig. 9f). With respect to interhemispheric differences, WACCM disagrees with the reanalyses: according to WACCM, M_y has a larger impact in the southern tropics than in the northern tropics, but M_y has a much larger impact in the northern tropics according to the reanalyses (Fig. 9e, f). These interhemispheric differences in the M_y contributions can be partly attributed to different forcings from the resolved waves between northern and southern tropics. The EPFD presents a stronger seasonality in the northern tropics than in the southern tropics in all of the datasets (Fig. S7 top row), which could partly explain the differences in the seasonality of M_y in the reanalyses, but it does not impact the M_y simulated by WACCM.

Figure 10 Monthly standard deviation over the 2005–2014 period at 15 hPa for (a, b) the N₂O volume mixing ratio (ppb), (c, d) the horizontal mixing term (ppbv d⁻¹) and (e, f)  the vertical residual advection term (ppbv d⁻¹) in (a, c, e) the Antarctic region (60–80^∘ S) and (b, d, f) the Arctic region (60–80^∘ N). The color code is shown in the legend. The yellow envelope shows the three realizations of the WACCM-CCMI simulation.

5 Interannual variability of the seasonal cycles

To analyze the interannual variability of the annual cycle, we compute for each month the 1σ standard deviations of the N₂O mixing ratio, M_y and A_z across the 10 simulated years. Figure 10 shows the annual cycles of these standard deviations for each dataset in the polar regions (60–80^∘ S and N) at 15 hPa for both hemispheres. First, we consider the variabilities of the N₂O mixing ratio. In the Antarctic, during austral winter and early spring the year-to-year change in the N₂O abundances is very small (Fig. 10a) because the duration, extension and strength of the polar vortex are very stable in a climatological sense, isolating air masses in the vortex from the highly variable midlatitudes (Waugh and Randel, 1999). The variability of the N₂O mixing ratio increases in October, i.e., during the breaking vortex period that is highly variable in time (Strahan et al., 2015). Furthermore, the midlatitude air masses, which have more variable composition, become free to reach the higher latitudes during this period. In the Arctic region, the interannual variability of the N₂O mixing ratio is also largest during springtime but this is spread over a longer period, i.e., from February to June (Fig. 10b), reflecting the large interannual variability in the duration, extension and zonal asymmetry of the Arctic polar vortex (Waugh and Randel, 1999). In both polar regions, WACCM-CCMI agrees well with BRAM2, whereas the CTM experiments simulate a smaller variability.

We now look at the interannual variability of A_z and M_y in the polar regions. Above the Antarctic, the interannual variabilities of A_z and M_y are maximum during spring (Fig. 10c, e), due to the large interannual variability in vortex breakup dates (Strahan et al., 2015). While the maximum variability of M_y is consistently reached in October in all of the reanalyses, WACCM-CCMI simulates an earlier maximum (September) that does not correspond to the maximum in its mean values. The lower wintertime variability of both A_z and M_y would increase if a longer period was considered that included the exceptional Antarctic vortices of 2002 (Newman and Nash, 2005) and 2019 (Yamazaki et al., 2019). Above the Arctic, M_y and A_z are most variable during winter, reflecting the frequent disruptions of the northern polar vortex by sudden stratospheric warmings (SSWs; Butler et al., 2017). A case study of the effect of a SSW on the N₂O TEM budget showed that A_z and M_y contribute more to this budget during the SSW event than in the seasonal mean (Randel et al., 1994). Thus, the large wintertime variabilities of A_z and M_y are explained by the occurrence of seven major SSWs detected in the reanalyses for the 2005–2014 period (Butler et al., 2017).

Figure 11 As for Fig. 10 but for the middle latitudes: (a, c, e) southern midlatitudes (40–60^∘ S); (b, d, f) northern midlatitudes (40–60^∘ N).

In Fig. 11, we show the interannual variabilities of the N₂O mixing ratio, M_y and A_z for each dataset in the surf zones (40–60^∘ S and N) at 15 hPa for both hemispheres. Regarding the N₂O mixing ratio, the interannual variability in the southern middle latitudes reaches the lowest values during austral winter. The datasets deliver very diverse values, with WACCM showing the largest variability and JRA-55 showing the least variability across the climatological year (Fig. 11a). In the northern midlatitudes, the interannual variability of the N₂O mixing ratio increases in late winter across all of the datasets as a response to the increased wintertime variability of the surf zone (Fig. 11b). The variability of WACCM-CCMI largely depends on the realization considered, except in October and November. Strong differences between ensemble members with respect to interannual variability indicate that the period considered is not long enough to explore the interannual variability in the northern midlatitudes and that the mean variability from this ensemble (with only three members) would not be representative of the internal variability of WACCM. The interannual variabilities of A_z and M_y in the southern midlatitudes are shown in Fig. 11c and e, respectively. As their mean value, A_z and M_y are most variable during austral spring and late summer in the reanalyses, whereas WACCM simulates an earlier peak during winter in the interannual variabilities of A_z and M_y compared with the reanalyses. In the northern midlatitudes, the interannual variabilities of A_z and M_y peak in winter, as do their mean values, and WACCM simulates smaller variabilities compared with the reanalyses (Fig. 11d, f).

Figure 12 As for Fig. 10 but for the tropics: (a, c, e) southern tropics (0–20^∘ S); (b, d, f) northern tropics (0–20^∘ N).

Figure 12 shows the annual cycles of the standard deviations of the N₂O mixing ratio, M_y and A_z for each dataset in the tropical regions (0–20^∘ S and N) at 15 hPa for both hemispheres. The interannual variability of the N₂O mixing ratio in both the southern and northern tropics depends considerably on the dataset (Fig.12a, b). WACCM-CCMI and the BASCOE reanalysis of the Aura MLS show very similar variabilities, especially in the southern tropics. As the QBO is the major source of variability in the tropical stratosphere (Baldwin et al., 2001), this confirms an earlier comparison that showed a good agreement between the WACCM model and MLS observations in the middle stratosphere in terms of the interannual variability of N₂O due to the QBO (Park et al., 2017). Among the CTM simulations, ERAI succeeds at delivering $\mathit{\sigma }\left(\stackrel{\mathrm{‾}}{X}\right)$ values that are as large as BRAM2 and WACCM-CCMI in the southern tropics, although this is not the case in the northern tropics.

The interannual variability of A_z in both hemispheres can be related to the impact of the QBO on the tropical upwelling (Flury et al., 2013). Among MERRA, ERAI and JRA-55, the fraction of variance in deseasonalized tropical upwelling ${\stackrel{\mathrm{‾}}{w}}^{\ast }$ that is associated with the QBO is the largest with ERAI (Abalos et al., 2015). Our findings support this conclusion, as the largest $\mathit{\sigma }\left(\stackrel{\mathrm{‾}}{A_z}\right)$ among the reanalyses is again found with ERAI (Fig. 12c, d). However, a detailed analysis of the impact of the QBO on the BDC as illustrated here goes beyond the scope of this study. The variability of M_y in the tropical regions is small compared with the extra-tropical regions (Fig. 12e, f), which is in agreement with calculations of effective diffusivity that show small variabilities within the tropical pipe (Abalos et al., 2016a). The reanalyses deliver a larger interannual variability in the northern tropics during boreal winter, whereas the variability of M_y presents a much weaker annual cycle in the southern tropics. WACCM-CCMI does not reproduce this hemispheric asymmetry, with a rather flat profile in both hemispheres and a clear underestimation in the northern tropics, as is also shown for its mean values. In the tropical regions, both the variabilities of M_y and A_z fail to explain the good agreement in the variability of N₂O between WACCM and BRAM2, as well as their disagreement with the dynamical reanalyses, because M_y and A_z directly contribute to the N₂O tendency rather than its mixing ratio.

6 Summary and conclusions

We have evaluated the climatological (2005–2014) N₂O transport processes in the stratosphere using the tracer continuity equation in the TEM formalism. In particular, we emphasized the horizontal mixing and the vertical advection terms (M_y and A_z, respectively). The upwelling term A_z reduces the N₂O concentrations in the tropics and increases them in the extra-tropics, whereas M_y tends to reduce the meridional gradients of N₂O and presents large hemispheric differences. As M_y or A_z contribute to the local rates of change of N₂O, this analysis is complementary to time-integrated diagnostics such as mAoA. The comparison investigates a variety of datasets, from a free-running chemistry climate model to a reanalysis where dynamics and chemistry are both constrained. The former comprises three realizations of the CCMI REF-C1 experiment by WACCM, and the latter is the chemical reanalysis of Aura MLS driven by ERA-Interim: BRAM2. The intercomparison also includes the BASCOE CTM driven by four dynamical reanalyses – ERAI, JRA-55, MERRA and MERRA-2 – in order to contribute to the S-RIP.

Considering the N₂O mixing ratio in the middle stratosphere, all of the datasets agree with respect to the annual cycle, with the large spread in the N₂O abundances of the CTM experiments (∼20 %) reflecting the diversity of mAoA obtained with the same model (Chabrillat et al., 2018). The upwelling term A_z also agrees among the datasets, especially in the NH where WACCM closely follows the reanalyses. The horizontal mixing term M_y in the NH is weaker in WACCM compared with the reanalyses. In the extra-tropics, this could be attributed to the weaker forcing from the planetary waves in WACCM compared with the reanalyses. The differences in M_y become striking in the wintertime Antarctic, where the polar vortex plays a major role. According to the reanalyses, horizontal mixing plays an important role in that region, but that is not found by WACCM. However, this large wintertime M_y in the reanalyses is challenged by an almost as large residual term. It should be noted that the residual term also includes effects from mixing by diffusion. An additional WACCM run with different gravity waves in the SH is used as a sensitivity test. Over the Antarctic, this test has small impact on M_y, but it significantly modifies A_z in the austral fall and winter due to the enhanced downwelling, and the N₂O mixing ratio during spring as a consequence of a more realistic timing of the vortex breakup.

The interannual variability of the mid-stratospheric horizontal mixing term M_y is largest in the polar regions. It is related to the variability in the vortex breakup dates during spring in the Antarctic, whereas it is related to the highly variable polar vortex in winter in the Arctic. The interannual variability of A_z is characterized by a large spread in the mid-stratospheric tropical regions where WACCM-CCMI and JRA-55 deliver a smaller contribution than the other reanalyses. This variability reflects the impact of the QBO on the tropical upwelling (Abalos et al., 2015).

The application of the TEM framework to tracer transport with reanalyses suffers from a poor closure of the budget in the polar regions. We chose to analyze these regions nonetheless, because the differences in M_y between WACCM and the reanalyses are larger than the residual term, but it remains important to better understand the causes of these large uncertainties. To this end, detailed studies of transport in the polar stratosphere are needed, e.g., comparing the residual circulations with indirect estimates derived from momentum and thermodynamic balances as well as evaluating the effective diffusivity in each dataset (Abalos et al., 2015, 2016a).

The next step of this research consists of the analysis of the interannual variations of the BDC, including the impact of the QBO and the El Niño–Southern Oscillation. Further extensions of this work would include the addition of new reanalysis products such as ERA5 and an intercomparison of several CCMs as already done for the residual circulation itself (Chrysanthou et al., 2019).