Effects of denitrification on the distributions of trace gas abundances in the polar regions: a comparison of WACCM with observations

Polar stratospheric clouds (PSCs) play a key role in the polar chemistry of the stratosphere. Nitric acid trihydrate (NAT) particles have been shown to lead to denitrification of the lower stratosphere. While the existence of large NAT particles (NAT “rocks”) has been verified by many measurements, especially in the Northern Hemisphere (NH), most current chemistry–climate models use simplified parameterizations, often based on evaluations in the Southern Hemisphere where the polar vortex is stable enough that accounting for NAT rocks is not as important as in the NH. Here, we evaluate the probability density functions of various gaseous species in the polar vortex using one such model, the Whole Atmosphere Community Climate Model (WACCM), and compare these with measurements by the Michelson Interferometer for Passive Atmospheric Sounding onboard the Environmental Satellite (MIPAS/Envisat) and two ozonesonde stations for a range of years and in both hemispheres. Using the maximum difference between the distributions of MIPAS and WACCM as a measure of coherence, we find better agreement for HNO₃ when reducing the NAT number density from the standard value of 10⁻² used in this model to $\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ for almost all spring seasons during the MIPAS period in both hemispheres. The distributions of ClONO₂ and O₃ are not greatly affected by the NAT density. The average difference between WACCM and ozonesondes supports the need to reduce the NAT number density in the model. Therefore, this study suggests using a NAT number density of $\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ for future simulations with WACCM.

1 Introduction

For decades, polar stratospheric clouds (PSCs) have been known to play a key role in explaining the stratospheric ozone hole (Solomon et al., 1986; Solomon, 1999; Tritscher et al., 2021). Reactions on their surfaces lead to the activation of chlorine reservoirs, known as chlorine activation. In addition, sedimentation of PSCs results in dehydration and denitrification of the lower stratosphere. Three types of PSCs with different composition and roles have been found to be important: (1) liquid supercooled ternary solution (STS) droplets (Carslaw et al., 1994) that are major contributors to chlorine activation (e.g., Peter, 1997); (2) ice PSCs that lead to dehydration of the stratosphere (e.g., Kelly et al., 1989); and (3) nitric acid trihydrate (NAT) particles that dominate the irreversible removal of nitric acid (HNO₃) from the lower stratosphere via sedimentation, known as denitrification, thereby potentially reducing the reformation of chlorine reservoir species which can affect ozone loss (e.g., Waibel et al., 1999).

The existence of NAT particles in the stratosphere has been verified by a variety of airborne measurements (Fahey et al., 1989; Voigt et al., 2000; Fahey, 2001; Molleker et al., 2014; Woiwode et al., 2016). As the polar vortex in the Southern Hemisphere (SH) is generally more stable than in the Northern Hemisphere (NH) (e.g., Schoeberl et al., 1992), the time period when denitrification can occur is much longer in the SH polar vortex than in the NH. On the other hand, denitrification in the NH has been found to occur locally, and the role of low-number-density large-size NAT particles, so-called NAT “rocks”, has been discussed and verified by measurements (Fahey, 2001; Fueglistaler et al., 2002; Drdla and Müller, 2012; Adriani et al., 2004; Woiwode et al., 2014, 2016). These large particles can lead to significant sedimentation, even on the shorter timescales needed to explain the occurrence of low HNO₃ volume mixing ratios (VMRs) in the NH spring.

Satellite instruments are able to provide measurements of PSCs and gaseous HNO₃ with daily near-global coverage (Höpfner et al., 2018; Spang et al., 2018; Pitts et al., 2018; Santee et al., 2007; Wespes et al., 2022). Höpfner et al. (2006) found the first evidence of a NAT belt around the Antarctic continent using the Michelson Interferometer for Passive Atmospheric Sounding onboard the Environmental Satellite (MIPAS/Envisat). Therefore, satellite measurements provide a unique opportunity for comparisons of HNO₃ and denitrification with current chemistry–climate models in both hemispheres. They also the allow study of the range of the probability distributions in the data for different years and conditions.

Accounting for denitrification is an important process in chemistry–climate and in chemistry–transport models, which is why many parameterizations with different levels of detail have been developed to account for the microphysics and sedimentation of NAT particles in these models (e.g., Considine et al., 2000; Carslaw et al., 2002; Grooß et al., 2005; Wohltmann et al., 2010; Wegner et al., 2013; Zhu et al., 2015; Kirner et al., 2011; Weimer et al., 2021). Some models account for the comprehensive microphysics of NAT particles (e.g., Zhu et al., 2015) to determine radii and concentrations which then redistribute gaseous HNO₃. Others use an intermediate approach, not considering the full microphysics but accounting for the amount of NAT with a measured NAT distribution (e.g., Kirner et al., 2011; Weimer et al., 2021). In addition, there are models with diagnostic approaches, as suggested by studies such as Considine et al. (2000), where NAT is formed based on the available gaseous HNO₃, sedimentation is computed and the HNO₃ contained in NAT particles is released back to the gaseous phase within the same model time step.

The standard version of the Whole Atmosphere component of the Community Earth System Model includes a diagnostic parameterization of NAT with a NAT number density set to a global value of 10⁻² cm⁻³ (Wegner et al., 2013). However, Wilka et al. (2021) found that this value leads to an overestimation of gaseous HNO₃ in the NH for the extraordinarily cold NH winter of 2019/2020 in comparison to measurements of the Microwave Limb Sounder. They found a better agreement with the measurements using an adopted NAT number density of $N_{\mathrm{NAT}}={\mathrm{10}}^{-\mathrm{5}}$ cm⁻³. Here, we investigate this further by applying an adapted version of the method by Zambri et al. (2021) to compare the distributions of various gaseous species within the polar vortex with MIPAS, which provided measurements during 2002–2012. MIPAS data are well suited to this task, as this instrument provides high signal-to-noise local measurements over multiple years in both hemispheres, with good dynamic range; other datasets may also be appropriate, but we use MIPAS here for this initial test. Further, using MIPAS and varying the NAT number density of the model, we can also investigate its associated impact on HNO₃, ClONO₂ and O₃ for many spring seasons in both hemispheres.

Details about the model and the observation as well as the methods to compare them can be found in Sect. 2. Section 3 shows the comparison between the datasets. Finally, some concluding remarks are given in Sect. 4.

2 Datasets

2.1 The Whole Atmosphere Community Climate Model (WACCM)

In this study, we use the Whole Atmosphere component (WACCM6) of the Community Earth System Model (CESM2.1) in Specified-Dynamics (SD) mode (Gettelman et al., 2019; Danabasoglu et al., 2020) to compare distributions of the trace gases with the satellite measurements. The model is relaxed towards the Modern-Era Retrospective Analysis for Research and Applications, version 2 (MERRA-2; Gelaro et al., 2017). In this study, we use a horizontal resolution of 1.9^∘ latitude × 2.5^∘ longitude and 88 vertical levels up to about 140 km. In the component set “FWmaSD”, a comprehensive chemistry for the middle atmosphere is included, as also used by studies such as Zambri et al. (2021).

Polar stratospheric clouds in WACCM are calculated using a diagnostic parameterization described by Considine et al. (2000), Kinnison et al. (2007) and Wegner et al. (2013). NAT particles are formed in thermodynamical equilibrium with the gaseous H₂O and HNO₃ below the NAT threshold temperature (Hanson and Mauersberger, 1988). The coexistence of NAT and STS is accounted for by allowing 20 % of HNO₃ to form NAT whereas the rest is available for STS (Wegner et al., 2013; Solomon et al., 2015). Sedimentation of the NAT particles, i.e., the vertical redistribution of gaseous HNO₃, is calculated using a simple upwind scheme (Considine et al., 2000). The radius of the particles is determined in this scheme by using the amount of condensed HNO₃ and the NAT number density, which is a global parameter in this scheme. Larger NAT density leads to smaller particles and vice versa, assuming the mass of HNO₃ condensed in NAT is constant in the grid box. Therefore, the NAT density is a parameter in this scheme that can be used to tune the denitrification in the model. Examination of observations of NAT particle number densities (Pitts et al., 2009, 2011) and an emphasis on HNO₃ data for the SH in 2005 led Wegner et al. (2013) to adopt a NAT value of 0.01 cm⁻³. Fahey (2001) measured NAT number densities between $\mathrm{2.3}\times {\mathrm{10}}^{-\mathrm{4}}$ and $\mathrm{2}\times {\mathrm{10}}^{-\mathrm{3}}$ cm⁻³. Pitts et al. (2009) and Pitts et al. (2011) used a wide range of N_NAT between 10⁻⁴ and 10⁻¹ cm⁻³ based on in situ measurements. In mountain waves, high NAT densities larger than 10⁻² cm⁻³ have been measured (e.g., Carslaw et al., 1998). Voigt et al. (2005) observed a NAT PSC with 10⁻⁴ cm⁻³. Here, we have had the benefit of the unusual NH year of 2020 and have placed more emphasis on both hemispheres to derive a parameterization that better represents the remaining gas-phase HNO₃ in the maximum number of years and for both hemispheres. However, it is important to emphasize that the model's NAT parameterization is subject to multiple simplifications of complex microphysics. For instance, NAT particles in a grid box with a horizontal extent of about 100×100 km² do not necessarily have the same radius but usually follow a multi-modal size distribution (e.g., Fahey, 2001). In addition, NAT particles are not allowed to change their size over time nor be transported while interacting with the atmosphere by nucleation and (re-)sublimation, and not all of the supersaturated amount of gaseous HNO₃ will become NAT. Therefore, observed NAT particle abundances may not be the best guide for this parameter choice.

We performed sensitivity simulations varying the NAT number density within the observed range from 10⁻² to 10⁻⁵ cm⁻³ in the current version 6 of WACCM. We also performed a simulation excluding stratospheric heterogeneous reactions, apart from the N₂O₅ + H₂O reaction, to evaluate the impact of heterogeneous processes on the gas-phase species. The N₂O₅ + H₂O reaction rate is nearly independent of temperature, also happens on the background aerosols, does not directly affect the halogen chemistry and is important for the partitioning of reactive nitrogen in the atmosphere; this is why we kept this reaction in that simulation. Without this reaction, the chemistry in the model would be completely unrealistic. All simulations are summarized in Table 1. They cover the satellite period starting from 1979 to the present.

Table 1 Simulations in this study. In the “noHetAll” simulation, the heterogeneous reactions are removed apart from the N₂O₅ +H₂O reaction.

2.2 The Michelson Interferometer for Passive Atmospheric Sounding (MIPAS)

The Michelson Interferometer for Passive Atmospheric Sounding (MIPAS) operated in limb geometry aboard the Environmental Satellite (Envisat) between July 2002 and April 2012 (Fischer et al., 2008). Envisat was placed in a sun-synchronous polar orbit at an altitude of around 800 km with more than 14 orbits per day. MIPAS measured a variety of trace gases including HNO₃, ClONO₂ and O₃ using a Fourier transform spectrometer in the infrared spectral range between 4.15 and 14.6 µm at tangent altitudes from 7 to 72 km (Fischer et al., 2008). The spatial resolution was approximately 3 km in the vertical and 30 km in the horizontal. The spectral resolution was 0.05 cm⁻¹ between 2002 and March 2004. Due to technical issues with the satellite and a corresponding gap in the measurements, MIPAS started measuring again in January 2005 with a reduced spectral resolution of 0.12 cm⁻¹.

Figure 1 Daily vortex-averaged HNO₃ VMR profiles in the analyzed pressure range between 30 and 150 hPa during local spring for MIPAS (a, f, k, p) and the WACCM simulations (b–e, g–j, l–o, q–t) in the NH (a–j) and SH (k–t) and for the years 2003 and 2011 (the first and second columns in the hemispheres). The number of profiles used to average each day depends on the availability of MIPAS data and the size of the polar vortex.

For this study, we use the V8 Level-2 MIPAS retrievals from the IMK-IAA processor (Kiefer et al., 2021, 2023; von Clarmann et al., 2009) of HNO₃, ClONO₂ and O₃ to compare their distributions with the WACCM simulations. The reported precision in the region between 30 and 150 hPa is 5 %–15 % (around 100 pptv) for HNO₃ (Sheese et al., 2016, 2017), 10 % (around 150 ppbv) for O₃ (Sheese et al., 2017; Kiefer et al., 2023) and 7 %–32 % (24–89 pptv) for ClONO₂ (Sheese et al., 2016; Höpfner et al., 2007).

The method used to compare the MIPAS data with WACCM is based on the approach by Zambri et al. (2021). It takes advantage of WACCM's ability to provide output at the locations and times closest to the MIPAS profiles. In combination with interpolation of the MIPAS vertical levels to the WACCM altitudes, this enables a direct comparison of the two datasets. As a result, probability density functions (PDFs) can be evaluated in the spatiotemporal range of interest and compared to the observations.

As NAT formation and denitrification are strongly temperature dependent and most efficient at the lowest temperatures, we compute PDFs for profiles inside the polar vortex, determined by MERRA-2 using the Nash criterion (Nash et al., 1996). As the sedimentation of the NAT particles takes several weeks (Tabazadeh et al., 2001), the largest effects of denitrification can be expected to be most amplified after the local winter, which is why we restrict our analysis to the early local spring, i.e., 1 February to 15 March in the Northern Hemisphere and 1 September to 15 October in the Southern Hemisphere, and to the pressure range from 30 to 150 hPa. Sensitivity studies that involved changing these two ranges did not affect the main message of the results shown in the next section. In order to make the datasets comparable, we remove the profiles from both WACCM and MIPAS data where negative values occur in the measurements. This does not have a significant effect on the results. The number of profiles depends on the area of the polar vortex and differs from year to year. Generally, the number of profiles is larger in the SH due to the larger size of the SH polar vortex compared with the NH.

Figure 2 Scatterplots of MIPAS and WACCM HNO₃ in (a) the NH and (b) the SH for the HetAll simulation using the three NAT number densities described in the text. The dots include all measurement and model data between 30 and 150 hPa. The dashed line corresponds to the one-to-one line.

Zambri et al. (2021) used the Kolmogorov–Smirnov test to evaluate the difference in the PDFs between the model and observations, which uses the number of data points and the maximum difference between the cumulative density functions (CDFs) to check whether the distributions are distinguishable in a statistical sense related to a significance level α. Here, we use the maximum difference between the CDFs F of MIPAS and WACCM “max(d)” to directly evaluate the difference in the distributions (Zambri et al., 2021):

$$\begin{array}{}\text{(1)}& \max \left(\mathrm{d}\right)={\max }_x\left|F_{\mathrm{MIPAS}}\left(x\right)-F_{\mathrm{WACCM}}\left(x\right)\right|,\end{array}$$

where x is the binned trace gas VMR. A more detailed discussion can be found in Zambri et al. (2021).

2.3 Ozonesondes

We also compare the WACCM simulations to balloon-borne in situ measurements of ozonesondes, made available by the World Ozone and Ultraviolet Data Centre. They use an electrochemical concentration cell to measure ozone profiles with a precision of 3 % to 5 % and an uncertainty of about ±10 % in the pressure range of interest here (Smit et al., 2007). We use ozonesondes of two stations: Eureka in the NH (80.04^∘ N, 86.1^∘ W) and Syowa in the SH (69.00^∘ S, 39.58^∘ E).

The method to compare the WACCM simulations with the ozonesonde measurements is similar to Wilka et al. (2021): the data are compared to the daily averaged WACCM ozone concentration at the grid point closest to the respective station for profiles within the polar vortex. Averages and spread of the differences between WACCM and the observations are evaluated.

3 Impact on the distribution of gas-phase species

In this section, we will analyze the influence of changing denitrification in WACCM on the trace gases HNO₃, ClONO₂ and O₃ during early spring in both hemispheres. Generally, these species are expected to be influenced by a changed NAT number density because it leads to a redistribution of HNO₃ (lower HNO₃ at high altitudes and higher HNO₃ at lower altitudes). During early spring, HNO₃ is photolyzed and forms NO₂ which combines with ClO to form ClONO₂. ClO is responsible for part of the catalytic ozone depletion and deactivated by the reaction into ClONO₂.

As a starting point, Fig. 1 shows vortex-mean profiles of the HNO₃ VMR during the local spring months for both hemispheres and the years 2003 and 2011 for MIPAS and all WACCM simulations of this study. Analysis of the other years of the MIPAS period showed similar results, which is why we focus on the two example years 2003 and 2011 here. The signature of denitrification is indicated in the MIPAS measurements (first row of Fig. 1, panels a, f, k and p) by the HNO₃ VMR decreasing over time at around 30 to 60 hPa, with corresponding increases due to renitrification at the end of the shown time series at lower levels, like in 2011 in the NH (Fig. 1f). In the simulation with the current operational NAT number density of 10⁻² cm⁻³, HetAll.1e-2 in the second row of Fig. 1 (panels b, g, l and q), the NAT particles do not get large enough for the simulated HNO₃ to compare well with the measured denitrification of the shown years. This means that larger HNO₃ VMRs than observed remain at higher altitudes at the end of the time series, whereas the VMR is underestimated at lower altitudes. In contrast to this, with the smallest tested number density of 10⁻⁵ cm⁻³, HetAll.1e-5 in the last row of Fig. 1 (panels e, j, o and t), the HNO₃ VMRs become too low, which indicates that the NAT particles become too large. By using a number density of $\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ (HetAll.5e-4 in the third row of Fig. 1, panels d, i, n and s), the differences compared to MIPAS HNO₃ are smallest and the altitude regions of denitrification and renitrification are comparable to the measurements.

In order to generalize this comparison, Fig. 2 shows correlation plots of the HetAll simulations using different NAT number densities compared with the MIPAS measurements for the whole pressure range (30–150 hPa) and all MIPAS years in the NH and SH. The HetAll.5e-4 simulation (green) has a relatively compact correlation with MIPAS in both hemispheres, whereas the other simulations show a larger scatter. The HetAll.1e-2 simulation shows HNO₃ VMRs that are too large compared with the observations because of NAT particles that are too small. As already suggested by Fig. 1, the HNO₃ VMRs in HetAll.1e-5 are too small compared with MIPAS as a result of overly large NAT particles which lead to sedimentation below the shown pressure range. As this figure includes the data of all MIPAS years and the whole pressure range, we will investigate the PDFs of this pressure range in the following, using the approach based on Zambri et al. (2021).

PDFs of HNO₃ for the whole MIPAS period from 2002 to 2012 and for the years 2003 and 2011 are shown in Fig. 3 within the polar vortex in both hemispheres. All individual years can be found in the Supplement (Fig. S1). The blue line corresponds to the simulation with the largest NAT density, the red line shows the simulation with smallest NAT density, and the green and yellow lines show the simulations using $N_{\mathrm{NAT}}=\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³, with the latter probing turning off all heterogeneous chemistry apart from N₂O₅ + H₂O. The black line corresponds to the MIPAS HNO₃ PDF.

Figure 3 HNO₃ PDFs within the polar vortex for MIPAS and the sensitivity simulations of WACCM for the NH (a–c) and SH (d–f) spring months in (a, d) the MIPAS period from 2002 to 2012, (b, e) 2003 and (c, f) 2011. Vertical dashed lines show the average VMR of the respective simulations. The number of profiles used to derive the PDFs in each panel is denoted by “nprof”. The maximum differences between the CDFs of MIPAS and the respective WACCM simulation are denoted by “max(d)” in the panels.

Figure 3 demonstrates that the NAT density has a large impact on the PDFs in the model between 30 and 150 hPa. In HetAll.1e-2, larger HNO₃ values are more frequent compared with the other simulations in all of the panels. In contrast to this, HetAll.1e-5 underestimates HNO₃, as manifested by a higher frequency at smaller values. The simulation with peak values closest to the MIPAS observations in all panels of Fig. 3 is HetAll.5e-4, indicating that heterogeneous chemistry is relevant and that $N_{\mathrm{NAT}}=\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ leads to PDFs comparable to MIPAS for all years (see also Fig. S1). However, the NAT parameterization in WACCM is not able to capture the minimum and maximum values; this is especially notable in the SH where VMRs as high as 12.5 ppbv have been measured by MIPAS, whereas the largest values in HetAll.5e-4 are around 10 ppbv (e.g., in Fig. 3d).

The overall agreement is also shown by the mean values illustrated by the dashed vertical lines in the figure. Moreover, the mean HNO₃ value for HetAll.5e-4 is closest to the MIPAS mean value in all of the panels. In addition, the maximum difference, max(d), of the HetAll.5e-4 simulations is smallest in the shown years compared with the other simulations.

In almost all years, heterogeneous chemistry leads to about 0.5 ppbv smaller HNO₃ mixing ratios in the NH (see also Fig. S1). In the SH, this effect is not as large as in the NH. In total, the effect of halogen heterogeneous chemistry on HNO₃ is small. The mean and shape of the distributions are similar for noHetAll.5e-4 and HetAll.5e-4.

As the number of profiles differs from year to year depending on the development of the polar vortex, more weight in the distributions of the whole period from 2002 to 2012 (Fig. 3a and d) is given to years with a larger polar vortex. This is why the distributions of all individual years are shown in the Supplement (Figs. S1, S2, S3). As a summary, the max(d) values of the individual years are shown in Tables 2 and 3 for the NH and SH, respectively. The simulation with the minimum difference is highlighted in bold for each year. Apart from some exceptions in the NH, the minimum difference occurs for the simulations using $N_{\mathrm{NAT}}=\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³. This is another indication that this value improves denitrification and the distribution of HNO₃ in the polar vortex in both hemispheres.

Table 2 Maximum differences between WACCM and MIPAS HNO₃ CDFs (“max(d)” in the figures) for all years (1 February to 15 March in the NH). Bold font is used to highlight the minimum values in each row.

Table 3 The same as Table 2 but for the SH spring (1 September to 15 October).

Exceptions are the years 2004 and 2012 in the NH, where the simulation with the smallest NAT density (HetAll.1e-5) shows the minimum differences compared with the MIPAS distribution (see Table 2). In 2004, relatively warm temperatures in the NH stratosphere lead to low denitrification (Manney et al., 2005). As shown in Fig. 4a, this is reflected by the WACCM HNO₃ distributions which do not show large differences when changing the NAT density during that year. However, there seems to be a systematic difference between WACCM and MIPAS in this year, maybe due to a lower accuracy of the instrument or, for instance, a mismatch in the location of the polar vortex in the reanalysis data used for the model and the temperature history during that winter. The 2011/2012 Arctic winter was characterized by an early breakdown of the polar vortex at the end of December (Roy and Kuttippurath, 2022), but there is some evidence of large NAT particles during the cold period during December 2011 (Woiwode et al., 2019). This is why the HNO₃ distributions change in WACCM and large NAT particles probably play a larger role in determining the distributions for that year (see Fig. 4b).

Figure 4 The same as Fig. 3 but for the years (a) 2004 and (b) 2012 in the NH.

Moderate changes due to the different NAT densities can also be found in other species. Figure 5 shows the distributions of chlorine nitrate (ClONO₂) for the same years and simulations as in Fig. 3. Both HetAll.1e-2 and HetAll.5e-4 have similar distributions and max(d) values. Depending on the year, minimum max(d) values occur in both simulations. HetAll.1e-5 seems to slightly underestimate the ClONO₂ mixing ratios in both hemispheres, illustrated by higher frequencies of low mixing ratios.

Figure 5 The same as Fig. 3 but for ClONO₂.

As a result of increased denitrification, impacts on the maximum values in the ClONO₂ distributions can be seen in both hemispheres. While the maximum ClONO₂ mixing ratios in the NH are closer to the MIPAS distributions using $N_{\mathrm{NAT}}={\mathrm{10}}^{-\mathrm{5}}$ cm⁻³, the HetAll.5e-4 simulation captures the maximum values in the SH distributions.

There is also a clear signature of heterogeneous chemistry in the ClONO₂ distributions. The deactivation reaction forming ClONO₂,

$$\begin{array}{}\text{(R1)}& \mathrm{ClO}+{\mathrm{NO}}_{\mathrm{2}}+\mathrm{M}\to {\mathrm{ClONO}}_{\mathrm{2}}+\mathrm{M},\end{array}$$

is usually faster than the deactivation reaction forming HCl,

$$\begin{array}{}\text{(R2)}& \mathrm{Cl}+{\mathrm{CH}}_{\mathrm{4}}\to \mathrm{HCl}+{\mathrm{CH}}_{\mathrm{3}},\end{array}$$

as long as NO₂ is available and the ozone VMRs are large enough (e.g., Kawa et al., 1992; Prather and Jaffe, 1990). At the edge regions of the polar vortex, some of the chlorine originally in HCl at the start of the winter shows up in ClONO₂ at the end of the winter, known as the ClONO₂ “collar” (e.g., Toon et al., 1989; Douglass et al., 1995; Solomon et al., 2016). Therefore, while the simulation without heterogeneous chemistry, noHetAll.5e-4, generally shows ClONO₂ VMRs that are too low compared with MIPAS, it is notable that the maximum values are increased by about 1 ppbv and are then comparable to MIPAS in all simulations with heterogeneous chemistry.

Similar impacts can be seen in the distributions of ozone (O₃), shown in Fig. 6. Ozone VMRs are about 1 ppmv too large in the SH for noHetAll.5e-4 compared with MIPAS. In the NH, the difference is smaller, but ozone is still overestimated in the simulation without heterogeneous chemistry.

Figure 6 The same as Fig. 3 but for O₃.

Wilka et al. (2021) showed that the agreement of WACCM model calculations of Arctic spring ozone losses in 2020 with ozonesonde data was heavily dependent on the denitrification parameterization used in the model, providing an important check. Figure 7 displays comparisons between WACCM and ozonesonde data at Syowa and Eureka for their respective spring seasons for all years between 2002 and 2012. A figure showing the comparison to Syowa without noHetAll.5e-4 can be found in the Supplement (Fig. S4). Similar comparisons are provided for the South Pole in Fig. S5. Because local Antarctic ozone abundances can drop to values lower than 1 ppmv between about 15 and 20 km in years with large ozone losses (see, e.g., Fig. 7a), care must be taken in interpreting percentage differences; we also provide absolute differences for that reason. The Syowa sondes agree considerably better with the model when the NAT number density of $\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ or less is used, compared with larger values of this parameter. At Eureka, results are far less sensitive due to smaller ozone losses, except in very cold years (like the year 2020, as shown in Wilka et al., 2021).

Figure 7 Mean and interquartile (0.25 and 0.75) range of (a, d) the volume mixing ratio as well as (b, e) relative differences and (c, f) absolute differences between the respective WACCM simulations and two ozonesonde stations for profiles within the polar vortex using the same definition as in the previous analysis. Panels (a)–(c) show Syowa Station in the Antarctic, and panels  (d)–(f) show Eureka Station in the Arctic. For WACCM, the grid point closest to the station is chosen. The average and distributions are taken for all profiles within the MIPAS period and the same months as for the satellite comparison. The first column (a, d) shows the O₃ VMR, the second column (b, e) shows the relative difference between WACCM and the ozonesonde measurements, and the third column (c, f) shows the absolute difference between WACCM and the ozonesonde measurements.

In summary, the results in this section showed that the diagnostic parameterization of NAT particles in WACCM is able to reflect the main features in the three key gas-phase species, HNO₃, ClONO₂ and ozone. When using a reduced NAT particle number density of $N_{\mathrm{NAT}}=\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ compared with the standard WACCM setup, the difference in the CDFs is minimized for almost all spring months of the MIPAS period from 2002 to 2012, and the bias between 30 and 150 hPa in comparison to the ozonesondes is reduced. This indicates that $N_{\mathrm{NAT}}=\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ improves the distribution of these species in the polar regions and that it should be used in future WACCM simulations.

4 Discussion and conclusions

In this study, we investigated the impact of the NAT number density (N_NAT) on the distributions of key gaseous species in the WACCM model. We compared PDFs within the polar vortex with the MIPAS satellite instrument which operated between 2002 and 2012. As a measure of the difference between the distributions of WACCM and MIPAS, we used the maximum difference, max(d), between the CDFs. To identify the impact of NAT on the distributions of the species, we performed sensitivity simulations varying N_NAT between 10⁻² and 10⁻⁵ cm⁻³. In the method applied here to compare the satellite measurements with the model, negative values in the satellite retrievals are cut off, which is why the method cannot be used for noisy data (such as N₂O₅ and ClO in MIPAS) and we restricted the analysis to the gaseous species HNO₃, ClONO₂ and O₃.

The diagnostic parameterization of NAT in the WACCM model was shown to be able to reproduce the general shape of the MIPAS distributions for almost all years of the MIPAS period in the polar vortices of both hemispheres. However, there was a general overestimation of gaseous HNO₃ when using the standard N_NAT of 10⁻² cm⁻³. Reducing this number concentration to $\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ was shown to also reduce the differences in the HNO₃ distributions between MIPAS and WACCM for almost all spring seasons in both hemispheres during the MIPAS period. Changes in the distributions of ClONO₂ and O₃ due to the new value of N_NAT were either negligible or the differences were reduced. Mean differences between the ozonesonde profiles and the WACCM simulations at the grid point closest to the Syowa and Eureka sonde stations were also shown to be decreased when reducing N_NAT. Therefore, this study suggests the use of $N_{\mathrm{NAT}}=\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ for future simulations. This new value is within the range of measured NAT number densities.

A simplified parameterization as applied in WACCM cannot capture all features of the distributions of all years. This was demonstrated by the years 2004 and 2012 where differences between the MIPAS and WACCM distributions were increased. Further, there are many uncertainties in the detailed chemistry from one year to another that may also be important for denitrification and ozone loss, including, for example, volcanic particle inputs, gravity and planetary wave amplitudes, and changes in circulation, to name only a few. It can also be seen that the smallest number density of 10⁻⁵ cm⁻³ is comparable to measurements during the exceptionally cold northern winter of 2019/2020 (Wilka et al., 2021). In addition, although overall differences in the PDFs are reduced with $N_{\mathrm{NAT}}=\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³, the Kolmogorov–Smirnov test, used by studies such as Zambri et al. (2021) with a similar approach to that used in this study, will fail if it does not employ a significance level α that is set to a very small value (O(10⁻⁷⁰)). Therefore, although the differences are decreased, the distributions of MIPAS and WACCM are still different in a statistical sense, probably due to the simplifications in the NAT parameterization of WACCM.

Future studies could apply this methodology to other long-term satellite measurements to evaluate the usage of this new NAT number density for further years, although their precision has to be low enough to allow for comparisons with individual profiles. Nevertheless, as a follow-up of the publication by Wilka et al. (2021), we recommend using $N_{\mathrm{NAT}}=\mathrm{5}\times {\mathrm{10}}^{-\mathrm{4}}$ cm⁻³ for future simulations with WACCM, indicating that NAT rocks play an important role, especially in the NH where we saw the largest impact on changing the NAT number density.