Theoretical analysis of the OH-initiated atmospheric oxidation reactions of imidazole

Abstract

Imidazoles are present in Earth's atmosphere in both the gas-phase and in aerosol particles, and have been implicated in the formation of brown carbon aerosols. The gas-phase oxidation of imidazole (C₃N₂H₄) by hydroxyl radicals has been shown to be preferentially initiated via OH-addition to position C5, producing the 5-hydroxyimidazolyl radical adduct. However, the fate of this adduct upon reaction with O₂ in the atmospheric gas-phase is currently unknown. We employed an automated approach to investigate the reaction mechanism and kinetics of imidazole's OH-initiated gas-phase oxidation, in the presence of O₂ and NO_x. The explored mechanism included reactions available to first-generation RO₂ radicals, as well as alkoxyl radicals produced from RO₂ + NO reactions. Product distributions were obtained by assembling and solving a master equation, under conditions relevant to the Earth's atmosphere. Our calculations show a complex, branched reaction mechanism, which nevertheless converges to yield two major closed-shell products: 4H-imidazol-4-ol (4H-4ol) and N,N′-diformylformamidine (FMF). At 298 K and 1 atm, we estimate the yields of 4H-4ol and FMF from imidazole oxidation initiated via OH-addition to position C5 to be 34 : 66, 12 : 85 and 2 : 95 under 10 ppt, 100 ppt and 1 ppb of NO respectively. This work also revealed O₂-migration pathways between the α-N-imino peroxyl radical isomers. This reaction channel is fast for the first-generation RO₂ radicals, and may be important during the atmospheric oxidation of other unsaturated organic nitrogen compounds as well.

---

1 Introduction

Imidazole is a five-membered aromatic heterocycle with two nitrogen atoms in meta positions: one that is pyrrole-like (>NH) and one that is pyridine-like (=N–). Alongside other organic nitrogen compounds,¹ imidazole and its derivatives (hence referred to as imidazoles) have been observed in the Earth's atmosphere, both in the gas-phase^2,3 and in aerosol particles.^4–7 Recently, imidazole's presence in the atmosphere has garnered attention due to its potential influence on the planet's climate.^4,8–10 More generally, imidazole is ubiquitous in nature due to its role in biochemistry,¹¹ where it can act as a versatile catalyst in enzyme active sites.¹²

The majority of atmospheric imidazoles are thought to be secondary in origin, i.e., formed in the atmosphere from reaction of precursor compounds. Numerous experimental and theoretical studies report that imidazoles can be formed from α-dicarbonyl compounds, such as glyoxal, reacting with ammonia or primary amines in aqueous aerosol particles.^13–20 Other possible sources of imidazoles to the atmosphere include anthropogenic emissions from industrial activities,⁴ biomass burning,²¹ and riverine water–air interface reactions.²² This latter source, involving oxidation of organic compounds by ozone at the river surface microlayer, may be a possible source of imidazole emissions directly into the atmospheric gas-phase. Although most studies investigating atmospheric imidazoles focus on condensed-phase chemistry, recent field experiments in the Hyytiälä boreal forest research station in Finland reported the detection of these compounds also in the gas-phase,^2,3 with concentrations of up to ∼183 ppt, which is comparable to the typical concentration of abundant amines.²³

Concerning their impact on the atmosphere, imidazoles have received much attention due to their role as components of brown carbon, in which they strongly affect the aerosol optical properties.^4,24,25 Moreover, studies show that imidazoles may act as photosensitizers, promoting accelerated aerosol particle growth via gas-particle interface reactions.^8,26,27 Aerosols still constitute one of the largest sources of uncertainty in climate model predictions,²⁸ largely due to insufficient understanding of the molecular-level processes controlling the formation, growth, ageing, composition, and physical properties of atmospheric aerosol particles.^9,29 Among these processes is the degradation of aerosol precursor species by atmospheric oxidation, which may lead to several different products, oftentimes with widely different properties.

Either in the gas-phase or dissolved in aqueous aerosol particles, imidazoles can react with OH radicals, the main atmospheric oxidant,^30–35 yet knowledge of their fate upon oxidation is limited. The mechanism for imidazole reaction with OH radical in aqueous solution has been previously investigated with experiments^32,33 and theory,³⁵ revealing that oxidation is preferentially initiated via electrophilic addition to position C5 leading to the 5-hydroxyimidazolyl radical adduct (5-OH) (Scheme 1).

Scheme 1 Imidazole's reaction with OH radical via attack on position C5.

There are only few studies investigating imidazole's gas-phase atmospheric chemistry.^31,36,37 A rate coefficient for the reaction of imidazole with OH radical was reported to be k(298 K) = 3.8 × 10⁻¹¹ cm³ molecule⁻¹ s⁻¹,³¹ as measured with flash photolysis-resonance fluorescence. Subsequent theoretical studies predicted that just like in the aqueous-phase, OH-addition to position C5 is favored (∼88% yield), while addition to position C2 is a minor channel.^36,37 However, no study has explored the atmospherically important reactions of the 5-OH product, an alkyl (C-centered) radical, with O₂.

Alkyl radicals formed in the atmosphere are prone to react fast with O₂, yielding peroxyl radical (RO₂) products:

R˙ + O₂ → RO₂˙ (1)

which may further react via a myriad of unimolecular or bimolecular reaction pathways.^38–40

RO₂ unimolecular reactions such as intramolecular H-shifts,

RO₂˙ → R′˙OOH (2)

in which the peroxyl group abstracts a hydrogen from another moiety in the molecule, or endo-cyclization,in which the peroxyl group attacks a C=C double-bond,⁴¹ may lead to propagation of the oxidative chain by reforming a C-centered radical, allowing for addition of a second O₂ molecule, in what are called autoxidation channels.^40,42,43 These channels may eventually produce highly functionalized, low-volatility products, contributing to enhanced aerosol particle growth.⁴⁴

If, however, a N-centered radical is formed (e.g., via H-shift from an amino group), subsequent reactions typically lead to termination of the oxidative chain. Their reactions with O₂ leads to the formation of a closed-shell imine plus HO₂,⁴⁵

RCH(–N˙R′)R′′ + O₂ → RC(=NR′)R′′ + HO₂˙ (4)

either via a step-wise addition/elimination mechanism or via a direct H-abstraction mechanism.^46–48 Since N-centered radicals are much less reactive towards O₂ than C-centered radicals, reactions with NO or NO₂ are also of importance,⁴⁹ yielding nitrosamines (RN(–NO)R′) or nitramines (RN(−NO₂)R′) respectively, known carcinogens.⁵⁰

RN˙R′ + NO_x˙ → RN(–NO_x)R′ (5)

Other RO₂ unimolecular channels may also lead to an early termination of the oxidative chain, such as HO₂-elimination reactions, in which a concerted H-shift from a β-peroxyl C, N or O atom (1,4 shift) and C–O(O) bond scission produces a C=C, C=N or C=O π-bond plus HO₂. In general, HO₂-elimination involving concerted H-shift from a carbon is negligible in the atmosphere,^48,51 whereas it becomes competitive, or even predominant, when it involves H-shift from a nitrogen^48,52–54

RO₂˙ → R′C(=NR′′′)R′′ + HO₂˙ (6)

or an oxygen:⁵⁵

RO₂˙ → R′C(=O)R′′ + HO₂˙ (7)

Relevant bimolecular channels of RO₂ include reaction with NO, NO₂, HO₂, and self- or cross-reaction with another R′O₂ radical.³⁸ These channels are often barrierless (or nearly so), since they involve radical recombination, and therefore their competitiveness is largely dictated by the favorability of unimolecular channels, as well as the concentration of reaction partners. For instance, typical NO_x concentrations in pristine environments are on the order of 1–10 ppt, where unimolecular, RO₂ + HO₂ and RO₂ + R′O₂ reactions can compete, while in more polluted conditions NO_x concentrations may be 3 orders of magnitude higher (1–10 ppb), and thus RO₂ + NO reactions tend to dominate.⁵⁶ Reaction of a RO₂ radical with NO leads predominantly to the formation of NO₂ plus the corresponding alkoxyl radical (RO),

RO₂˙ + NO˙ → RO˙ + NO₂˙ (8)

another key reaction intermediate in the atmospheric oxidation of organic compounds that is able to react via multiple channels. These channels include H-abstraction by O₂, leading to a carbonyl plus HO₂, intramolecular H-shifts, and β-scissions,

RO˙ → R′C(=O)R′′ + R′′′˙ (9)

in which cleavage of a C_α–X_β bond produces a carbonyl and a X-centered radical.⁵⁷ β-Scission reactions can be extremely fast, especially in RO radicals derived from cyclic compounds, whose ring-strain is released if cleavage occurs at an endocyclic bond. Moreover, if cleavage occurs at a C_α–N_β bond, a N-centered radical is formed, with potential to produce nitrosamines and nitramines.

Since peroxyl radicals derived from cyclic compounds tend to be, in general, less prone to react via unimolecular pathways than acyclic analogues,⁵⁸ RO₂ + NO reactions could be important during imidazole oxidation even under lower NO concentrations. The rigid structure of a cyclic RO₂ imposes constraints on the orientation of functional groups, thus hindering or completely preventing access to some sites for intramolecular –OO˙ radical attack. On the other hand, the presence of a nitrogen in the ring allows for the possibility of HO₂-elimination with imine formation, reaction which can be fast even for cyclic reactants.^48,54,59 Furthermore, aromatic homocyclic compounds display some mechanistic peculiarities during atmospheric oxidation, such as a prevalence of RO₂ endo-cyclization pathways (Reaction (3)) producing bicyclic endoperoxides.^39,60,61 The mechanism of heterocyclic aromatic atmospheric oxidation is less understood, and so far it is uncertain whether the homocyclic reactivity trends apply also for these compounds, or if they follow entirely different pathways.

To gain fundamental understanding of imidazoles in the atmosphere, we investigate the mechanism of imidazole gas-phase oxidation initiated by OH radical addition to position C5, and follow its chemistry upon the formation of RO₂ in the presence of NO (Reactions (1)–(8)). The explored oxidation mechanism also includes the further reactions available to RO radicals produced viareaction (8). We employ automated high-level quantum chemistry methods to explore the reaction pathways involved in the oxidation mechanism systematically, calculate reaction rate coefficients, and explore the overall reaction kinetics to predict product distributions under a range of atmospherically relevant conditions.

2 Methods

2.1 Potential energy surface exploration

Unimolecular isomerization and decomposition reactions were explored using KinBot,^59,62,63 a workflow code for automated gas-phase reaction search and kinetics. KinBot automatically creates and submits calculations to explore and characterize the relevant regions of a potential energy surface (PES) as described below. Overall, starting from an initial reactant structure well (e.g. RO₂), KinBot proposes the set of unimolecular isomerization and dissociation reactions that this species can undergo, and generates a guess for the related transition state geometries based on reaction templates. After filtering out the invalid reactions and the ones with barriers above a user-defined energy cutoff, the process is repeated for each new well, until no more new stationary points are found.

Stationary point geometries were first optimized at the L1 = B3LYP/6-31+G* level of theory.^64–69 We discarded all reactions with an L1 barrier above a threshold, for which we found that 20 kcal mol⁻¹ relative to the initial reactant was a meaningful choice in all cases. Additionally, intrinsic reaction coordinate (IRC) calculations were used for each of the remaining saddle points to identify the connected minima. Saddle-points are discarded if neither or both of the connected minima corresponded to the reactant. KinBot searched the conformational space of each stationary point at the L1 level. Subsequently, all unique conformers with a Boltzmann factor larger than 0.05 at 300 K were re-optimized at the L2 = ωB97X-D/aug-cc-pVTZ level of theory.^70–73 Finally, single-point L3 = ROHF-ROCCSD(T)-F12a/cc-pVDZ-F12^74,75 calculations were done for the lowest energy conformers of each stationary point. As is recommended for this type of calculation,⁷⁵ we changed the value of the geminal Slater exponent β to 0.9, using the keyword “gem_beta = 0.9”. In order to understand the uncertainties associated with our choice of method for L3, we performed ROHF-ROCCSD(T)-F12a/cc-pVTZ-F12⁷⁵ and ROHF-UCCSD(T)-F12a/cc-pVDZ-F12 single-point calculations on a subset of important species. Comparison of the results obtained is shown in Table S11 (ESI†).

We selected the lowest energy structures in two ways, based on the zero-point energy (ZPE) corrected energies and based on the 300-K free-energy. The former definition was used in the multi-conformer transition state theory (MC-TST) calculations, while the latter was used to select the structure in the master-equation, where we used a single conformer for each stationary point. Note that in the MC-TST calculations we applied the same L3//L2 energy correction to the higher energy conformers as to the lowest one.

The L1 and L2 calculations were done via KinBot's atomic simulation environment (ASE)⁷⁶ interface to Gaussian,⁷⁷ while Molpro v2022.2^78–80 was used for the L3 energies offline with automatically generated inputs.

For a few reactions, such as bimolecular reactions connecting successive parts of the oxidation mechanism, the PES was explored manually. In those cases, conformer search was performed with the Spartan v20⁸¹ software, employing the Merck-molecular force field (MMFF94) method.⁸²

2.2 Reaction kinetics

We estimated thermal rate coefficients for each unimolecular reaction using transition state theory (TST),⁸³ both using only the lowest-energy conformer and a multi-conformer approach (MC-TST),⁷³ assuming rigid-rotor-harmonic-oscillator (RRHO) partition functions for each structure. Product yields were estimated by assembling and solving a master equation, using the Master Equation Solver for Multi-Energy well Reactions (MESMER) software version 7.0⁸⁴ using the lowest-free energy conformers only.

2.2.1 Reaction rate coefficients. Thermal rate coefficients (k(T)) were estimated using the Eyring equation⁸³where k_B is the Boltzmann constant, h is the Planck constant, and T is the temperature. Δ^‡E_v=0 is the reaction energy barrier, defined as the difference in ZPE-corrected energies between transition state and reactant, calculated at L3//L2. Q_TS and Q_R are the total molecular partition functions of the transition state and reactant respectively, calculated under the RRHO approximation. The transmission coefficient κ was estimated using the Eckart tunnelling correction.⁸⁵

In order to assess the effect of multiple conformations on reaction rates, we also used a MC-TST approach⁷³ to estimate k(T) values using a sum of Boltzmann-weighed conformer contributions.

where Q_R,i and Q_TS,i are the total molecular partition functions of conformer i of the reactant and of the transition state, respectively. E_R,i − E_R,0 and E_TS,i − E_TS,0 are the ZPE-corrected energy differences between local minima and the global minimum (subscript 0) for the reactant and the transition state respectively, calculated at the L2 level.

2.2.2 Product yields. Oxidation product yields were obtained by solving a master equation, and simulating the time-evolution of species concentrations under atmospheric conditions. Zero-point corrected energies, vibrational frequencies, and rotational constants calculated for the lowest-free energy conformation of each stationary point, under the RRHO approximation, were used as input to the master equation. Initial association reactions, such as OH radical addition to imidazole and O₂ addition to the subsequently produced 5-OH adduct, were assumed to be barrierless reactions and were treated using the inverse Laplace transform (ILT) method.^86,87 A temperature-dependent rate coefficient of 1.7 × 10⁻¹² exp(+930 K/T) cm³ molecule⁻¹ s⁻¹ was assigned to the imidazole + OH reaction, corresponding to experimentally measured values.³¹ For all reaction steps involving O₂ addition to alkyl radicals, whose rate coefficients are unknown, a temperature-dependent rate coefficient value typical of such reactions, 2.0 × 10⁻¹² × (T/300 K)⁻¹ cm³ molecule⁻¹ s⁻¹, was used.

Unimolecular reactions (isomerizations and dissociations) were modelled with standard Rice–Ramsperger–Kassel–Marcus (RRKM) theory, using the Eckart tunnelling correction.⁸⁵ Each dissociation transition state was assumed to be connected by isomerization to a post-reactive complex (PC), whose decomposition into bimolecular products was modelled in a separate, irreversible and barrierless reaction step, using ILT with a rate coefficient of 10⁻¹⁰ cm³ molecule⁻¹ s⁻¹ corresponding to the reverse association reaction. Finally, RO₂ + NO reactions were modelled as simple bimolecular sinks, with a rate coefficient value of 2.7 × 10⁻¹² exp(+360 K/T) cm³ molecule⁻¹ s⁻¹.⁸⁸

In the simulations a high concentration of OH radicals (10¹⁵ molecules cm⁻³cf ∼ 10⁶ molecules cm⁻³ ambient OH concentration) was used with the intent of starting the oxidative chain within ∼1 μs. This enabled us to see reactive intermediates on short time-scales and resolve low-yield intermediates better in the numerical simulations, but it did not affect any other reaction step in the master equation. For O₂, a mixing ratio of 0.20946 (5.16 × 10¹⁸ molecules cm⁻³ at 1 atm, 298 K) was used. The RRKM-ME calculations were done in the 250–350 K temperature range and using three NO concentrations (10 ppt, 100 ppt and 1 ppb). We extrapolated the range of investigated temperatures beyond realistic atmospheric conditions, especially towards the hotter end, with the intention of gaining further insight into the kinetics of oxidation.

We employed the single exponential down model to treat collisional energy transfer, using a 〈ΔE〉_down parameter value of 225 cm⁻¹. Lennard-Jones parameters for modelled wells were calculated from pure-compound critical properties,^89,90 which were in turn estimated with the Joback group-additivity method.⁹¹ N₂ was set as the bath gas, with Lennard-Jones parameters ε/k_B = 91.85 K and σ = 3.919. We employed 50 cm⁻¹ energy grains, which spanned a maximum energy of 20 k_BT above the highest-energy stationary point. Test calculations done with a smaller energy grain size (30 cm⁻¹) resulted in negligible differences in product yields (maximum difference of ∼0.04%). The time evolution of species distributions was obtained by summing over populations in energy grains at each time step. We used the quad-double (qd) numerical precision in all ME calculations.

Once the product distributions from the reaction of first-generation RO₂ radicals were obtained, we explored the reaction pathways available to the open-shell products with significant yields, using KinBot and repeating the computational routine described above (PES exploration and RRKM-ME simulation), with a few differences. The fate of each major alkoxyl radical product was investigated with a separate RRKM-ME calculation, setting the parent RO₂ radical as the initial reactant. Formation of the alkoxyl radical occurs with two steps, both modelled with standard ILT: first, RO₂ + NO reaction produces a ROONO adduct, which then dissociates into RO + NO₂ in a second step. We set the NO and NO₂ concentrations to 100 ppt. An O₂ addition reaction step is added for each possible alkyl radical intermediate, but this time modelled as simple bimolecular sinks. Once a major, second-generation RO₂ product is identified, the RRKM-ME calculation is run again, but this time treating the corresponding O₂-addition channel with standard ILT, and including subsequent unimolecular reaction steps.

Finally, to assess the uncertainty of calculated product yields due to the uncertainty in RRKM-ME model parameters and the calculated PES, we employed the automated procedure present in KinBot.⁵⁹ Here the approach was adapted for the MESMER code. A detailed description for this procedure is given in the ESI.†

3 Results and discussion

3.1 First generation peroxyl radicals

We start by exploring the fate of imidazole after initial OH radical addition to position C5, producing the C-centered alkyl radical, 5-OH. Adduct 5-OH has no unimolecular reaction channel that is able to compete with O₂-addition, as the only pathway found during reaction search (within the 20 kcal mol⁻¹ energy barrier cutoff) was dissociation back to reactants (imidazole + OH). Direct H-abstraction from C5 by O₂, leading to the aromatic species 1H-imidazol-5-ol (5-ol), is also unable to compete with O₂-addition. The calculated energy barrier for this reaction is 9.93 kcal mol⁻¹, and the associated rate coefficient value estimated with TST is 1.14 × 10⁻²¹ cm³ molecule⁻¹ s⁻¹.

Reaction of O₂ with 5-OH may initially lead to the formation of one of four different RO₂ isomers. The odd-electron density in 5-OH is delocalized over two carbon atoms, C2 and C4, due to resonance stabilization, thus the addition of O₂ may occur at either of these positions. Moreover, depending on the side of the ring at which the O₂ attack occurs, addition to positions C2 and C4 may lead to two RO₂ stereoisomers each. Anti-addition (relative to the OH moiety) may lead to peroxyl radicals 2a-RO₂ or 4a-RO₂, while syn-addition may lead to peroxyl radicals 2s-RO₂ or 4s-RO₂. The stationary points on the PES, including O₂-addition reactions, as well as the unimolecular reaction pathways available to each of the RO₂ species, are shown in Fig. 1.

Fig. 1 PES of the studied reaction pathways available to the first-generation RO₂ from imidazole atmospheric oxidation, initiated via OH radical-addition to position C5. Relative ZPE-corrected energies shown in kcal mol⁻¹, calculated at the ROHF-ROCCSD(T)-F12a/cc-pVDZ-F12//ωB97X-D/aug-cc-pVTZ level of theory. Solid black lines indicate unimolecular steps, solid gray lines indicate barrierless association and dissociation steps, and dotted lines indicate minor channels. The dashed circle indicates the initial reactants (imidazole + OH + O₂), and dashed boxes indicate the major products from major unimolecular channels (Z-FMF + OH, and 4H-4ol + HO₂).

The 4-RO₂ adducts are 3–4 kcal mol⁻¹ more stable than their 2-RO₂ counterparts. This difference can be explained by the hyperconjugative interaction between the amine nitrogen's lone-pair and the imine's π-antibonding orbital (n_N → π_CN*), which is retained in the C4-addition products (in the amidine group), but lost upon O₂-addition to position C2. The energy differences between syn- and anti-adducts are smaller, and the reasons are less obvious. Unlike in anti-isomers, the peroxyl and hydroxyl groups in syn-isomers are able to interact via an intramolecular H-bond, albeit a weak one, since ring-strain prevents an optimal O–H–O orientation. The 2s-RO₂ isomer is therefore 0.8 kcal mol⁻¹ lower in energy compared to 2a-RO₂. However, the opposite trend is observed for the 4-isomers: 4s-RO₂ is 0.4 kcal mol⁻¹ higher in energy than 4a-RO₂. Here, destabilization brought by steric repulsion between the neighboring peroxyl and hydroxyl groups, which is significant with β-OH syn-substitution, overcomes the stabilizing effect of the H-bond.

In the 2-RO₂ isomers, the peroxyl group is bonded to an α-amino C, and thus both of these species may undergo HO₂-elimination (1,4 H-shift from the nitrogen concerted with cleavage of the C–O(O) bond), resulting in C=N double-bond formation via transition states TS5a or TS5s. These two HO₂-elimination pathways have similar reaction barrier heights, differing by about 0.4 kcal mol⁻¹, and lead to the same products, the cyclic diimine 4H-imidazol-4-ol (4H-4ol) plus HO₂ radical. In contrast, the reaction channels which are available to the 4-RO₂ isomers depend on the reactant's stereochemistry that determines which moieties are accessible for attack of the peroxyl group. Unlike the syn isomer 4s-RO₂, the anti isomer 4a-RO₂ may also react via HO₂-elimination (via transition state TS6a), but in this case with concerted H-transfer from the α-hydroxyl carbon, producing the aromatic species 1H-imidazol-5-ol (5-ol) plus HO₂ radical. Even though this channel recovers the aromaticity of the heterocycle, leading to a product that is more stable than 4H-4ol by 16.7 kcal mol⁻¹, its reaction barrier is the highest among the ones studied here for first-generation RO₂ species.

The peroxyl group in the syn isomer 4s-RO₂ has access to the alcohol group, and this species may react via an 1,5 H-shift (via TS7s), producing the hydroperoxy-alkoxyl radical (HOOQO) M7s. This channel has a slightly lower reaction barrier height than the HO₂-elimination pathways, and it initiates a cascade of exothermic unimolecular reaction steps, ultimately leading to the formation of another closed-shell product, (Z)-N,N′-diformylformamidine (Z-FMF), plus OH radical. Most substituted RO radicals may undergo fast β-scission, producing a carbonyl and a C-centered radical. Species M7s is especially prone to react in this manner (via TS8s) because the reaction barrier associated with β-scission is only 0.9 kcal mol⁻¹. Release of ring-strain and the development of a high degree of electron delocalization in TS8s may explain the increased reactivity observed. The (HOO)C–C(O) distance in M7s is also abnormally large (∼2 Å), and for this reason, this part of the PES was calculated manually. The product of this β-scission, species M8s, is an α-OOH alkyl radical, a class of intermediates which are known to undergo spontaneous OH radical elimination,⁹² producing a closed-shell carbonyl compound

RC˙(–OOH)R′ → RC(=O)R′ + ˙OH (12)

Since most conformations of the reactant spontaneously decompose into the products (Z-FMF and OH radical) during geometry optimization, we assumed that dissociation of M8s occurs without a barrier.

Finally, O₂-migration reaction steps connect the two anti-isomers (2a-RO₂ ⇄ 4a-RO₂) and the two syn-isomers (2s-RO₂ ⇄ 4s-RO₂), via TS4a and TS4s respectively, with the lowest energy barriers among the RO₂ unimolecular channels. These reactions resemble endo-cyclization steps (Reaction (3)), but involving a concerted attack of the peroxyl group on the C-terminus of the C=N double-bond and C–O(O) bond cleavage. As a result, the peroxyl group and the C=N double-bond exchange positions.

We note that T1 diagnostics and the percentage of triples (T) contribution to the total atomization energy %TAE[(T)]^93,94 suggest that transition states TS4a, TS4s, TS5a, TS5s, TS6a, TS7s, and TS8s have a mild degree of multi-configurational character (see Section S9 in the ESI†). Test calculations with a unrestricted coupled-cluster wavefunction (ROHF-UCCSD(T)-F12a/cc-pVDZ-F12) show a substantially lower (∼2.1 kcal mol⁻¹) relative energy for TS7s, in comparison to the restricted open-shell method used for L3 (see Table S11, ESI†). However, since this difference does not change the order of importance among exit reaction channels, uncertainties associated with the use of a ROCCSD(T) method (as opposed to UCCSD(T)⁹⁵) are not expected to significantly affect branching ratios.

While the formation of electronically-excited RO₂ intermediates directly from the association of an alkyl radical with O₂ may be possible,⁹⁶ a thorough exploration of excited-state reaction mechanisms is beyond the scope of this study. We argue, however, that if excited RO₂ species are formed from association of 5-OH and O₂, the likely major outcome is fast relaxation to the electronic ground state via internal conversion.⁹⁷ Further discussion is presented in Section S8 of the ESI† alongside test calculations.

The calculated reaction barrier heights and associated thermal rate coefficients, estimated with TST, for each unimolecular reaction described above are shown in Table 1. To capture the chemical evolution of the first-generation RO₂ system, we invoked an RRKM-ME model, in which we also included the reaction of each RO₂ radical with NO (Reaction (8)) as loss terms as described in the Methods section, to determine the competition between unimolecular and bimolecular channels. The alkoxyl radicals produced this way are 4a-RO, 4s-RO, 2a-RO and 2s-RO. The simulation results (298 K and 10 ppt NO) are summarized in Fig. 2.

Table 1 Reaction barrier heights (Δ^‡E_v=0) calculated as the difference in zero-point corrected energy between TS and reactant, calculated at L3//L2 level, for unimolecular reaction steps available to first-generation RO₂ intermediates. Thermal rate coefficients estimated with lowest-conformer TST (k_LC-TST) and multi-conformer TST (k_MC-TST) at 298 K, ratio of lowest-conformer partition functions (Q_TS/Q_R), and Eckart tunnelling factors (κ). *: non-elementary reaction steps. Rate coefficient values estimated assuming a direct connection between 4a-RO₂ and TS5a, and between 4s-RO₂ and TS5s. **: Eckart tunnelling factor lower than unity arrises because of quantum mechanical reflection⁹⁸

Reaction step	Description	Δ^‡Ev=0 (kcal mol^−1)	k LC-TST (s^−1)	k MC-TST (s^−1)	k MC/kLC	Q TS/QR	κ
2a-RO2 → 4a-RO2	O2-Migration	10.05	1.69 × 10^5	8.24 × 10^4	0.49	0.61	1.04
2a-RO2 → PC5a	HO2-Elimination	17.57	4.14	3.05	0.74	1.67	3.02
4a-RO2 → 2a-RO2	O2-Migration	14.04	2.01 × 10^2	30.94	0.49	0.47	1.04
4a-RO2 → PC6a	HO2-Elimination	24.29	1.78 × 10^−5	5.54 × 10^−6	0.31	1.04	1.77
4a-RO2 → PC5a*	HO2-Elimination*	21.56	3.82 × 10^−3	1.46 × 10^−3	0.38	1.31	3.02
2s-RO2 → 4s-RO2	O2-Migration	11.89	2.49 × 10^4	1.05 × 10^4	0.42	2.05	1.02
2s-RO2 → PC5s	HO2-Elimination	19.17	0.328	0.145	0.44	2.40	2.46
4s-RO2 → 2s-RO2	O2-Migration	14.39	3.65 × 10^2	93.6	0.42	1.54	1.02
4s-RO2 → M7s	1,5 H-shift	20.50	8.37 × 10^−3	2.83 × 10^−3	0.34	0.38	3.75
4s-RO2 → PC5s*	HO2-Elimination*	21.67	3.62 × 10^−3	1.30 × 10^−3	0.36	1.81	2.46
M7s → M8s	β-Scission	0.93	1.28 × 10^12	9.88 × 10^11	0.77	1.03	0.96**

---

Fig. 2 Time evolution of species distributions from RRKM-ME calculation for the first-generation RO₂ radicals, at 298 K, 1 atm and 10 ppt NO, with imidazole + OH as starting reactants. Vertical dashed lines indicate relevant time-scales: (A) onset of oxidation; (B) 4s-RO₂ ⇄ 2s-RO₂ and 4a-RO₂ ⇄ 2a-RO₂ equilibria are achieved; (C) 4a-RO₂ ⇄ 4s-RO₂ equilibrium is approached; (D) onset of unimolecular reactions; (E) onset of bimolecular RO₂ + NO → RO + NO₂ reactions; (F) end of reaction.

Significantly faster than other competing channels, the O₂-migration reactions serve to rapidly establish equilibrium between the 2- and 4-RO₂ isomer pairs, highly favoring the 4-RO₂ isomers. Past point B in Fig. 2, i.e., after ∼1 ms most of the produced 2a-RO₂ and 2s-RO₂ radicals have been converted into their respective more stable 4a-RO₂ and 4s-RO₂ isomers. Due to these fast O₂-migration steps, the 4-RO₂ isomers also have prompt, indirect access to the HO₂-elimination pathways available to the 2-RO₂ isomers via TS5a and TS5s.

Equilibrium between the syn and anti 4-RO₂ isomers is approached after around 1 s (Fig. 2(C)), with a portion of the 4s-RO₂ species being converted into the slightly more stable (∼0.3 kcal mol⁻¹) 4a-RO₂, indicating that dissociation of the RO₂ adducts back to the reactants (5-OH + O₂) is significant. Some studies investigating the reaction of O₂ with other resonance stabilized radicals reported a similar behaviour.⁹⁹ Due to partial or complete loss of resonance stabilization, the RO₂ radicals produced in these reactions correspond to relatively shallow wells, leading to faster back-reaction rates. In the absence of other fast RO₂ sinks, which seems to be the case during the OH-initiated oxidation of imidazole, the R˙ + O₂ ⇄ RO₂˙ reaction may approach equilibrium on relatively short timescales. As a consequence, the calculated product yields have a low dependence on the uncertainty in the rate coefficients of the four O₂-addition reactions, which were estimated as described in the Methods section.

As seen in Table 2 and Fig. 3, results from RRKM-ME calculations show that unimolecular channels via TS5a, TS5s and TS7s, as well as bimolecular RO₂ + NO → RO + NO₂ channels involving the 4a-RO₂ and 4s-RO₂ isomers are all competitive at 298 K and 1 atm. The branching ratios of competitive channels are, however, strongly dependent on the NO concentration. We illustrate this effect in Fig. 3, where we show the branching fraction predictions and their uncertainty at various NO concentrations for the main channels. The only negligible channel under most conditions considered here is HO₂-elimination via TS6a to produce 5-ol. This finding is in line with other studies indicating that this type of reaction (HO₂-elimination involving H-shift from a C atom) is a negligible pathway in the Earth's atmosphere.^48,51,59 Production of 5-ol via direct H-abstraction from 5-OH by O₂ is also negligible, as confirmed by a test RRKM-ME calculation where this channel was included alongside the RO₂ channels.

Table 2 Product yields obtained from ME calculations for reaction of first-generation RO₂ radicals, at 298.15 K and 1 atm. Standard deviations (±σ) calculated with the logit-transformed distribution of yields obtained with uncertainty analysis

Product	Yield (%) (±σ)
	[NO] = 10 ppt	[NO] = 100 ppt	[NO] = 1 ppb
5-ol	0.1 (^+0.1−0.0)	0.0 (^+0.1−0.0)	0.0 (^+ 0.1−0.0)
4H-4ol	33.7 (^+15.6−9.4)	12.4 (^+6.4−4.7)	1.8 (^+1.1−0.6)
Z-FMF	47.0 (^+11.3−16.6)	17.3 (^+9.2−7.2)	2.4 (^+1.8−1.0)
4s-RO	4.9 (^+2.2−2.3)	18.1 (^+4.4−6.6)	25.3 (^+9.8−6.0)
4a-RO	14.2 (^+5.0−5.9)	51.9 (^+9.9−11.4)	70.1 (^+6.0−11.1)
2-RO	0.1 (^+0.0−0.0)	0.3 (^+0.1−0.1)	0.4 (^+0.1−0.1)

---

Fig. 3 Product branching ratios obtained from ME calculations for reaction of adduct 5-OH with O₂, at a range of temperatures (250 K, 262.5 K, 275 K, 287.5, 298 K, 312.5 K, 325 K, 337.5 K, 350 K) and NO concentrations (10 ppt, 100 ppt and 1 ppb). All calculations were done at 1 atm. Gray circles represent branching ratio values obtained from uncertainty analysis runs. Gray shaded areas represent standard deviations, calculated with the logit-transformed distribution of branching ratios from uncertainty analysis. Dark gray: ±σ; light gray: ±2σ.

Under more pristine conditions ([NO] = 10 ppt) unimolecular channels predominate, products 4H-4ol + HO₂ and Z-FMF + OH showing yields of 33.7 (^+15.6_−9.4)% and 47.0 (^+11.3_−16.6)% respectively at 298 K, whereas bimolecular channels are minor, with products 4a-RO + NO₂ and 4s-RO + NO₂ showing yields of 14.2 (^+5.0_−5.9)% and 4.9 (^+2.2_−2.3)% respectively. Rising NO concentrations are accompanied by a sharp rise in the yield of alkoxyl radicals, which become the major products at intermediate (100 ppt) and high (1 ppb) NO concentration levels. Unlike with unimolecular channels, competition among the four bimolecular channels depends on the relative concentration of the RO₂ reactants. Thus, the equilibrium established between the RO₂ isomers is reflected on the relative yield of RO products: the yield of 4a-RO is about 3 times larger than that of 4s-RO, while 2a-RO and 2s-RO have a combined yield of only 0.4 (±0.1)% under high NO concentration (1 ppb).

Product yields from first-generation RO₂ also display a strong temperature dependence. As expected, unimolecular channels are favored at higher temperatures, while the barrierless bimolecular RO₂ + NO → RO + NO₂ channels, whose rate coefficients have a weak negative temperature dependence, are favored at lower temperatures. Thus, rising temperatures counteract the effect of rising NO concentrations, as seen in Fig. 3. Also, upwards of 312 K at [NO] = 100 ppt, 4H-4ol becomes the major product, with a yield of 45.5%, surpassing that of Z-FMF (42.4%). While only one reaction channel leads to Z-FMF, two HO₂-elimination channels can produce 4H-4ol, both with very similar energy barrier heights and loose transition states (see ratio of TS and reactant partition functions Q_TS/Q_R in Table 1). So, as the temperature rises, entropy favors formation of 4H-4ol. The equilibrium with the four RO₂ isomers discussed previously is more clearly observed when comparing the RRKM-ME time evolution of concentrations calculated at different temperatures, shown in Fig. S5 (ESI†). Lower temperatures see a delay in the establishment of these equilibria, which is reflected in the ratio of RO stereoisomer yields. At 250 K, RO₂ back-reactions are minor and only occur prior to thermalization of the RO₂ adducts, seen by the slightly higher concentration of 4a-RO₂ relative to 4s-RO₂.

The dependence of product yields on pressure is weak, as seen from Table S7 (ESI†), with yields changing by a maximum of ∼0.5% in the 0.1–10 atm range, at 298 K. This possibly indicates that the kinetics for this part of the oxidation mechanism is close to the high-pressure limit, the regime at which TST is applicable.

3.2 Alkoxyl radical products

Given that the alkoxyl radicals 4a-RO and 4s-RO appear to be major products of imidazole oxidation by OH radicals, being competitive even under low NO concentrations, we also investigated their fate in the atmosphere, as illustrated by the PES in Fig. 4. The RO species were assumed to be produced via the following:

4s-RO₂˙ + NO˙ → 4s-ROONO ΔE = −24.8 kcal mol⁻¹ (13)

4s-ROONO → 4s-RO˙ + NO₂˙ ΔE = +12.7 kcal mol⁻¹ (14)

4a-RO₂˙ + NO˙ → 4a-ROONO ΔE = −24.5 kcal mol⁻¹ (15)

4a-ROONO → 4a-RO˙ + NO₂˙ ΔE = +13.4 kcal mol⁻¹ (16)

Fig. 4 Potential energy surface of studied reaction pathways available to the major alkoxyl radical products. Relative zero-point energy values shown in kcal mol⁻¹, calculated at the ROHF-ROCCSD(T)-F12a/cc-pVDZ-F12//ωB97X-D/aug-cc-pVTZ level of theory. Dotted lines indicate minor channels. Dashed circles indicate the initial reactants 4s-RO and 4a-RO (“s” and “a” stand for syn- and anti-respectively), and dashed boxes indicate major products.

Three initial unimolecular reaction pathways are available to both the syn isomer (via TS12s, TS13s and TS14s) and the anti isomer (via TS12a, TS13a and TS14a) of the major alkoxyl radical intermediate 4-RO. Each of these channels correspond to the cleavage of one of the RO˙ C_α–X_β bonds (β-scissions), leading to the formation of a carbonyl group and a X-centered radical.

The RO channels involving cleavage of the C–C bond, via TS12a and TS12s, have the lowest associated reaction energy barriers (4.17 and 5.04 kcal mol⁻¹ respectively), and initially lead to different products. In the syn RO isomer, the oxyl and alcohol groups interact with a weak H-bond, which is preserved (and perhaps strengthened) as the C–C bond is broken during reaction via TS12s, as revealed by IRC calculations. Further along the post-TS reaction coordinate, the alcohol hydrogen is spontaneously transferred to the carbonyl oxygen atom in a 1,7 H-shift, producing the C-centered radical M12s. On the other hand, no such intramolecular H-bond is present in the anti RO isomer, and C–C bond scission proceeds with no concerted rearrangements, leading to species M12a. However, torsion about the C˙–N bond in M12a (through transition state TS15a) brings the alcohol group in proximity of the carbonyl, allowing for the 1,7 H-shift to spontaneously occur, thus connecting M12a to M12s via a low torsional barrier (∼3.6 kcal mol⁻¹).

Resonance brought by N-imino and N-amido substitution to the radical center makes species M12s especially stable, being almost 15 kcal mol⁻¹ lower in energy than M12a. M12s is also completely planar and assumes a 8-membered ring-like structure, held together by an intramolecular H-bond. Planarity favors conjugation of the functional groups, but the ring-like conformation requires large bond angles, introducing strain to the system. The H-bond is lost with inversion at the conjugated C–N–C moiety (TS17s), producing species M17s, but it also leads to strain relief while maintaining the planar geometry. As a consequence, M17s is just slightly less stable than M12s, by about 0.6 kcal mol⁻¹.

The channels involving cleavage of the C–N bond, via TS13a and TS13s, have slightly higher reaction energy barriers (9.57 and 7.37 kcal mol⁻¹ respectively), both leading to the same product, the N-centered (amidinyl) radical M13. Similarly, cleavage of the C–H bond (H-elimination) in either RO stereoisomer via TS14a or TS14s leads to the same outcome, the closed-shell product imidazol-4-one (P4) plus an H atom, but crossing much higher reaction barriers (16.02 and 15.38 kcal mol⁻¹ respectively).

All of the C–C and C–N bond scission transition states (TS12a, TS12s, TS13a, and TS13s) are stabilized by release of ring strain, which can explain the low energy barriers observed overall (<10 kcal mol⁻¹). Additionally, TS12a and TS12s are further stabilized by a hyperconjugative interaction of the developing, singly-occupied p_C-orbital with the hydroxyl O lone-pair (n_O → p_C), and to a lesser extent, with the amidine N(H) lone-pair (n_N → p_C). These interactions can explain the lower C–C bond scission energy barriers relative to C–N bond scissions. The effects underlying the observed difference in reactivity between stereoisomers is less clear, but a tentative explanation could be the presence of an intramolecular H-bond in 4s-RO and the maintenance of such bond in TS12s and TS13s. Similar to what is observed with 4s-RO₂, an H-bond between the oxyl and hydroxyl groups in 4s-RO is possible but sub-optimal, due to a constraint imposed by the cyclic structure. This constraint is relaxed as the C–N bond is lengthened in TS13s, allowing for a stronger H-bond, which grants extra stabilization to the transition state. Thus, the energy barrier to C–N bond scission is ∼3 kcal mol⁻¹ lower for 4s-RO than it is for 4a-RO, where no such H-bond is possible.

Subsequent reactions available to the two C-centered radicals (M12a and M12s) and the N-centered radical (M13) produced from the RO β-scission channels were also explored. Several of these reaction pathways were found during reaction search for M12a and M12s. However, most have high (>20 kcal mol⁻¹) associated energy barriers, and are not expected to compete with O₂-addition. Thus, for species M12a and M12s, we limited our analysis to the pathways whose reaction barriers are lower in energy than the original RO reactants 4a-RO and 4s-RO. These reactions are included in the PES in Fig. 4.

We also note here that, given the high degree of conjugation present in species M12s and M12a, it is not immediately obvious whether equilibrium between their respective conformations (or E–Z isomers) can be achieved, especially considering the short lifetime of C-centered radicals. For this reason, we treated M12a and M12s as single-conformer species, corresponding to the structures obtained from IRC end-points, connected to TS12a and TS12s respectively. An exhaustive search for reaction paths connecting E–Z isomers would go beyond the scope of this work, but a few critical torsional and inversion transition states were found, either accidentally by KinBot's algorithm (TS15a, TS17a, TS18a) or by manual search (TS17s). Thus, torsional isomerism in species M12a and M12s was treated explicitly, but including only the steps necessary to properly connect reaction channels to either M12a or M12s. For example, one of the reaction channels found is the one leading from species M12a to species M19a, involving a 1,5 H-shift. In order to react in this manner, species M12a must first undergo torsion about the imine C=N bond (TS17a), followed by torsion about the C˙–N bond (TS18a), leading to species M18a, which may then finally react via H-shift (TS19a) to yield product M19a. The N-inversion transition state TS17s was included in the analysis because, as shall be discussed further in Section 3.3, O₂-addition to species M12s or M17s leads to different RO₂ radicals which are not necessarily interconnected via low-barrier channels.

Species M12a and M12s may also react by cyclization (via TS16a or TS16s), where the carbonyl O attacks the radical center, forming a 6-membered ring. A more detailed discussion on the cyclization and H-shift channels available to M12a or M12s is provided in Section S2 of the ESI.†

Species M13, the product of RO C–N bond scission, may react further via a 1,4 H-shift (TS20), where the hydrogen is transferred from the α-OH carbon to the imine N, or via a 1,5 H-shift (TS21), where the hydrogen is transferred from the aldehyde group instead.

The calculated reaction barrier heights and associated thermal rate coefficients for the unimolecular channels discussed in this section are shown in Table 3. Product yields obtained from ME calculations at 298.15 K and 1 atm are shown in Table 4. For the sake of limiting computational costs, a separate ME calculation was done for each RO stereoisomer, produced from their respective parent RO₂ radical viareactions (13)–(16). Reaction of O₂ with C-centered radicals M12a, M12s, M16a, M16s, M17s, M19a, M20, and M21 was treated as a simple bimolecular sink, with a rate coefficient value of k_O2 (298 K) = 2 × 10⁻¹² cm³ molecule⁻¹ s⁻¹. Typically, N-centered radicals formed in the atmosphere may also react with O₂, NO or NO₂ (Reactions (4) and (5)). Reaction of M13 with O₂ is expected to lead to the formation of a nitrile plus HO₂,¹⁰⁰ whereas reaction with NO or NO₂ may produce a nitrosamine or a nitramine. Only reaction with O₂ was included in the ME modelling for bimolecular loss terms involving M13, but treated as a simple bimolecular sink leading to a single product with a rate coefficient value of 10⁻¹⁶ cm³ molecule⁻¹ s⁻¹, corresponding to an upper limit to rate coefficients typically observed for the reaction of N-centered radicals with O₂.⁴⁷

Table 3 Reaction barrier heights (Δ^‡E_v=0) calculated as the difference in zero-point corrected energy between TS and reactant, calculated at the L3//L2 level, for unimolecular reaction steps available to the major RO intermediates and its bond scission products. Thermal rate coefficients estimated with lowest-conformer TST (k_LC-TST) and multi-conformer TST (k_MC-TST) at 298 K, ratio of lowest-conformer partition functions (Q_TS/Q_R), and Eckart tunnelling factors (κ)

Reaction step	Description	Δ^‡Ev=0 (kcal mol^−1)	k LC-TST (s^−1)	k MC-TST (s^−1)	k MC/kLC	Q TS/QR	κ
4a-RO → M12a	β-Scission	4.17	5.78 × 10^9	9.01 × 10^9	1.56	0.91	1.17
4a-RO → M13	β-Scission	9.57	6.24 × 10^5	3.90 × 10^5	0.63	0.73	1.40
4a-RO → P4 + H	H-Elimination	16.02	6.93	5.61	0.81	0.62	—
4s-RO → M12s	β-Scission	5.04	1.47 × 10^9	1.00 × 10^9	0.68	1.06	1.10
4s-RO → M13	β-Scission	7.37	2.32 × 10^7	1.58 × 10^7	0.68	0.70	1.35
4s-RO → P4 + H	H-Elimination	15.38	48.8	33.3	0.68	1.47	—
M12a → M16a	Cyclization	7.05	6.66 × 10^6	7.43 × 10^6	1.12	0.14	1.10
M12a → M17a	N-Inversion	9.65	6.37 × 10^5	—	1	1.18	1.04
M12a → M12s	C–N torsion	3.58	1.35 × 10^10	—	1	0.87	1.05
M17a → M18a	C–N torsion	6.13	3.52 × 10^8	—	1	1.55	1.14
M12s → M16s	Cyclization	17.30	0.212	—	1	0.13	1.29
M12s → M17s	N-Inversion	11.98	7.92 × 10^3	—	1	0.73	1.05
M13 → M20	1,4 H-shift	17.73	56.4	39.6	0.70	0.26	351.24
M13 → M21	1,5 H-shift	16.22	3.04 × 10^2	2.14 × 10^2	0.70	0.51	74.58

---

Table 4 Product yields obtained from ME calculations for reaction of major RO radical products, at 298.15 K and 1 atm. Species R-O₂ correspond to the products of reaction of R with O₂. Standard deviations (±σ) calculated with the logit-transformed distribution of yields obtained with uncertainty analysis. *: combined yields of M12s-O₂ and M17s-O₂

Product	Yield (%) (±σ)
	Initial reactant: 4a-RO2	Initial reactant: 4s-RO2
M12s-O2/M17s-O2*	97.8 (^+1.0−1.4)	95.8 (^+1.6−3.3)
M13-O2	∼ 0	1.5 (^+1.5−0.8)
M16a-O2	0.4 (^+0.3−0.1)	∼ 0
M19a-O2	1.9 (^+1.3−0.8)	∼ 0
M20-O2	∼0	0.5 (^+1.1−0.3)
M21-O2	∼0	2.4 (^+1.8−1.4)

---

Results from ME calculations revealed that both 4a-RO and 4s-RO react almost solely via their respective C–C bond scission channels, producing species M12a and M12s respectively. For 4s-RO, C–N bond scission (via TS13s) is slightly competitive, with a yield of ∼4.4% at 298 K, whereas it is a negligible channel for the anti isomer (via TS13a). The channels involving C–H bond scission via TS14a and TS14s, expected to be negligible due to their much higher energy barriers, were not included in the ME treatment.

As for the fate of the subsequently produced C-centered radicals, the major loss channel of species M12a is conversion to M12s via bond torsion/H-shift, as described above, with a yield of about 97.4% at 298 K. Thus, the great majority of the RO radicals (either syn or anti) produced from imidazole + OH reacts further leading to the same C-centered radical, M12s. Once formed, M12s can undergo fast N-inversion via TS17s to produce M17s. The only competitive loss channel for M12s and M17s is reaction with O₂, yielding second-generation RO₂ radicals. These RO₂ are the major products from both RO isomers, with summed yields of about 95.7% and 97.4% from 4s-RO and 4a-RO respectively. About 85% of the second-generation RO₂ is produced from O₂-addition to M17s, whereas 15% is produced from O₂-addition to M12s, a distribution which indicates that equilibrium between M12s and M17s is approached, favoring the latter. Due to a less rigid structure, M17s is ∼0.2 kcal mol⁻¹ lower in free energy relative to M12s.

Cyclization constitutes a minor channel for M12a, with a yield of 0.4% at 298 K, and it is negligible for M12s. The H-shift channel initiated via the torsional transition state TS17a is also minor for M12a, but with a slightly higher yield (1.9% at 298 K), despite involving a reaction energy barrier ∼2.6 kcal mol⁻¹ higher than cyclization. This observation, which may be explained by vibrationally hot molecules reacting at non-thermal rates, is further discussed in Section S2 of the ESI.†

Concerning the fate of the N-centered radical M13, results indicate that its reaction with O₂ is a minor channel, with a maximum yield of only 1.5% at 298 K from the total 4s-RO, the 1,5 H-shift channel via TS21 being slightly more competitive, with a yield of 2.4%. Reaction of N-centered radicals with NO or NO₂ can lead to nitrosamines (RN–NO) or nitramines (RN–NO₂) respectively, compounds which are known to be toxic. However, these channels are expected to be negligible for 4s-RO, given that reaction of M13 with O₂ is a barely competitive sink, even when assuming an upper limit to the rate coefficient. Under this assumption, reaction with O₂ has a pseudo-first order rate coefficient of ∼500 s⁻¹, which is 20 times faster than the upper limit for the rate of reaction with NO or NO₂ (∼25 s⁻¹, assuming collision-limit k_NOx = 10⁻¹⁰ cm³ molecule⁻¹ s⁻¹ and [NO_x] = 10 ppb).

In contrast to what was seen with first-generation RO₂, the fate of alkoxyl radicals 4s-RO and 4a-RO shows a weak temperature dependence in the range analysed here (see Tables S3–S5 in the ESI†). Variation in product yields is, at most, just over 1% for 4a-RO, and ∼3% for 4s-RO.

3.3 Second generation peroxyl radicals

The reaction mechanisms which determine the fate of the second-generation RO₂ radicals are seemingly rather complex; however, they all converge to the same outcome. The investigated second-generation RO₂, as well as the unimolecular reactions available to these radicals are shown in Scheme 2. Reaction of M12s or M17s (circled in Scheme 2) with O₂ may initially lead to two different RO₂ products each, depending on which C atom of the conjugated C–N–C moiety O₂-addition occurs. O₂-addition to M12s may lead to Z2-RO₂ or Z4-RO₂, while O₂-addition to M17s may lead to E2-RO₂ or E4-RO₂. Analogous to what is observed with the first-generation RO₂, an O₂-migration step connects the Z-RO₂ isomer pair via transition state TS25 and another connects the E-RO₂ isomer pair via TS26. Since the α-peroxyl carbon in these RO₂ is sp³ hybridized, torsion about the adjacent C–N bond is possible, which allows for two additional connections by O₂-migration: TS27 connects E2-RO₂ to Z4-RO₂ and TS28 connects Z2-RO₂ to E4-RO₂. Direct connection between E4-RO₂ and Z4-RO₂ or between E2-RO₂ and Z2-RO₂ is possible only via imine N-inversion channels (via TS36 and TS37), which have energy barriers too high (>22 kcal mol⁻¹) to be important. KinBot's conformer sampling algorithm accounts for E–Z isomerism in alkenes, but the version used here does not separate those isomers in imines, and for this reason we chose to employ only the reaction search module of the code for second-generation RO₂. Conformer sampling and high-level quantum chemistry calculations were done manually.

Scheme 2 Studied reaction pathways available to the second-generation RO₂ from imidazole + OH reaction, beginning with M12s and M17s, circled. Relative zero-point-corrected energies given in kcal mol⁻¹, calculated at the ROHF-ROCCSD(T)-F12a/cc-pVDZ-F12//ωB97X-D/aug-cc-pVTZ level of theory. Major products are shown in dashed squares.

The second-generation RO₂ species correspond to significantly deeper wells than the first-generation ones, and most have access to decomposition channels with lower reaction barriers than their O₂-migration steps, so rapid equilibrium between RO₂ isomers is not achieved. Species Z2-RO₂, E2-RO₂ and Z4-RO₂ can react via HO₂-elimination with C=O double-bond formation, all yielding FMF, the first RO₂ producing the (Z) isomer and the latter two producing the (E) isomer. While HO₂-elimination in Z2-RO₂ and E2-RO₂ involves the typical concerted 1,4 H-shift/C–O bond scission, in Z4-RO₂ this reaction occurs with a 1,6 span instead, which is possible due to the presence of an imine group between the reacting OH and OO groups. Species E4-RO₂ cannot react in this manner, since the OH group is not accessible to the peroxyl group. However, E4-RO₂ has access to a low-barrier (10.2 kcal mol⁻¹, via TS26) O₂-migration, connecting it to E2-RO₂, and granting the former indirect access to the latter's HO₂-elimination channel. Reaction channels involving HO₂-elimination with C=N double-bond formation were also found for Z2-RO₂, E4-RO₂ and Z4-RO₂, but these have significantly higher reaction barriers (see Fig. S4, ESI†) and are therefore expected to be negligible.

The reaction pathways available to the second-generation RO₂ are included in the RRKM-ME modelling for alkoxyl radicals 4s-RO and 4a-RO, connecting the PESs in Fig. 4 and Scheme 2. Results from ME calculations, shown in Table 5, indicate that virtually all of the second-generation RO₂ radicals formed in the atmosphere undergo HO₂-elimination to yield FMF, about 91% of which formed as the (E) isomer and about 9% as the (Z) isomer. This ratio between E–Z isomers is largely insensitive to the temperature or to the stereochemistry of the precursor RO radical. Nevertheless, the energy barrier of imine N-inversion in FMF is ∼15 kcal mol⁻¹, and inter-conversion between these isomers is relatively fast, with a rate coefficient of ∼100 s⁻¹ at 298 K. Thus, under these conditions, an equilibrium favoring the more stable (Z)-FMF is expected to be established on the time-scales relevant to the fate of this product. A (Z)/(E) ratio of 57/43 is expected at equilibrium, at 298.15 K.

Table 5 Product yields obtained from ME calculations for reaction of second-generation RO₂ radical products, at 298.15 K and 1 atm. Standard deviations (±σ) calculated with the logit-transformed distribution of yields obtained with uncertainty analysis

Product	Yield (%) (±σ)
	Initial reactant: 4a-RO2	Initial reactant: 4s-RO2
E-FMF	91.2 (^+3.4−7.0)	91.2 (^+4.3−7.4)
Z-FMF	8.8 (^+6.4−3.3)	8.8 (^+5.4−4.2)

---

3.4 Multiple conformer effects

In comparison to the standard lowest-conformer TST, the inclusion of multiple-conformer effects with MC-TST lowers the estimated rate coefficients of most of the studied reactions by a similar factor, as seen from Tables 1, 3, and Table S6 (ESI†). Thus, relying on error cancellation, we expect that the neglect of multiple conformations in the RRKM-ME calculations introduces small errors, which are well within the variation ranges assumed during uncertainty analysis of product distributions. The exceptions to this trend are 4a-RO C–C bond scission (4a-RO → M12a), M12a cyclization (M12a → M16a), O₂-migration via TS25 (Z2-RO₂ → E2-RO₂), and HO₂-elimination via TS31. However, apart from cyclization, these reactions are either the only competitive channels for their respective reactants, or much slower than competing channels, so that a slight increase in rate coefficients is not expected to affect product yields. Treatment of multiple-conformations leads to a small increase in the rate coefficient for M12a cyclization (a factor of 1.12) and this reaction competes only with single-conformer channels.

3.5 Uncertainties in product yields

Uncertainties in product yields are given as standard deviation ranges in Fig. 3 and Table 2 for first-generation RO₂, in Table 4 for the alkoxyl radical products, and in Table 5 for second-generation RO₂. The first stage of the oxidation mechanism, involving reactions with the first-generation RO₂, is by far the most uncertain when it comes to predicting product yields. The largest sources of uncertainty are likely to be the rates of unimolecular reaction pathways, more specifically, the HO₂-elimination steps via TS5a or TS5s and the 1,5 H-shift step via TS7s. The calculated energy barrier heights for these three reactions are very similar (0.8, 1.2 and 0.1 kcal mol⁻¹, relative to 5-OH + O₂, for TS5a, TS5s, and TS7s respectively), and therefore a variation of ±0.5 kcal mol⁻¹ to the energy of their transition states can completely alter their order of importance. This is clearly observed from the broad product yield uncertainty ranges under conditions in which unimolecular channels dominate, such as low [NO] and high temperatures.

Smaller, but still present, the uncertainty associated with bimolecular RO₂ channels is mostly due to the ±0.25 kcal mol⁻¹ variation in the energy of RO₂ wells, altering their relative concentrations at equilibrium (or near-equilibrium), than it is to the variation in k_NO rate coefficients. The rate coefficients for RO₂ + NO reactions are not very sensitive to the identity of the R group,⁸⁸ hence a small variation range was applied to k_NO during the analysis (see Methods section). Surprisingly, variation in the rate coefficients of O₂-addition reactions (k_O2), which can have values spanning different orders of magnitude, and are unknown for the systems investigated here, had a minimal impact in the overall product yield uncertainty. As discussed in Section 3.1, the equilibrium established (or approached) between the four RO₂ isomers nullifies the uncertainties associated with their formation rates. A test RRKM-ME uncertainty analysis at 298 K and [NO] = 100 ppt, applying random variations only to 5-OH + O2 k_O2 values, produces negligible uncertainties in product yields, with a maximum standard deviation of +0.03%/−0.05%.

Concerning the latter stage of the oxidation mechanism, the fate of the alkoxyl radical products is much less uncertain, and largely independent of the ratio of stereoisomers produced – 97.8 (^+1.0_−1.4)% of the 4a-RO and 95.8 (^+1.6_−3.3)% of the 4s-RO formed react leading to M12s (and M17s), ∼100% of which produces FMF upon O₂-addition/HO₂-elimination. The ratio of E–Z isomers of FMF produced immediately after oxidative chain termination is somewhat uncertain, however, as discussed in Section 3.3, equilibrium between these isomers is expected to be rapidly established, favoring Z-FMF.

3.6 Comparison to oxidation of similar compounds

Similar to aliphatic N-containing heterocycles whose RO₂ chemistry has been studied, such as pyrrolidine,⁵⁹ piperazine,⁵⁴ piperidine,¹⁰¹ and 1,2,3,6-tetrahydropyrazine (THP),⁴⁸ HO₂-elimination with imine formation is also an important reaction channel during imidazole OH-initiated gas-phase oxidation. However, unlike with these aliphatic species, this reaction is not the major channel for imidazole, due to a few structural and mechanistic differences, as outlined below.

Relative to their respective α-amino RO₂ reactants, the reaction barrier heights of HO₂-elimination channels in the oxidation of imidazole, pyrrolidine, piperazine, piperidine and THP are all relatively low, with calculated Δ^‡E values ranging from 14.4 to 19.2 kcal mol⁻¹. As it seems, a more important factor controlling differences in reactivity here is the well depth of those RO₂. The RO₂ derivatives of saturated species (pyrrolidine, piperazine, and piperidine), and the β-N-imino RO₂ derivative of THP correspond to deep wells (−33.9, −35.4, −35.1, and −33.3 kcal mol⁻¹ respectively), and their respective HO₂-elimination barriers are submerged far below the entrance energy level (R˙ + O₂), meaning that excess energy imparted to the system during the preceding O₂-addition step may propel these reactions to happen at non-thermal rates. On the other hand, since imidazole is an aromatic species whose oxidation is initiated preferentially via OH-addition, rather than via H-abstraction, subsequent O₂-addition reactions are bound to be accompanied by loss of resonance stabilization, producing RO₂ radicals corresponding to shallower wells. Thus, reaction barriers to HO₂-elimination are not submerged for imidazole, albeit nearly so, and non-thermal reactions are less likely. THP is not completely saturated, bearing an endocyclic imine group, and its oxidation shares some characteristics with that of imidazole. The α-C-imino RO₂ derivative of THP, whose parent alkyl radical is also resonance stabilized, corresponds to a shallow well (−14.8 kcal mol⁻¹), with its available HO₂-elimination channel having a slightly submerged barrier (−0.4 kcal mol⁻¹).⁴⁸

Another important factor is that O₂ attack at the other available site in a parent resonance stabilized alkyl radical can compete with the α-amino substitution required for HO₂-elimination. In imidazole, the effect of this competition is partly counteracted by the availability of fast O₂-migration channels, which promote an equilibrium between the γ-amino RO₂ (4-RO₂) and the α-amino RO₂ (2-RO₂). This means that the γ-amino RO₂ species also have access to the HO₂-elimination channel, but surpassing a higher reaction energy barrier, since these RO₂ correspond to deeper wells. In the oxidation of THP, however, no such competition is possible. In the α-C-imino alkyl radical, resonance occurs at a conjugated C–C=N moiety, and O₂ attack on the nitrogen cannot compete with C-addition.

Finally, the reactivity of imidazole's RO₂ derivatives is significantly affected by OH-substitution, which is absent from the aliphatic heterocycles considered in this comparison. The OH group not only grants access to a slightly faster reaction channel, the 1,5 H-shift via TS7s, but it can alter the energetics of competing pathways, via intramolecular H-bonds and steric effects, as discussed in Section 3.1. Moreover, OH-substitution adds further complexity to the oxidation mechanism by introducing stereoisomeric effects.

Many of the characteristics of imidazole's oxidation mechanism are also observed for homocyclic aromatics, such as a OH-addition-initiated oxidation, shallow RO₂ wells, and stereoisomeric effects.^102,103 The major differences displayed by oxidation of homocyclic aromatics, such as toluene¹⁰² and benzene,¹⁰³ are the absence of competitive HO₂-elimination channels, which are in general slower if involving H-shift from a C atom, and the predominance of endo-cyclization reactions leading to bicyclic endoperoxides. Endo-cyclization transition states were found for imidazole's RO₂ derivatives during reaction search, but the corresponding energy barriers exceed the 20 kcal mol⁻¹ cutoff at L1. Manual calculations at the L3//L2 level revealed that these reactions indeed cannot be competitive for imidazole's products (see Fig. S1 in the ESI†), with barriers >12 kcal mol⁻¹ relative to 5-OH + O₂. The high barriers may be due to the extra ring-strain brought to an already strained five-membered cycle, and the placement the radical center on the less-accommodating nitrogen bridge atom. That being said, a similarly high barrier was found for analogous reactions during the OH-initiated oxidation of THP, a less strained six-membered heterocycle,⁴⁸ possibly indicating that the channel is unfavorable for N-containing RO₂ mainly due to destabilization brought by formation of the N-centered radical.

Perhaps a unique mechanistic feature of the oxidation of imidazole is the availability of O₂-migration channels. The molecular structure of the O₂-migration transition state TS4a is depicted in Fig. 5. To the best of our knowledge, no such reaction channel has been reported before in mechanistic studies of atmospheric oxidation reactions. As of yet, it is unclear whether these reactions are possible for other N-containing aromatics, but the availability of O₂-migrations for the second-generation RO₂ studied here indicates that these reactions should only require the formation of α-N-imino RO₂ radicals. In the OH-initiated oxidation of oxazole,¹⁰⁴ for example, O₂-migration is presumably possible and competitive. In the oxidation of pyrrole,¹⁰⁵ however, a α-N-imino RO₂ cannot be formed, and O₂-migration may not be available. Despite this, drawing from knowledge of imidazole's oxidation mechanism, we would expect pyrrole's RO₂ derivatives to react preferentially via similar unimolecular channels: HO₂-elimination with imine formation and 1,5 H-shift from the alcohol group. A previous computational study¹⁰⁵ has investigated the fate of the first-generation RO₂ from pyrrole's OH-initiated oxidation, however, neither of these channels were reported. We speculate that these reactions do exist for pyrrole, and could very well be the major loss channels, alongside bimolecular RO₂ + NO reactions.

Fig. 5 Structure of O₂-migration transition state TS4a. C atoms in gray, H atoms in white, O atoms in red, and N atoms in blue.

3.7 Atmospheric implications

Our calculations predict that over 95% of imidazole oxidized by OH radicals in the atmospheric gas-phase, via addition to position C5, leads to oxidative chain termination within two RO₂ generations, producing two closed-shell products, the cyclic diimine 4H-4ol and the acyclic amidine FMF. A scheme depicting the overall mechanism of imidazole's OH-initiated oxidation via adduct 5-OH, simplified to include only important reaction steps, is shown in Scheme 3. Relative yields of 4H-4ol and FMF are relatively uncertain and sensitive to NO_x concentrations, but the latter species is consistently predicted to be the major product, with ∼66%, ∼85%, and ∼95% yield under 10 ppt, 100 ppt and 1 ppb of NO respectively, at 298 K. The minor product 4H-4ol is predicted to be formed with ∼34%, ∼12%, and ∼2% yield under 10 ppt, 100 ppt and 1 ppb of NO respectively, at 298 K. The reaction channels leading to N-centered radical M13 are minor, with a maximum yield of ∼2% at 298 K and 1 ppb NO, considering the combined uncertainties in product yields. Moreover, subsequent unimolecular reactions available to M13 are fast, leading to C-centered radicals M20 and M21. Thus, we expect that formation of nitrosamines or nitramines via reaction of M13 with NO or NO₂ to be negligible. Yield of autoxidation channels, leading to highly functionalized products, is also expected to be negligible. The C-centered radical with the highest yield (M19a), predicted to be formed with a maximum of ∼2% yield at 298 K and 1 ppb NO, has access to a potentially predominant O₂-addition/HO₂-elimination channel leading to FMF. It is apparent that no matter how complex imidazole's oxidation may be, reaction pathways seem to converge to the formation of the dicarbonyl amidine FMF. In fact, this compound was also proposed as a bio-degradation product of ornidazole,¹⁰⁶ an imidazole derivative used as antibiotic.

Scheme 3 Overall mechanism of OH-initiated oxidation of imidazole via adduct 5-OH. Dashed boxes indicate major products. Minor products and less important reaction steps omitted for clarity.

Concerning the atmospheric implications of the major product FMF, we speculate on its fate in the atmosphere. Upon reaction with another OH radical in the gas-phase, FMF could produce an isocyanate (R–N=C=O) via three different initial H-abstraction pathways. Initial H-abstraction by OH radicals from the amide (–NH–CH=O) carbon^107,108 or nitrogen, followed by O₂-addition/HO₂-elimination, could both lead to N-formylformimidoyl isocyanate (O=CH–N=CH–NCO). Similarly, initial H-abstraction/O₂-addition at the =N–CH=O group could lead to the exact same isocyanate, but via a HO₂-elimination transition state with a larger 1,6 span, like that observed for Z4-RO₂ (TS29z). If, however, the amidine hydrogen is abstracted, and assuming HO₂-elimination is both available and competitive for a R–N=C(OO)–NH–R type radical, it could lead to a carbodiimide (R–N=C=N–R′). Nothing is known about carbodiimides in the atmosphere, however, formation of isocyanates has indeed been reported to be possible from the OH-oxidation of amides^107,108 and imines.⁴⁸ Given that isocyanates are typically toxic, their potential formation from FMF oxidation warrants further investigation. Nevertheless, FMF is likely to be highly water soluble due to a presence multiple H-bond acceptors in its structure, and therefore it may partition into aerosol particles efficiently. Assuming a condensation sink rate typical of more pristine environments (10⁻⁴ s⁻¹), and an imidazole concentration of 183 ppt, we estimate a steady-state gas-phase concentration for FMF of up ∼46 ppt. In aerosols, FMF may participate in aqueous-phase reactions with potentially many outcomes. One such reaction, which is perhaps the most important, is amidine hydrolysis leading to formamide plus N-formylformamide. Speculating further, another possibility is the formation of 1,3,5-triazine, a six-membered aromatic compound with three alternated nitrogens, via a multi-step reaction analogous to imidazole's own proposed formation route in aerosols, from reaction of glyoxal with ammonia.¹⁹

4 Conclusions

In this work we used quantum chemistry and theoretical reaction kinetics to investigate the atmospheric gas-phase chemistry of imidazole. We employed an automated kinetics workflow code to explore and characterize the stationary points on the potential energy surface relevant for imidazole's atmospheric oxidation initated by OH-radical addition to position C5, and propagated onwards by RO₂ chemistry. We estimated product distributions by assembling and solving a master equation under a range of different atmospheric conditions. The explored mechanism included pathways available to the alkoxyl radical products from RO₂ + NO reactions, and the subsequently formed second-generation RO₂ radicals. Our calculations predict a complex, branched oxidation mechanism, but with converging pathways leading to two major closed shell species, the cyclic diimine 4H-imidazol-4-ol (4H-4ol) and the acyclic N,N′-diformylformamidine (FMF). Relative yields of 4H-4ol and FMF are strongly dependent on NO concentrations, but the latter species consistently appears as the main product.

Overall, we estimate the branching percentage for FMF and 4H-4ol formation from imidazole's oxidation initiated via OH-addition to position C5 to be 66 : 34, 85 : 12, and 95 : 2 under 10 ppt, 100 ppt and 1 ppb of NO respectively, at 298 K. Under more pristine atmospheric conditions ([NO] = 10 ppt) and 298 K, major reaction channels for first-generation peroxyl radicals are HO₂-elimination leading to 4H-4ol (∼34%), a 1,5 H-shift channel followed by β-scission and OH-elimination leading to FMF (∼47%), and RO₂ + NO reactions leading to alkoxyl radicals (∼19%). Under more polluted conditions ([NO] = 1 ppb), however, bimolecular reactions with NO dominate (∼96%). The major alkoxyl radical products may further react, the majority (∼96–98%) of which undergoing C–C bond scission, followed by O₂-addition/HO₂-elimination to yield FMF.

Our calculations also revealed the availability of O₂-migration reaction pathways for α-N-imino peroxyl radicals. This type of reaction involves a concerted attack of the peroxyl group on the C-terminus of the imine and C–O(O) bond scission, exchanging the positions of the two groups:

RR′C=N–C(OO)R′′R′′′ ⇄ RR′(OO)C–N=CR′′R′′′ (17)

For imidazole's first-generation peroxyl radicals, these reactions are fast (10²–10⁵ s⁻¹), and serve to establish equilibrium between the connected RO₂ isomers. If available, this mechanism may be important in the atmospheric oxidation of other unsaturated organic nitrogen compounds.

This work provides a detailed description of imidazole chemistry in the atmospheric gas-phase, and may serve as a blueprint to understand the fate of important substituted imidazole compounds, as well as other N-containing aromatics.

Author contributions

Thomas Golin Almeida: conceptualization, methodology, software, validation, formal analysis, investigation, writing – original draft, visualization; Carles Martí: software, writing – review & editing; Theo Kurtén: funding acquisition, supervision, writing – review & editing; Judit Zádor: conceptualization, methodology, software, supervision, writing – review & editing; Sommer L. Johansen: conceptualization, supervision, funding acquisition, writing – review & editing.

Data availability

Research data files are available in the CERN Zenodo repository: https://doi.org/10.5281/zenodo.13143561. The raw data (energies, cartesian coordinates, and vibrational frequencies), required to reproduce the results of this work, is available in the ESI.†

Conflicts of interest

There are no conflicts to declare.

Acknowledgements

The authors thank Arkke Eskola, Timo Pekkanen, Siddharth Iyer, and Rashid Valiev for insightful discussions on the topic of this work. The work of JZ and CM was supported by the U.S. Department of Energy, Office of Science, Office of Basic Energy Sciences, Division of Chemical Sciences, Geosciences and Biosciences under the Gas Phase Chemical Physics program. SLJ was supported by the Laboratory Directed Research and Development program at Sandia National Laboratories. The work of TGA and TK was funded by the University of Helsinki (CHEMS) and the Research Council of Finland (center of excellence VILMA, grant no. 346369). Computational resources were provided by CSC – IT Center for Science, in Finland. This article has been authored by employees of National Technology & Engineering Solutions of Sandia, LLC under contract no. DE-NA0003525 with the U.S. Department of Energy (DOE). The employees co-own right, title and interest in and to the article and are responsible for its contents. The United States Government retains and the publisher, by accepting the article for publication, acknowledges that the United States Government retains a non-exclusive, paid-up, irrevocable, world-wide license to publish or reproduce the published form of this article or allow others to do so, for United States Government purposes. The DOE will provide public access to these results of federally sponsored research in accordance with the DOE Public Access Plan https://www.energy.gov/downloads/doe-public-access-plan.

Notes and references

1. X. Ge, A. S. Wexler and S. L. Clegg, Atmos. Environ., 2011, 45, 524–546.
2. E. D. Pusfitasari, J. Ruiz-Jimenez, I. Heiskanen, M. Jussila, K. Hartonen and M.-L. Riekkola, Sci. Total Environ., 2022, 808, 152011.
3. E. D. Pusfitasari, J. Ruiz-Jimenez, A. Tiusanen, M. Suuronen, J. Haataja, Y. Wu, J. Kangasluoma, K. Luoma, T. Petäjä, M. Jussila, K. Hartonen and M.-L. Riekkola, Atmos. Chem. Phys., 2023, 23, 5885–5904.
4. M. Teich, D. van Pinxteren, S. Kecorius, Z. Wang and H. Herrmann, Environ. Sci. Technol., 2016, 50, 1166–1173.
5. Y. Takao, T. Atarashi, T. Kubo, M. Nagae and T. Nakayama, Asian J. Atmos. Environ., 2019, 13, 259–265.
6. M. Teich, M. Schmidtpott, D. van Pinxteren, J. Chen and H. Herrmann, J. Sep. Sci., 2020, 43, 577–589.
7. C. He, H. Wang, D. Gong, S. Lv, G. Wu, R. Wang, Y. Chen, Y. Ding, Y. Li and B. Wang, Sci. Total Environ., 2022, 806, 150804.
8. M. T. Martins-Costa, J. M. Anglada, J. S. Francisco and M. F. Ruiz-López, Chem. Sci., 2022, 13, 2624–2631.
9. M. Hallquist, J. C. Wenger, U. Baltensperger, Y. Rudich, D. Simpson, M. Claeys, J. Dommen, N. Donahue, C. George, A. Goldstein, J. F. Hamilton, H. Herrmann, T. Hoffmann, Y. Iinuma, M. Jang, M. E. Jenkin, J. L. Jimenez, A. Kiendler-Scharr, W. Maenhaut, G. McFiggans, Th. F. Mentel, A. Monod, A. S. H. Prévôt, J. H. Seinfeld, J. D. Surratt, R. Szmigielski and J. Wildt, Atmos. Chem. Phys., 2009, 9, 5155–5236.
10. T. Drugé, P. Nabat, M. Mallet, M. Michou, S. Rémy and O. Dubovik, Atmos. Chem. Phys., 2022, 22, 12167–12205.
11. D. L. Nelson and M. M. Cox, Lehninger Principles of Biochemistry, W.H. Freeman, New York, NY, 4th edn, 2004.
12. V. B. Silva and E. S. Orth, Quím. Nova, 2021, 44, 318–333.
13. M. Galloway, P. Chhabra, A. Chan, J. Surratt, R. Flagan, J. Seinfeld and F. Keutsch, Atmos. Chem. Phys., 2009, 9, 3331–3345.
14. D. O. D. Haan, A. L. Corrigan, K. W. Smith, D. R. Stroik, J. J. Turley, F. E. Lee, M. A. Tolbert, J. L. Jimenez, K. E. Cordova and G. R. Ferrell, Environ. Sci. Technol., 2009, 43, 2818–2824.
15. N. Sedehi, H. Takano, V. A. Blasic, K. A. Sullivan and D. O. De Haan, Atmos. Environ., 2013, 77, 656–663.
16. X. Lian, G. Zhang, Y. Yang, Q. Lin, Y. Fu, F. Jiang, L. Peng, X. Hu, D. Chen, X. Wang, P. Peng, G. Sheng and X. Bi, Environ. Sci. Technol. Lett., 2020, 8, 9–15.
17. G. T. Drozd, K. S.-M. Brown and H. Q. Karp, ACS Earth Space Chem., 2022, 6, 1772–1781.
18. J. Kua, H. E. Krizner and D. O. De Haan, J. Phys. Chem. A, 2011, 115, 1667–1675.
19. V. P. Tuguldurova, A. V. Fateev, O. K. Poleshchuk and O. V. Vodyankina, Phys. Chem. Chem. Phys., 2019, 21, 9326–9334.
20. H. Debus, Justus Liebigs Ann. Chem., 1858, 107, 199–208.
21. A. Laskin, J. S. Smith and J. Laskin, Environ. Sci. Technol., 2009, 43, 3764–3771.
22. Y. Wang, H. Deng, P. Li, J. Xu, G. Loisel, H. Pang, X. Xu, X. Li and S. Gligorovski, Environ. Sci. Technol. Lett., 2022, 9, 493–500.
23. M. van Pinxteren, K. W. Fomba, D. van Pinxteren, N. Triesch, E. H. Hoffmann, C. H. Cree, M. F. Fitzsimons, W. von Tümpling and H. Herrmann, Atmos. Environ., 2019, 203, 183–195.
24. T. Moise, J. M. Flores and Y. Rudich, Chem. Rev., 2015, 115, 4400–4439.
25. A. Laskin, J. Laskin and S. A. Nizkorodov, Chem. Rev., 2015, 115, 4335–4382.
26. K. Z. Aregahegn, B. Nozière and C. George, Faraday Discuss., 2013, 165, 123–134.
27. W.-Y. Li, X. Li, S. Jockusch, H. Wang, B. Xu, Y. Wu, W. G. Tsui, H.-L. Dai, V. F. McNeill and Y. Rao, J. Phys. Chem. A, 2016, 120, 9042–9048.
28. S. Szopa, V. Naik, B. Adhikary, P. Artaxo, T. Berntsen, W. Collins, S. Fuzzi, L. Gallardo, A. Kiendler-Scharr, Z. Klimont, H. Liao, N. Unger and P. Zanis, in Short-Lived Climate Forcers, ed. V. Masson-Delmotte, P. Zhai, A. Pirani, S. Connors, C. Péan, S. Berger, N. Caud, Y. Chen, L. Goldfarb, M. Gomis, M. Huang, K. Leitzell, E. Lonnoy, J. Matthews, T. Maycock, T. Waterfield, O. Yelekçi, R. Yu and B. Zhou, Cambridge University Press, Cambridge, United Kingdom and New York, NY, USA, 2021, pp. 817–922.
29. D. Spracklen, J. Jimenez, K. Carslaw, D. Worsnop, M. Evans, G. Mann, Q. Zhang, M. Canagaratna, J. Allan, H. Coe, G. McFiggans, A. Rap and P. Forster, Atmos. Chem. Phys., 2011, 11, 12109–12136.
30. T. Felber, T. Schaefer and H. Herrmann, J. Phys. Chem. A, 2019, 123, 1505–1513.
31. F. Witte and C. Zetzsch, International Symposium on Gas Kinetics (9th) Held in Bordeaux, France on 20-25 July 1986. Abstracts, 1986.
32. A. Samunl and P. Neta, J. Phys. Chem., 1973, 77, 1629–1635.
33. K. Bansal and R. Sellers, J. Phys. Chem., 1975, 79, 1775–1780.
34. G. Lassmann, L. A. Eriksson, F. Himo, F. Lendzian and W. Lubitz, J. Phys. Chem. A, 1999, 103, 1283–1290.
35. J. Llano and L. A. Eriksson, J. Phys. Chem. B, 1999, 103, 5598–5607.
36. Z. Safaei, A. Shiroudi, E. Zahedi and M. Sillanpää, Phys. Chem. Chem. Phys., 2019, 21, 8445–8456.
37. T. V.-T. Mai and L. K. Huynh, Phys. Chem. Chem. Phys., 2019, 21, 21162–21165.
38. J. J. Orlando and G. S. Tyndall, Chem. Soc. Rev., 2012, 41, 6294–6317.
39. D. R. Glowacki and M. J. Pilling, ChemPhysChem, 2010, 11, 3836–3843.
40. J. Zádor, C. A. Taatjes and R. Fernandes, Prog. Energy Combust. Sci., 2011, 37, 371–421.
41. S. Richters, M. Pfeifle, M. Olzmann and T. Berndt, Chem. Commun., 2017, 53, 4132–4135.
42. J. D. Crounse, L. B. Nielsen, S. Jørgensen, H. G. Kjaergaard and P. O. Wennberg, J. Phys. Chem. Lett., 2013, 4, 3513–3520.
43. J. D. Savee, E. Papajak, B. Rotavera, H. Huang, A. J. Eskola, O. Welz, L. Sheps, C. A. Taatjes, J. Zádor and D. L. Osborn, Science, 2015, 347, 643–646.
44. F. Bianchi, T. Kurtén, M. Riva, C. Mohr, M. P. Rissanen, P. Roldin, T. Berndt, J. D. Crounse, P. O. Wennberg, T. F. Mentel, J. Wildt, H. Junninen, T. Jokinen, M. Kulmala, D. Worsnop, J. A. Thornton, N. Donahue, H. G. Kjaergaard and M. Ehn, Chem. Rev., 2019, 119, 3472–3509.
45. C. R. Lindley, J. G. Calvert and J. H. Shaw, Chem. Phys. Lett., 1979, 67, 57–62.
46. Y. Tang and C. J. Nielsen, Atmos. Environ., 2012, 55, 185–189.
47. C. Liu, F. Ma, J. Elm, Z. Fu, W. Tang, J. Chen and H.-B. Xie, Chemosphere, 2019, 237, 124411.
48. T. Golin Almeida and T. Kurtén, ACS Earth Space Chem., 2022, 6, 2453–2464.
49. Y. G. Lazarou, K. G. Kambanis and P. Papagiannakopoulos, J. Phys. Chem., 1994, 98, 2110–2115.
50. S. Ravnum, E. Rundén-Pran, L. Fjellsbø and M. Dusinska, Regul. Toxicol. Pharmacol., 2014, 69, 250–255.
51. N. Hyttinen, H. C. Knap, M. P. Rissanen, S. Jørgensen, H. G. Kjaergaard and T. Kurten, J. Phys. Chem. A, 2016, 120, 3588–3595.
52. G. da Silva, B. B. Kirk, C. Lloyd, A. J. Trevitt and S. J. Blanksby, J. Phys. Chem. Lett., 2012, 3, 805–811.
53. M. P. Rissanen, A. J. Eskola, T. L. Nguyen, J. R. Barker, J. Liu, J. Liu, E. Halme and R. S. Timonen, J. Phys. Chem. A, 2014, 118, 2176–2186.
54. Z. Ren and G. da Silva, ACS Earth Space Chem., 2019, 3, 2510–2516.
55. C. T. Cheng, M. N. Chan and K. R. Wilson, J. Phys. Chem. A, 2016, 120, 5887–5896.
56. M. E. Jenkin, R. Valorso, B. Aumont and A. R. Rickard, Atmos. Chem. Phys., 2019, 19, 7691–7717.
57. J. J. Orlando, G. S. Tyndall and T. J. Wallington, Chem. Rev., 2003, 103, 4657–4690.
58. M. E. Jenkin, R. Valorso, B. Aumont, A. R. Rickard and T. J. Wallington, Atmos. Chem. Phys., 2018, 18, 9297–9328.
59. J. Zádor, C. Martí, R. Van de Vijver, S. L. Johansen, Y. Yang, H. A. Michelsen and H. N. Najm, J. Phys. Chem. A, 2023, 127, 565–588.
60. M. E. Jenkin, R. Valorso, B. Aumont, A. R. Rickard and T. J. Wallington, Atmos. Chem. Phys., 2018, 18, 9329–9349.
61. A. W. Birdsall and M. J. Elrod, J. Phys. Chem. A, 2011, 115, 5397–5407.
62. R. Van de Vijver and J. Zádor, Comput. Phys. Commun., 2020, 248, 106947.
63. J. Zádor, R. Van de Vijver, A. Dewyer and C. Martí, KinBot, 2023, https://github.com/zadorlab/KinBot.git.
64. A. D. Becke, J. Chem. Phys., 1993, 98, 1372–1377.
65. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785.
66. R. Ditchfield, W. J. Hehre and J. A. Pople, J. Chem. Phys., 1971, 54, 724–728.
67. W. J. Hehre, R. Ditchfield and J. A. Pople, J. Chem. Phys., 1972, 56, 2257–2261.
68. P. C. Hariharan and J. A. Pople, Theor. Chim. Acta, 1973, 28, 213–222.
69. T. Clark, J. Chandrasekhar, G. W. Spitznagel and P. V. R. Schleyer, J. Comput. Chem., 1983, 4, 294–301.
70. J.-D. Chai and M. Head-Gordon, Phys. Chem. Chem. Phys., 2008, 10, 6615–6620.
71. T. H. Dunning Jr, J. Chem. Phys., 1989, 90, 1007–1023.
72. R. A. Kendall, T. H. Dunning Jr and R. J. Harrison, J. Chem. Phys., 1992, 96, 6796–6806.
73. K. H. Møller, R. V. Otkjær, N. Hyttinen, T. Kurtén and H. G. Kjaergaard, J. Phys. Chem. A, 2016, 120, 10072–10087.
74. G. Knizia, T. B. Adler and H.-J. Werner, J. Chem. Phys., 2009, 130, 054104.
75. K. A. Peterson, T. B. Adler and H.-J. Werner, J. Chem. Phys., 2008, 128, 084102.
76. A. H. Larsen, J. J. Mortensen, J. Blomqvist, I. E. Castelli, R. Christensen, M. Dułak, J. Friis, M. N. Groves, B. Hammer, C. Hargus, E. D. Hermes, P. C. Jennings, P. B. Jensen, J. Kermode, J. R. Kitchin, E. L. Kolsbjerg, J. Kubal, K. Kaasbjerg, S. Lysgaard, J. B. Maronsson, T. Maxson, T. Olsen, L. Pastewka, A. Peterson, C. Rostgaard, J. Schiøtz, O. Schütt, M. Strange, K. S. Thygesen, T. Vegge, L. Vilhelmsen, M. Walter, Z. Zeng and K. W. Jacobsen, J. Phys.: Condens. Matter, 2017, 29, 273002.
77. M. J. Frisch, G. W. Trucks, H. B. Schlegel, G. E. Scuseria, M. A. Robb, J. R. Cheeseman, G. Scalmani, V. Barone, G. A. Petersson, H. Nakatsuji, X. Li, M. Caricato, A. V. Marenich, J. Bloino, B. G. Janesko, R. Gomperts, B. Mennucci, H. P. Hratchian, J. V. Ortiz, A. F. Izmaylov, J. L. Sonnenberg, D. Williams-Young, F. Ding, F. Lipparini, F. Egidi, J. Goings, B. Peng, A. Petrone, T. Henderson, D. Ranasinghe, V. G. Zakrzewski, J. Gao, N. Rega, G. Zheng, W. Liang, M. Hada, M. Ehara, K. Toyota, R. Fukuda, J. Hasegawa, M. Ishida, T. Nakajima, Y. Honda, O. Kitao, H. Nakai, T. Vreven, K. Throssell, J. A. Montgomery, Jr., J. E. Peralta, F. Ogliaro, M. J. Bearpark, J. J. Heyd, E. N. Brothers, K. N. Kudin, V. N. Staroverov, T. A. Keith, R. Kobayashi, J. Normand, K. Raghavachari, A. P. Rendell, J. C. Burant, S. S. Iyengar, J. Tomasi, M. Cossi, J. M. Millam, M. Klene, C. Adamo, R. Cammi, J. W. Ochterski, R. L. Martin, K. Morokuma, O. Farkas, J. B. Foresman and D. J. Fox, Gaussian 16, Revision C.02, Gaussian, Inc., Wallingford, CT, 2016.
78. H.-J. Werner and P. J. Knowles, et al., MOLPRO, version 2022.2, a package of ab initio programs.
79. H.-J. Werner, P. J. Knowles, G. Knizia, F. R. Manby and M. Schütz, Wiley Interdiscip. Rev.: Comput. Mol. Sci., 2012, 2, 242–253.
80. H.-J. Werner, P. J. Knowles, F. R. Manby, J. A. Black, K. Doll, A. Heßelmann, D. Kats, A. Köhn, T. Korona and D. A. Kreplin, et al. , J. Chem. Phys., 2020, 152, 144107.
81. Spartan'20, Wavefunction, Inc., Irvine, CA, 2020.
82. T. A. Halgren, J. Comput. Chem., 1999, 20, 730–748.
83. H. Eyring, J. Chem. Phys., 1935, 3, 107–115.
84. D. R. Glowacki, C.-H. Liang, C. Morley, M. J. Pilling and S. H. Robertson, J. Phys. Chem. A, 2012, 116, 9545–9560.
85. C. Eckart, Phys. Rev., 1930, 35, 1303.
86. J. W. Davies, N. J. B. Green and M. J. Pilling, Chem. Phys. Lett., 1986, 126, 373–379.
87. S. H. Robertson, M. J. Pilling, D. L. Baulch and N. J. B. Green, J. Phys. Chem., 1995, 99, 13452–13460.
88. M. E. Jenkin, R. Valorso, B. Aumont and A. R. Rickard, Atmos. Chem. Phys., 2019, 19, 7691–7717.
89. C. W. Gao, J. W. Allen, W. H. Green and R. H. West, Comput. Phys. Commun., 2016, 203, 212–225.
90. L. S. Tee, S. Gotoh and W. E. Stewart, Ind. Eng. Chem. Fundam., 1966, 5, 356–363.
91. K. G. Joback and R. C. Reid, Chem. Eng. Commun., 1987, 57, 233–243.
92. L. Vereecken, T. L. Nguyen, I. Hermans and J. Peeters, Chem. Phys. Lett., 2004, 393, 432–436.
93. A. Karton, E. Rabinovich, J. M. Martin and B. Ruscic, J. Chem. Phys., 2006, 125, 114108.
94. A. Karton, S. Daon and J. M. Martin, Chem. Phys. Lett., 2011, 510, 165–178.
95. K. Madajczyk, P. S. Żuchowski, F. Brzek, Ł. Rajchel, D. Kdziera, M. Modrzejewski and M. Hapka, J. Chem. Phys., 2021, 154, 134106.
96. A. V. Copan, H. F. Schaefer III and J. Agarwal, Mol. Phys., 2015, 113, 2992–2998.
97. R. R. Valiev and T. Kurten, R. Soc. Open Sci., 2020, 7, 200521.
98. H. S. Johnston and J. Heicklen, J. Phys. Chem., 1962, 66, 532–533.
99. T. T. Pekkanen, G. Lendvay, M. Döntgen, R. S. Timonen and A. J. Eskola, Proc. Combust. Inst., 2023, 39, 189–197.
100. A. J. C. Bunkan, Y. Tang, S. R. Sellevag and C. J. Nielsen, J. Phys. Chem. A, 2014, 118, 5279–5288.
101. W. Tan, L. Zhu, T. Mikoviny, C. J. Nielsen, A. Wisthaler, B. D’Anna, S. Antonsen, Y. Stenstrøm, N. J. Farren, J. F. Hamilton, G. A. Boustead, T. Ingham and D. E. Heard, J. Phys. Chem. A, 2024, 128, 2789–2814.
102. L. Xu, K. H. Møller, J. D. Crounse, H. G. Kjaergaard and P. O. Wennberg, Environ. Sci. Technol., 2020, 54, 13467–13477.
103. R. Wu, S. Pan, Y. Li and L. Wang, J. Phys. Chem. A, 2014, 118, 4533–4547.
104. A. Shiroudi, M. A. Abdel-Rahman, A. M. El-Nahas and M. Altarawneh, New J. Chem., 2021, 45, 2237–2248.
105. T. V.-T. Mai, H. T. Nguyen and L. K. Huynh, Chemosphere, 2021, 263, 127850.
106. J. Zeng, S. Xu, K. Lin, S. Yao, B. Yang, Z. Peng, T. Hao, X. Yu, T. Zhu, F. Jiang and J. Sun, Water Res., 2024, 249, 120940.
107. N. Borduas, G. da Silva, J. G. Murphy and J. P. Abbatt, J. Phys. Chem. A, 2015, 119, 4298–4308.
108. N. Borduas, J. P. Abbatt, J. G. Murphy, S. So and G. da Silva, Environ. Sci. Technol., 2016, 50, 11723–11734.

---

Footnote

† Electronic supplementary information (ESI) available: Additional computational details and cartesian coordinates for optimized molecular geometries. See DOI: https://doi.org/10.1039/d4cp02103g

---