Comparison of saturation vapor pressures of α-pinene + O₃ oxidation products derived from COSMO-RS computations and thermal desorption experiments

Accurate information on gas-to-particle partitioning is needed to model secondary organic aerosol formation. However, determining reliable saturation vapor pressures of atmospherically relevant multifunctional organic compounds is extremely difficult. We estimated saturation vapor pressures of α-pinene-ozonolysis-derived secondary organic aerosol constituents using Filter Inlet for Gases and AEROsols (FIGAERO)–chemical ionization mass spectrometer (CIMS) experiments and conductor-like screening model for real solvents (COSMO-RS). We found a good agreement between experimental and computational saturation vapor pressures for molecules with molar masses around 190 g mol⁻¹ and higher, most within a factor of 3 comparing the average of the experimental vapor pressures and the COSMO-RS estimate of the isomer closest to the experiments. Smaller molecules likely have saturation vapor pressures that are too high to be measured using our experimental setup. The molecules with molar masses below 190 g mol⁻¹ that have differences of several orders of magnitude between the computational and experimental saturation vapor pressures observed in our experiments are likely products of thermal decomposition occurring during thermal desorption. For example, dehydration and decarboxylation reactions are able to explain some of the discrepancies between experimental and computational saturation vapor pressures. Based on our estimates, FIGAERO–CIMS can best be used to determine saturation vapor pressures of compounds with low and extremely low volatilities at least down to 10⁻¹⁰ Pa in saturation vapor pressure.

1 Introduction

Secondary organic aerosol (SOA) is formed in the gas phase by the condensing of organic molecules with low volatilities. In atmospheric science, organic compounds are often grouped based on their saturation vapor pressures into volatile organic compounds (VOCs), intermediate-volatility organic compounds (IVOCs), semi-volatile organic compounds (SVOCs), low-volatility organic compounds (LVOCs), extremely low volatility organic compounds (ELVOCs) and ultra-low-volatility organic compounds (ULVOCs) (Donahue et al., 2012; Schervish and Donahue, 2020). In the ambient air, ULVOCs can nucleate to initiate SOA formation (Kirkby et al., 2016; Bianchi et al., 2016), while ELVOCs, LVOCs and SVOCs can condense on existing particles to contribute to the growth of SOA (Ehn et al., 2014). A large source of organic compounds in the atmosphere comprises biogenic volatile organic compounds (BVOCs) emitted by plants (Jimenez et al., 2009; Hallquist et al., 2009). These BVOCs are oxidized in the gas phase by oxidants, such as OH, O₃ and NO₃, to form less volatile compounds through addition of oxygen-containing functional groups. In order to determine the role of different oxidation products in SOA formation, it is essential to have reliable methods to estimate the volatility of complex organic molecules formed in the atmosphere.

Monoterpenes (C₁₀H₁₆) are an abundant class of BVOCs emitted by various plants (Guenther et al., 1995). The oxidation mechanisms of monoterpene reactants vary significantly, leading to products with different SOA formation capabilities (Thomsen et al., 2021). Additionally, the initial oxidant (i.e., OH, O₃ and NO₃) affects the SOA formation rates of the oxidation products (Kurtén et al., 2017). For example, α-pinene, a very abundant monoterpene in the atmosphere, has been widely studied in both laboratory and field experiments (Docherty et al., 2005; Hall and Johnston, 2011; Hao et al., 2011; Ehn et al., 2012, 2014; Kristensen et al., 2014; Lopez-Hilfiker et al., 2015; McVay et al., 2016; D'Ambro et al., 2018; Huang et al., 2018; Claflin et al., 2018; Ye et al., 2019), and α-pinene oxidation, with both O₃ and OH, is efficient at producing oxygenated organic molecules and SOA. Consequently, the oxidation mechanism and potential structures of α-pinene + O₃ products have been extensively studied both experimentally and computationally (Rissanen et al., 2015; Berndt et al., 2018; Kurtén et al., 2015; Iyer et al., 2021; Lignell et al., 2013; Aljawhary et al., 2016; Mutzel et al., 2016; Kristensen et al., 2020; Thomsen et al., 2021). Various methods have been used for estimating the saturation vapor pressures or saturation mass concentrations of α-pinene ozonolysis products (Ehn et al., 2014; Kurtén et al., 2016; D'Ambro et al., 2018; Buchholz et al., 2019; Peräkylä et al., 2020; Räty et al., 2021). However, measuring these saturation vapor pressures accurately is extremely difficult (Seinfeld and Pankow, 2003). Additionally, different experimental and theoretical methods are known to produce very different saturation vapor pressures (Bilde et al., 2015; Kurtén et al., 2016; Bannan et al., 2017; Wania et al., 2017; Ylisirniö et al., 2021). For example, the agreement between different experiments is better measuring at a subcooled state compared to a solid state, perhaps due to ambiguity of the physical state of the solid-state samples (Bilde et al., 2015).

During recent years, saturation vapor pressures of atmospherically relevant multifunctional organics have been derived from their desorption temperatures using mass spectrometers equipped with the Filter Inlet for Gases and AEROsols (FIGAERO; Lopez-Hilfiker et al., 2014). Using this method, the saturation vapor pressures can be estimated from the desorption temperatures of the molecules. However, the measurements need to be calibrated using compounds with known saturation vapor pressures in order to find the correlation between saturation vapor pressure and desorption temperature. For example, a recent experimental study highlighted how different sample preparation methods affect measured desorption temperatures of FIGAERO calibration experiments (Ylisirniö et al., 2021), and aerosol particle size and operational parameters generally affect the measurement results as well (Schobesberger et al., 2018; Thornton et al., 2020). Additionally, when the desorbed molecules are detected using a chemical ionization mass spectrometer (CIMS), only the elemental compositions are obtained without any information on the chemical structures, based on which saturation vapor pressures could be further constrained.

Another possible complication in thermal desorption experiments is thermal decomposition reactions during the heating of the sample. For example, Stark et al. (2017) studied the effect of thermal decomposition on the determination of volatility distributions from FIGAERO–CIMS measurements. They concluded that most of the condensed-phase species decompose during thermal desorption experiments, in agreement with several other studies of laboratory and ambient FIGAERO measurements (Lopez-Hilfiker et al., 2015, 2016; Schobesberger et al., 2018). Most recently, Yang et al. (2021) found that decarboxylation and dehydration reactions are significant in FIGAERO measurements for multifunctional carboxylic acids that have more than four oxygen atoms, a degree of unsaturation between 2 and 4, and a maximum desorption temperature (T_max) higher than 345 K.

Among theoretical models, the conductor-like screening model for real solvents (COSMO-RS; Klamt, 1995; Klamt et al., 1998; Eckert and Klamt, 2002) has been seen as the most promising method for calculating partitioning properties because it does not require calibration, unlike group-contribution methods (Wania et al., 2014). During recent years, this quantum-chemistry-based method has been used to estimate saturation vapor pressures of atmospherically relevant multifunctional compounds (Wania et al., 2014, 2015; Kurtén et al., 2016; Wang et al., 2017; Krieger et al., 2018; Kurtén et al., 2018; D'Ambro et al., 2019; Hyttinen et al., 2020, 2021b). For example, Kurtén et al. (2016) compared COSMO-RS-derived saturation vapor pressures of 16 α-pinene ozonolysis products with those estimated with various group-contribution methods. They found that COSMO-RS (parametrization BP_TZVPD_FINE_C30_1501 implemented in the COSMOtherm program; COSMOtherm, 2015) predicts saturation vapor pressures of up to 8 orders of magnitude higher than group-contribution methods, such as EVAPORATION (Compernolle et al., 2011) and SIMPOL.1 (Pankow and Asher, 2008). The COSMOtherm15-estimated saturation vapor pressures indicated that the studied highly oxidized monomers derived from the ozonolysis of α-pinene were likely classified as SVOCs with saturation vapor pressures higher than 10⁻⁵ Pa (Kurtén et al., 2016). However, the parametrization in COSMOtherm has a large effect on the calculated properties since the model is parametrized using a set of well-known compounds with experimental properties available. There have been significant improvements since the BP_TZVPD_FINE_C30_1501 parametrization used by Kurtén et al. (2016), especially with better description of the effect of hydrogen bonding on thermodynamic properties. This is an important factor in calculating properties of multifunctional compounds that are able to form intramolecular hydrogen bonds (H-bonds). For example, Hyttinen et al. (2021b) found that with an improved conformer sampling method (recommended by Kurtén et al., 2018) and a newer parametrization (BP_TZVPD_FINE_19), COSMOtherm-estimated saturation vapor pressures of the two most highly oxygenated α-pinene ozonolysis monomer products studied by Kurtén et al. (2016) are up to 2 orders of magnitude lower than SIMPOL.1 estimates, while COSMOtherm predicted higher saturation vapor pressures than SIMPOL.1 for 15 different α-pinene + OH-derived dimers.

In this study, we investigate the saturation vapor pressures of SOA constituents formed in α-pinene ozonolysis, using both FIGAERO–CIMS experiments and the COSMO-RS theory. We are especially interested in whether the calibration done using compounds with saturation vapor pressures limited to the LVOC–SVOC range is valid for estimating saturation vapor pressures of ELVOCs and ULVOCs. We compare saturation vapor pressures derived from both experiments and calculations (different isomers) in order to evaluate the experimental method. Additionally, we investigate the prevalence of thermal decomposition in our experiment.

2 Methods

2.1 Chamber experiments

The experiments were conducted at a 9 m³ Teflon environmental reaction chamber. The chamber is located at the University of Eastern Finland (Kuopio, Finland). During the experiment, the chamber was operated as a batch reactor; i.e., the experimental conditions were set at the start of the experiment, and after the chemistry was initiated, the proceeding changes in the gas and particle phase in the closed system were sampled. The chamber is set on a foldable frame which allows the chamber to collapse when deflated, maintaining a constant pressure. The chamber and the instruments were situated inside a temperature-controlled environment (temperature set to 295.15 K). Before the experiment, the chamber was flushed overnight with dry, clean air to reduce the impact of evaporation of residues from preceding experiments from the walls.

To prepare the chamber for the experiment, it was first filled with clean air, which was sampled by a proton-transfer-reaction time-of-flight mass spectrometer (PTR-ToF-MS, Ionicon, Inc.), and a Filter Inlet for Gases and AEROsols (FIGAERO) coupled with a time-of-flight chemical ionization mass spectrometer (ToF CIMS) to determine the chamber background. The next section will provide a more thorough description of the instruments. After the chamber was filled close to operational capacity (9 m³), α-pinene was introduced into the chamber. This was done by flushing dry purified air through an α-pinene diffusion source and into the chamber until target concentration (11 ppb) was reached. α-Pinene levels were monitored with an online PTR-ToF-MS. Polydisperse ammonium sulfate seed aerosol (∼ 10 000 cm⁻³, maximum number concentration at ∼ 80 nm) was added to provide condensation nuclei and to prevent possible nucleation during the experiments. Lastly, 30 ppb of externally generated ozone (using an ozone generator with a UV lamp of wavelength 185 nm) was introduced into the reaction chamber to start the chemistry. The experiment duration was 8 h from when the chemistry started (ozone was added). There was practically no change in the chamber size during the experiment due to the low sampling flows compared to the total chamber volume.

2.1.1 Instrumentation

In this study, we analyzed particle-phase composition measurements performed with a Filter Inlet for Gases and AEROsols (FIGAERO) inlet system coupled with a time-of-flight chemical ionization mass spectrometer with iodide ionization (I-CIMS, Aerodyne Research Inc.), a system that allows for measurement of both gas-phase and particle-phase compounds with a single instrument (Lopez-Hilfiker et al., 2014, 2015; Ylisirniö et al., 2021). In the FIGAERO inlet, the aerosol particles are collected on a Teflon filter (Zefluor 2 µm PTFE membrane filter, Pall Corp.) while simultaneously analyzing the gas phase. After a predetermined collection time (here 45 min) is finished, the sampled particle matter is evaporated using a gradually heated nitrogen flow with a heating rate of 11.7 K min⁻¹ and the evaporated molecules are carried into the detector instrument I-CIMS. Integrating over the heating time will give the total signal of a particular compound in the sample being processed. The working principle of the I-CIMS has been introduced elsewhere (Lee et al., 2014; Iyer et al., 2017), but in short, oxidized gas-phase constituents are detected by clustering negatively charged iodide anions (I⁻) with suitable organic compounds. Clustering of the organic molecules and I⁻ happens in an ion molecule reaction chamber (IMR), which is actively controlled to be at 10⁴ Pa pressure.

The particle sampling period was set to 45 min, and the particle analysis period consisted of a 15 min ramping time (when the filter was heated linearly from room temperature to 473.15 K) and a 15 min soak period (where the filter temperature was kept at 473.15 K). Thus, there is 45 min of gas-phase measurements followed by a 30 min gap while particle chemical composition is being analyzed. Seven particle samples were collected during the 8 h SOA experiment.

2.1.2 Data analysis

All FIGAERO–CIMS data were preprocessed with tofTools (version 611) running in MATLAB R2019b (MATLAB, 2019) and further processed with custom MATLAB scripts. Saturation vapor pressures of the oxidized organics were estimated based on their thermograms, i.e., signal as a function of temperature along the heating of the particle sample in FIGAERO. We used 20 s averaging in the thermograms. The temperature axis calibration sample was made as described by Ylisirniö et al. (2021), by using an atomizer to produce a particle population with a similar size distribution to the one present in the chamber experiment. Polyethylene glycol (PEGn with n equal to 6, 7 and 8) with known saturation vapor pressures (see Fig. S1 and Table S1 in the Supplement) was used to produce the calibration particle population. Following the calibration fit, the saturation vapor pressure (p_sat in Pa) of a molecule can be calculated from the temperature of the highest signal (T_max in K):

$$\begin{array}{}\text{(1)}& p_{\text{sat}}=e^{-\mathrm{0.1594}\cdot T_{max}+\mathrm{40.13}}.\end{array}$$

Ylisirniö et al. (2021) found a good exponential correlation between the temperature of the highest signal and saturation vapor pressure ranging up to p_sat = 5 × 10⁻⁴ Pa (PEG5). However, like theirs, our calibration only reached down to 9 × 10⁻⁸ Pa in saturation vapor pressure, which introduces an additional source of uncertainty to the saturation vapor pressures estimated from the experiments. In addition to the linear correlation between T_max and log ₁₀p_sat, Ylisirniö et al. (2021) proposed a polynomial calibration curve, which leads to lower saturation vapor pressure estimates at higher desorption temperatures (T_max > 350 K). With our three calibration points, it is impossible to find a reliable polynomial fit to extrapolate to a higher T_max. Instead, we assume a similar difference between the two calibration curves to what was estimated by Ylisirniö et al. (2021). For example, using our linear fit, 392 K corresponds to 2 × 10⁻¹⁰ Pa, but in the polynomial fit, the same T_max corresponds to about 10⁻¹¹ Pa (see Table S2 and Fig. S2 for more values).

The variation in T_max values between three calibration runs varies from 0.5 K for the smallest, PEG6 (282.3 g mol⁻¹), to 7.6 K for the largest, PEG8 (370.4 g mol⁻¹). With our calibration curve, these differences correspond to a factor of 1.1 and 3.3 variation in the saturation vapor pressures, respectively. Saturation vapor pressures were calculated for multiple α-pinene-derived SOA constituents from six different samples, i.e., six different subsequent thermal desorptions, during the one 8 h experiment. The first sample of our experiment was omitted because the signals were much lower in the first sample than in the other samples. This was likely caused by lower concentrations of oxidation products in the chamber at the beginning of the experiment. In our experiment, the variation in T_max values between the different thermal desorption cycles ranged from 2.0 to 11.1 K. The variation in T_max values increases with the increasing molar mass (see Fig. S3). The 11.1 K variation corresponds to a factor of 5.8 variation in p_sat. Most of the studied compounds have saturation vapor pressures within a factor of 4 from the six measurement cycles.

We used desorption temperatures to estimate saturation vapor pressures even though the particle-to-gas partitioning in our experiment is also affected by the activity coefficient of the compound in the sample. For example, Ylisirniö et al. (2021) found a 5–7 K difference in the temperatures of maximum desorption signal between pure PEG and PEG-400 mixture (average molecular mass ∼ 400 g mol⁻¹), which they attributed to the additional compounds in the mixture. In the case of similar multifunctional compounds, the activity coefficients of individual compounds in the mixture (estimated using COSMOtherm) are likely to be close to unity, with respect to the pure compound reference state (the compound has similar chemical potentials in a pure state and in the mixture, which leads to activity coefficients close to 1 in COSMOtherm calculations). We therefore assume that the mixture is ideal and estimate saturation vapor pressures from desorption temperatures.

2.2 COSMOtherm calculations

Our experiments provided us with elemental compositions of compounds in our SOA sample and saturation vapor pressures corresponding to each composition. To compare with the experiments, we computed saturation vapor pressures of potential ozonolysis product structures corresponding to the measured elemental compositions using the COSMO-RS theory with the newest BP_TZVPD_FINE_21 parametrization, implemented in the COSMOtherm program Release 2021 (BIOVIA COSMOtherm, 2021). COSMO-RS uses statistical thermodynamics to predict properties of molecules in both condensed and gas phases. The interactions between molecules in the condensed phase are described using the partial charge surfaces of the molecules derived from quantum chemical calculations.

For our COSMO-RS calculations, we selected conformers containing no intramolecular H-bonds, detailed previously by Kurtén et al. (2018), Hyttinen and Prisle (2020), and Hyttinen et al. (2021b). This method has been shown to provide more reliable saturation vapor pressure estimates for multifunctional oxygenated organic compounds even if they are able to form intramolecular H-bonds (Kurtén et al., 2018). Additionally, Hyttinen and Prisle (2020) found that in COSMOtherm, conformers containing multiple intramolecular H-bonds are given high weights in the conformer distribution due to their low COSMO energies, even if conformers containing no intramolecular H-bonds would be more stable in the condensed phase. The conformer search was performed using the Merck molecular force field (MMFF94; Halgren, 1996) and the systematic algorithm (sparse systematic algorithm for isomers that have more than 100 000 possible conformers) of Spartan'14 (Wavefunction Inc., 2014). Instead of omitting conformers containing intramolecular H-bonds after running the quantum chemical calculations, as recommended by Kurtén et al. (2018), we removed conformers containing intramolecular H-bonds already after the initial conformer search step in order to decrease the number of density functional theory (DFT) calculations needed for the input file generation (see Sect. S1 in the Supplement for the details).

The quantum chemical single-point calculations and geometry optimizations were performed using the COSMOconf program (BIOVIA COSMOconf, 2021), which utilizes the TURBOMOLE program package (TURBOMOLE V7.4.1, 2019). First, single-point calculations at a low level of theory (BP/SV(P)-COSMO) were used to remove conformers with similar chemical potentials. After a geometry optimization at the same level, duplicate conformers with similar geometries and chemical potentials were omitted. Duplicates were also removed after the final optimization at the higher level of theory (BP/def-TZVP-COSMO), and final single-point energies were calculated for the remaining conformers at the highest level of theory available in the current COSMOtherm version (BP/def2-TZVPD-FINE-COSMO, currently only available for single-point calculations). Finally, the intramolecular hydrogen bonding of each remaining conformer was checked using the pr_steric keyword in COSMOtherm and up to 40 conformers containing no intramolecular H-bonds were selected for our saturation vapor pressure calculations (see Kurtén et al., 2018, for more details). The gas-phase energies of the selected conformers were obtained by optimizing the condensed-phase geometries and calculating single-point energies at the levels of theory corresponding to the COSMO calculations (BP/def-TZVP-GAS and BP/def2-TZVPD-GAS, respectively). Gas-phase conformers containing intramolecular H-bonds (formed in the gas-phase geometry optimization) were omitted, and the gas-phase single-point energy of the corresponding condensed-phase geometry was used instead. When all conformers found contained intramolecular H-bonds, conformers containing a single intramolecular H-bond were selected for the COSMO-RS calculation (see Sect. S1).

In COSMOtherm, the saturation vapor pressure (p_sat,i in mbar) of a compound is estimated using the free-energy difference of the compound in the pure condensed phase ($G_i^{\left(l\right)}$ in kcal mol⁻¹) and in the gas phase ($G_i^{\left(g\right)}$ in kcal mol⁻¹):

$$\begin{array}{}\text{(2)}& p_{\text{sat},i}=e^{-\left(G_i^{\left(l\right)}-G_i^{\left(g\right)}\right)/RT}.\end{array}$$

Here R is the gas constant (in $\mathrm{kcal}{\mathrm{K}}^{-\mathrm{1}}{\mathrm{mol}}^{-\mathrm{1}}$) and T is the temperature (in K).

3 Results and discussion

3.1 Saturation vapor pressures

We selected 26 elemental compositions (20 monomers and 6 dimers) from our FIGAERO–CIMS measurements for the comparison with COSMOtherm-estimated saturation vapor pressures. All elemental compositions that contain up to 10 carbon atoms are assumed to be monomers (containing carbon atoms only from the original reactant α-pinene), while compounds with 11–20 carbon atoms are assumed to be dimers (covalently bound accretion products of two monomers). For COSMOtherm analysis, we selected one to seven isomer structures that can be formed from α-pinene ozonolysis for each elemental composition. The degree of unsaturation of the studied monomers and dimers is 1–4 and 4–5, respectively, determining how many double bonds or ring structures each isomer must contain.

The structures of the studied monomers were formed based on structures suggested by previous experimental and computational studies (Lignell et al., 2013; Aljawhary et al., 2016; Mutzel et al., 2016; Kristensen et al., 2014; Kurtén et al., 2015; Iyer et al., 2021). These structures are shown in Figs. S4–S8. For dimer calculations, we selected elemental compositions that can be formed using the studied monomer structures, assuming a loss of H₂O₂, H₂O or no atoms from the original monomers. With a loss of H₂O₂ or H₂O, a dimer can be formed by a recombination of two hydroperoxy or hydroxy groups to form a peroxide or an ether. Additionally, if one or both of the monomers are carboxylic acids, the dimer contains an ester or an acid anhydride ($\mathrm{RC}(=\mathrm{O})\mathrm{OC}(=\mathrm{O}){\mathrm{R}}^{\prime }$) group, respectively. A dimer can also be formed in a condensed-phase reaction between a hydroxide and an aldehyde to form a hemiacetal (ROR^′OH). In hemiacetal formation, no atoms are lost from the reactant monomers. In order to reduce the number of computationally heavy dimer calculations, we selected only one pair of monomer isomers with the same elemental composition for each dimerization reaction. For most of the monomers used to form the studied dimers, the best agreement between experimental and computational saturation vapor pressures was found with the isomer that had the lowest COSMOtherm-estimated p_sat. We therefore mainly chose the monomer isomers with the lowest p_sat to form the studied dimer isomers. Table S4 shows which monomers were used to form each of the studied dimer isomers.

Figure 1 shows COSMOtherm-estimated saturation vapor pressures of the studied isomers, as well as vapor pressures derived from the experimental T_max values. The agreement between COSMOtherm-estimated and experimentally determined saturation vapor pressures is good for molar masses higher than 190 g mol⁻¹. Even with a limited selection of dimer structures, the agreement between COSMOtherm and FIGAERO–CIMS is very good, and even better agreement could likely be found by selecting additional dimer isomers for COSMOtherm calculations. Using a polynomial correlation between T_max and log ₁₀p_sat, the p_sat estimates of the studied monomers (highest T_max at 378 K) would likely decrease by 1 order of magnitude or less (see Fig. S2). With such a small decrease, all of the studied monomers would still be classified as LVOCs, with the exception of C₉H₁₈O₁₀, which would be classified as an ELVOC. The experimental saturation vapor pressures of the studied dimers (excluding C₁₈H₂₆O₆) would decrease by 1 to 2 orders of magnitude, which would improve the agreement between the experimental and calculated saturation vapor pressures. Until more accurate calibration of the FIGAERO–CIMS instrument becomes available, the experimental p_sat from the linear and polynomial fits can be used as upper and rough lower-limit estimates, respectively.

Figure 1 Saturation vapor pressures of the studied α-pinene-derived ozonolysis products as a function of their molar mass at 298.15 K. The colors represent the number of oxygen atoms in each isomer. The isomers shown with filled markers have COSMOtherm-estimated p_sat values closest to those of the experiments. The bars of the experimental values show the range of saturation vapor pressures from the six samples used, instead of error estimates of the measurements. Suspected thermal decomposition products are shown with the lighter gray color. The dashed lines indicate different volatilities using the classification of Donahue et al. (2012) and Schervish and Donahue (2020) and assuming ideality (γ = 1) in the conversion from mass concentration to vapor pressure.

The large discrepancy between the experimental and calculated p_sat of the lowest-molar-mass molecules (light gray bars in Fig. 1) suggests that the measured T_max values are related to the thermal decomposition temperatures of larger compounds, rather than to saturation vapor pressures of the measured elemental compositions. It is unlikely that COSMOtherm would overestimate saturation vapor pressures by several orders of magnitude using the newest parametrization and improved conformer selection (Kurtén et al., 2018). Additionally, if the low-molar-mass compounds are IVOCs (p_sat > 10⁻² Pa), as predicted by COSMOtherm, they are not likely to contribute to the SOA formation. Conversely, the calibration curve sets a practical upper limit to experimentally derivable p_sat based on the experiment temperature and premature evaporation. For example, the upper limit p_sat corresponding to the initial temperature of the experiment (T_max = 294.15 K) is 1.2 × 10⁻³ Pa. However, the highest experimental saturation vapor pressure among the studied molecules is 8.5 × 10⁻⁶ Pa, which corresponds to $T_{max}=\mathrm{325}$ K. This may indicate that the SOA constituents selected for our analysis do not contain SVOCs and the selected elemental compositions corresponding to SVOCs in the experiments were in fact thermal decomposition products rather than oxidation products of α-pinene ozonolysis.

Both the computational and experimental p_sat values correlate with molar mass, the O:C ratio having a smaller effect on p_sat. In addition to molar mass, saturation vapor pressure is known to depend on the functional groups of the molecule, as is seen in the difference of several orders of magnitude in the COSMOtherm estimates of different isomers at the same elemental compositions. We have previously noted that the addition of a CH₂ (∼ 14 g mol⁻¹) to a multifunctional molecule has an effect of 0.5 orders of magnitude on saturation vapor pressure, and the addition of an oxygen atom (∼ 16 g mol⁻¹) similarly has an effect of 0.5–1 orders of magnitude on saturation vapor pressure depending on the functional group (Hyttinen et al., 2021b). This means that addition of an oxygen atom may decrease the saturation vapor pressure (in Pa per g mol⁻¹) either less or more than the addition of a CH₂ group depending on the oxygen-containing functional group. The COSMOtherm-estimated saturation vapor pressures can also vary by more than an order of magnitude for different stereoisomers with identical functional groups (Kurtén et al., 2018).

Based on our FIGAERO–CIMS measurements, the studied monomer products derived from α-pinene ozonolysis present in the SOA are LVOCs or ELVOCs, while the studied dimers are mainly ELVOCs. We would like to note that this does not reflect the composition of α-pinene ozonolysis SOA but simply represents the set of elemental compositions selected for the analysis. The lowest experimental p_sat among the studied elemental compositions is at 5.4 × 10⁻¹¹ Pa. The saturation vapor pressure corresponding to the upper limit temperature of our experiment ($T_{max}=\mathrm{473.15}$ K) is 4.7 × 10⁻¹⁶ Pa (linear calibration curve), which means that saturation vapor pressures below 4.7 × 10⁻¹⁶ Pa cannot be estimated in our experiments. It is also possible that saturation vapor pressures of dimers with the lowest volatilities (p_sat < 10⁻¹¹ Pa) cannot be estimated using thermal desorption, as the molecules would thermally decompose before evaporating from the sample (Yang et al., 2021).

3.2 Correlation between monomer and dimer vapor pressures

The COSMOtherm calculations of dimers are computationally more demanding than those of monomers, due to a larger size and higher number of possible conformers. In group-contribution methods, such as SIMPOL.1, the saturation vapor pressure of a compound is estimated as the sum of contributions of each of the functional groups in the molecule:

$$\begin{array}{}\text{(3)}& {\log }_{\mathrm{10}}{p}_{\text{sat},i}={\mathrm{\Sigma }}_k{\mathit{\nu }}_{k,i}{b}_k,\end{array}$$

where ν_k,i is the number of functional groups of type k in compound i.

We used the same approach to estimate the saturation vapor pressures of dimers and compared those values to saturation vapor pressures estimated using COSMOtherm. However, instead of using the functional groups of the dimer, we used the contributions of the two monomers that formed the dimer. This way, the group-contribution term b_k was replaced by COSMOtherm-estimated saturation vapor pressures of the monomers multiplied with a scaling factor (S_n) to account for the changing functional groups and loss of atoms in dimerization reaction n.

$$\begin{array}{}\text{(4)}& p_{\text{sat},\text{dimer}}=S_n{p}_{\text{sat},\text{monomer1}}{p}_{\text{sat},\text{monomer2}}\end{array}$$

Many of the monomer isomers had to be altered slightly to accommodate the chosen dimerization reactions, which makes a direct comparison between monomers and dimers impossible. We therefore only investigate the acid anhydride formation, for which we have the most dimers to compare. For this comparison we formed additional C₁₃H₁₈O₉ dimers from all carboxylic acid isomers of C₅H₈O₆ and C₈H₁₂O₄ in order to better test the effect of molar mass on S.

Figure 2 shows the correlation between the COSMOtherm-estimated saturation vapor pressures of dimers and those of the monomers that were used to form the dimers. We see that the product of monomer vapor pressures is 1–3 orders of magnitude higher than the dimer vapor pressure. There is also a size dependence in the scaling factor; the values of S as a function of dimer size are shown in Fig. S14. Of the studied acid anhydride dimers, C₁₃H₁₈O₉ (the smallest dimer) isomers have scaling factors of 10⁻³–10⁻² and C₁₇H₂₄O₁₀ (the largest dimer) isomers of 10⁻²–10⁻¹. As a comparison, SIMPOL.1 predicts S = 1.1 × 10⁻² for the acid anhydride (ketone and ester) formation from two carboxylic acid monomers, with no size dependence.

Figure 2 Correlation between COSMOtherm-estimated saturation vapor pressures of the studied dimers (p_sat) and the product monomer vapor pressures (p_sat,monomer) at 298.15 K. The dimers are ordered from the smallest molar mass to the highest. The deviation from the 1:1 line represents the S value of Eq. (4).

The correlation between COSMOtherm-estimated saturation vapor pressures of monomers and dimers can be used to obtain rough saturation vapor pressures estimates of a larger number of dimer compounds by computing only the saturation vapor pressures of their constituent monomers. This reduces the computational cost of dimer calculations since the number of conformers and the calculation times increase exponentially with the size of the molecule. For example, the COSMOtherm-estimated saturation vapor pressures of the studied dimer with the highest molar mass (C₁₇H₂₄O₁₀) are much lower than the experimental ones, with a difference of around a factor of 2 between the highest COSMOtherm estimate and the lowest experimental value. Assuming that C₁₇H₂₄O₁₀ is an acid anhydride formed from C₈H₁₂O₄ and C₉H₁₄O₇, trying all combinations of the studied carboxylic acid isomers gives a p_sat range of 5.6 × 10⁻¹⁴–7.0 × 10⁻⁹ Pa using Eq. (4) and S = (1.0–9.5) × 10⁻². This range overlaps with the experimental range of 8.1 × 10⁻¹¹–2.7 × 10⁻¹⁰ Pa (see Fig. 1).

3.3 Thermal decomposition

A recent study by Yang et al. (2021) proposed two major thermal decomposition pathways for multifunctional carboxylic acids occurring in FIGAERO–CIMS: dehydration (Reaction R1) and decarboxylation (Reaction R2).

$$\begin{array}{}\text{(R1)}& {\displaystyle }\mathrm{HOC}(=\mathrm{O})\mathrm{RC}(=\mathrm{O})\mathrm{OH}\to \mathrm{R}\mathrm{1}\mathrm{C}(=\mathrm{O})\mathrm{OC}(=\mathrm{O})\mathrm{R}\mathrm{1}+{\mathrm{H}}_{\mathrm{2}}\mathrm{O}\text{(R2)}& {\displaystyle }\mathrm{RC}(=\mathrm{O}){\mathrm{CH}}_{\mathrm{2}}\mathrm{C}(=\mathrm{O})\mathrm{OH}\to \mathrm{RC}(=\mathrm{O}){\mathrm{CH}}_{\mathrm{3}}+{\mathrm{CO}}_{\mathrm{2}}\end{array}$$

We selected four elemental compositions (C₇H₁₀O₄, C₇H₁₀O₆, C₈H₁₂O₄ and C₈H₁₂O₆) to investigate the two possible thermal decomposition reactions. The different isomers of C₇H₁₀O₄, C₇H₁₀O₆, C₈H₁₂O₄ and C₈H₁₂O₆ and their thermal decomposition reactants are shown in Figs. S11 and S12. These elemental compositions were selected because their experimental and COSMOtherm-estimated saturation vapor pressures had large differences. Additionally, their elemental compositions are possible products of both of the studied reactions. Both reactions are possible only if the product contains fewer than 10 carbon atoms (the monomer reactant of the dehydration reaction can contain up to 10 carbon atoms) and the degree of unsaturation is at least 3 (the product of dehydration contains at least one ring structure and two double bonds). The thermal decomposition reactants also fulfill the number-of-oxygen and degree-of-unsaturation criteria given by Yang et al. (2021). Other likely decomposition products among the studied monomers are C₄H₄O₅, C₄H₄O₆, C₉H₁₂O₃, C₁₀H₁₆O₃ and C₉H₁₄O₄ (see Fig. 1).

Figure 3 Saturation vapor pressures of potential products (red markers) and the reactants (blue markers) of thermal decomposition reactions (CO₂ loss: × or +; H₂O loss: ∘). The structures of the studied thermal decomposition reactants and products of C₇H₁₀O₆ are shown here as an example; the structures for the other three elemental compositions are shown in Figs. S11 and S12. Note that the reactant molecules are plotted with the same molar mass as the corresponding product molecule, instead of the molar mass of the reactant molecule.

Figure 3 shows the COSMOtherm-estimated p_sat of the studied thermal decomposition product isomers of C₇H₁₀O₄, C₇H₁₀O₆, C₈H₁₂O₄ and C₈H₁₂O₆ in red markers. The corresponding reactants are shown in blue markers at the product molar mass, and experimental p_sat values of the product elemental compositions are given as a range of the six measurement points. The studied thermal decomposition reaction is possible if the COSMOtherm-estimated saturation vapor pressure of the reactant molecule is lower than the experimental saturation vapor pressure of the product elemental composition. Otherwise, the reactant would desorb from the sample before the thermal decomposition reaction has taken place. For example, the elemental composition of C₇H₁₀O₄ at 158.15 g mol⁻¹ and its reactants C₈H₁₀O₆ (decarboxylation) and C₇H₁₂O₅ (dehydration) all have higher estimated saturation vapor pressures than the one derived experimentally (though Reactant-1 vapor pressure is close to the experimental one; see Fig. 3). This indicates that the measured C₇H₁₀O₄ is likely not a product of dehydration. The decarboxylation reaction is a possible source of the measured C₇H₁₀O₄, assuming under- or overestimation of the saturation vapor pressure by our experiments or COSMOtherm, respectively. Another possibility is that the measured C₇H₁₀O₄ is a fragmentation product of some other thermal decomposition reaction, where the reactant has an even lower saturation vapor pressure. For C₈H₁₂O₄, C₇H₁₀O₆ and C₈H₁₂O₆, some of the studied thermal decomposition reactants have saturation vapor pressures lower than the experimental p_sat. This means that the reactant molecules would remain in the sample at the measured T_max (i.e., potential thermal decomposition temperature). The measured C₈H₁₂O₄ is more likely a product of decarboxylation than dehydration because the proposed dehydration reactant has a higher COSMOtherm-estimated saturation vapor pressure than the experimental p_sat of the product C₈H₁₂O₄. C₇H₁₀O₆ and C₈H₁₂O₆ have similar estimated and experimental saturation vapor pressures, and the measured molecules can therefore be either thermal decomposition products or simply relatively low volatility isomers.

The saturation vapor pressures of the thermal decomposition reactant molecules are 3.3–6.5 (on average 4.7) orders of magnitude lower than the saturation vapor pressures of the corresponding product molecules. The differences between SIMPOL.1-estimated saturation vapor pressures of the reactants and products are 4.9 and 3.9 orders of magnitude for the dehydration and decarboxylation reactions, respectively. Based on this, it is unlikely that the detected molecule is formed in either of these specific thermal decomposition reactions if the COSMOtherm-estimated p_sat of the detected molecule is more than 7 orders of magnitude higher than the p_sat derived from FIGAERO–CIMS experiments. In those cases, the reactant is likely a larger monomer or even a dimer that decomposes to form two large fragments.

It is also possible that other molecules detected in our FIGAERO–CIMS experiments are thermal decomposition products formed during the heating of the sample, though it is impossible to determine if this is true based only on information available from our measurements and calculations. If the decomposition temperature is lower than the T_max of the decomposition product molecule, the measured T_max values can correspond to the saturation vapor pressures of the decomposition products. However, this possibility was not taken into account when we selected the isomers for the COSMOtherm calculations.

Table 1 Saturation vapor pressures of carboxylic acids in pascals. COSMOtherm-, FIGAERO–CIMS- and SIMPOL.1-derived saturation vapor pressures are given at 298.15 K.

Experiments: ^a 296 K, Bilde and Pandis (2001). ^b 298.15 K, Lienhard et al. (2015). ^c 298.15 K, Babar et al. (2020). ^d 298 K, Kostenidou et al. (2018).
COSMOtherm: ^e 298.15 K, Kurtén et al. (2016) (different isomers).
^f The isomers detected in the FIGAERO–CIMS experiments may be different or products of thermal decomposition.

3.4 Comparison with previous studies

Recently, Thomsen et al. (2021) identified multiple carboxylic acids in SOA formed in α-pinene ozonolysis experiments using an ultra-high-performance liquid chromatograph (UHPLC). Out of the compounds included in Thomsen et al. (2021), elemental compositions corresponding to diaterpenylic acid acetate (DTAA, C₁₀H₁₆O₆), 3-methyl-1,2,3-butanetricarboxylic acid (MBTCA, C₈H₁₂O₆), OH-pinonic acid (C₁₀H₁₆O₄), oxo-pinonic acid (C₁₀H₁₄O₄), pinonic acid (C₁₀H₁₄O₃) and terpenylic acid (C₈H₁₂O₄) were measured in our FIGAERO–CIMS experiments. In addition, an elemental composition corresponding to pinic acid (C₉H₁₄O₄) was seen in our experiments, but we were not able to determine T_max values from the thermogram. Previously, Kurtén et al. (2016) calculated saturation vapor pressures of C₁₀H₁₆O₄, C₁₀H₁₆O₆ and C₁₀H₁₆O₈ with several isomers not included in our calculations using COSMOtherm15. In Table 1, we have summarized the saturation vapor pressures of the carboxylic acids identified by Thomsen et al. (2021), as well as all studied isomers of C₁₀H₁₆O₄, C₁₀H₁₆O₆ and C₁₀H₁₆O₈, from our FIGAERO–CIMS measurements, COSMOtherm and SIMPOL.1 calculations, and previous studies. The experimental saturation vapor pressures from previous studies (Bilde and Pandis, 2001; Lienhard et al., 2015; Babar et al., 2020) are given for the specific isomer, while COSMOtherm15 values (Kurtén et al., 2016) are for various other isomers.

Based on COSMOtherm-estimated saturation vapor pressures, it is unlikely that the T_max values of pinic acid, pinonic acid, terebic acid and terpenylic acid could be determined in our FIGAERO–CIMS experiments due to their high volatilities. The FIGAERO–CIMS-derived saturation vapor pressures of C₈H₁₂O₆ and C₁₀H₁₆O₆ agree with previous measurements of MBTCA and DTAA by Babar et al. (2020) and Kostenidou et al. (2018), respectively. The difference of 3 orders of magnitude in experimentally determined saturation vapor pressures of MBTCA (Lienhard et al., 2015; Babar et al., 2020; Kostenidou et al., 2018) demonstrates how different experimental methods give widely different values.

Kurtén et al. (2016) computed saturation vapor pressures of isomers that are potential products of gas-phase autoxidation, rather than products of condensed-phase reactions. The isomers in Kurtén et al. (2016) therefore contain mainly carbonyl, hydroperoxide and peroxy acid groups. Using a combination of group-contribution methods and COSMOtherm15, they concluded that molecules with high oxygen content are likely LVOCs. However, systematic conformer sampling and newer COSMOtherm parametrizations can lead to p_sat estimates in COSMOtherm that are several orders of magnitude lower (Hyttinen et al., 2021b). Two of the C₁₀H₁₆O₆ isomers studied here were taken from Kurtén et al. (2016). Our saturation vapor pressure estimates are 2–4 orders of magnitude lower than those estimated by Kurtén et al. (2016) (see Table S5). Our calculations and experiments show that most of the studied dimers (C₁₃ and higher carbon numbers) are likely ELVOCs (around p_sat < 10⁻⁹ Pa) and the studied monomers with high molar masses (i.e., C₉–C₁₀ and O₁₀) may be ELVOCs, while the studied monomers with lower molar masses (around 190 < M_w < 275 g mol⁻¹) are likely LVOCs (around 10⁻⁹ < p_sat < 10⁻⁵ Pa), with the exception of some higher p_sat isomers at lower molar masses (M_w < 235 g mol⁻¹).

We additionally compared our COSMOtherm vapor pressures with those calculated with SIMPOL.1. The comparison is shown in Fig. S13. We can see that with the molecules in this study, SIMPOL.1 is more likely to overestimate than underestimate COSMOtherm-estimated saturation vapor pressures. COSMOtherm-estimated saturation vapor pressures are up to a factor of 430 higher and up to a factor of 3.5 × 10⁴ lower than those estimated using SIMPOL.1.

4 Conclusions

We have shown that COSMOtherm-estimated saturation vapor pressures agree (for M_w > 190 g mol⁻¹) with those derived from particle-phase thermal desorption measurements of the α-pinene ozonolysis SOA system, taking into account the possibility of thermal decomposition. With our limited set of compounds, we cannot determine the lower-limit saturation vapor pressure for which our experimental method is valid. Additionally, our limited set of calibration compounds further restricts our ability to reliably estimate saturation vapor pressures of the lowest-volatility compounds. However, molecules with ultra-low volatilities likely do not evaporate from the sample during the experiments without fragmenting and are therefore not detected by FIGAERO–CIMS. The measured α-pinene ozonolysis monomer products selected from our SOA sample are mainly LVOCs, and dimers are mainly ELVOCs. The smaller monomers (M_w < 190 g mol⁻¹) with the highest saturation vapor pressures (IVOCs) were likely not present in the sample aerosol collected from the chamber; instead, they are likely products of thermal decomposition formed from larger compounds during the experiment.

Comparison between estimated and experimental p_sat can provide insight about the possible chemical structures of SOA constituents. Based on our results, the commonly used FIGAERO–CIMS instrument is best suited for measuring saturation vapor pressures of monoterpene-derived highly oxygenated monomers in the LVOC and ELVOC range with M_w > 190 g mol⁻¹. Hence, it is reliable for estimating saturation vapor pressures of oxidation products of monoterpenes, such as α-pinene, keeping in mind that the smallest measured molecules are likely products of thermal decomposition. COSMOtherm can be used to estimate saturation vapor pressures of compounds for which p_sat is outside the applicable range of FIGAERO–CIMS experiments, i.e., IVOCs, SVOCs and ULVOCs, if the exact structures of the molecules are known.

In conclusion, this study gives us useful information for studying saturation vapor pressures of multifunctional compounds and further information on the gas-to-particle partitioning of the compounds, which is key when the SOA formation is investigated. Recently, it has been shown that SOA formation has a clear effect on both direct and indirect radiative forcing (Yli-Juuti et al., 2021), highlighting the atmospheric relevance of our study.