Effects of organic coating on the nitrate formation by suppressing the N₂O₅ heterogeneous hydrolysis: a case study during wintertime in Beijing–Tianjin–Hebei (BTH)

Although stringent emission mitigation strategies have been carried out since 2013 in Beijing–Tianjin–Hebei (BTH), China, heavy haze with high levels of fine particulate matter (PM_2.5) still frequently engulfs the region during wintertime and the nitrate contribution to PM_2.5 mass has progressively increased. N₂O₅ heterogeneous hydrolysis is the most important pathway of nitrate formation at nighttime. In the present study, the WRF-Chem model is applied to simulate a heavy haze episode from 10 to 27 February 2014 in BTH to evaluate contributions of N₂O₅ heterogeneous hydrolysis to nitrate formation and effects of organic coating. The model generally performs reasonably well in simulating meteorological parameters, air pollutants, and aerosol species against observations in BTH. N₂O₅ heterogeneous hydrolysis with all the secondary organic aerosol assumed to be involved in coating considerably improves the nitrate simulations compared to the measurements in Beijing. On average, organic coating decreases nitrate concentrations by 8.4 % in BTH during an episode, and N₂O₅ heterogeneous hydrolysis with organic coating contributes about 30.1 % of nitrate concentrations. Additionally, the reaction also plays a considerable role in the heavy haze formation, with a PM_2.5 contribution of about 11.6 % in BTH. Sensitivity studies also reveal that future studies need to be conducted to predict the organic aerosol hygroscopicity for accurately representing the organic coating effect on N₂O₅ heterogeneous hydrolysis.

1 Introduction

Within recent decades, China has been suffering from pervasive and persistent haze pollution caused by elevated levels of fine particulate matters (PM_2.5), particularly in Beijing–Tianjin–Hebei (BTH) (Guo et al., 2014; Gao et al., 2016; Wang et al., 2016). Numerous studies have revealed that the inorganic aerosols, including nitrate, sulfate, and ammonium, are the most abundant component of PM_2.5 during haze pollution episodes in BTH, and that the evolution of the haze pollution is characterized by the formation of substantial amounts of sulfate and nitrate (Sun et al., 2013, 2015; Zhang et al., 2013; Zhao et al., 2013). Since 2013, several aggressive emission control strategies have been implemented in China, including desulfurization and dedusting for coal combustion, vehicle restriction and execution of stringent emission standards in key industries (Tao et al., 2017). However, the control of emissions of nitrate gaseous precursors does not seem to be effective since many observations have shown that the nitrate aerosol concentration has progressively increased in recent several years (Zhang et al., 2012, 2015; Sun et al., 2015; Tao et al., 2017).

In the atmosphere, nitrate aerosol is formed via nitrous acid (HNO₃) to balance the inorganic cations in the aerosol phase. HNO₃ is produced through four pathways (Kim et al., 2014): (1) the reaction of OH and NO₂ (main gas-phase pathway and usually considered the daytime pathway because the OH radical is severely limited at night due to lack of O₃ and peroxide photolysis), (2) NO₃ radical reaction with hydrocarbons, (3) aqueous reaction of NO₃ radical to form HNO₃, and (4) NO₃ conversion to N₂O₅ with subsequently heterogeneous chemical conversion to form HNO₃. The last pathway is referred to as the most important pathway during nighttime since both NO₃ and N₂O₅ are photolytically liable, or even under a heavy haze situation with weak sunlight and high relative humidity (RH) (Brown et al., 2016).

The heterogeneous hydrolysis of N₂O₅ on the surface of deliquescent aerosols to form HNO₃ is quantified by the reaction probability (${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$) (Bertram and Thornton, 2009; Chen et al., 2018; Davis et al., 2008; Riemer et al., 2003). ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ has been measured by previous laboratory experiments, dependent on particulate chemical composition, RH, temperature, aerosol surface area, and water content (Chang et al., 2011), and on the order of 10⁻² (Zheng et al., 2015). In modeling studies, various parameterizations of ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ are used to simulate the nitrate formation. Dentener and Crutzen (1993) first used 0.1 as the representative ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ in a three-dimensional global model. Riemer et al. (2003) developed a ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ parameterization on the surface of aerosols containing sulfate and nitrate (hereafter referred to as Riemer03), which has widely been used and further improved in air quality models. Davis et al. (2008) implemented a ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ parameterization on the surface of particles containing sulfate, nitrate, and ammonium as a function of RH, temperature, and phase state to improve simulations of N₂O₅ hydrolysis. Bertram and Thornton (2009) developed a parameterization to consider the influence of chloride salts on ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ as a function of RH.

The coating of particles by organic materials has been reported to inhibit N₂O₅ uptake (Anttila et al., 2006) and suggested as a possible explanation for field observations of suppressed N₂O₅ uptake (Brown et al., 2006). Evans and Jacob (2005) have incorporated a ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ parameterization on surfaces of sulfate particles as a function of RH and temperature into the GEOS-CHEM model, including the effects of dust, sea salt, sulfate, elemental carbon, and organic carbon but ignoring nitrogen-containing species. Riemer et al. (2009) developed a N₂O₅ uptake parameterization (Riemer09) based on the laboratory results of Anttila et al. (2006), which combines nitrogen-containing and organic effects on N₂O₅ hydrolysis. The parameterization has been used to estimate the maximum effect of organic coating by assuming that all available secondary organic compounds (SOCs) contribute to the coating in a 3-D model. The results show that SOC could suppress N₂O₅ uptake significantly, reducing particulate nitrate concentrations by up to 90 %. Lowe et al. (2015) further combined the organic coating and chloride salt effects on ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ in the WRF-Chem model. Most recently, Chen et al. (2018) developed a new ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ parameterization with respect to RH, temperature, and aerosol composition, showing that organic coating effect on ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ is not as important as expected over western and central Europe. However, there is still a lack of modeling studies focused on the effect of organic coating on ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ and particulate nitrate formation in China. Wang et al. (2017) evaluated the potential particulate nitrate formation through the N₂O₅ hydrolysis reaction without considering the organic coating effect during a haze pollution episode in Beijing, and found that the observed nitrate concentration (20.6 µg m⁻³ on average) is lower than the assessment (57.0 µg m⁻³ on average). Considering the high organic aerosol concentration and increasing trend of particulate nitrate during haze days in BTH, it is imperative to assess the effect of organic coating on N₂O₅ hydrolysis and its consequent contribution to nitrate formation.

In the present study, based on Riemer09 parameterization, the contribution of the organic coating effect on N₂O₅ hydrolysis to nitrate formation is investigated using the WRF-Chem model. The model configuration and methodology are described in Sect. 2. Results and sensitivity studies are presented in Sect. 3. Discussion and summary are given in Sect. 4.

2 Model and methodology

2.1 WRF-Chem model and configuration

A modified version of the WRF-Chem model (Grell et al., 2005) is used in this study, which is developed by Li et al. (2010, 2011a, b, 2012) at the Molina Center for Energy and the Environment. A new flexible gas-phase chemical module has been developed and implemented into the version of the WRF-Chem model, which can be utilized with different chemical mechanisms, including CBIV, RADM2, and SAPRC. The gas-phase chemistry is solved by an Eulerian backward Gauss–Seidel iterative technique with a number of iterations, inherited from NCAR-HANK (Hess et al., 2000). In the study, the SAPRC99 chemical mechanism is used based on the available emission inventory. For the aerosol simulations, the CMAQ/Models-3 aerosol module (AERO5) developed by U.S. EPA has incorporated into the model (Binkowski and Roselle, 2003). Briefly, wet deposition uses the method in the CMAQ module and dry deposition of chemical species is parameterized following Wesely (1989). The photolysis rates are calculated using the Fast Tropospheric Ultraviolet and Visible (FTUV) radiation model (Tie et al., 2003; Li et al., 2005), with the aerosol and cloud effects on the photochemistry (Li et al., 2011a).

ISORROPIA (version 1.7) is used to predict the thermodynamic equilibrium between the ammonia–sulfate–nitrate–chloride–water aerosols and their gas-phase precursors of H₂SO₄–HNO₃–NH₃–HCl–water vapor (Nenes et al., 1998). It is worth noting that the most recent extension of ISORROPIA, known as ISORROPIA II, has incorporated a larger number of aerosol species (Ca, Mn, K salts) and is designed to be a superset of ISORROPIA (Fountoukis and Nenes, 2007). Considering that crustal species are not considered in the study, ISORROPIA (version 1.7) is still used to calculate inorganic components and ISORROPIA II is imperative to be incorporated into the WRF-Chem model in future studies. In addition, a parameterization of sulfate heterogeneous formation involving aerosol liquid water (ALW) has been developed and implemented into the model, which has successfully reproduced the observed rapid sulfate formation during haze days (Li et al., 2017). The sulfate heterogeneous formation from SO₂ is parameterized as a first-order irreversible uptake by ALW surfaces, with a reactive uptake coefficient of $\mathrm{0.5}\times {\mathrm{10}}^{-\mathrm{4}}$ assuming that there is enough alkalinity to maintain the high iron-catalyzed reaction rate.

The organic aerosol (OA) module is based on the volatility basis set (VBS) approach with aging and detailed information can be found in Li et al. (2011b). The POA components from traffic-related combustion and biomass burning are represented by nine surrogate species with saturation concentrations (C^∗) ranging from 10⁻² to 10⁶ µg m⁻³ at room temperature (Shrivastava et al., 2008), and assumed to be semi-volatile and photochemically reactive (Robinson et al., 2007). The secondary organic aerosol (SOA) formation from each anthropogenic or biogenic precursor is calculated using four semi-volatile organic compounds (VOCs) with effective saturation concentrations of 1, 10, 100, and 1000 µg m⁻³ at 298 K. The SOA formation via the heterogeneous reaction of glyoxal and methylglyoxal is parameterized as a first-order irreversible uptake by aerosol particles with an uptake coefficient of $\mathrm{3.7}\times {\mathrm{10}}^{-\mathrm{3}}$ (Liggio et al., 2005; Zhao et al., 2006; Volkamer et al., 2007). The OA module has reasonably reproduced the POA and SOA concentration against measurements, and detailed model performance can be found in Li et al. (2011b), Feng et al. (2016), and Xing et al. (2019).

The anthropogenic emission inventory with a horizontal resolution of 6 km is developed by Zhang et al. (2009), with the base year of 2013, including industry, transportation, power plant, residential, and agriculture sources. The Model of Emissions of Gases and Aerosols from Nature (MEGAN) is used to calculate the biogenic emissions online (Guenther et al., 2006).

Table 1 WRF-Chem model configurations.

A heavy haze episode from 10 to 27 February 2014 in BTH is simulated in association with the field observation of air pollutants and secondary inorganic aerosols. Detailed model configuration can be found in Table 1 and the simulation domain is presented in Fig. 1.

Figure 1 WRF-Chem simulation domain with topography. The red filled circles show the locations of the cities with ambient air quality monitoring sites, and the size of the circles represents the number of sites in each city. The blue filled rectangle denotes the CRAES observation site in Beijing.

2.2 Parameterization of the heterogeneous hydrolysis of N₂O₅

The reaction of N₂O₅ heterogeneous hydrolysis on the surface of deliquescent aerosols to form HNO₃ can be represented as

$$\begin{array}{}\text{(1)}& {\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}+{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\to \mathrm{2}\cdot {\mathrm{HNO}}_{\mathrm{3}}.\end{array}$$

This reaction is usually implemented into air quality transport models as a first-order loss:

$$\begin{array}{}\text{(2)}& {\displaystyle \frac{\partial \left[{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}\right]}{\partial t}}=-k_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}\cdot \left[{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}\right].\end{array}$$

[N₂O₅] represents the N₂O₅ concentration in the atmosphere. The loss rate constant, $k_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$, is parameterized in the following way:

$$\begin{array}{}\text{(3)}& k_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}={\displaystyle \frac{\mathrm{1}}{\mathrm{4}}}\cdot c_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}\cdot S\cdot {\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}},\end{array}$$

where $c_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ is the average molecular velocity of N₂O₅, and S is the available aerosol surface area density.

In this study, the parameterization of ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ follows Riemer03 and Riemer09. In the parameterization, primary emission compounds such as elemental carbon, insoluble organic matter (mostly part of POA), insoluble inorganic matter, and mineral dust particles are assumed to serve as a nucleus of aerosols. Condensation of soluble chemical components and further water vapor on the surface of the nucleus forms an aqueous layer. The nucleus and the aqueous layer are assumed to be a unified “core” (aqueous core) in the Riemer03 and Riemer09 parameterizations. In the Riemer03 parameterization, soluble inorganic components including sulfate and nitrate are taken into consideration for suppressing the N₂O₅ heterogeneous hydrolysis uptake in the aqueous core, and the parameterization of ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ is defined as

$$\begin{array}{}\text{(4)}& {\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}=f\cdot {\mathit{\gamma }}_{\mathrm{1}}+\left(\mathrm{1}-f\right)\cdot {\mathit{\gamma }}_{\mathrm{2}},\end{array}$$

with γ₁=0.02 and γ₂=0.002, and f is defined as

$$\begin{array}{}\text{(5)}& f={\displaystyle \frac{m_{{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}}}{m_{{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}}+m_{{\mathrm{NO}}_{\mathrm{3}}^{-}}}}.\end{array}$$

$m_{{\mathrm{SO}}_{\mathrm{4}}^{\mathrm{2}-}}$ and $m_{{\mathrm{NO}}_{\mathrm{3}}^{-}}$ are the aerosol mass concentrations of soluble sulfate and nitrate.

In the Riemer09 parameterization, unreactive organic layers are further considered for the suppression of N₂O₅ hydrolysis by covering the aqueous core. Organic layers may be formed by secondary organic aerosols, and such layers may consist of a single layer of molecules (monolayered coatings) or of several molecule layers (multilayered coatings) on the surface of the aqueous core. These organic layers are assumed to be an organic “coating” (shell) in the Riemer09 parameterization. The resistor scheme to calculate ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ in the Riemer09 parameterization is parameterized as follows:

$$\begin{array}{}\text{(6)}& {\displaystyle \frac{\mathrm{1}}{{\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}}}={\displaystyle \frac{\mathrm{1}}{{\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}},\mathrm{core}}}}+{\displaystyle \frac{\mathrm{1}}{{\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}},\mathrm{coat}}}},\end{array}$$

where ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}},\mathrm{core}}$ is the reaction probability of the aqueous core, which is calculated using Eq. (4), and ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}},\mathrm{coat}}$ is the pseudo-reaction probability of the organic coating calculated by the following formulation:

$$\begin{array}{}\text{(7)}& {\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}},\mathrm{coat}}={\displaystyle \frac{\mathrm{4}\cdot R\cdot T\cdot H_{\mathrm{org}}\cdot D_{\mathrm{org}}\cdot R_{\mathrm{c}}}{c_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}\cdot l\cdot R_{\mathrm{p}}}},\end{array}$$

where R is the universal gas constant, T is temperature, H_org is the Henry's law constant for N₂O₅ in the organic coating, and D_org is the diffusion coefficient for N₂O₅ in the organic coating. H_org and D_org depend on the physicochemical properties of the compounds comprising the organic coating. In the Riemer09 scheme, H_org⋅D_org is defined as $\mathrm{0.03}\cdot H_{\mathrm{aq}}\cdot D_{\mathrm{aq}}$. H_aq is the Henry's law constant of N₂O₅ for the aqueous phase (H_aq=5000 M atm⁻¹) and D_aq is the diffusion coefficient of N₂O₅ in the aqueous phase ($D_{\mathrm{aq}}={\mathrm{10}}^{-\mathrm{9}}$ m² s⁻¹). R_p, R_c, and l are the radius of the particle, radius of the inorganic core, and thickness of the coating, respectively. R_p, R_c, and l are calculated as follows:

$$\begin{array}{}\text{(8)}& {\displaystyle }{R}_{\mathrm{p}}=R_{\mathrm{c}}+l,\text{(9)}& {\displaystyle }l=R_{\mathrm{p}}\cdot \left(\mathrm{1}-{\mathit{\beta }}^{\frac{\mathrm{1}}{\mathrm{3}}}\right),\text{(10)}& {\displaystyle }\mathit{\beta }={\displaystyle \frac{V_{\mathrm{inorg}}}{V_{\mathrm{inorg}}+V_{\mathrm{org}}}},\end{array}$$

where V_inorg and V_org are the volume of inorganic and organic materials, respectively.

2.3 Statistical methods for model evaluation

In this study, the mean bias (MB), root-mean-square error (RMSE), the index of agreement (IOA), mean fractional bias (MFB), and mean fractional error (MFE) are used to evaluate the model performance in simulating air pollutants.

$$\begin{array}{}\text{(11)}& {\displaystyle }\mathrm{MB}={\displaystyle \frac{\mathrm{1}}{N}}\sum _{i=\mathrm{1}}^N\left(P_i-O_i\right)\text{(12)}& {\displaystyle }\mathrm{RMSE}={\left[{\displaystyle \frac{\mathrm{1}}{N}}\sum _{i=\mathrm{1}}^N{\left(P_i-O_i\right)}^{\mathrm{2}}\right]}^{\frac{\mathrm{1}}{\mathrm{2}}}\text{(13)}& {\displaystyle }\mathrm{IOA}=\mathrm{1}-{\displaystyle \frac{\sum _{i=\mathrm{1}}^N{\left(P_i-O_i\right)}^{\mathrm{2}}}{\sum _{i=\mathrm{1}}^N{\left(\left|P_i-\stackrel{\mathrm{‾}}{O}\right|+\left|O_i-\stackrel{\mathrm{‾}}{O}\right|\right)}^{\mathrm{2}}}}\text{(14)}& {\displaystyle }\mathrm{MFB}={\displaystyle \frac{\mathrm{1}}{N}}\sum _{i=\mathrm{1}}^N{\displaystyle \frac{P_i-O_i}{\left(P_i+O_i\right)/\mathrm{2}}}\text{(15)}& {\displaystyle }\mathrm{MFE}={\displaystyle \frac{\mathrm{1}}{N}}\sum _{i=\mathrm{1}}^N{\displaystyle \frac{\left|P_i-O_i\right|}{\left(P_i+O_i\right)/\mathrm{2}}}\end{array}$$

Here P_i and O_i are the simulated and observed variables, respectively. N is the total number of the simulations for comparisons, and $\stackrel{\mathrm{‾}}{O}$ donates the average of the observation. The IOA ranges from 0 to 1, with 1 showing a perfect agreement of the simulation with the observation.

2.4 Air pollutant observations

Simulations are compared to available meteorological and air pollutant observations to evaluate the model performance. The meteorological parameters including temperature, RH, wind speed, and direction are obtained from the website http://www.meteomanz.com (last access: 20 June 2019). The hourly observations of PM_2.5, O₃, SO₂, NO₂, and CO concentrations are released by the China National Environmental Monitoring Center and can be downloaded from the website http://106.37.208.233:20035 (last access: 20 June 2019).

Additionally, hourly organic carbon (OC) and elemental carbon (EC) concentrations are measured using a thermal–optical reflectance carbon analyzer (OCEC RT-4, Sunset Lab, USA) at the Chinese Research Academy of Environmental Sciences (CRAES, 40.04^∘ N, 116.40^∘ E) in Beijing (Wei et al., 2014; Liu et al., 2018). Hourly sulfate, nitrate, ammonium, and other inorganic ions are sampled and analyzed by ion chromatography (URG 9000S, Thermo Fisher Scientific, USA).

The OC ∕ EC ratio approach is used to derive the SOA mass concentrations from EC and OC filter measurements as follows (Strader, 1999; Cao et al., 2004):

$$\begin{array}{}\text{(16)}& {\displaystyle }\mathrm{POC}=\mathrm{EC}\times \left({\displaystyle \frac{\mathrm{POC}}{\mathrm{EC}}}\right),\text{(17)}& {\displaystyle }\mathrm{SOC}=\mathrm{OC}-\mathrm{POC},\text{(18)}& {\displaystyle }\mathrm{SOA}=\mathrm{SOC}\times \left({\displaystyle \frac{\mathrm{SOA}}{\mathrm{SOC}}}\right),\end{array}$$

where POC and SOC are the primary OC and secondary OC, respectively. In the present study, POC∕EC and SOA∕SOC are assumed to be 2.4 and 1.6, respectively, based on previous studies (Cao et al., 2007; Aiken et al., 2008; Yu et al., 2009), and detailed information about the approach can be found in Feng et al. (2016). It is worth noting that those assumed POC∕EC and SOA∕SOC could potentially affect the model–measurement comparisons.

3 Results and discussion

3.1 Synoptic conditions during the wintertime of 2014

Based on the NCEP FNL reanalysis data (https://rda.ucar.edu/datasets/ds083.2, last access: 20 June 2019), we have initially performed the analysis of synoptic conditions using the wind, temperature, relative humidity, and geopotential height fields at 500 and 850 hPa averaged from 10 to 27 February 2014 over China (Fig. 2). At 500 hPa, flat westerly winds prevail over BTH and its surrounding area, indicating stagnant atmospheric circulation conditions (Fig. 2a). Moreover, the flat isotherm distribution is similar to that of the isobar at 500 hPa, showing that there is no obvious exchange of cold and warm air masses, which together with the flat westerly leads to weak turbulent mixing in the vertical direction and stable weather conditions (Fig. 2b). At 850 hPa, the southeast coastal areas of China are controlled by the anticyclone, whose center is located over the South China Sea (Fig. 2c). In eastern China influenced by the anticyclone, the weak southerly wind prevails over the BTH and its surrounding regions, providing a favorable condition for stagnant weather conditions and the formation of air pollution. With the prevailing southerly wind, the warm and humid air flow and the polluted air mass are subject to being transported from south to north, aggravating the air pollution in BTH. In addition, high-relative-humidity conditions facilitate heterogeneous reactions for secondary aerosol formation (Fig. 2d).

Figure 2 Distributions of average winds (black flag vectors), geopotential heights (blue lines), temperature (red lines), and relative humidity (contour fill) at (a, b) 500 hPa and (c, d) 850 hPa from 10 to 27 February 2014, respectively.

3.2 Model performance

In order to quantify effects of the N₂O₅ heterogeneous hydrolysis and organic coating on the nitrate formation, three experiments have been performed in the study. In the base case, the Riemer09 parameterization is used to take into consideration the organic coating effect on the N₂O₅ heterogeneous hydrolysis by assuming that all the SOA is involved in coating (hereafter referred to as the B-case). In the first sensitivity case, the contribution of N₂O₅ heterogeneous hydrolysis to nitrate formation is not considered (hereafter referred to as the H0-case). In the second sensitivity case, the organic coating effect is not considered in the Riemer09 parameterization (hereafter referred to as the C0-case). The simulation results in the B-case are compared to observations in BTH.

3.2.1 Meteorological parameter simulations in Beijing

Considering that the meteorological conditions play a crucial role in air pollution simulations, which determine accumulation or dispersion of pollutants, verifications are first performed for the simulations of meteorological fields. Figure 3 presents the temporal profile of the simulated and observed temperature, RH, wind speed, and wind direction averaged over 12 meteorological sites in Beijing from 10 to 27 February 2014. The WRF-Chem model reproduces the temporal variation in the surface temperature during the whole episode well. The MB and RMSE are −0.2 and 1.7 ^∘ C, and the IOA reaches 0.94, indicating good agreement of the simulations with observations (Table 2). The simulated temporal RH variations are also consistent with observations, with a MB, RMSE, and IOA of 2.6 %, 10.9 %, and 0.89, respectively. In addition, the model tracks the temporal variations in the surface wind reasonably well, with IOAs of 0.73 and 0.66 for the wind speed and direction, respectively.

Figure 3 Temporal variations in simulated (red line) and observed (black dots) meteorological parameters of near-surface (a) temperature, (b) relative humidity, (c) wind speed, and (d) wind direction averaged at 12 meteorological sites in Beijing from 10 to 27 February 2014.

Table 2 Statistics for model performance.

a, b, and c represent the meteorological parameter averaged over 12 meteorological sites in Beijing, the air pollutant averaged over all ambient monitoring stations in BTH, and the aerosol component at the CRAES site in Beijing, respectively.

3.2.2 Air pollutant simulations in BTH

Figure 4 shows the relationship between observed and simulated mass concentrations of PM_2.5, O₃, SO₂, NO₂, and CO in Beijing, Tianjin, and Hebei from 10 to 27 February 2014. The correlation coefficient (R) of PM_2.5 mass concentrations between observations and simulations in Beijing, Tianjin, and Hebei is 0.83, 0.80, and 0.90, respectively, indicating a good performance of the WRF-Chem model in simulating the PM_2.5 concentration in BTH. The correlation of O₃ and NO₂ mass concentrations between observations and simulations is not as good as that of PM_2.5 concentrations in BTH, with an R between 0.6 and 0.8. Apparently, the R of SO₂ simulations with observations shows that the WRF-Chem model still has difficulties in simulating SO₂ concentrations in BTH well, particularly in Hebei. Except uncertainties from SO₂ emissions, such as source intensities and distributions and diurnal profiles, the bias of simulated wind fields also substantially influences the SO₂ simulation (Bei et al., 2017). In particular, SO₂ is principally emitted by the point source, including the power plants and agglomerated industrial zones, so the SO₂ simulations are more sensitive to the wind field simulation uncertainties. In terms of R, the SO₂ simulations in Beijing and Tianjin are better than those in Hebei, indicating that the SO₂ emissions in Beijing and Tianjin are generally determined by area sources, i.e., the residential living, but the point source dominates the SO₂ concentration in Hebei. Considering the long lifetime of CO in the atmosphere, the CO simulation is decided by its emission and the meteorological fields. The R of CO simulations with observations in BTH ranges from around 0.6 to 0.7, showing that the CO emissions used in the study and simulated meteorological fields are generally reasonable.

Figure 4 Relationships between observed and simulated mass concentrations of PM_2.5, O₃, NO₂, SO₂, and CO in Beijing, Tianjin, and Hebei from 10 to 27 February 2014. The red line is the linear regression between observations and simulations, and the black dashed line presents the 1 : 1 line.

Figure 5 Comparison of observed (black dots) and simulated (red line) diurnal profiles of near-surface hourly (a) PM_2.5, (b) O₃, (c) NO₂, (d) SO₂, and (e) CO averaged over all ambient monitoring stations in BTH from 10 to 27 February 2014.

Figure 5 presents the diurnal profiles of simulated and observed PM_2.5, O₃, NO₂, SO₂, and CO mass concentrations averaged over all ambient monitoring stations in BTH during the simulated episode. The WRF-Chem model reproduces the diurnal variations in the PM_2.5 mass concentrations against observations in BTH well. The MB and RMSE are −6.3 and 27.6 µg m⁻³, respectively, and the IOA is 0.96. The model generally replicates the haze-developing stage well, but fails to capture the observed spikes of PM_2.5 mass concentrations, which might be caused by the uncertainty of the simulated meteorological fields or irregular air pollutant emissions (Bei et al., 2017). The simulated O₃ diurnal variations are in good agreement with observations, with a MB and IOA of 1.4 µg m⁻³ and 0.91, respectively. The model tracks the observed diurnal variations in NO₂ mass concentrations with an IOA of 0.92 well, but it slightly overestimates NO₂ concentrations compared to observations with a MB of 6.6 µg m⁻³. However, during nighttime, the model overestimation is considerable, which is perhaps due to the model biases in modeling nighttime planetary boundary layer (PBL). Although the model reasonably yields the variation trend of the observed SO₂ concentration, with an IOA of 0.85, the dispersion of the simulated SO₂ concentration is rather large, with a RMSE of 27.8 µg m⁻³. In addition, the model overestimates the SO₂ concentration compared to observations, and the MB is 7.6 µg m⁻³, which might be mainly caused by the emission inventory that has undergone noticeable changes since implementation of emission control strategies in 2013 in BTH. The model performs well in simulating CO diurnal variations against observations, with a MB and IOA of 0.2 µg m⁻³ and 0.90, respectively.

Figure 6 Spatial distributions of average (a) PM_2.5, (b) O₃, (c) NO₂, and (d) SO₂ mass concentrations from 10 to 27 February 2014. Colored dots, colored contour, and black arrows are observations, simulations, and simulated surface winds, respectively.

Figure 6 presents the distributions of simulated and observed near-surface mass concentrations of PM_2.5, O₃, NO₂, and SO₂ along with the predicted wind fields averaged during the episode. Generally, the simulated wind in BTH is weak during the episode and the southerly wind prevails, corresponding to the synoptic situation at 850 and 500 hPa well, which is favorable for the accumulation of air pollutants. The observed PM_2.5 concentrations are more than 115 µg m⁻³ on average, showing that BTH suffers from heavy haze pollution (Fig. 6a). The model generally reproduces the spatial distribution of PM_2.5 concentrations against observations well, with the PM_2.5 concentration exceeding 150 µg m⁻³ in the plain area of BTH. The simulated and observed O₃ mass concentrations are less than 50 µg m⁻³ in the plain area of BTH, and in several megacities, including Beijing, Tianjin, Baoding, and Shijiazhuang, the O₃ concentrations are less than 30 µg m⁻³ (Fig. 6b). The low O₃ concentrations during the episode are generally caused by the weak insolation during wintertime, which is unfavorable for photochemical reactions, and the titration due to high NO_x emissions in BTH (Fig. 6c). The simulated NO₂ concentrations are generally more than 40 µg m⁻³, consistent with the observations at monitoring sites in BTH. The simulated and observed SO₂ mass concentrations in cities or their surrounding areas are still rather high, exceeding 50 µg m⁻³ (Fig. 6d). Elevated SO₂ concentrations in BTH during wintertime are to some degree contributed by the residential coal combustion (Li et al., 2018). High levels of NO₂ and SO₂ show that stringent emission mitigation strategies still need to be implemented in BTH.

3.2.3 Sulfate, ammonium, and SOA simulations in Beijing

The SOA and sulfate concentration directly influences the N₂O₅ heterogeneous hydrolysis in the Riemer09 parameterization, and the ammonium aerosol concentration substantially affects the nitrate aerosol formation. Therefore, Fig. 7 presents the temporal profiles of observed and calculated SOA, sulfate, and ammonium mass concentrations at the CRAES site in Beijing from 10 to 27 February 2014. The model reasonably tracks the diurnal variation in the SOA concentration compared to observations, with a MB and IOA of −1.2 µg m⁻³ and 0.83, respectively. The observed SOA concentration exhibits rather large fluctuations, which are not reproduced well by the model. The simulated sulfate trend is generally in agreement with observations with an IOA of 0.88, but there are considerable model biases. During the first pollution event, the model reasonably reproduces the sulfate increase during the haze-developing stage, but the early falloff of sulfate concentrations during the dissipation stage causes a substantial underestimation. However, during the second pollution event, the model considerably overestimates the sulfate concentration against the measurement from 22 to 26 February 2014. The ammonium simulation is slightly better than that of sulfate, with an IOA of 0.90.

Figure 7 Comparison of observed (black dots) and simulated (red line) diurnal profiles of hourly (a) SOA, (b) sulfate, and (c) ammonium concentrations at the CRAES site in Beijing from 10 to 27 February 2014.

Furthermore, the MFB and MFE between simulations and observations are also calculated to evaluate the model performance in simulating meteorological parameters and air pollutants (Table 2). Boylan and Russell (2006) proposed that MFB should be within ±60 % and MFE should be below 75 % for a satisfactory model performance. For the simulation in the B-case, MFB values are within 27 % and MFE values are below 55 %, indicating that the model performance is satisfactory.

In summary, the WRF-Chem model performs reasonably well in simulating air pollutants and aerosol species, showing that the simulated meteorological fields and emissions used in the study are generally reasonable.

3.3 Contributions of the N₂O₅ heterogeneous hydrolysis and organic coating to the nitrate formation

Figure 8 provides the nitrate temporal variations in the three cases against observations at CRAES site in Beijing from 10 to 27 February 2014. When the N₂O₅ heterogeneous hydrolysis is not considered in the H0-case, although the model tracks the observed nitrate variations well with an IOA of 0.91, it considerably underestimates nitrate concentration against the measurement, with a MB of −17.0 µg m⁻³. When the N₂O₅ heterogeneous hydrolysis is taken into consideration based on the Riemer09 parameterization without organic coating in the C0-case, the nitrate simulation is improved compared to that in the H0-case, with an IOA of 0.95. However, the model begins to overestimate the nitrate concentration compared to the measurement, with a MB of 5.4 µg m⁻³. In the B-case, when all the SOA is assumed to be involved in coating to suppress the N₂O₅ heterogeneous uptake on surfaces of deliquescent aerosols, the model performs best in simulating the nitrate variation compared to the measurement, with a MB and IOA of 0.1 µg m⁻³ and 0.96, respectively. The remarkable consistency of the simulated nitrate in the B-case with the measurement indicates that the organic coating plays an important role in improving the nitrate simulation. It is worth noting that the MB for nitrate aerosols at the CRAES site in the B-case is close to zero, but the RMSE is still rather large, reaching 19.0 µg m⁻³, showing considerable underestimation and overestimation, caused by uncertainties of meteorological fields and emissions. For example, the model overestimates nitrate concentrations on 11, 13, and 14 February and underestimates on 24 February against measurements. In addition, the early occurrence of intensified winds in the morning on 16 February in simulations causes rapid falloff of nitrate concentrations, leading to substantial model biases.

Figure 8 Temporal variations in observed (black dot) and simulated (green line: H0-case; blue line: C0-case; red line: B-case) nitrate concentrations at the CRAES site in Beijing from 10 to 27 February 2014.

Figure 9 Spatial distributions of average nitrate contributions of (a) the N₂O₅ heterogeneous hydrolysis and (b) organic coating in BTH from 10 to 27 February 2014.

Figure 9a presents the distribution of contributions of the N₂O₅ heterogeneous hydrolysis to the nitrate formation averaged during the episode by differentiating simulations in the B-case and H0-case. The contribution of the N₂O₅ heterogeneous hydrolysis to the nitrate formation is substantial in BTH, exceeding 15 µg m⁻³ in the plain area. Although the O₃ concentration is fairly low in BTH during the episode (Fig. 6b), particularly during nighttime (Fig. 5b), the elevated NO₂ level still facilitates the N₂O₅ formation to warrant occurrence of the N₂O₅ heterogeneous hydrolysis. Previous studies have revealed that N₂O₅ heterogeneous hydrolysis is vital in nitrate formation. For example, Wang et al. (2017) calculated the daily average nitrate formation potential from the N₂O₅ heterogeneous hydrolysis in Beijing, showing that the reaction accounts for 52 % of the total nitrate formation. Su et al. (2017) investigated the contribution of N₂O₅ heterogeneous hydrolysis to the nitrate formation in Beijing during autumn in 2015 and found that the reaction causes a 21.0 % enhancement of nitrate concentrations. In the present study, the nitrate contribution of the N₂O₅ heterogeneous hydrolysis is 29.4 % in Beijing during the episode on average, which is close to the result in Su et al. (2017) but much lower than that in Wang et al. (2017). The average nitrate contribution of the reaction in BTH is about 30.1 %, showing that the reaction constitutes an important nitrate source during the haze pollution episode. Additionally, the N₂O₅ heterogeneous hydrolysis contributes 11.6 % of the PM_2.5 concentration on average, playing a considerable role in the haze formation in BTH.

However, it is worth noting that the brute force method (BFM) is used to quantify the contribution of the N₂O₅ heterogeneous hydrolysis to the nitrate formation (Dunker et al., 1996). The BFM is generally used to assess the importance of some source, but it lacks consideration of interactions of the complicated physical and chemical processes in the atmosphere (Zhang and Ying, 2011). Therefore, in the study, the contribution of the N₂O₅ heterogeneous hydrolysis to the nitrate formation might be underestimated, considering the competition of inorganic cations from HNO₃ formed through gas-phase reactions and sulfate aerosols in the atmosphere. It is imperative to use the source-oriented base module to evaluate the nitrate contribution of the reaction.

Figure 10 Temporal variation in the simulated ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ in the B-case in Beijing from 10 to 27 February 2014.

Figure 9b shows the distribution of the average decrease in nitrate concentrations due to suppression of organic coating during the episode by differentiating simulations in the B-case and C0-case. The organic coating reduces the nitrate concentration by more than 5 µg m⁻³ in the plain area of BTH, and on average the decrease in nitrate aerosols is 4.7 µg m⁻³ or 8.4 % in BTH during the episode. Riemer et al. (2009) have shown that when the nitrate levels are high (above 15 µg m⁻³), the organic coating decreases nitrate concentrations by 10 %–15 % over Europe. However, Chen et al. (2018) have demonstrated that the suppression of organic coating is negligible over western and central Europe, with an influence on nitrate concentrations of less than 2 % on average and 20 % at the most significant moment. Apparently, except for N₂O₅ and water-soluble OA in the atmosphere, the effect of organic coating is also dependent on NH₃, RH, and temperature. Hence, the inconsistency between the model results about the organic coating effect can be attributed to the variation in simulation conditions. For example, in order to obtain substantial effects of the organic coating, N₂O₅, SOA, and NH₃ need to be present when RH is high and temperature is low. However, those conditions are rarely fulfilled simultaneously over western and central Europe, causing a negligible effect of organic coating (Chen et al., 2018). Additionally, Wang et al. (2017) have indicated that the evaluated nitrate level with the N₂O₅ heterogeneous hydrolysis in Beijing is much higher than the observation, which they have attributed to atmospheric dilution and deposition. It is worth noting that the organic coating effect might constitute one of the most possible reasons for the overestimation of nitrate concentrations, considering the elevated SOA level in Beijing, which suppresses the N₂O₅ heterogeneous hydrolysis and results in high observed N₂O₅ concentrations (Wu et al., 2017). Figure 10 presents the temporal variation in the simulated ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ in Beijing during the episode. The simulated ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ fluctuates between 0.009 and 0.02 when organic coating is included, with an average of 0.013. The estimated ${\mathit{\gamma }}_{{\mathrm{N}}_{\mathrm{2}}{\mathrm{O}}_{\mathrm{5}}}$ in Beijing by Wang et al. (2017) ranges from 0.025 to 0.072 without consideration of the suppression of organic coating, indicating that organic coating substantially hinders the N₂O₅ heterogeneous hydrolysis, likely causing the observed high level of N₂O₅ during nighttime.

It is worth noting that, in the study, the assumption of metastable aerosols is used or the water-soluble aerosol is assumed to be only in a liquid state in simulations. However, Wang et al. (2008) have highlighted the effect of the hysteresis of particle-phase transitions on the distribution of solid and aqueous aerosols. The aerosol phase is generally regulated by the hysteresis loop. Atmospheric particles containing inorganic salts remain solid until the RH reaches the DRH (deliquescence relative humidity). At the DRH, the solid particle spontaneously absorbs water to become a saturated aqueous solution. However, the liquid particle does not crystallize when the RH is below the DRH (Seinfeld and Pandis, 2006). Therefore, another possible pathway exists to suppress the N₂O₅ hydrolysis; i.e., the inorganic particles might be in solid phase without organic coating. Further studies need to be conducted to evaluate the hysteresis effect on the N₂O₅ hydrolysis and organic coating.

3.4 Studies of organic aerosol hygroscopicity sensitivity to nitrate formation

In Sect. 3.3, the WRF-Chem model considerably improves nitrate simulations when considering the N₂O₅ heterogeneous hydrolysis and organic coating effects. Organic aerosols (OAs) are broadly classified as primary OA (POA) directly emitted and SOA formed in the atmosphere, some of which are water soluble. In order to explore the effects of different OA coatings on the nitrate formation, an additional four sensitivity studies are conducted, in which half of SOA (C1-case), all SOA (C2-case), all SOA and half of POA (C3-case), and all SOA and POA (C4-case) are involved in coating.

Figure 11 shows the Taylor diagram (Taylor, 2001) to present the variance, bias, and correlation of the observed and simulated nitrate concentrations in the four sensitivity cases at the CRAES site during the episode. In the C1-case, when half of SOA is considered to be involved in coating, the simulated nitrate concentration is the most consistent with the observation, with a correlation coefficient of 0.96. In the C2-case with all SOA assumed to be engaged in coating, the correlation coefficient decreases to 0.95. The normalized standardized deviation (NSD) is 1.02 for the C1-case and C2-case, showing the model overestimation in these two cases. With half of POA involved in coating, the NSD is very close to 1.0 (0.99), indicating the simulated nitrate concentration in the C3-case is almost the same as the observation on average, but the correlation coefficient of 0.94 is less than those in the C1-case and C2-case. When all of OA is assumed as the coating, the bias between simulated and observed nitrate concentrations is the largest, and the effect of POA on suppressing the N₂O₅ heterogeneous hydrolysis might be overestimated in the C4-case.

Figure 11 Taylor diagram (Taylor, 2001) to present the variance, bias, and correlation of the observed and simulated nitrate concentrations at the CRAES site in Beijing from 10 to 27 February 2014.

Sensitivity results show that the effects of different organic compounds on suppressing the N₂O₅ heterogeneous hydrolysis to form nitrate vary, depending on the content of water-soluble OA. Laboratory and field measurements have revealed that OA becomes progressively oxidized and more hygroscopic during the aging process in the atmosphere (Jimenez et al., 2009). The OA hygroscopicity constitutes a necessary prerequisite for accurately representing the organic coating effect on the N₂O₅ heterogeneous hydrolysis. According to the simulations in the present study, in BTH, not all SOA can serve as the coating to suppress the nitrate formation, and the effect of POA on coating might be neglected. Xing et al. (2019) have shown that in BTH, the heterogeneous SOA formed by irreversible uptake of glyoxal and methylglyoxal on wet aerosol surfaces contributes about 30 % of the SOA mass during haze days. Considering the possible heterogeneous SOA contribution of other carbonyl compounds and the atmospheric aging of OA, about half of SOA should likely be hygroscopic and involved in coating.

4 Conclusion

Nitrate aerosol has constituted a main component of PM_2.5 with implementation of aggressive emission control strategies since 2013 in BTH. In the study, the Riemer09 parameterization is implemented into the WRF-Chem model to simulate the nitrate formation from the N₂O₅ heterogeneous hydrolysis referred to as the most important pathway of the nitrate formation at nighttime. A heavy haze episode from 10 to 27 February 2014 in BTH is simulated using the WRF-Chem model to verify the effect of organic coating on the N₂O₅ heterogeneous hydrolysis and its consequent contribution to nitrate formation. Analyses of synoptic fields show a stagnant weather condition with the prevailing southerly wind in the low-level atmosphere in BTH and surrounding areas during the episode, facilitating accumulation of air pollutants and heavy haze formation.

The WRF-Chem model performs reasonably in predicting the temporal variations in the meteorological parameters compared to observations in Beijing. The model generally reproduces the temporal variations in and spatial distributions of air pollutants well against observations at monitoring sites in BTH. In addition, the simulated diurnal profiles of sulfate, ammonium, and SOA are also in good agreement with the measurements at the CRAES site in Beijing.

The Riemer09 parameterization with all the SOA assumed to be involved in coating considerably improves the nitrate simulations compared to the measurements at the CRAES site in Beijing. When organic coating is not considered in the Riemer09 parameterization, the model overestimates the nitrate concentration against the measurements. On average, organic coating decreases nitrate concentrations by 4.7 µg m⁻³ or 8.4 % in BTH during the episode. Furthermore, the N₂O₅ heterogeneous hydrolysis with organic coating contributes about 30.1 % of nitrate concentrations and 11.6 % of the PM_2.5 concentration in BTH, playing a considerable role in the haze formation.

Sensitivity studies reveal that the OA hygroscopicity is a necessary prerequisite for accurately evaluating the organic coating effect on the N₂O₅ heterogeneous hydrolysis. In the present study, POA might not serve as coating and about half of SOA should be involved in coating to suppress the nitrate formation. Future studies still need to be conducted to further predict the OA hygroscopicity, in order to more precisely represent the organic coating effect on the N₂O₅ heterogeneous hydrolysis in chemical transport models.