Secondary organic aerosol formation from smoldering and flaming combustion of biomass: a box model parametrization based on volatility basis set

Residential wood combustion remains one of the most important sources of primary organic aerosols (POA) and secondary organic aerosol (SOA) precursors during winter. The overwhelming majority of these precursors have not been traditionally considered in regional models, and only recently were lignin pyrolysis products and polycyclic aromatics identified as the principal SOA precursors from flaming wood combustion. The SOA yields of these components in the complex matrix of biomass smoke remain unknown and may not be inferred from smog chamber data based on single-compound systems. Here, we studied the ageing of emissions from flaming and smoldering-dominated wood fires in three different residential stoves, across a wide range of ageing temperatures (−10, 2 and 15 ^∘C) and emission loads. Organic gases (OGs) acting as SOA precursors were monitored by a proton-transfer-reaction time-of-flight mass spectrometer (PTR-ToF-MS), while the evolution of the aerosol properties during ageing in the smog chamber was monitored by a high-resolution time-of-flight aerosol mass spectrometer (HR-ToF-AMS). We developed a novel box model based on the volatility basis set (VBS) to determine the volatility distributions of the oxidation products from different precursor classes found in the emissions, grouped according to their emission pathways and SOA production rates. We show for the first time that SOA yields in complex emissions are consistent with those reported in literature from single-compound systems. We identify the main SOA precursors in both flaming and smoldering wood combustion emissions at different temperatures. While single-ring and polycyclic aromatics are significant precursors in flaming emissions, furans generated from cellulose pyrolysis appear to be important for SOA production in the case of smoldering fires. This is especially the case at high loads and low temperatures, given the higher volatility of furan oxidation products predicted by the model. We show that the oxidation products of oxygenated aromatics from lignin pyrolysis are expected to dominate SOA formation, independent of the combustion or ageing conditions, and therefore can be used as promising markers to trace ageing of biomass smoke in the field. The model framework developed herein may be generalizable for other complex emission sources, allowing determination of the contributions of different precursor classes to SOA, at a level of complexity suitable for implementation in regional air quality models.

1 Introduction

Atmospheric aerosols impact visibility, human health, and climate on a global scale (Stocker et al., 2013; World Health Organization, 2013). A thorough understanding of their chemical composition, sources, and processes is a fundamental prerequisite to develop appropriate mitigation policies. Laboratory experiments using smog chambers enable the detailed examination of the gas-phase composition and ageing of different emissions such as biomass smoke (e.g. Bruns et al., 2016; Bian et al., 2017), car exhaust (Gordon et al., 2014a, b; Platt et al., 2017; Gentner et al., 2017; Pieber et al., 2018), aircraft exhaust (Miracolo et al., 2011; Kılıç et al., 2018) or cooking emissions (Klein et al., 2016). Results from these studies consistently show that the measured concentrations of secondary organic aerosol (SOA), formed upon oxidation and partitioning of the oxidized vapours, greatly exceed estimated concentrations based on the oxidation of volatile organic compounds (VOCs) traditionally assumed to be the dominant SOA precursors (Jathar et al., 2012). The SOA formed from these chemically speciated VOCs is defined as traditional SOA (T-SOA) and is explicitly accounted for in chemical transport models (CTMs). However, Robinson et al. (2007) suggested that a significant fraction of the unexplained SOA is due to the oxidation of lower-volatility organics, i.e. semi-volatile and intermediate-volatility organic compounds (SVOC and IVOC, respectively), collectively referred to as non-traditional SOA (NT-SOA) precursors (Donahue et al., 2009).

In spite of its importance, incorporating NT-SOA into current organic aerosol (OA) models remains challenging without the identification and the quantification of the most important precursors (Jathar et al., 2012). For simplification purposes several methods based on the volatility of the emissions and a volatility-based oxidation mechanism have been developed. Currently the volatility basis set (VBS) scheme is considered to be the most suitable approach to simulate the ageing processes of non-speciated organic vapours (Donahue et al., 2006). The VBS scheme represents OA as a discrete volatility-resolved mass distribution. Reactions are described by the transfer of OA mass between volatility bins, thereby accounting for the contribution of non-traditional vapours to SOA formation without the need to incorporate explicit chemical mechanisms. Robinson et al. (2007) proposed that SVOCs, IVOCs and their products react with hydroxyl radicals (OH) to form products that are an order of magnitude lower in volatility than their precursors. Pye and Seinfeld (2010) proposed a single-step mechanism for the non-speciated SVOCs, where the products of oxidation were 2 orders of magnitude lower in volatility than the precursors. They used SOA yield (defined as SOA mass formed divided by reacted precursor mass) data for naphthalene as a surrogate for all non-speciated IVOCs, even though these are thought to be mainly branched and cyclic alkanes (Robinson et al., 2007, 2010; Schauer et al., 1999). Both methods have been implemented in plume models as well as regional and global chemical transport models and have reduced discrepancies between measured and predicted SOA concentrations and properties (Shrivastava et al., 2008; Dzepina et al., 2009; Pye and Seinfeld, 2010; Jathar et al., 2011). However, considerable uncertainties remain in the relative contributions of non-traditional precursors to different emissions, their ability to form SOA and their reaction rate constants (Jathar et al., 2014a). Limitations in SOA modelling are also a direct consequence of limitations in measurements: namely undetected or unidentified precursors and a limited number of studies available investigating the influence of different parameters such as temperature, emission load and combustion regimes. For instance, the overwhelming majority of smog chamber studies have been conducted under summertime conditions (20–30 ^∘C), preventing the assessment of temperature effects on both SOA-producing reactions and the partitioning thermodynamics (Jathar et al., 2013).

Similar limitations apply to the consideration of emissions in models. Biomass combustion is a major source of gas- and particle-phase air pollution on urban, regional and global scales (Grieshop et al., 2009; Lanz et al., 2010; Crippa et al., 2013; Gobiet et al., 2014; Chen et al., 2017; Bozzetti et al., 2017). Globally, approximately 3 billion people burn biomass or coal for residential heating and cooking (World Energy Council, 2016), often using old and highly polluting appliances. Emissions from these devices are highly variable depending on fuel type and fuel moisture (McDonald et al., 2000; Schauer et al., 2001; Fine et al., 2002; Pettersson et al., 2011; Eriksson et al., 2014; Reda et al., 2015; Bertrand et al., 2017) and typically include a complex mixture of non-methane organic gases (NMOGs), primary organic aerosol (POA) and black carbon (BC). Once emitted into the atmosphere, organic compounds can react with oxidants such as OH radicals, ozone (O₃) and nitrate radicals (NO₃). These reactions remain poorly understood, which greatly hinders the quantification of wood combustion SOA in ambient air. Bruns et al. (2016) investigated the SOA formation from residential log wood combustion from a single type of stove under stable flaming conditions only. They reported that T-SOA precursors, included in models account for only 3 % to 27 % of the measured SOA whereas 84 % to 116 % was from NT-SOA precursors including in total 22 individual compounds and two lumped compound classes, mainly consisting of polycyclic aromatic hydrocarbons from incomplete combustion (e.g. naphthalene) and cellulose and lignin pyrolysis products (e.g. furans and phenols, respectively). The estimated SOA concentrations were based on the literature SOA yields of single precursors, obtained from smog chamber experiments, and a good agreement was observed between predicted and measured SOA. However, the method suffers from two drawbacks. First, the dependence of the yields on the organic aerosol loading and temperature was not considered. Second, although the relative contributions of different precursors to SOA were estimated, thermodynamic parameters for chemical transport models (CTMs) were not determined. Based on the same experiments, the lumped concentrations of the 22 non-traditional volatile organic compounds and two compound classes were constrained in a box model (Ciarelli et al., 2017a). Improved parameters were retrieved describing the volatility distributions and the production rates of oxidation products from the overall mixture of precursors present in biomass smoke. While this method is well suited for CTMs (Pandis et al., 2013; Ciarelli et al., 2017b), it does not provide any information about the contributions of the different chemical classes to the aerosol. Similar limitations are associated with the study of other emissions, e.g. fossil fuel combustion or evaporation (Jathar et al., 2013, 2014b). The development of models capable of simulating the contribution of the different chemical species to the aerosol at different conditions is especially important in the light of the current development of highly time resolved chemical ionization mass spectrometry, capable of quantifying these products. To realize the full potential of the data acquired by this instrumentation, a modelling framework capable of predicting the production rates and the partitioning between the gas and the particle phases of the oxidation products from complex emissions is required.

Here, we extend the past analysis investigating the most recent smog chamber data of residential wood combustion based on 14 experiments performed in 2014–2015 under various conditions. Different experimental temperatures of the smog chamber were investigated, namely −10, 2 and 15 ^∘C. Three different stove types were tested, including conventional and modern residential burners. Different emission load and different hydroxyl (OH) radical exposure were examined. Moreover distinct combustion regimes were sampled across the different experiments for the first time to investigate the secondary organic aerosol chemical composition and yields from flaming and smoldering emissions. Integrated VBS-based model and novel parameterization methods based on a genetic algorithm (GA) approach were developed to predict the contribution of the oxidation products of different chemical classes present in complex emissions and to better explain the SOA formation process, providing useful information to regional air quality models. Overall, this study presents a general framework which can be adapted to assess SOA closure for complex emissions from different sources.

2 Methods

2.1 Smog chamber set-up and procedure

Two smog chamber campaigns were conducted to investigate SOA production from multiple domestic wood combustion appliances as a function of combustion phase, initial fuel load and OH exposure. These experiments were previously described in detail (Bruns et al., 2016; Ciarelli et al., 2017a; Bertrand et al., 2017, 2018a) and are summarized here. Experiments from Bertrand et al. (2017, 2018a) will be referred to as Set 1 and experiments from Bruns et al. (2016) and Ciarelli et al. (2017a) as Set 2.

The emissions were generated by three different log wood stoves for residential wood combustion: stove 1 manufactured before 2002 (Cheminées Gaudin Ecochauff 625), stove 2 fabricated in 2010 (Invicta Remilly) and stove 3 (Avant, 2009, Attika). For each stove three to four replicate experiments were performed with a loading of 2–3 kg of beech wood having a total moisture content ranging between 2 % and 19 %. The fire was ignited with three starters made of wood wax, wood shavings, paraffin and five types of natural resin. The starting phase was not studied. In total, 14 experiments were performed, consisting of two experiments at −10 ^∘C, seven experiments at 2 ^∘C and five experiments at 15 ^∘C. These experiments cover the typical range of European winter temperatures and are summarized in Table 1.

Table 1 Experimental parameters for the 14 smog chamber experiments used in this study, including smog chamber temperature, stove type, modified combustion efficiency (MCE), and the initial concentrations of POA and OGs.

NA: not available.

Ward and Hardy (1991) define the flaming and smoldering conditions according to the modified combustion efficiency, $\text{MCE}={\mathrm{CO}}_{\mathrm{2}}/\left(\mathrm{CO}+{\mathrm{CO}}_{\mathrm{2}}\right)$. Specifically, MCE > 0.9 is identified as the flaming condition, while MCE < 0.85 is identified as the smoldering condition. MCE values for the different experiments are reported in Table 1. According to this parameter, Set 1 and Set 2 experiments were dominated by smoldering and flaming, respectively. Practically, we achieved the different burning conditions by varying the amount of air in the stoves, therefore changing the combustion temperature. For Set 1, closing the air window decreased the flame temperature, resulting in a transition from a flaming to a smoldering fire. This could be visibly identified, together with the development of a thick white smoke from the chimney. We note that this conduct is very common in residential stoves, to keep the fire running for longer. Meanwhile, for Set 2, after lighting the fire, we kept a high air input to maintain a flaming fire. At the same time, we monitored the MCE in real time and only injected the emissions into the chamber when the MCE increased above 0.9.

The experiments were performed in a flexible Teflon bag of nominally 7 but typically about 5.5 m³ equipped with UV lamps (40 lights, 90–100 W, Cleo Performance, Philips, wavelength λ<400 nm) enabling photo-oxidation of the emissions (Platt et al., 2013; Bruns et al., 2015). The chamber is located inside a temperature-controlled housing. Relative humidity was maintained at 50 % and three different temperatures were investigated. Emissions from the stoves were sampled from the chimney into the chamber through heated (140 ^∘C) stainless-steel lines to reduce the loss of semi-volatile compounds. An ejector dilutor was installed (Dekati Ltd, DI-1000) to dilute emissions in the chamber by a factor of 10 before sampling. The sample injection lasted for approximately 30 min for each experiment, and was followed by an injection of 1 µL d9-butanol (98 %, Cambridge Isotope Laboratories), which was used to estimate the OH exposure as described in Sect. 3.1.1 (Barmet et al., 2012). The chamber was then allowed to equilibrate for 30 min to ensure stabilization and homogeneity and to fully characterize the primary emissions before ageing. OH radicals were produced by UV irradiation of nitrous acid (HONO) injected after chamber equilibration, generated as described in Taira and Kanda (1990) by reaction of diluted sulfuric acid (H₂SO₄) and sodium nitrate (NaNO₂) in a gas flask followed by transportation into the chamber with a carrier gas with a flow rate of 1 L min⁻¹. The smog chamber was then irradiated with UV lights for approximately 4 h to simulate atmospheric ageing.

Before and after each experiment, the smog chamber was cleaned with humidified pure air (100 % RH) and O₃ (1000 ppm) under irradiation with UV lights for at least 1 h, followed by flushing with pure dry air for at least 10 h. The background particle- and gas-phase concentrations were then measured in the clean chamber.

The total amount and composition of the emissions depend on the oxygen supply, temperature, fuel elemental composition and combustion conditions, which can be broadly classified as flaming or smoldering (Koppmann et al., 2005; Sekimoto et al., 2018). Flaming combustion occurs at high temperature and consists of volatilization of hydrocarbons from the thermal decomposition of biomass leading to rapid oxidation and efficient combustion, producing CO₂, water and black carbon (BC). Instead, smoldering combustion is flameless and can be initiated by weak sources of heat and results in less efficient combustion of fuel, leading to gas-phase products (mainly CO, CH₄ and volatile organic compounds).

2.2 Instrumentation

We characterized the emissions with a suite of gas- and particle-phase instrumentation. Organic gases (OGs) were measured by a proton-transfer-reaction time-of-flight mass spectrometer (PTR-ToF-MS 8000, Ionicon Analytik). A detailed description of the instrument can be found in Jordan et al. (2009). The PTR-ToF-MS was operated under standard conditions (ion drift pressure of 2.2 mbar and drift intensity of 125 Td) in H₃O⁺ mode, allowing the detection of OGs with a proton affinity higher than that of water. Quantification was performed by using known individual reaction rate constants (Cappellin et al., 2012), otherwise a value of $\mathrm{2}\times {\mathrm{10}}^{-\mathrm{9}}$ cm³ s⁻¹ was assumed. The effective rate constants applied to both Set 1 and Set 2 can be found in Bruns et al. (2017). Data were analysed using the Tofware software 2.4.2 (PTR module as distributed by Ionicon Analytik GmbH, Innsbruck, Austria) running in Igor Pro 6.3.

Non-refractory primary and aged particle composition was monitored by a high-resolution time-of-flight aerosol mass spectrometer (HR-ToF-AMS, Aerodyne Research Inc.) (DeCarlo et al., 2006). The HR-ToF-AMS is described in detail elsewhere (Bruns et al., 2016; Bertrand et al., 2017) and summarized here. The instrument was operated under standard conditions (temperature of vaporizer 600 ^∘C, electronic ionization (EI) at 70 eV, V mode) with a temporal resolution of 10 s. Data analysis was performed in Igor Pro 6.3 (Wave Metrics) using SQuirrel 1.57 and Pika 1.15Z assuming a collection efficiency of 1. The O : C ratio was determined according to Aiken et al. (2008).

The concentration of equivalent black carbon (eBC) was determined from the absorption coefficient measured with a seven-wavelength Aethalometer (Magee Scientific Aethalometer model AE33) with a time resolution of 1 min at a wavelength of 880 nm (Drinovec et al., 2015).

The particle number concentration and size distribution (16 to 914 nm) were provided by a scanning mobility particle sizer (SMPS, consisting of a custom-built differential mobility analyser (DMA) and a condensation particle counter (CPC 3022, TSI)) with a time resolution of 5 min. Supporting gas measurements included a CO₂ analyser (LI-COR), a CH₄ monitor, a total hydrocarbon (THC) monitor (flame ionization detector, THC monitor Horiba APHA-370), and NO and NO₂ (NO_x analyser, Thermo Environmental) monitors.

3 Data analysis

The data analysis entails three steps detailed in this section. The first subsection describes the determination of the amount of oxidized OGs in the chamber. The second subsection details the determination of the amount of SOA formed in the chamber. The last subsection describes the box model used for the parameterization of SOA formation from the OGs.

3.1 Organic gas loss in the chamber

In the chamber, OGs were oxidized to several oxidation products, referred to as oxidized OGs (CG, condensable gases) in the following analysis. According to their volatility, these products may remain in the gas phase or partition to the particle phase, thereby contributing to SOA formation.

We described the change in any OG concentration over time as a combination of its loss and production as follows:

$$\begin{array}{}\text{(1)}& {\displaystyle }\begin{array}{rl}& {\displaystyle \frac{\text{d}\left[\mathrm{OG}\right]}{\text{d}t}}=P-\\ & \left(\sum k_{\text{dil}}\times \left[\mathrm{OG}\right]+k_{\mathrm{OH}}\times \left[\text{OH}\right]\times \left[\mathrm{OG}\right]+k_{\text{other}}\times \left[\mathrm{OG}\right]\right)\end{array}.\end{array}$$

Here, P corresponds to the production of an OG in the chamber, e.g. from the oxidation of other primary OGs. k_dil is the dilution rate constant in reciprocal seconds. k_OH [OH] [OG] in molecules per cubic centimetre per second represents the consumption rate due to oxidation by OH, where k_OH is the reaction rate constant and [OH] is the OH concentration. k_other[OG] in molecules per cubic centimetre per second is the loss rate of OG by other processes, where k_other is the reaction rate constant in reciprocal seconds. The loss of some OGs could not be explained by their reaction with OH and dilution alone for Set 1, so we added this additional term, which is discussed after the first two processes are constrained. We considered primary OGs that exhibited a clear decay with time to be strictly of primary origin, and hence neglected their production from other OGs (i.e. P=0). This assumption signifies that the yields estimated under our conditions are upper limits. In reality, the detection of aromatic hydrocarbons (e.g. single-ring aromatic hydrocarbons, SAHs, and polycyclic aromatic hydrocarbons, PAHs) by the PTR-ToF-MS may be affected by the interference due to fragmentation during ionization of their oxidation products (Gueneron et al., 2015). On the other hand, directly emitted oxygenated aromatics could themselves be the oxidation products of aromatic hydrocarbons and their production may continue during the experiment. However, the assumption of P=0 does not introduce a significant error for most OGs with significant primary emissions because the observed OG decay was consistent with their OH reaction rate constant for Set 2 as demonstrated by Bruns et al. (2017) for the +15 ^∘C conditions. In the following, we describe the processes governing the changes in the OG concentrations in the chamber and the approaches adopted for the determination of the different parameters in Eq. (1).

3.1.1 Reaction with OH radical

The OH exposure, which is the integrated OH concentration over time, was estimated based on the differential reactivity of two OGs. Specifically, we used d9-butanol (fragment at mass-to-charge ratio m∕z 66.126, [C₄D₉]⁺) and naphthalene (fragment at m∕z 129.070, [C₁₀H₈]H⁺). These compounds are selected because they can be unambiguously detected (no isomers or interferences, high signal-to-noise ratio), are not produced during the experiment and have OH reaction rate constants that are precisely measured and significantly different from each other. The OH exposure can be expressed as follows:

$$\begin{array}{}\text{(2)}& \text{OH exposure}=\left({\displaystyle \frac{\ln {\left(\frac{\mathrm{d}\mathrm{9}-\mathrm{butanol}}{\mathrm{naphthalene}}\right)}_{\mathrm{0}}-\ln {\left(\frac{\mathrm{d}\mathrm{9}-\mathrm{butanol}}{\mathrm{naphthalene}}\right)}_t}{k_{\mathrm{OH},\mathrm{but}}-k_{\mathrm{OH},\mathrm{naph}}}}\right),\end{array}$$

where (d9-butanol ∕ naphthalene)₀ is the ratio between these compounds at t=0 (before lights were turned on), (d9-butanol ∕ naphthalene)_t is the ratio measured at time t, and k_OH,but and k_OH,naph are the OH reaction rate constants of d9-butanol and naphthalene, respectively ($k_{\mathrm{OH},\mathrm{but}}=\mathrm{3.14}\times {\mathrm{10}}^{-\mathrm{12}}$ cm³ molec.⁻¹ s⁻¹ and $k_{\mathrm{OH},\mathrm{naph}}=\mathrm{2.30}\times {\mathrm{10}}^{-\mathrm{11}}$ cm³ molec.⁻¹ s⁻¹) (Atkinson and Arey, 2003).

For Set 1, the OH exposure at the end of each experiment ranged between 5 and 8×10⁶ molec. cm⁻³ h, corresponding to approximately 5–8 h in the atmosphere (given global average and typical wintertime OH concentrations of 1×10⁶ molec. cm⁻³). For Set 2, higher OH exposures were reached (3 to 7×10⁷ molec. cm⁻³ h at the end of each experiment, corresponding to 2–3 d in the atmosphere). This is likely because both sets of experiments utilized a similar HONO molar flow (and thus similar OH production rate), but higher OG concentrations were reached in Set 1, which could possibly have resulted in a higher OH sink. We calculated the OH concentration, [OH] in Eq. (1), numerically as d(OH exposure)∕dt.

3.1.2 Smog chamber dilution

For Set 2, dilution was dominated by the constant injection of HONO to the chamber and accounted for as described in Bruns et al. (2016). For Set 1 the dilution rate of primary OGs in the chamber was calculated as follows. The integrated dilution over time, K_dil, was determined as the ratio between the d9-butanol concentration corrected for the reaction of d9-butanol with OH and the d9-butanol concentration at $t=\mathrm{0}\left(\right[\mathrm{d}\mathrm{9}-\mathrm{butanol}{]}_{\mathrm{0}}$):

$$\begin{array}{}\text{(3)}& K_{\text{dil}}={\displaystyle \frac{\left[\mathrm{d}\mathrm{9}-\mathrm{butanol}\right]\times e^{k_{\mathrm{OH},\mathrm{but}}\times \text{OH exposure}}}{{[\mathrm{d}\mathrm{9}-\mathrm{butanol}]}_{\mathrm{0}}}}.\end{array}$$

Figure S1 shows calculated values for K_dil as a function of time. Increasingly high dilution at the end of these experiments (Fig. S1), between 20 % and 35 % (i.e. a dilution ratio of 0.8 to 0.65), is a result of constant injection of HONO and excessive sampling at high rates, resulting in inputs from laboratory air (likely through leaks in the Teflon bag or the connections of chamber inlets and outlets). The inputs from the HONO pure carrier gas and of laboratory air are roughly comparable (about 12 % vs. 8 %–23 % at the end of the experiment). Other than dilution, the effect of both inputs on the gas-phase chemical composition is negligible. Dilution rates are non-linear, increasing as the experiment progresses due to continuous dilution within a decreasing chamber volume. The dilution rate constant, k_dil, in Eq. (1), is the differential of K_dil over time.

3.1.3 Losses by other processes

To corroborate the estimated OH exposures and dilution rates, we examined the loss of prominent OGs with known reaction rate constants against OH. For Set 2 and as demonstrated by Bruns et al. (2017), the OG decay is consistent with their estimated loss based on their dilution and reaction with OH. By contrast, for Set 1, the decay of some OGs could not be solely explained by their reaction with OH and dilution, suggesting additional reactions with oxidants other than OH as discussed below. The OG total consumption by this process (∫k_other[OG]) was estimated as the difference between the total measured decay of the OG of interest and the fraction consumed by both dilution and oxidation by OH radicals. The reaction rate constants for several precursors towards OH (k_OH) are not available in the literature. In addition, many fragments may have several isomers, each of which is associated with different rate constants. Effective rate constants for all precursors considered were estimated from their decay in Set 2, where the combination of OH reaction and dilution fully explained the decay of OGs with known OH reaction rate constants.

3.1.4 Precursor classification

A common set of 263 ions was extracted from the PTR-ToF-MS. Among these ions, 86 showed a clear decay with time and were thus identified and selected as potential SOA precursors. Previous work based on Set 2 experiments showed that the PTR-ToF-MS measures the most important SOA precursors, which explained the measured SOA mass within 40 % uncertainty and without systematic bias (Bruns et al., 2016). Therefore, these compounds are expected to capture the dominant fraction of SOA mass, although we cannot rule out losses in the PTR-ToF-MS inlet or small contributions from other precursors such as alkanes. The compound identification was supported by previous publications (McDonald et al., 2000; Fine et al., 2001; Nolte et al., 2001; Schauer et al., 2001; Stockwell et al., 2015), including gas chromatography–mass spectrometry (GC–MS) analysis when available.

The size of our dataset does not allow us to retrieve the volatility distribution for single precursors, which would entail the determination of more than 86 free parameters. This is especially the case as the time series of precursors, decaying with oxidation, are typically strongly correlated, which prevents resolving systematic differences between the yields of the different single precursors. Therefore, lumping is needed to decrease the model degree of freedom. Accordingly, precursors are grouped in six chemical classes: furans, single-ring aromatic hydrocarbons (SAHs), polycyclic aromatic hydrocarbons (PAHs), oxygenated aromatics (OxyAH) and organic compounds containing six and more or fewer than six carbon atoms (OVOC_c ≥ 6 and OVOC_c < 6, respectively) (Table S1). The selected precursors in each class are the same for each experiment and each dataset. This lumping approach is based on the two main objectives of the study:

1. Compare the SOA yields of specific precursors determined in complex emission experiments with those determined in single-compound systems.

2. Identify the main SOA precursors in flaming and smoldering wood emissions, at different temperatures.

To be able to compare to literature yields, we lumped species that have similar yields in the same chemical class: e.g. at an organic aerosol concentration of 10 µg m⁻³ the yield of PAHs is ∼ 20 % (objective 1). In addition, we classified the precursors based on the pathway by which they are emitted, which will allow us to determine which compounds dominate SOA formation in flaming and smoldering emissions (objective 2). We differentiated between oxygenated aromatics, mainly emitted through lignin pyrolysis, furans emitted through cellulose pyrolysis, and single-ring aromatics and PAHs generated from incomplete combustion, especially from flaming wood. The remaining SOA precursors are all oxygenated OGs; therefore we separated them according to their carbon number knowing that larger precursors will have higher yields than smaller precursors. Our ability to precisely extract yields specific to a precursor class heavily relies on differences in the oxidation rates or emission patterns of the precursors. Therefore, the classification approach adopted here, where classes are expected to have different contributions during different experiments (Bhattu et al., 2019), facilitates the extraction of yields of the different classes.

3.2 Calculation of the OA mass in the chamber

The total organic aerosol measured by the HR-ToF-AMS was corrected for particle losses in the chamber due to gravitational and diffusional deposition. To assess the total wall losses due to both processes, we assumed that the condensable vapours partition only to the suspended aerosols but not to the wall.

Assuming that black carbon is inert in the chamber, it was possible to use its decay to estimate the particle loss to the walls. The aerosol attenuation measured at 880 nm (at the end of each experiment, ∼ 4 h) with an Aethalometer was used to estimate the particle loss rate to the wall. This attenuation is proportional to the eBC mass concentration and within uncertainties independent of the ageing extent, as demonstrated in Kumar et al. (2018). Using eBC as a tracer, we inherently assumed that eBC and OA were internally mixed and homogeneously distributed over the aerosol size range. The decay of eBC due to both dilution and deposition onto the chamber walls was parametrized as follows:

$$\begin{array}{}\text{(4)}& {\displaystyle \frac{\text{d}\left[\mathrm{eBC}\right]}{\text{d}t}}=-k_{\text{dil}}\left[\mathrm{eBC}\right]-k_{\text{wall}}\left[\mathrm{eBC}\right],\end{array}$$

where k_wall is the first-order wall loss rate constant used to correct the measured OA concentration for wall losses, ranging between 4 and $\mathrm{8}\times {\mathrm{10}}^{-\mathrm{5}}$ s⁻¹, and k_dil is the dilution rate constant determined above.

The wall-loss-corrected organic aerosol, OA_WLC, was calculated using Eq. (5):

$$\begin{array}{}\text{(5)}& {\mathrm{OA}}_{\mathrm{WLC}}\left(t\right)=\mathrm{OA}\left(t\right)+\underset{\mathrm{0}}{\overset{t}{\int }}{k}_{\text{wall}}\times \mathrm{OA}\left(t\right)\times \text{d}t,\end{array}$$

where OA(t) is the measured organic aerosol concentration in microgram per cubic metre. The total OA present in the chamber was estimated as the suspended OA concentration measured by the HR-ToF-AMS plus the estimated OA lost to the wall. This concentration was directly compared to the condensable gas (CG) concentration estimated according to Eq. (A8) presented in Appendix A. As mentioned, our approach did not take into account the losses of precursor vapours or their oxidation products in the gas phase onto the walls. We note that these processes were unlikely to have a substantial effect on the precursors considered, which were largely highly volatile species, even at lower temperature. Based on the calculation of the equilibrium constant of semi-volatile species on the walls by Bertrand et al. (2018b), we estimated that at 293 K the fraction of these compounds absorbed on the walls is < 5 %. Meanwhile, the walls could indeed act as a sink for the semi-volatile oxidation products. This effect was not taken into account in the current study, but we expect that it was minimized under our conditions, by the high OA concentration in the chamber and high production rates (Zhang et al., 2014; Nah et al., 2017).

3.3 Modelling SOA formation

The general aim of the model is the determination of the parameters describing the volatility distributions of the oxidation products from different precursor classes and their temperature dependence. A simplified schematic of the modelling framework is described in Fig. 1. It consists of (1) a box model that describes the partitioning of the condensable gases generated through oxidation, (2) the model input parameters obtained from the smog chamber, (3) the model output parameters and (4) the model optimization based on a genetic algorithm (GA). Each of these parts is described in the following sections.

Figure 1 Schematic of the modelling framework. The box model simulates the formation of SOA from each precursor class j in volatility bins i. The best solution of the initialized input parameter volatility distribution (described as the mean value of logC* for each precursor class μ_j and standard deviation σ) and enthalpy of vaporization (ΔH_vap) parameters are optimized by a genetic algorithm, using minimum mean bias and root-mean-square error (RMSE) between modelled and measured OA concentration as the fitness function. Green boxes represent measured data from the chamber experiments.

3.3.1 Box model

We assume the partitioning of CG between the gas and the particle phases to obey Raoult's law (Strader et al., 1999), where the aerosol can be described as a pseudo-ideal organic solution, of SOA and POA species. The volatility basis set (VBS, implemented by Koo et al., 2014, in the Comprehensive Air quality Model with eXtensions, CAMx) was used to classify the oxidation products of the different precursors into surrogates with different volatility, distributed into discrete logarithmically spaced bins (Donahue et al., 2006).

We considered the most basic mechanism by which SOA may form. That is, the oxidation products from the different precursor classes described above instantaneously partition into the condensed phase depending on their volatility. No additional reactions in the gas or particle phase were considered (e.g. reaction with oxidants, photolysis or oligomerization). In addition, we neglected the contribution of primary oxidation products of the gas-phase semi-volatile species (co-emitted with POA) compared to the OGs detected by the PTR-ToF-MS, based on the findings of Bruns et al. (2016) and Ciarelli et al. (2017a). Finally, we considered the species in the gas and the particle phases to be permanently at equilibrium, as condensation is expected to be faster than oxidation (timescales for oxidized vapour condensation < 1 min assuming no particle-phase diffusion limitations; Bertrand et al., 2018b). While including additional processes in the model is feasible, this would result in a significantly higher-dimensional parameter space, which cannot be unambiguously inferred from the present data. We consider that without supportive data, e.g. chemically resolved characterization of the particle-phase species, such reactions could not be well constrained or even deduced from structure activity relationships, given the many unknowns in complex emissions. Therefore, such a simplified scheme of SOA formation from complex emissions may be compared in the future with chemically resolved data to help the identification of additional mechanisms that were not considered here.

The derivation of the thermodynamic equations governing SOA formation from precursors, implemented in the box model, is detailed in Appendix A and only a brief description of the model principles is given here. In the following, let i and j be the indices for the different volatility bins and precursor classes, respectively. The model determines the molar distribution of the oxidation products from different precursor classes in the different volatility bins, Υ_i,j, together with the compounds' enthalpy of evaporation, $\mathrm{\Delta }{H}_{\text{vap},i,j}$. The latter describes the temperature dependence of the oxidation products' effective molar saturation concentration, $x_{i,j}^{*}$. For this, the model iteratively solves Eqs. (A6) and (A7) at every experimental time step (time resolution of 10 s) for all experiments, to retrieve the surrogate molar concentrations in the particle phase, x_i,j|_p, and the total surrogates' molar concentration in the condensed organic phase, x_OA (see Appendix A).

3.3.2 Model inputs

The model uses as main inputs the molar concentrations of the condensable gases from different precursors (in total n= six precursor classes) in both phases, x_j|_g + p. The concentration of condensable gases in both the gas and particle phases is derived from the consumption rates of VOC_j determined by the PTR-ToF-MS, by numerically integrating Eq. (A8). x_j|_g + p is related to the concentrations of the different surrogates from a precursor class j in different volatility bins, $x_{i,j}{\mathrm{|}}_{\text{g}+\text{p}}$ (Eq. A6), through their yields, Υ_i,j, according to Eq. (A9). The number of volatility bins, m, is set to six, approximately corresponding to the following mass saturation concentrations: $C_{i,j}^{*}\left(\mathrm{\mu }\mathrm{g}{\mathrm{m}}^{-\mathrm{3}}\right)=\left\{{\mathrm{10}}^{-\mathrm{1}};{\mathrm{10}}^{\mathrm{0}};{\mathrm{10}}^{\mathrm{1}};{\mathrm{10}}^{\mathrm{2}};{\mathrm{10}}^{\mathrm{3}};{\mathrm{10}}^{\mathrm{4}}\right\}$.

In addition to x_j|_g + p, the model needs as input $x_i^{{\text{OM}}_{\text{p}}}{\mathrm{|}}_{\text{p}+\text{g}}$, the molar concentration of primary organic matter from a volatility bin i in both gas and particle phases (Eq. A7). $x_i^{{\text{OM}}_{\text{p}}}{\mathrm{|}}_{\text{g}+\text{p}}$ is inferred from the measured POA concentrations injected in the chamber at the beginning of the experiment and using the volatility distribution function of wood combustion emissions in May et al. (2013). It is assumed constant with ageing. The computation of the fraction of POA_i in the condensed phase is similar to that for SOA species in Eq. (A6).

The secondary surrogates' elemental composition (C#_i,j, O#_i,j and H#_i,j) is also used as model input to compute the surrogates' molecular weight, MW_i,j, required for C_OA calculations (see Sect. 3.3.3). A single C#_i,j value is calculated per chemical class, based on the average $C{\mathit{\#}}_j^{\text{VOC}}$ of the respective precursor class, and considering $C{\mathit{\#}}_{i,j}=C{\mathit{\#}}_j^{\text{VOC}}-\mathrm{\Delta }C$, where ΔC is the average loss in carbon due to fragmentation during the oxidation of precursors from all classes. ΔC is determined by systematically changing its value in multiple model runs and selecting the value that best explains the observed O : C ratios (see Fig. 8b). Likewise, a single H#_i,j value is assumed per chemical class, considering that H#_j∕C#_j equals $H{\mathit{\#}}_j^{\text{VOC}}/C{\mathit{\#}}_j^{\text{VOC}}$. Finally, O#_i,j is constrained by the C#_j and the surrogate volatility ($C_{i,j}^{*}$) based on the simplification of the SIMPOL model (Pankow and Asher, 2008), provided by Eq. (3) in Donahue et al. (2011). Based on this relationship, O#_i,j increases with decreasing C#_j and $C_{i,j}^{*}$. The O : C ratio of primary emissions is constrained in the model to the measured O : C in the beginning of each experiment, by setting $C{\mathit{\#}}^{{\text{OM}}_{\text{p}}}$ and calculating $O{\mathit{\#}}_i^{{\text{OM}}_{\text{p}}}$ using the same methodology as for O#_i,j. The resulting $C{\mathit{\#}}^{{\text{OM}}_{\text{p}}}$ and $O{\mathit{\#}}_i^{{\text{OM}}_{\text{p}}}$ and the corresponding primary organic matter molecular weight, MW${}_i^{{\text{OM}}_{\text{p}}}$, as well as C#_j, O#_i,j, H#_j and MW_i,j, is reported in Table S4.

3.3.3 Model outputs

The model provides the Υ_i,j and $\mathrm{\Delta }{H}_{\text{vap},i,j}$ parameters. To reduce the model's degree of freedom we consider a single ΔH_vap for all surrogates from different chemical classes in different volatility bins. Υ_i,j is considered to follow a kernel normal distribution as a function of log C^*, ${\mathrm{{\rm Y}}}_{i,j}\sim N\left({\mathit{\mu }}_j,\mathit{\sigma }\right)$, where μ_j is the median value of log C^* and σ is the standard deviation. This step (1) ensures positive Υ_i,j parameters, (2) significantly reduces the model's degree of freedom and (3) allows constraint of the total concentration of surrogates from a certain chemical class: ${\sum }_i^m{\mathrm{{\rm Y}}}_{i,j}=\mathrm{1}$.

The set of Υ_i,j and $\mathrm{\Delta }{H}_{\text{vap},i,j}$ parameters is determined by minimizing the sum of mean bias (MB) and the root-mean-square error (RMSE) between modelled mass concentrations of the particulate organic phase, C_OA, calculated using Eq. (6), and concentrations measured by the AMS.

$$\begin{array}{}\text{(6)}& C_{\text{OA}}={\left.\sum _j^n{\sum }_i^m{\text{MW}}_{i,j}{x}_{i,j}\right|}_{\text{p}}+{\left.{\text{MW}}_i^{{\text{OM}}_{\text{p}}}{x}_i^{{\text{OM}}_{\text{p}}}\right|}_{\text{p}}\end{array}$$

The model fitted to the measured C_OA was also validated by external AMS measurements of the O : C ratio determined through high-resolution analysis. The modelled O : C ratio was calculated at every experimental time step as follows:

$$\begin{array}{}\text{(7)}& \mathrm{O}:\mathrm{C}={\displaystyle \frac{{\left.\sum _j^n{\sum }_i^{m}O{\mathit{\#}}_{i,j}{x}_{i,j}\right|}_{\text{p}}+{\left.O{\mathit{\#}}_i^{{\text{OM}}_{\text{p}}}{x}_i^{{\text{OM}}_{\text{p}}}\right|}_{\text{p}}}{{\left.\sum _j^n{\sum }_i^{m}C{\mathit{\#}}_{i,j}{x}_{i,j}\right|}_{\text{p}}+{\left.C{\mathit{\#}}_i^{{\text{OM}}_{\text{p}}}{x}_i^{{\text{OM}}_{\text{p}}}\right|}_{\text{p}}}}.\end{array}$$

3.3.4 Model optimization

The model is optimized to determine the volatility distributions of the oxidation products from different precursor classes described by μ_j, σ and their temperature dependence described by ΔH_vap to best fit the observed OA concentrations. For the model optimization, we used a genetic algorithm (GA), a metaheuristic procedure inspired by the theory of natural selection in biology, including selection, crossover and mutation processes, to efficiently generate high-quality solutions to optimize problems (Goldberg et al., 2007; Mitchell, 1996). The GA is initiated with a population of randomly selected individual solutions. The performance of each of these solutions is evaluated by a fitness function, and the fitness values are used to select more optimized solutions, referred to as parents. The new generation of solutions (denoted children) is produced by either randomly changing a single parent (as mutation) or combing the vector entries of a pair of parents (as crossover). The evolution process will be repeated until the termination criterion is reached, here maximum iteration time. In this study, a population of 50 different sets of model parameters (μ_j, σ and ΔH_vap) was considered for each GA generation. The sum of mean bias and RMSE between measured and modelled C_OA of the 14 experiments was used as a fitness function to evaluate the solutions. We assume the termination criterion is reached if no improvement in the fitness occurs after 50 generations, with a maximum of 500 total iterations allowed. The GA calculations were performed using the package “GA” for R (Scrucca et al., 2017). A bootstrap method was then adopted to quantify the uncertainty in the constrained parameters.

Figure 2 Comparison of the OG emissions between the two sets of experiments. Panels (a) and (b) display average primary OG mass spectra from Set 1 (representing the smoldering phase) and Set 2 (representing mainly the flaming phase), respectively. Spectra are normalized to the initial total OG concentration in micrograms per cubic metre. Compounds are colour coded by chemical classes. (c) The p value vs. fold change comparing the fingerprints of primary OGs between the two sets of experiments. The fold change was calculated as the ratio of the intensities of each ion normalized to the total signal, between Set 2 and Set 1 averaged across experiments. Data points above p=0.05 have significantly different contributions to the total OGs between the two sets of experiments. Blue coloured data points on the right-hand side designate compounds enriched in the emissions from Set 2 while green coloured data points on the left-hand side designate compounds enriched in the emissions from Set 1. (d) Spearman correlation matrix for the primary OG mass spectra between all experiments highlighting the variability in the composition of the primary emissions. Each experiment is identified by an index (see Table 1) and the experiments from Set 2 were reordered according to similarity with Set 1.

4 Results and discussion

4.1 Comparison of primary emissions across experiments

A larger amount of primary OG was emitted into the chamber in Set 1, with concentrations ranging from 950 to 7860 µg m⁻³, while in Set 2 the primary OG concentrations ranged from 300 to 1360 µg m⁻³. In the same way, the measured OA at the beginning of each test (POA) ranged from 10 to 180 µg m⁻³ and from 9 to 22 µg m⁻³ for Set 1 and Set 2, respectively (Table 1). The OG composition for Set 1 and Set 2 is summarized in Fig. 2 showing the mean PTR-ToF-MS mass spectra (Fig. 2a, b), the relative contributions of the different compounds for different datasets (Fig. 2c) and the variability in composition among all experiments (Fig. 2d). Set 1 shows higher relative contributions of furans and OVOC_c < 6 (Fig. 2a), while the contributions of PAH, SAH and OxyAH are higher in Set 2 (Fig. 2b). The OxyAH compounds, mainly methyl- and methoxyphenols, are produced by lignin pyrolysis (Fine et al., 2001) while furans are formed from cellulose pyrolysis (Mettler et al., 2012). The chemical class referred to as “Others” comprises compounds that do not show a clear decay upon oxidation and are therefore not considered SOA precursors in the following analysis. The category others is dominated by acetic acid, previously reported as a major species in residential wood burning emissions (Bhattu et al., 2019). The majority of the compounds differ among datasets and the most significant difference estimated through the p value (probability associated with a Student's t test) occurs for the OVOC_c < 6, which is about a factor of 6 higher than the significance threshold (p=0.05) for Set 1. In order to investigate the similarity between all experiments a Spearman correlation matrix was calculated. Experiments from Set 1 appear to be consistently similar to each other while the experiments from Set 2 are significantly different among each other in terms of composition of the primary emissions. A possible reason for such a discrepancy was the difficulty in injecting flaming emissions without any significant smoldering contribution for Set 2. This hypothesis is supported by the strong similarity between precursor compounds measured in some experiments supposed to represent flaming-phase emissions only but apparently included some smoldering as well (9, 12, 13) and experiments from Set 1 (1–7). Figure S2 reports the relative contributions of primary OGs (before photo-oxidation) for all compound classes for the 14 experiments. Note that these trends are not correlated with the modified combustion efficiency (MCE), reported in Table 1, defined as ${\mathrm{CO}}_{\mathrm{2}}/(\mathrm{CO}+{\mathrm{CO}}_{\mathrm{2}})$, which was constant at 0.97 g g⁻¹, for Set 2, but ranged from 0.8 g g⁻¹ to 0.91 g g⁻¹ for Set 1.

Figure 3 (a) Box plot for the relative contributions to the total primary emissions of the different precursor classes and of OM_p averaged between experiments. OM_p is the total semi-volatile and low-volatility organic matter calculated by means of the VBS model assuming the volatility distribution from May et al. (2013) and using as input the measured organic aerosol (OA) mass. The top and bottom whiskers represent the 90th and 10th percentiles, respectively, while the top, middle and bottom lines of the boxes show the 75th, 50th and 25th percentiles, respectively. The circles represent each single experiment from the two datasets investigated. (b) Average contributions of different precursor families and OM_p for Set 1 and Set 2.

Figure 3 shows the contribution of each class of precursor compounds and of the primary semi-volatile organic matter (OM_P) to the total primary emissions. OM_P is the total organic matter in the semi-volatile and low-volatility range in the particle and the gas phases (saturation concentration < 1000 µg m⁻³ at 298 K) (see Sect. 3.3.2). Overall, the highest average relative contributions are related to OVOC_c < 6 followed by furans and OxyAH, but we also note an average large contribution by SAH for Set 2. The two sets of experiments investigated clearly show different primary composition of emissions in terms of dominant contributions; in Set 1 OVOC_c < 6 and OxyAH dominate by far the total primary emissions while in Set 2 the main species influencing the total primary emissions are OxyAH, OVOC_c < 6 and SAH with roughly similar contributions (see Fig. 3b). Moreover, the calculated averaged OM_p ∕ OGs ratios are around 0.05 and 0.03 for Set 2 and Set 1, respectively.

OGs undergo oxidation during atmospheric ageing to form a complex mixture of products, some of which remain in the gas phase while others have sufficiently low volatility to partition to the particle phase. The consumption of the different OG classes over time is shown in Fig. 4a, b for Set 1 and Set 2, respectively. The general trends manifest that PAH and OxyAH are the most reactive classes, exhibiting an average consumption of up to 80 % at the end of the experiments (after ∼ 4 h of ageing) while for both datasets SAH appears to be the least reactive class (with an average consumption between 10 % and 20 % at the end of the experiments). Relevant compounds in the latter class are benzene (C₆H₆), toluene (C₇H₈) and xylene (C₈H₁₀); their slow reactivity is consistent with literature reaction rate constants toward OH, from the NIST database (NIST chemistry WebBook, Linstrom and Mallard, 2018), of $\mathrm{1.22}\times {\mathrm{10}}^{-\mathrm{12}}$, $\mathrm{6.13}\times {\mathrm{10}}^{-\mathrm{12}}$ and $\mathrm{7.51}\times {\mathrm{10}}^{-\mathrm{12}}$ (cm³ molec.⁻¹ s⁻¹), respectively.

The two datasets differ for the OH dose; we observe in Set 2 (Expt. 8–14) an overall higher consumption of all precursor classes due to the higher OH dose (representing a longer ageing time in an ambient atmosphere). PAH shows the highest reactivity followed by OxyAH while for Set 1 (Expt. 1–7) the fastest class of compounds to react is the OxyAH followed by PAH and OVOC_c < 6. Moreover, despite the lower OH exposure reached for Set 1, the consumption of OxyAH and furans is substantially higher at comparable exposure levels.

Table 2 Reported average reaction rate constants (10⁻¹¹ cm³ molec.⁻¹ s⁻¹) towards OH per family at the beginning of each experiment, including average reactivity (AVG) and standard deviation (SD).

Figure 4 Average consumption of SOA precursor classes against average OH exposure for Set 1 (a) and Set 2 (b). The observed decay of SOA precursor families (OGs) as described in Eq. (1) is due to both oxidation processes and dilution in the chamber. (c) SOA-to-POA ratio for each experiment against average OH exposure coloured according to experimental temperatures (2, −10, 15 ^∘C).

Figure 5 Box plot of mass-weighted average OH reaction rate constants (k_OH) determined for each precursor class from Bruns et al. (2016) for Set 2 only (see Table S3). The individual k_OH values for all compounds are also shown for all experiments, colour coded according to the experimental temperatures. The top and bottom whiskers represent the 90th and 10th percentiles, while the top, middle and bottom lines of the boxes show the 75th, 50th and 25th percentiles, respectively.

In the same way, we also observe a higher SOA production for Set 1 compared to Set 2 at comparable OH exposure. The SOA-to-POA ratio ranges between 2 and 6, similar to ratios observed in previous studies (Heringa et al., 2011; Bruns et al., 2015; Grieshop et al., 2009; Tiitta et al., 2016).

Overall the same chemical classes appear to behave differently across the different sets of experiments. Such an inconsistency in behaviour is due to either differences in the chemical composition within the same class or due to additional reactivity occurring in Set 1.

To investigate the chemical differences within the same class of compounds across different experiments, Table 2 reports the average reaction rate constants against OH of the different chemical classes calculated at the beginning of each experiment following Eq. (8).

$$\begin{array}{}\text{(8)}& {\stackrel{\mathrm{‾}}{k}}_{{\mathrm{OH}}_{j,k}}=\sum _i{k}_{{\mathrm{OH}}_{c,j}}\times {\displaystyle \frac{{\mathrm{OG}}_{c,j,k}}{{\mathrm{OG}}_{j,k}}}\end{array}$$

Here, c represents the single compound, j the family and k the experiment. k_OH is the reaction rate constant toward OH and OG refers to the primary OGs.

OH reaction rate constants ($k_{{\mathrm{OH}}_{c,j}}$) for each compound were calculated from Set 2 only (Fig. 5), whereas the decay of the precursors contributing most to SOA formation during ageing was compared with the expected decay based on literature. The good agreement indicates that for these experiments the consumption of the precursors was dominated by OH (Bruns et al. 2017). The average OH reaction rate constants ($k_{{\mathrm{OH}}_{c,j}}$) are reported in Table S3. They are determined for each precursor class and calculated with a first-order exponential fitting on the precursors' decay curves previously corrected for dilution.

Figure 6 Fraction of consumed precursor compounds for Set 1 by OH oxidation (a), dilution (b) and other reactivity (c) at the end of the experiments. Each point corresponds to a single compound averaged among experiments and normalized to the initial concentration. The colour legend represents the statistically significant deviation from zero reactivity with the investigated reactant (p value).

The average reaction rate constants per family (${\stackrel{\mathrm{‾}}{k}}_{{\mathrm{OH}}_{j,k}})$ are similar among the same families for different experiments, suggesting that the variable behaviour of the chemical classes across different experiments was due to differences in the reactive environment rather than a different chemical composition within a given class.

As introduced in Sect. 3.1.3, the total measured decay of OGs in Set 1 could not be fully explained by dilution and reactivity against OH, suggesting the presence of an additional loss process. To assess the remaining oxidation processes, k_OH values were used to estimate the missing loss process for Set 1 according to Eq. (1). In this way the consumed fraction due to OH chemistry, dilution in the chamber and the additional reactivity was calculated for each OG compound family and is shown in Fig. 6. The additional reactivity appears to contribute to the total precursor consumption for most of the classes, with particular relevance for the OxyAH and PAH classes. This is in contrast with the calculated loss for Set 2 (Fig. S3) where the dominant consumption is due to OH.

One possible hypothesis is that the remaining loss process might be due to reaction with the nitrate radical (NO₃), which absorbs in the visible region (∼ 500–650 nm) and thus is not efficiently photolysed by the black lights used here (Reed et al., 2016). Figure S4 shows that for Set 1 the NO₃ reaction rate constants ($k_{{\mathrm{NO}}_{\mathrm{3}}}$) for compounds found in the NIST database (Linstrom and Mallard, 2018) (see Table S1) are well correlated with the amount reacted, making nitrate chemistry a likely loss process (Schwantes et al., 2017). We note that while the OH production rate in the chamber is similar for the two sets of experiments, given the same injection rate of HONO, the OH total reactivity is significantly higher for Set 1 because of the injection of ∼ 5 times higher OG concentrations. As a result, the OH concentration is around 1.5×10⁶ vs. 10×10⁶ molec. cm⁻³ for Set 1 and Set 2, respectively, increasing the availability of OGs for consumption by other processes during Set 1.

4.2 Model evaluation

As previously mentioned, the reaction of the OGs produces oxidation products (CG, condensable gases), which can consequently partition between the gas and particle phases. Their concentration, estimated by accounting for the production minus the loss in the chamber, as described in Eq. (1), was used as input for the box model with the volatility basis set (VBS) scheme. Figure S5 shows the volatility distributions of different precursor classes for Set 1 and Set 2 to assess the influence of different processes occurring in different experiments.

We note that for Set 1 the lower-volatility bins exhibit higher contributions compared to Set 2 but still within 2 standard deviations, such that it is difficult to distinguish statistically different yields for most of the cases. Hence, the total volatility distribution was used to calculate the mass yield instead of a specific one for each dataset (Fig. 9). Mass yields calculated for the specific datasets Set 1 and Set 2 are reported in Fig. S6.

Figure 7 (a) Comparison between the sum of simulated organic aerosol (OA) concentrations from different precursor families and the measured OA concentration as a function of the OH exposure for each smog chamber experiment. Set 1 (Expt. 1–7, 2 ^∘C) and Set 2 (Expt. 8–9, −10 ^∘C and Expt. 10–14, 15 ^∘C). (b) Probability distribution of mean bias weighted by the average measured OA concentration. The resulting mean bias is −7.2 µg m⁻³ (∼ 15 %) and the root-mean-square error (RMSE) is 37.4 µg m⁻³. The model tends to underestimate the measured OA for Set 1 with a mean bias of 7 % while it tends to overestimate the measured OA for Set 2 with a mean bias of 8 %.

Figure 7 shows the modelled and measured OA mass for all 14 experiments, where Set 1 accounts for both OH and NO₃ chemistry while Set 2 includes OH chemistry only. The modelled OA is divided into POA and six SOA classes attributed to the respective precursor classes. Overall, the model performance is satisfactory, although in general the final OA concentration is slightly overpredicted while the initial production rate is underpredicted. Most of the SOA is attributed to furans (30.8 %), OxyAH (19 %) and OVOC_c < 6 (12.5 %) for Set 1 while for Set 2 there is a generally lower contribution from OVOC_c ≥ 6 (7.8 %) and a higher contribution from PAH (12 %), especially for experiments 08, 13 and 14. The mean bias between measured and modelled OA averaged over all experiments is −7.2 µg m⁻³, which corresponds to ∼ 15  % on average.

Figure 8 (a) Comparison between modelled and measured O : C ratio as a function of the OH exposure for each experiment. Set 1 (Expt. 1–7, 2 ^∘C) and Set 2 (Expt. 8–9, −10 ^∘C and Expt. 10–14, 15 ^∘C). (b) Probability distribution of the relative bias (normalized by the averaged measured O : C ratios). The resulting mean relative bias is 0.006 and the root-mean-square error (RMSE) is 0.06.

4.3 Investigation of OA chemical and physical properties

Comparisons between measured and modelled O : C ratios are reported in Fig. 8. The oxidation products' elemental compositions from which the modelled O : C ratio is calculated are presented in Table S4. The oxidation products' carbon numbers, which best explained the observed O : C ratio, correspond to a set ΔC of 0.6 (Fig. S7). There is a general increase in the O : C ratio with time. Model and observations match in terms of average O : C ratio for each experiment but the temporal evolution of the ratio is not well predicted, suggesting that there are additional processes that are not taken into account in the model. As the model was initiated using the measured POA O : C ratio at OH exposure equal to zero, an agreement between model and measurements can be observed at this time of the experiment. We do not find any systematic correlation of the bias with chamber conditions except for lower concentrations where experiments exhibit higher O : C ratios at the end. Overall for experiments conducted at lower temperature (−10 ^∘C) the model tends to overestimate the O : C while for higher-temperature experiments (15 ^∘C) the model clearly under-predicts the ratio upon ageing, especially at the end of the experiments, indicating the presence of compounds with a higher number of oxygen molecules (lower number of carbon molecules) than predicted.

Figure 9 shows mass yield curves for each class of compounds, in comparison with mass yields of several single compounds from literature (Table S2). The single compounds were selected according to their presence in the current study and compared to the respective chemical family. The reported published yields were discriminated according to experimental NO_x regimes, as low-NO_x conditions generally lead to higher OA yields. Our study is more likely representative of a high-NO_x regime and thus assesses the SOA forming potential for this atmospherically relevant condition. For the model and measurement conditions (C_OA∼20–600 µg m⁻³) the following median mass yields ranges were found: 7 %–20 % for furans, 10 %–25 % for SAH, 14 %–32 % for PAH, 9 %–24 % for OxyAH, 23 %–46 % for OVOC_c ≥ 6 and 6 %–18 % for OVOC_c < 6. Mass yields were also calculated for the two separate datasets and reported in Fig. S6 in order to assess differences driven by specific combustion regimes. Considering ambient-relevant conditions of C_OA∼50 µg m⁻³, we note a generally good agreement between the two datasets except for OxyAH, which exhibit higher median yields (representing smoldering phase) of ∼ 25 % in Set 1 and only ∼ 10 % in Set 2 (representing flaming phase). Set 2 on the other hand exhibited higher median yields for the OVOC_c ≥ 6 family (∼ 20 % compared to ∼ 16 % for Set 1).

The effect of temperature and OH exposure on OA concentrations and yields is shown in Fig. 10 for different primary OM loads (total primary gaseous and particulate OM_p + OGs) of 6, 60 and 600 µg m⁻³. The range of temperatures investigated varies between 255 and 315 K. We find a general increase in total OA concentration with increasing OH exposure, decreasing experimental temperature and higher initial loads, as expected. The average increase in OA concentration is 0.001, 0.03 and 0.6 µg m⁻³ K⁻¹ for 6, 60 and 600 µg m⁻³, respectively. Concerning SOA yields the temperature effect is also a function of OH exposure and aerosol load; SOA yields increase by 0.0001, 0.0006 and 0.002 g g⁻¹ K⁻¹ on average for 6, 60 and 600 µg m⁻³, respectively, with a higher effect predicted at lower temperature. We note overall an average yield increase by a factor of 3–4 for a 10-fold increase in the primary OM loads at the highest OH exposure considered (8×10⁷ molec. cm⁻³ h). Set 2 exhibits higher yields because of lower contributions from the OVOC_c < 6 family, which does not produce significant amounts of SOA.

Figure 9 Mass yields for each class of compounds. The solid lines represent the median values while the lower and upper limits are the 25th and 75th percentiles, respectively. The different markers in each plot are yields published in the literature for different single compounds (see Table S2). The colour code denotes different NO_x regimes (red denoting high NO_x, blue low or no NO_x and grey not specified NO_x regimes). The shaded background represents the range of our experiments (20–600 µg m⁻³); outside this shaded area yields are extrapolated from the model.

Figure 10 Modelled OA concentrations (a) and yields (b) under different temperature, OH exposure and initial OA concentrations (6, 60 and 600 µg m⁻³). Upper and lower panels are based on Set 1 and Set 2, respectively. Temperature is provided in kelvin (K) to avoid confusion with the experimental data in degrees Celsius (^∘C).

Figure 11 Fractional contributions of the six precursor classes to SOA formation under different temperature, OH exposure and initial OA concentrations. Left and right panels are based on Set 1 and Set 2, respectively. Temperature is provided in kelvin (K) to avoid confusion with the experimental data in degrees Celsius (^∘C).

SOA yields increase with increasing C_OA due to additional partitioning but also due to changes in the chemical composition and volatility of SOA species since they age differently with different experimental temperature and concentrations. Compounds with different oxygen-to-carbon ratios lead to different functionality, polarity and vapour pressure upon ageing. Moreover, different temperatures result in different evaporation enthalpies which consequently influence the compounds' volatility and lifetime. The modelled SOA ΔH_vap for each family of precursors is 17.5 kJ mol⁻¹ after GA calculation. This value is within the ranges of values reported in literature, where values between 11 and 44 kJ mol⁻¹ were reported for biogenic and anthropogenic SOA precursors depending upon the reactant hydrocarbon mixture and NO_x concentration. For SOA formed from both α-pinene and toluene, a negative correlation between ΔH_vap and NO_x concentration was observed (Offenberg et al., 2006).

The modelled fractional contributions of the six different precursor classes to SOA are shown in Fig. 11 for Set 1 and Set 2. The most dominant contribution is from the OxyAH family apart from Set 1 at high OH exposure where the contributions from precursors with higher volatility (furans and OVOC_c < 6) are more strongly temperature-dependent. In detail, the OVOC_c < 6 family exhibits a higher contribution with higher initial load and higher OH exposure. The SAH and PAH families have relevant contributions for Set 2 because these compounds are strongly emitted during the flaming phase. Since the SAH family consists of compounds that are less reactive than other families, they become relevant just at high OH exposure, while the compounds in the PAH family react faster and show a decreasing contribution with increasing OH exposure.

5 Conclusions

We performed box model simulations, based on the volatility basis set (VBS) approach, of residential wood combustion smog chamber experiments conducted at different temperatures, different combustion conditions and using different residential stoves. Primary emissions of SOA precursor compounds (organic gases, OGs) and organic aerosol (OA) as well as their evolution during ageing in the smog chamber were simultaneously monitored by a PTR-ToF-MS and an HR-ToF-AMS, respectively. This enabled the identification of the nature of SOA precursors lumped into different classes according to their chemical composition.

The knowledge about the nature of SOA precursors was used to better constrain model parameters, in the oxidation products' production rates and elemental composition. Using the measured OA mass, we were able to determine the volatility distributions and ΔH_vap for the products formed from the oxidation of the dominant precursor compound classes. We estimated the contributions of different compound classes to SOA and evaluate how the variability in the emission composition under a wide range of conditions would influence the SOA yield predictions. Investigation of different experimental temperatures allowed the evaluation of the model evaporation enthalpies, which have a decisive influence on the volatility of the emissions and hence their atmospheric lifetimes. Upon ageing, compounds with lower atomic O : C ratios are converted through the oxidation pathway to products with higher functionality, higher polarity and lower saturation vapour pressure. As a result, a part of these products (re-)condenses to the particle phase with partial pressures determined by their volatility, ambient temperature and concentration of the particulate organic mixture. While the degree of oxygenation increases during ageing, organic species may also fragment into more volatile compounds, being eventually converted into CO₂. Understanding the balance between oxygenation and fragmentation, their effect on volatility of emissions and timescale of these processes is essential to predict the evolution of the OA concentration.

Overall we developed a framework useful to constrain complex emissions and suitable for sophisticated mass spectrometry analysis with the novelty and ability of identifying the contributions of different classes of OG precursors to SOA formation. The main focus of the study included the investigation of smoldering vs. flaming emissions, resulting in predominant contributions of different classes of compounds according to the combustion phase investigated. Smoldering-phase emissions were dominated by the OVOC_c < 6 compound family while the flaming phase exhibited higher contributions by the SAH and PAH families. For both phases, SOA formation is found to be dominated by OxyAH (e.g. phenols and cresols), emitted from lignin pyrolysis. These species were therefore predicted to be important markers to be monitored in air pollution studies in order to estimate the SOA forming potential from real emissions.

Appendix A

In this appendix, we present the derivation of the thermodynamic equations used in the model. We considered the species in the gas and the particle phases to be permanently at equilibrium, as condensation is expected to be faster than oxidation (timescales for oxidized vapour condensation < 1 min assuming no particle-phase diffusion limitations; Bertrand et al., 2018b). Accordingly, the relation between the aerosol and gas-phase activities satisfies the following expression:

$$\begin{array}{}\text{(A1)}& x_{i,j}{\mathrm{|}}_{\text{g}}={\mathit{\gamma }}_{i,j}{\mathit{\chi }}_{i,j}{x}_{i,j}^{\mathrm{0}}.\end{array}$$

Here, x_i,j|_g denotes the gas-phase molar concentration of a surrogate in a volatility bin i formed from a precursor class j. The product ${\mathit{\gamma }}_{i,j}{\mathit{\chi }}_{i,j}$ represents the activity of the same surrogate in the particle phase, where γ_i,j and χ_i,j are the activity coefficient and the fraction of the surrogate in the particle phase, respectively. $x_{i,j}^{\mathrm{0}}$ is the equilibrium molar saturation concentration of the pure surrogate, related to its equilibrium vapour pressure, $p_{i,j}^{\mathrm{0}}$, according to Eq. (A2).

$$\begin{array}{}\text{(A2)}& x_{i,j}^{\mathrm{0}}={\displaystyle \frac{p_{i,j}^{\mathrm{0}}}{RT}}\end{array}$$

Here, R and T are the ideal gas constant and the temperature, respectively. The effective molar saturation concentration ($x_{i,j}^{*}$), which takes into account the influence of nonideal mixing on the compounds' activity, can be defined as the product of γ_i,j and $x_{i,j}^{\mathrm{0}}$:

$$\begin{array}{}\text{(A3)}& x_{i,j}^{*}={\mathit{\gamma }}_{i,j}{x}_{i,j}^{\mathrm{0}}.\end{array}$$

Similar to $p_{i,j}^{\mathrm{0}}$, $x_{i,j}^{*}$ can be written as a function of temperature, according to the Clausius–Clapeyron relationship based on Eq. (A4).

$$\begin{array}{}\text{(A4)}& \ln \left({\displaystyle \frac{x_{i,j}^{*}\left(T\right)}{x_{i,j}^{*}\left(T_{\text{ref}}\right)}}\right)=-{\displaystyle \frac{\mathrm{\Delta }{H}_{\text{vap},i,j}}{R}}\left({\displaystyle \frac{\mathrm{1}}{T}}-{\displaystyle \frac{\mathrm{1}}{T_{\text{ref}}}}\right)\end{array}$$

Here, T and T_ref are the experimental and reference (T_ref=298 K) temperatures, respectively. $\mathrm{\Delta }{H}_{\text{vap},i,j}$ is the effective enthalpy of evaporation. It includes the effects of temperature on (1) the pure compound saturation vapour pressure ($p_{i,j}^{\mathrm{0}}$), (2) the compounds' mixing properties in the condensed phase, i.e. γ_i,j, and (3) the radical chemistry reaction rate constants and branching ratios in the gas phase (Stolzenburg et al., 2018).

$x_{i,j}^{*}$ is related to the Donahue effective mass saturation concentration (Donahue et al., 2012), $C_{i,j}^{*}$, which is the inverse of the Pankow equilibrium constant (Pankow, 1987), through the compounds' molecular weight, MW_i,j, as indicated in Eq. (A5).

$$\begin{array}{}\text{(A5)}& C_{i,j}^{*}={\text{MW}}_{i,j}{x}_{i,j}^{*}\end{array}$$

To constrain the model to the measurements, it is of convenience to rearrange Eqs. (A1) and (A3) as a function of the surrogate total concentration, $x_{i,j}{\mathrm{|}}_{\text{g}+\text{p}}=x_{i,j}{\mathrm{|}}_{\text{g}}+x_{i,j}{\mathrm{|}}_{\text{p}}$, and the total molar concentration of species in the particle phase (x_OA), which yields the following expression:

$$\begin{array}{}\text{(A6)}& {\left.x_{i,j}\right|}_{\text{p}}={\left.{\left(\mathrm{1}+{\displaystyle \frac{x_{i,j}^{*}}{x_{\text{OA}}}}\right)}^{-\mathrm{1}}{x}_{i,j}\right|}_{\text{g}+\text{p}}.\end{array}$$

The parameters in Eq. (A6) were obtained as follows.

• The modelled molar concentration of the particulate organic phase was expressed as the sum of the concentration of all surrogates in the particle phase Eq. (A7).

  $$\begin{array}{}\text{(A7)}& x_{\text{OA}}=\sum _j^n{\sum }_i^m{x}_{i,j}{\left|{}_{\text{p}}+x_i^{{\text{OM}}_{\text{p}}}\right|}_{\text{p}}\end{array}$$

  Here, m and n are the total number of volatility bins and precursor chemical classes, respectively. $x_i^{{\text{OM}}_{\text{p}}}{\mathrm{|}}_{\text{p}}$ is the particle-phase molar concentration of primary organic matter in volatility bin i. The $x_i^{{\text{OM}}_{\text{p}}}{\mathrm{|}}_{\text{p}}$ is calculated based on $x_i^{{\text{OM}}_{\text{p}}}{\mathrm{|}}_{\text{g}+\text{p}}$ in both phases, following a similar computation as for SOA (Eq. A6). $x_i^{{\text{OM}}_{\text{p}}}{\mathrm{|}}_{\text{g}+\text{p}}$ was inferred from the measured POA concentrations injected into the chamber at the beginning of the experiment and using the volatility distribution function of wood combustion emissions in May et al. (2013).

• $x_{i,j}{\mathrm{|}}_{\text{g}+\text{p}}$ was derived from the precursor oxidation rates measured by the PTR-ToF-MS. The change in the total concentration of oxidation products from a precursor class j in both gas and particle phases, x_j|_g + p, was expressed as follows:

  $$\begin{array}{}\text{(A8)}& {\displaystyle }\begin{array}{rl}{\displaystyle \frac{\text{d}{x}_j{\mathrm{|}}_{\text{g}+\text{p}}}{\text{d}t}}& =k_{\mathrm{OH}}\times \left[\mathrm{OH}\right]\times \left[{\mathrm{OG}}_j\right]\\ & +k_{\text{other}}\times \left[{\mathrm{OG}}_j\right]-k_{\text{dil}}\times x_j{\mathrm{|}}_{\text{g}+\text{p}}\end{array}.\end{array}$$

  Here, [OG_j] is the molar concentration of total OG precursors in class j. For the definition of the other parameters, the reader is referred to Eq. (1). x_j|_g + p is calculated by numerically integrating Eq. (A8). x_j|_g + p is related to the concentration of the different surrogates from a precursor class j in different volatility bins, $x_{i,j}{\mathrm{|}}_{\text{g}+\text{p}}$ (Eq. A6), through their yields, Υ_i,j:

  $$\begin{array}{}\text{(A9)}& {\left.x_{i,j}\right|}_{\text{g}+\text{p}}={\left.{\mathrm{{\rm Y}}}_{i,j}{x}_j\right|}_{\text{g}+\text{p}}.\end{array}$$

These yields, which represent the surrogate volatility distributions, were determined by the model.