Influence of biomass burning vapor wall loss correction on modeling organic aerosols in Europe by CAMx v6.50

Increasing evidence from experimental studies suggests that the losses of semi-volatile vapors to chamber walls could be responsible for the underestimation of organic aerosol (OA) in air quality models that use parameters obtained from chamber experiments. In this study, a box model with a volatility basis set (VBS) scheme was developed, and the secondary organic aerosol (SOA) yields with vapor wall loss correction were optimized by a genetic algorithm based on advanced chamber experimental data for biomass burning. The vapor wall loss correction increases the SOA yields by a factor of 1.9–4.9 and leads to better agreement with measured OA for 14 chamber experiments under different temperatures and emission loads. To investigate the influence of vapor wall loss correction on regional OA simulations, the optimized parameterizations (SOA yields, emissions of intermediate-volatility organic compounds from biomass burning, and enthalpy of vaporization) were implemented in the regional air quality model CAMx (Comprehensive Air Quality Model with extensions). The model results from the VBS schemes with standard (VBS_BASE) and vapor-wall-loss-corrected parameters (VBS_WLS), as well as the traditional two-product approach, were compared and evaluated by OA measurements from five Aerodyne aerosol chemical speciation monitor (ACSM) or aerosol mass spectrometer (AMS) stations in the winter of 2011. An additional reference scenario, VBS_noWLS, was also developed using the same parameterization as VBS_WLS except for the SOA yields, which were optimized by assuming there is no vapor wall loss. The VBS_WLS generally shows the best performance for predicting OA among all OA schemes and reduces the mean fractional bias from −72.9 % (VBS_BASE) to −1.6 % for the winter OA. In Europe, the VBS_WLS produces the highest domain average OA in winter (2.3 µg m⁻³), which is 106.6 % and 26.2 % higher than VBS_BASE and VBS_noWLS, respectively. Compared to VBS_noWLS, VBS_WLS leads to an increase in SOA by up to ∼80 % (in the Balkans). VBS_WLS also leads to better agreement between the modeled SOA fraction in OA (fSOA) and the estimated values in the literature. The substantial influence of vapor wall loss correction on modeled OA in Europe highlights the importance of further improvements in parameterizations based on laboratory studies for a wider range of chamber conditions and field observations with higher spatial and temporal coverage.

1 Introduction

Organic aerosol (OA) accounts for a substantial fraction of atmospheric particulate matter (Jimenez et al., 2009), which is closely associated with human health impacts and climate change (Cohen et al., 2017; Kanakidou et al., 2005; Lelieveld et al., 2015). Organic aerosol originates from a variety of natural and anthropogenic sources (Hallquist et al., 2009), among which residential biomass burning has been recognized as the dominant source for both primary (POA) and secondary (SOA) organic aerosols in Europe during wintertime (Butt et al., 2016; Jiang et al., 2019b; Qi et al., 2019). Despite its substantial contribution to OA, biomass burning OA is largely underestimated by chemical transport models (CTMs) (Ciarelli et al., 2017a; Hallquist et al., 2009; Robinson et al., 2007; Theodoritsi and Pandis, 2019; Woody et al., 2016).

Many efforts have been devoted to understanding and diminishing the gap between modeled and observed OA from biomass burning. One of the major reasons for underestimated OA is the absence of semi-volatile organic compounds (SVOCs) from residential biomass burning in the current emission inventories (Denier van der Gon et al., 2015). A smog chamber study showed that the precursors traditionally included in CTMs account for only ∼3 %–27 % of the observed SOA from residential biomass burning (Bruns et al., 2016). In order to compensate for the effects from missing precursors, various modeling studies have treated POA as semi-volatile and adopted different scaling approaches to calculate SVOC and IVOC emissions. The most commonly used method is to increase the POA emissions by a factor of 3 (Ciarelli et al., 2017a; Fountoukis et al., 2014; Jiang et al., 2019b; Tsimpidi et al., 2010), while recent studies have also developed new profiles based on nonmethane organic compounds (NMOCs) (Cai et al., 2019; Lu et al., 2018). However, an increasing number of laboratory experimental studies have found that SVOC and IVOC emissions have high variability depending on different burning conditions and fuel types (Hatch et al., 2015, 2017; Jen et al., 2019; Koss et al., 2018; Sekimoto et al., 2018), and the estimation of SVOCs and IVOCs in modeling studies remains to be improved. Meanwhile, increasing evidence from chamber experiments has demonstrated that losses of semi-volatile vapors to chamber walls could lead to a substantial underestimation of OA (Akherati et al., 2020; Bertrand et al., 2018; Bian et al., 2015; Krechmer et al., 2016; Loza et al., 2010; Matsunaga and Ziemann, 2010; Zhang et al., 2014). Unlike particle wall losses – which are routinely corrected in chamber studies – the effects of vapor wall losses are rarely investigated and considered in modeling practices.

Zhang et al. (2014) reported that vapor wall losses may lead to an underestimation of SOA by a factor of 1.1–4.2, depending on different NO_x conditions. This factor has been adopted by several CTM studies to scale up the yields of SOA. For instance, Baker et al. (2015) tested the sensitivity of CMAQ to vapor wall loss by increasing the yields of semi-volatile gases by a factor of 4, with the traditional two-product approach used for OA simulations. This factor was also implemented in a box model with a volatility basis set (VBS) scheme (Hayes et al., 2015), which distributed organic species into logarithmically spaced volatility bins, and was shown to improve the model performance for predicting SOA (Donahue et al., 2011, 2006; Hodzic et al., 2010; Robinson et al., 2007). Nevertheless, recent studies have shown that vapor wall losses lead to even larger variability in SOA yields according to different chamber conditions and precursor species (Akherati et al., 2019; Cappa et al., 2016). On the other hand, some studies took vapor wall loss corrections into account in CTMs using the SOA yields generated from the statistical oxidation model (SOM). Cappa et al. (2016) and Akherati et al. (2019) used the traditional two-product model to fit vapor-wall-loss-corrected SOA yields and applied the yields in the regional CTM UCD/CIT (University of California at Davis/California Institute of Technology). They reported, however, that the two-product fits might not be sufficiently robust. Furthermore, Hodzic et al. (2016) used a vapor-wall-loss-corrected VBS parameterization in the global model GEOS-Chem based on chamber experiments conducted on individual precursors, which are highly dependent on the experimental conditions. Each of these studies clearly called for a better assessment of the uncertainties across the entire range of precursor compounds as well as under different chamber conditions.

Here, we (1) developed a VBS-based box model and fit the vapor-wall-loss-corrected SOA yields of biomass burning IVOCs based on 14 chamber experiments under different temperature and emission loads, (2) implemented the vapor-wall-loss-corrected VBS parameters in the regional chemical transport model Comprehensive Air Quality Model with extensions (CAMx), and (3) investigated the role of vapor wall loss correction in model performance by comparing modeled organic aerosols from traditional and modified VBS OA schemes with ambient observations at multiple European sites. Biomass burning in this study refers to residential biomass burning, while wildfires and prescribed burning are not included.

2 Parameterization method

2.1 Chamber experimental data

The parameterization of the VBS scheme was based on experimental data from two smog chamber campaigns in 2014–2015. It includes 14 experiments conducted under various temperature conditions (−10, 2, 15 ^∘C) and covered a wide range of emission loads (from 19 to 284 µg m⁻³). Emissions were generated by combustion of beechwood in three different woodstoves, including conventional and modern burners, manufactured in 2002–2010. Beechwood is selected as it is one of the major forest types in Europe, and beechwood is widely used for residential heating and cooking in Europe. Although different biomass fuel types may largely affect the emitted organic gas species and affect SOA formation, a recent study showed that the effect of biomass fuel type on SOA formation is much smaller than the effects of initial organic mass (OM) load and hydroxyl (OH) radical exposure (Lim et al., 2019). The organic gases covering 86 intermediate-volatility and semi-volatile organic compounds (IVOCs and SVOCs), which are SOA precursors, were measured by a proton-transfer-reaction time-of-flight mass spectrometer (PTR-ToF-MS). The PTR-ToF-MS was operated under standard conditions in H₃O⁺ mode, as introduced in Stefenelli et al. (2019). A common set of 263 ions was extracted from the measurements, and among these ions, 86 showed clear decay with time and were identified as potential SOA precursors. These are listed in Table S1 of Stefenelli et al. (2019). The aerosol evolution was monitored by a high-resolution time-of-flight aerosol mass spectrometer (HR-ToF-AMS). The particle wall loss has been corrected as described in Stefenelli et al. (2019). The conditions of each chamber experiment are shown in Table S1. A more detailed description of the experiments can be found in Stefenelli et al. (2019), Bertrand et al. (2017), and Bruns et al. (2016).

2.2 VBS box model

A VBS box model was developed to simulate the formation and evolution of primary and secondary OA in the chamber. In the model, we assumed that the condensable gases generated from oxidation of the precursors could (1) partition to the particle phase, (2) be lost on the chamber wall, and/or (3) be diluted by other gases injected into the smog chamber. CAMx includes four types of precursors from anthropogenic sources, i.e., toluene, xylene, benzene, and IVOCs, which includes all the other unspeciated organic gases. According to our measurements, the traditional anthropogenic precursors toluene, xylene, and benzene only account for ∼15 % of the total organic gases. To facilitate the implementation of the optimized parameters in CAMx, all the measured SOA precursors including the traditional ones were lumped into one surrogate as IVOCs with the same reaction rate and volatility distribution. In comparison, Stefenelli et al. (2019) assigned the same set of compounds to six different classes according to their properties (reaction rates, expected SOA yields, etc.) as well as their origins and occurrence in the emissions. These included furans and methoxy-phenols from the pyrolysis of cellulose and lignin, respectively, single-ring and polyaromatic hydrocarbons from flaming combustion, and oxygenated non-aromatic compounds with fewer and more than six carbon atoms. The current lumping approach of all these species into one surrogate, despite variations in their properties, is more adapted for implementation into CAMx and for assessing vapor wall losses with additional parameters included in the box model. The organic compounds were distributed into six logarithmically spaced volatility bins corresponding to saturation concentrations of 10⁻¹, 10⁰, 10¹, 10², 10³, and 10⁴ µg m⁻³. The change in the organic gas concentration (C) for a constituent within the volatility bin i (C_i) can be described by Eq. (1), where P is the production of organic gas (OG) in the chamber due to oxidization of precursors, k_cs is the condensation sink (s⁻¹) describing the speed of condensable gases condensing on existing aerosol particles, k_w is the rate constant of vapor lost to the wall, k_dil is the dilution rate, and $C_{{\mathrm{eq}}_{i,\mathrm{p}}}$ and $C_{{\mathrm{eq}}_{i,\mathrm{w}}}$ represent the gas-phase equilibrium concentrations to the aerosol particles and chamber wall surface, respectively:

$$\begin{array}{}\text{(1)}& \begin{array}{rl}{\displaystyle \frac{\mathrm{d}{C}_i}{\mathrm{d}t}}& =P\cdot {\mathit{\zeta }}_i-k_{\mathrm{cs}}\left(C_i-C_{{\mathrm{eq}}_{i,\mathrm{p}}}\right)\\ & -k_{\mathrm{w}}\left(C_i-C_{{\mathrm{eq}}_{i,\mathrm{w}}}\right)-k_{\mathrm{dil}}\cdot C_i.\end{array}\end{array}$$

The production rates of oxidized organic gases (P) are used as inputs of the box model. They are determined by the consumption rates of precursors measured by PTR taking into account their dilution. ζ_i represents the mass fraction of primary and oxidation products in a volatility bin i. ζ_i of POA from biomass burning is obtained from May et al. (2013), with values of 0.2, 0.1, 0.1, 0.2, 0.1, and 0.3 for compounds in each volatility bin. ζ_i of oxidation products is assumed to follow a kernel normal distribution as a function of log C^∗, $\mathit{\zeta }\sim N(\mathit{\mu },{\mathit{\sigma }}^{\mathrm{2}}$), where μ is the median value of log C^∗ and σ is the standard deviation, which will be optimized as described in Sect. 2.3. The assumption of a normal distribution could ensure positive ζ_i values, allow constraining the total mass fraction of the certain surrogate to equal 1, and reduce the model's degree of freedom significantly, as reported in Stefenelli et al. (2019). The time series of k_dil is obtained from Stefenelli et al. (2019). The k_cs of each experiment is obtained from Bertrand et al. (2018). The k_w varies significantly depending on the chamber conditions such as the chamber size and relative humidity. Zhang et al. (2014) reported k_w values of $\mathrm{2.5}\times {\mathrm{10}}^{-\mathrm{4}}$ s⁻¹ and $\mathrm{1}\times {\mathrm{10}}^{-\mathrm{4}}$ s⁻¹ for toluene and other VOCs, respectively, while it is much higher in recent studies, such as $\mathrm{1.2}\times {\mathrm{10}}^{-\mathrm{3}}$ to $\mathrm{2.4}\times {\mathrm{10}}^{-\mathrm{3}}$ s⁻¹ in Krechmer et al. (2016), $\mathrm{1.28}\times {\mathrm{10}}^{-\mathrm{3}}$ s⁻¹ in Akherati et al. (2020), and $\sim \mathrm{1}\times {\mathrm{10}}^{-\mathrm{3}}$ to $\mathrm{3.3}\times {\mathrm{10}}^{-\mathrm{3}}$ s⁻¹ in Bertrand et al. (2018). To cover the wide range of vapor wall loss, we tested three k_w values of 0.0020, 0.0033, and 0.0040 s⁻¹ based on the conditions of our chamber. A base case was also developed by assuming there is no vapor wall loss in the chamber (k_w=0). The condensation of a species in the particle phase (C_p) can then be described by Eq. (2).

$$\begin{array}{}\text{(2)}& {\displaystyle \frac{\mathrm{d}{C}_{i,\mathrm{p}}}{\mathrm{d}t}}=k_{\mathrm{cs}}\left(C_i-C_{{\mathrm{eq}}_{i,\mathrm{p}}}\right)-k_{\mathrm{dil}}\cdot C_{i,\mathrm{p}}\end{array}$$

Following the partitioning model of Pankow (1994), the gas-phase concentrations at equilibrium with respect to the particle phase ($C_{{\mathrm{eq}}_{i,\mathrm{p}}}$) and to the chamber wall ($C_{{\mathrm{eq}}_{i,\mathrm{w}}}$) are determined by their partitioning coefficients ξ_i and ξ_i,w (Donahue et al., 2009), as shown in Eqs. (3) and (4).

$$\begin{array}{}\text{(3)}& {\displaystyle }{C}_{{\mathrm{eq}}_{i,\mathrm{p}}}=(C_{i,\mathrm{g}}+C_{i,\mathrm{p}})\cdot [\mathrm{1}-{\mathit{\xi }}_i],{\mathit{\xi }}_i={\left(\mathrm{1}+{\displaystyle \frac{C_i^{\ast }}{C_{\mathrm{OA}}}}\right)}^{-\mathrm{1}}\text{(4)}& {\displaystyle }\begin{array}{rl}{C}_{{\mathrm{eq}}_{i,\mathrm{w}}}& =(C_{i,\mathrm{g}}+C_{i,\mathrm{w}})\cdot [\mathrm{1}-{\mathit{\xi }}_{i,\mathrm{w}}],\\ {\mathit{\xi }}_{i,\mathrm{w}}& ={\left(\mathrm{1}+{\displaystyle \frac{C_i^{\ast }}{C_{\mathrm{wall}}}}\right)}^{-\mathrm{1}}\end{array}\end{array}$$

Here, C^∗ represents the saturation concentration, C_OA is the wall-loss-corrected OA concentration measured by the AMS, and C_wall is the equivalent organic mass concentration at the wall determined in Bertrand et al. (2018). The Clausius–Clapeyron equation (Eq. 5) was applied to take into account the effects of temperature on C^∗:

$$\begin{array}{}\text{(5)}& C^{\ast }=C_{T\mathrm{0}}^{\ast }\cdot {\displaystyle \frac{T_{\mathrm{0}}}{T}}\cdot \exp \left({\displaystyle \frac{\mathrm{\Delta }{H}_{\mathrm{vap}}/\mathrm{8.314}}{\mathrm{1}/T_{\mathrm{0}}-\mathrm{1}/T}}\right),\end{array}$$

where $C_{T\mathrm{0}}^{\ast }$ is the mass saturation concentration under the reference temperature (T₀). T is the temperature of each experiment, while T₀ equals 298 K. ΔH_vap (J) is the enthalpy of vaporization at the reference temperature, and 8.314 is the universal gas constant (J mol⁻¹ K⁻¹). $\mathrm{\Delta }{H}_{\mathrm{vap}}=\left\{\mathrm{70}\mathrm{000}-\mathrm{11}\mathrm{000}\times \log C^{\ast }\right\}$ is adopted for the primary set (May et al., 2013), while ΔH_vap of the oxidized products is determined during model optimization. The C_wall was determined in previous studies to be of the order of a few milligrams per cubic meter (mg m⁻³; Bertrand et al., 2018). In this study, we run the box model for three different C_wall values (1, 5, 25 mg m⁻³) with a reference temperature of 2 ^∘C (275.15 K) according to Bertrand et al. (2018).

2.3 Model optimization

The model is optimized to constrain the volatility distribution (as a function of log C^∗, $\mathit{\zeta }\sim N(\mathit{\mu },{\mathit{\sigma }}^{\mathrm{2}})$) and ΔH_vap of the oxidized products. A genetic algorithm (GA) is used to find the best-fit parameters leading to the lowest average root mean square error (RMSE) and mean bias (MB) between modeled and measured OA concentrations for all 14 experiments. The genetic algorithm is a metaheuristic algorithm inspired by the natural selection process to generate optimized solutions (Mitchell, 1996). It begins by creating an initial population of individual solutions (20 different combinations of μ, σ, and ΔH_vap here) within certain upper and lower bounds, called parents. The performance of each solution is evaluated by a fitness function, which is the sum of RMSE and MB between the modeled and measured OA concentrations of 14 experiments in this study. A new generation of solutions is then formed either by making random changes to a single parent (called mutation) or by combining the vector entries of a pair of parents (called crossover). The process will be repeated until reaching the stopping conditions, which involves either the number of iterations reaching 50 or the stall generations (generation with no significant change in the fitness function) reaching 20. The GA is conducted using the genetic algorithm solver in the MATLAB R2019a Global Optimization Toolbox (The MathWorks, Inc).

3 Modeling approach

3.1 Regional chemical transport model CAMx

The regional model CAMx version 6.50 (Ramboll, 2018) was used to model organic aerosol in Europe (15^∘ W–35^∘ E, 35–70^∘ N) for the whole year of 2011, with a horizontal resolution of 0.25^∘ × 0.125^∘ and 14 terrain-following vertical layers from ∼20 m above the ground up to 460 hPa. The Carbon Bond 6 Revision 2 (CB6r2) gas-phase mechanism (Hildebrandt Ruiz and Yarwood, 2013) was selected. The gas–aerosol partitioning of inorganic aerosols was simulated by the ISORROPIA thermodynamic model (Nenes et al., 1998). For organic aerosols, several OA schemes, including both the traditional two-product approach (SOA chemistry partitioning scheme, SOAP) and the VBS scheme with different parameterizations, were applied (see Sect. 3.2).

The meteorological parameters were prepared with the Weather Research and Forecasting model (WRF version 3.7.1; Skamarock et al., 2008) based on the 6 h European Centre for Medium-Range Weather Forecasts (ECMWF) reanalysis global data (Dee et al., 2011). The meteorological parameters were evaluated and reported in a previous study (Jiang et al., 2019a), which showed that most of the meteorological parameters met the criteria for meteorological model performance by Emery (2001). The initial and boundary conditions were obtained from the global model MOZART-4/GEOS-5 (Horowitz et al., 2003). Inputs of ozone column densities were produced based on the Total Ozone Mapping Spectrometer (TOMS) data by the National Aeronautics and Space Administration (NASA, https://acd-ext.gsfc.nasa.gov/anonftp/toms/omi/data/Level3e/ozone/, last access: 17 March 2021), and the photolysis rates were then calculated by the Tropospheric Ultraviolet and Visible (TUV) Radiation Model version 4.8 (NCAR, 2011). The source-specific anthropogenic emissions were based on the European emission inventory TNO-MACC-III (Monitoring Atmospheric Composition and Climate) (Kuenen et al., 2014). The biogenic emissions (isoprene, monoterpenes, sesquiterpenes, soil NO) were simulated by the PSI model developed at the Laboratory of Atmospheric Chemistry at the Paul Scherrer Institute (Andreani-Aksoyoglu and Keller, 1995; Jiang et al., 2019a; Oderbolz et al., 2013). More details about the model inputs can be found in our previous studies performed using the same input data (Jiang et al., 2019a, b).

Table 1 Description of the different OA schemes.

^a IVOB is the abbreviation for “IVOCs from biomass burning” in CAMx.
^b The yield values correspond to volatility bins with saturation concentrations of 10⁻¹, 10⁰, 10¹, 10², and 10³ µg m⁻³.
^c SOAP does not separate IVOCs from biomass burning and other anthropogenic sectors, and therefore it is not comparable with the SOA yields for IVOBs.

3.2 Parameterization of OA schemes

To investigate the effects of vapor-wall-loss-corrected yields and to compare to other modifications and/or parameterizations that are currently strongly debated in the community, five simulations with different OA schemes were conducted in this study (Table 1). Besides VBS_WLS, which uses the optimized parameterization with vapor wall loss correction for the biomass burning sector, SOAP and VBS_BASE represent the two standard parameterizations in CAMx; VBS_3POA represents a common approach to offset the missing SVOC emissions in recent modeling studies without vapor wall loss, and VBS_noWLS is another reference case without vapor wall loss, which uses exactly the same parameters as VBS_WLS except for the SOA yields from IVOCs. Details about each OA scheme are introduced below.

• SOAP. The SOAP (SOA chemistry partitioning) module is a semi-volatile equilibrium scheme based on the traditional two-product approach. POA emissions are assumed to be inert in SOAP. The updated parameterization of SOAP2.1 in CAMx v6.50 used the aerosol yield data that correct for vapor wall losses in smog chamber experiments based on Zhang et al. (2014).

• VBS_BASE. The VBS_BASE used the standard VBS parameterization in CAMx v6.50. The IVOC emissions from different sources were calculated based on the literature. IVOCs from gasoline and diesel vehicles were calculated as 25 % and 20 % of NMVOC emissions from gasoline and diesel vehicles, respectively (Jathar et al., 2014). IVOC emissions from residential biomass burning were estimated as 4.5 times POA emissions based on Ciarelli et al. (2017a). The IVOC emissions from other anthropogenic sources were calculated as 1.5 times POA as proposed by Robinson et al. (2007).

• VBS_3POA. An increasing number of experimental and modeling studies have reported a considerable contribution of semi-volatile organic compounds (SVOCs) to SOA formation (Bruns et al., 2016; Ciarelli et al., 2017b; Denier van der Gon et al., 2015; Hatch et al., 2017; Woody et al., 2015), while SVOCs are absent in the current emission inventories. Despite considerable variability of SVOC emissions from biomass burning according to recent studies, the VBS_3POA is supposed to be a reference case representing the commonly used approach without vapor wall loss, and therefore we adopted the routine approach of multiplying the POA emissions by a factor of 3 to offset the influence of missing SVOC emissions. This approach has been widely used in modeling studies (Ciarelli et al., 2016, 2017a; Shrivastava et al., 2011; Tsimpidi et al., 2010). All the other parameters were kept the same as the standard VBS parameterization in CAMx v6.50. The VBS_BASE IVOC emissions were adopted here.

• VBS_WLS. The VBS_WLS used the optimized parameters by the VBS box model, including the emissions and vapor-wall-loss-corrected yields for IVOCs from residential biomass burning, and the ΔH_vap of the oxidized products. The modified parameters for volatility-bin-specific yields and ΔH_vap of the oxidized products from IVOCs can be found at https://doi.org/10.5281/zenodo.3998342 (Jiang, 2020). The optimized mass yields in the box model were converted to molar yields using the default molecular weights in CAMx (Table 1). Both the optimized and default molar yields have a sum larger than 1 as the VBS scheme accounts for both oxygenation and fragmentation (Koo et al., 2014). The reaction rate with OH (k_OH) was calculated based on the measurements following Stefenelli et al. (2019). Based on the chamber measurements, the IVOC emissions from residential biomass burning are ∼13.7 times the primary OM load (Fig. S1), among which the traditional precursors in CAMx from biomass burning (toluene, xylene, and benzene) account for ∼15 % of the total emissions. To avoid double counting of these traditional precursors, which are already included in the emission inventory, we applied a factor of 12 to calculate the IVOC emissions from biomass burning. The IVOC emissions from other sources were estimated by using the same approach as in VBS_3POA.

• VBS_noWLS. The VBS_noWLS was designed as a reference for VBS_WLS, which adopted the same parameters as VBS_WLS except for the yields. The VBS_noWLS used the fitted yields from the box model assuming that there is no vapor wall loss (k_w=0).

3.3 Model evaluation

The general model performance for the major air pollutants (SO₂, NO₂, O₃, PM_2.5) was reported in our previous study (Jiang et al., 2019b), which was comparable to other modeling studies in Europe. OA measurements and source apportionment studies using positive matrix factorization (PMF) analysis from five Aerodyne aerosol chemical speciation monitor (ACSM) or aerosol mass spectrometer (AMS) stations in the winter of 2011 were used to evaluate modeled primary and secondary organic aerosol by different OA schemes: Zurich (Canonaco et al., 2013), Marseille (Bozzetti et al., 2017), the SIRTA (Site Instrumental de Recherche par Télédétection Atmosphérique) facility located in the Paris region (Zhang et al., 2019), and Bologna and San Pietro Capofiume (SPC) (Paglione et al., 2020). For Zurich and SIRTA, only data collected from late autumn to early spring (January, February, March, November, and December) – when emissions from biomass burning are relatively high – were used for the statistical analysis, although the observations covered longer time periods. The spatial distribution and observation periods of each station are shown in Fig. S2. Statistical metrics, including mean bias (MB), mean error (ME), root mean square error (RMSE), mean fractional bias (MFB), and mean fractional error (MFE), between modeled and observed primary and secondary OA were calculated.

Figure 1 Comparison between measured and modeled OA with an optimized parameterization under k_w=0.0033 s⁻¹ and C_wall=5 mg m⁻³ (a, b). Relation between the end-point wall loss factor R_wall of each experiment and initial OM loads under different temperature (c). The gray dashed lines in (a) represent modeled OA with the same parameterization but k_w=0.

4 Results and discussion

4.1 Modeled and measured OA from chamber experiments

The optimized parameters were then applied to the box model to simulate OA production for 14 chamber experiments. Figure 1 shows the comparison between measured OA and modeled primary and secondary OA under the median chamber conditions (k_w=0.0033 s⁻¹, C_wall=5 mg m⁻³) for each experiment. The model reproduces the process of OA formation for most of the experiments well, except for experiments no. 9 and no. 14, which have relatively lower OM loads (26 and 48 µg m⁻³ for Exp9 and Exp14, respectively). This can be partially explained by the different weighting impact for experiments with high or low OM loads. The experiments with higher OM loads normally have larger MB and RMSE at the beginning of optimization and therefore have a higher impact during the model optimization. A direct consequence is that the optimized parameters would work better for experiments with higher OM loads. However, the model performance in each experiment could also be influenced by a series of other factors such as temperature and chamber conditions. While the model simulation without vapor wall loss correction largely overestimates OA at the initial time point and underestimates the final OA (Fig. S3), the agreement between the modeled and measured trends improved when the vapor wall loss was taken into account. The mean bias (MB) and root mean square error (RMSE) between the modeled and observed OA in 14 experiments are 6.7 and 42.2 µg m⁻³ for the case under median chamber conditions, which are 48 % and 12 % lower than in the case without vapor loss correction (MB $=-\mathrm{12.8}$, 47.8 µg m⁻³). To investigate the role of vapor wall loss in modeled OA, another set of simulations was performed in which we used the same optimized parameterization under the median chamber conditions but set k_w=0. In these cases, the modeled OA concentrations (dashed line in Fig. 1) were based on the assumption that there is no vapor wall loss. The wall loss ratio R_wall, which is defined as the ratio between the modeled OA concentration without (k_w=0) and with (k_w=0.0033 s⁻¹, C_wall=5 mg m⁻³) vapor wall loss, was calculated for the end point of each experiment (Fig. 1c). The R_wall values varied from 1.5 (Exp2) to 3.2 (Exp11) among the 14 experiments and showed a clear dependence on the initial OA loads.

Figure 2 Dependence of the wall loss factor $R_{\mathrm{wall}}(C_{\mathrm{OA},k_{\mathrm{w}}=\mathrm{0}}/C_{\mathrm{OA},\mathrm{optimal}{}_{k_{\mathrm{w}}}}$) on initial organic mass load, temperature, and condensation sink.

To further understand the factors influencing R_wall, we conducted a series of model simulations with and without vapor wall loss under different initial organic mass loads, temperature, and condensation sink inputs (Fig. 2). Higher k_w and C_wall lead to higher R_wall values for all the cases, and different chamber conditions (k_w, C_w) could result in a different R_wall by a factor of 1.2–1.6, depending on different temperature, OM loads, and condensation sinks. The R_wall values generally decrease with increasing initial OM loads, which is consistent with the fact that R_wall values for Exp8–14 are higher than Exp1–7. The increased R_wall with the increasing temperature explains why Exp10 through Exp14 (T=15 ^∘C) have higher R_wall than Exp8 and Exp9 ($T=-\mathrm{10}$ ^∘C), while they have similar OM load levels. The condensation sink is inversely correlated with R_wall, indicating that the higher the rate of condensable gases condensing on the existing particles, the lower the vapor loss to the chamber wall, and therefore the lower the effect of vapor wall loss on modeled OA.

Figure 3 Optimized yield factors (a) and the mass yield of SOA from biomass burning at 298 K (b) with and without vapor wall loss correction. The blue bars (a) and line (b) with vapor wall loss refer to median chamber conditions with k_w=0.0033 s⁻¹ and C_wall=5 mg m⁻³.

The optimized volatility distribution for the secondary condensable gases from biomass burning (ppm per ppm IVOC) based on different wall loss assumptions (k_w>0 or k_w=0) is displayed in Fig. 3a. The optimized yields considering vapor wall loss lead to a 3.3 times higher mass in the low-volatility bins ($\log C^{\ast }\le \mathrm{0}$) compared to that assuming k_w=0, indicating significant effects of vapor wall loss correction on predicting SOA production. To give a more direct view of the effects of vapor wall loss on the SOA yield, we integrated the mass of SOA for all the volatility bins at 298 K (Fig. 3b). The mass yield under the median chamber conditions for vapor wall loss (k_w=0.0033 s⁻¹, C_wall=5 mg m⁻³) is higher than the base case without considering the vapor wall loss about by factors of 4.9 (when C_OA= 0.1 µg m⁻³) to 1.9 (when C_OA=1000 µg m⁻³). The influence of vapor wall loss on mass yield decreases with decreasing temperature. At 0 ^∘C, the mass yield with vapor wall loss correction is higher than the base case by factors of 4.3 (when C_OA=0.1 µg m⁻³) to 1.7 (when C_OA=1000 µg m⁻³).

Table 2 Statistical results for model performance in simulating OA, SOA, and POA. The number of daily average observations from five ACSM/AMS stations is 216.

Figure 4 Concentrations of measured and modeled OA, POA, and SOA at five ACSM or AMS stations in winter (a) and mean bias for different OA schemes (b). The lines inside boxes represent median values, and the yellow triangles represent mean values.

4.2 Performance of CAMx with different OA schemes

The modeled OA concentrations with different OA schemes were compared with measurements from five ACSM/AMS stations in winter. The statistical results are shown in Table 2, and the distributions of OA concentrations and the mean bias between modeled and measured primary and secondary OA are displayed in Fig. 4. OA is underestimated overall with all OA schemes. The VBS schemes lead to a better model performance than the two-product approach SOAP, except for VBS_BASE with the default VBS parameterization. These results are consistent with a previous study using CAMx (Meroni et al., 2017), in which the better performance of SOAP compared to the default VBS was reported as a result of error compensation. The improved performance of modified VBS (3POA, noWLS, WLS) for OA mainly comes from the contribution of SOA (Table 2). The modeled SOA by 3POA and noWLS is very similar, and therefore the analysis below will focus on the comparison between noWLS and WLS, for which the only difference is that WLS uses vapor-wall-loss-corrected yields for IVOCs from biomass burning, while noWLS uses the fitted yields assuming no vapor wall loss (k_w=0). WLS reduces the MFB between the modeled and measured SOA from 52.5 % in noWLS to 20.0 %. WLS shows a better average MB than noWLS; however, it also increases the upper whisker of the MB (Fig. 4b), largely affected by overestimated SOA in Bologna and SPC.

Figure 5 Measured and modeled daily average OA using different OA schemes in winter. ZRH: Zurich, BLQ: Bologna, MRS: Marseille, SIRTA: Paris SIRTA, SPC: San Pietro Capofiume.

Limited by the availability of OA measurements, the effects of vapor wall loss correction on model performance present a clear site dependence in this study. The modeled and measured daily average OA concentrations at each site are shown in Fig. 5. The temporal variations of primary and secondary OA at these sites can be found in Fig. S4. VBS_WLS leads to the best performance for both OA and SOA in Marseille and SIRTA, in spite of an overall underestimation (Fig. S4b, c). In Zurich, the vapor-wall-loss-corrected yields for biomass burning improve the model performance in February and March, while there is an overestimation of OA and SOA for all the OA schemes in November (Fig. S4a). The largest contribution to OA during this period was found to be from biogenic SOA, which was potentially overestimated due to overestimated temperatures during the same time period (Jiang et al., 2019b). Bologna and SPC are located in the Po Valley where biomass burning contributes most to winter OA (Jiang et al., 2019b), and therefore higher effects from vapor wall loss correction on SOA are observed compared to other sites. At SPC, fog scavenging processes played an important role in OA during the measurements (Gilardoni et al., 2014); however, the meteorological model failed to reproduce the fog events due to the coarse resolution in this study (Jiang et al., 2019b). Consequently, both VBS_WLS and noWLS lead to an overestimation of OA and SOA, while SOAP and VBS_BASE show better performance, probably due to compensation for errors (Fig. S4e). In Bologna, a significant overestimation of temperature was found on 2 to 6 December (Jiang et al., 2019b), leading to a significant underestimation of SOA for all the OA schemes (Fig. S4d). Excluding this period, the modeled SOA by VBS_WLS is 89 % higher than the measurements, while the modeled SOA concentrations by the other schemes are closer to the measurements, with relative differences of −64 % for SOAP, −10 % for VBS_BASE, and 4 % for VBS_noWLS.

The distinct performance of vapor-wall-loss-corrected VBS at different sites could arise from various factors. It might come from the high uncertainties of SVOC and IVOC emissions from biomass burning, which were estimated by the same factor for the whole domain but were reported to have substantial inter-country variations (Denier van der Gon et al., 2015). Missing formation and removal processes such as photolytic and heterogeneous oxidation in the model could also result in different model performance for specific sites. In addition, in spite of the advanced chamber measurements we used to optimize the yield parameters covering a wide range of precursor species and multiple temperature and chamber conditions, the fitted vapor-wall-loss-corrected parameterization is still highly uncertain. To achieve a more robust parameterization and to further improve the model performance for OA, more studies on SVOC and IVOC emissions, as well as the formation and removal mechanisms of SOA based on extensive laboratory studies and field observations with higher spatial and temporal coverage, are needed.

Figure 6 Modeled OA, SOA, and POA in winter (DJF, December–January–February) by different OA schemes.

4.3 Effects of vapor wall loss correction on modeled OA in Europe

4.3.1 OA

The modeled OA results in Europe for the whole year of 2011 with different OA schemes were compared to investigate the effects of OA schemes and the vapor wall loss correction. Among all the sources, residential biomass burning contributed 16.3 %–52.6 % POA and 5.9 %–28.9 % SOA in winter (Jiang et al., 2019b), indicating the potential roles of vapor wall loss for the biomass burning sector. Figure 6 shows the modeled OA, SOA, and POA in winter (December–January–February). VBS_WLS leads to the highest domain average OA (2.3 µg m⁻³), which is 16.4 %, 26.2 %, 38.7 %, and 106.6 % higher than VBS_3POA, VBS_noWLS, SOAP, and VBS_BASE, respectively. The VBS schemes generally produce higher OA than SOAP, except for the default parameterization (VBS_BASE) in which the lack of SVOC emissions is not considered. However, SOAP leads to the second highest SOA after VBS_WLS, especially in northern Europe where monoterpene emissions from coniferous forests are relatively high. This is mostly because of the high terpene SOA yields in SOAP2.1, which were reduced in the later version of the CAMx model (CAMx v7.0, http://www.camx.com, last access: 17 March 2021). The vapor-wall-loss-corrected yields lead to increased SOA in large areas of central and southern Europe (Fig. 7). The largest difference is predicted for the Po Valley and Romania regions with a high residential biomass burning impact (Fig. S5). The overall relative differences between VBS_WLS and VBS_noWLS are more than 80 %, and the highest grid-scale increment reaches 5.6 µg m⁻³ in the region of the Balkans. The modeled POA concentrations are similar to those in the VBS case with correction for SVOC (3POA, noWLS, WLS), with domain average concentrations ranging from 0.9 (noWLS) to 1.1 (3POA) µg m⁻³, and therefore no significant effects were observed from vapor wall loss correction (Fig. 6). The POA simulated by VBS_BASE (0.3 µg m⁻³) is even lower than SOAP (0.7 µg m⁻³), as POA is semi-volatile and could evaporate and react with oxidants to form secondary products in VBS, while SOAP assumes POA to be inert.

Figure 7 Differences in modeled OA, SOA, and POA in winter (DJF, December–January–February) by VBS schemes with (VBS_WLS) and without (VBS_noWLS) vapor wall corrections.

The effects of different VBS schemes on OA are much smaller in summer (Fig. S6). Despite a slight increase from the VBS_BASE (1.2 µg m⁻³), the modeled OA by the three modified VBS schemes is quite similar (1.4–1.5 µg m⁻³). The effects of vapor-wall-loss-corrected yields for biomass burning emissions are negligible due to low emissions in summer (Fig. S7). SOAP produced the highest OA (2.1 µg m⁻³) in summer due to the high SOA yields from monoterpenes as explained before.

Figure 8 Modeled fractions of annual mean SOA to total OA (fSOA) using different OA schemes. Modeled results for VBS_3POA are very similar to VBS_noWLS and are therefore not shown here.

4.3.2 Fraction of SOA in OA

The effects of the updated VBS schemes on the fraction of annual average SOA in total OA (fSOA = SOA $/$ OA) are shown in Fig. 8. The VBS schemes lead to a higher fSOA (domain average 71.4 %–87.3 %) compared to SOAP (domain average 69.9 %) in most of the domain except for northern Europe, where SOAP produces high biogenic SOA. The increased POA emissions to offset the missing SVOC emissions (3POA, noWLS, WLS) decrease the fSOA compared to the default VBS parameterization (BASE), while the vapor wall loss correction yields (WLS) result in ∼5.8 % higher fSOA than noWLS for the domain average, and the largest grid-scale increase reaches 43.4 % in the Balkans. The absolute differences between fSOA for WLS and noWLS are relatively higher in rural areas than urban areas, where fSOA is lower due to high primary emissions.

Figure 9 Comparison between modeled and measured fSOA from the literature over the year (see data and sources in Table S2). The shading indicates the confidence intervals of the regression lines.

The modeled fSOA values were compared with measurements from previous studies in Europe (Crippa et al., 2014; Jiang et al., 2019b). The measured fSOA from the literature covered 18 sites and different seasons between 2008 and 2011 (Table S2). SOAP tends to underestimate the fSOA, while VBS_BASE significantly overpredicts the fSOA (Fig. 9). Both WLS and noWLS tend to underestimate the high fSOA and overestimate the low fSOA. VBS_WLS has 5 % higher fSOA than VBS_noWLS and shows the highest agreement on the range of fSOA with the measurements and the average fSOA values (measured: 69.6 %; VBS_WLS: 69.1 %). The largest improvements occur in winter when the vapor-wall-loss-corrected yields of biomass burning emissions largely increase the SOA production.

5 Conclusions

In this study, we optimized the SOA yields for a VBS-based box model using 14 chamber experiments with biomass burning and implemented the fitted VBS parameters (SOA yields, IVOC emissions from biomass burning, and enthalpy of vaporization) in the regional air quality model CAMx v6.5. The influence of the vapor wall loss correction on the model performance was investigated by comparing modeled primary and secondary OA with traditional and modified OA schemes, including the two-product approach (SOAP), the standard VBS (VBS_BASE), VBS with 3 times the POA to compensate for the missing SVOCs (VBS_3POA), VBS with vapor wall loss correction (VBS_WLS), and an additional reference scenario with the same parameterizations as in VBS_WLS except for using the default SOA yields from biomass burning IVOCs (VBS_noWLS).

The vapor wall loss correction increases the mass distributed in the low-volatility bins ($\log C^{\ast }\le \mathrm{0}$) by a factor of 4.3 and increases the SOA yields by a factor of 1.9–4.9 (at 298 K). Comparison of the modeled results with different OA schemes to field measurements from five ACSM/AMS stations in Europe in winter suggests that VBS_WLS generally has the best performance to predict OA, which lowers the highest mean fractional bias from −72.9 % (VBS_BASE) to −1.6 % for OA and from −77.8 % (SOAP) to 20.0 % for SOA. In Europe, the VBS_WLS produces the highest domain average OA in winter (2.3 µg m⁻³), which is 106.6 % and 26.2 % higher than VBS_BASE and VBS_noWLS, respectively. The largest influence of vapor wall loss correction was predicted in Romania where the VBS_WLS increases SOA by ∼80 % compared to VBS_noWLS due to high emissions from residential biomass burning. VBS_WLS also leads to the highest agreement with measurements for the SOA fraction in OA (fSOA) from the literature.

The optimized parameterization with vapor wall loss correction in this study is expected to provide some insight to improve SOA underestimation in CTMs. Despite the overall improvement of model performance for predicting SOA, the VBS_WLS was found to increase the mean bias at specific sites compared to noWLS. To achieve a more robust parameterization and to further improve the model performance, complementary studies on SVOC and IVOC emissions, as well as on the formation and removal mechanisms of SOA based on extensive laboratory studies and field observations with higher spatial and temporal coverage, are still needed.