Tam V.-T.
Mai
*abc and
Lam K.
Huynh
*cd
aUniversity of Science, Ho Chi Minh City, 227 Nguyen Van Cu, Ward 4, District 5, Ho Chi Minh City, Vietnam
bMolecular Science and Nano-Materials Lab, Institute for Computational Science and Technology, SBI Building, Quang Trung Software City, Tan Chanh Hiep Ward, District 12, Ho Chi Minh City, Vietnam
cVietnam National University HoChiMinh City, Quarter 6, Linh Trung Ward, Thu Duc City, Ho Chi Minh City, Vietnam
dInternational University, Block 6, Linh Trung Ward, Thu Duc City, Ho Chi Minh City, Vietnam. E-mail: hklam@hcmiu.edu.vn
First published on 21st July 2021
In this short communication, we resolve the large discrepancy (∼750 times) in the rate constants between computational and experimental studies for the atmospheric oxidation reaction of oxazole initiated by OH radicals. This achievement is due to the use of the rigorous stochastic Rice–Ramsperger–Kassel–Marcus based master equation (RRKM-ME) rate model, which allows mimicking the dynamics of the title reaction without introducing any assumption (e.g., fast pre-equilibrium between the reactants and van der Waals pre-reactive complex). The model shows a large tunneling contribution at low temperatures and the minor role of the hindered internal rotation treatment. Detailed mechanistic insights (e.g., the NTC behaviors) are also revealed to advance related applications of the title system.
Fig. 1 schematically presents the detailed energy profile for the oxazole + OH → products reaction, consisting of three H-abstraction and three OH-addition channels, explored at the M06-2X/aug-cc-pVTZ level. Due to the somewhat similar systems of imidazole4/pyrrole5 + OH, the interaction between an oxazole molecule and an OH radical can form two pre-active van der Waals complexes or reactant-complexes (namely, RC1 and RC2) whose structures are presented in Fig. 2, located at −5.3 and −3.5 kcal mol−1, respectively, relative to the entrance channel. Energetically, the former complex, RC1, is well consistent with those found in imidazole4/pyrrole5 + OH (e.g., −6.8 and −4.9 kcal mol−1, respectively). Using the same electronic structure methods, Shiroudi and co-workers1 reported two complexes RC1′ and RC2′, lying at higher energies (i.e., −3.3 and −3.5 kcal mol−1, respectively). It is worth noting that different RC1/RC1′ structures and similar RC2/RC2′ structures, via different transition states, were found in Shiroudi et al.'s and our studies (cf.Fig. 2). In particular, in our study, RC1 is found to be about 2.0 kcal mol−1 lower (i.e., −5.3 vs. −3.3 kcal mol−1) due to the intramolecular hydrogen bond formed between the H atom of the OH radical and the N atom of oxazole. Also, RC1 and RC2 are found to go through transition states (TS1, TS2 & TS4) and (TS3, TS5 & TS6) while RC1′ and RC2′ proceed via (TS1 and TS4) and (TS2, TS3, TS5 & TS6), respectively.
![]() | ||
Fig. 1 Reaction energy profile (0 K), calculated at the M06-2X/aug-cc-pVTZ level, for the reaction of oxazole + OH → products. ZPE correction is included, and units are in kcal mol−1. Numbers in parentheses, “(![]() |
![]() | ||
Fig. 2 Optimized structures of reactant-complexes RC1, RC1′, RC2 and RC2′ obtained at the M06-2X/aug-cc-pVTZ level of theory. Bond lengths and dihedral angles are in Å and degree (°), respectively. Values in “(![]() |
From RC1, two H-abstraction (viaTS1 and TS2) and one OH-addition (viaTS4) channels can occur to form (P1 + H2O), (P2 + H2O) and I1, respectively; while the channels via transition states TS3 (direct loss-H pathway) and TS5/TS6 (OH-addition) can proceed from RC2. Note that the unimolecular decomposition of the formed adducts (i.e., I1, I2, and I3) via the β-scission of the C–H bonds is further studied to provide more insights into the reaction mechanism of the title reaction in this work. The optimized structures (reactants, complexes, TSs, intermediates, and products) are shown in Fig. S1 (ESI†). In general, our calculated numbers consistently match with those obtained from Shiroudi and co-workers1 with the mean absolute deviation (MAD) values of 0.25 and 0.35 kcal mol−1, calculated using the same methods, M06-2X/aug-cc-pVTZ and ROCBS-QB3//M06-2X/aug-cc-pVTZ, respectively (cf. Table S3, ESI†).
With the calculated equilibrium concentration ratio Keq(298 K) = [RC1]eq/[RC2]eq = 55.7 and the low intrinsic barrier (comparable with the RC2 energy) between the two complexes, RC1 is expected to dominate from both thermodynamic and kinetic points of view. In fact, two kinetic models (namely, the so-called “RC1&RC2 model” using explicitly two complexes RC1 & RC2 as described by the PES in Fig. 1 and the so-called “RC1 model” using the lowest complex RC1 as described by the PES in Fig. S2, ESI†) predict very close rate constant values with each other (i.e., the largest difference of 40.3% at high temperatures and the MAD value of ∼18.3%, cf.Fig. 4). Therefore, model RC1, in addition to model RC1&RC2, was used in some analyses to simplify the calculations (e.g., species profile simulations).
Based on the PES characterized at the M06-2X/aug-cc-pVTZ level, the predicted kinetic behaviors of the oxazole + OH → products reaction, including the time-resolved species profiles and rate coefficients, are reported thoroughly. As seen in Fig. 3 for the condition of 300 K and 760 torr, the species profiles significantly change with time; thus, at least, those species appearing in the figure (i.e., OH, RC1, I1, I3, P4 + H, P2 + H2O, P1 + H2O, P5 + H and P6 + H) should be included to characterize the dynamic/kinetic behaviors correctly. The computed total rate coefficients with the HIR and tunneling corrections, ktot, are in good agreement with those measured by Witte et al.2,3 (e.g., 8.57 × 10−12 (model RC1) vs. 8.29 × 10−12 cm3 per molecule per s at T = 324 K, cf.Fig. 4) with the MAD of less than 10%, at T = 299–468 K & P = 100 torr. In contrast, Shiroudi's numbers are much lower than the experimental values2,3 with the average factor of ∼750 (e.g., 1.04 × 10−14vs. 8.29 × 10−12 cm3 per molecule per s at T = 324 and P = 100 torr).
![]() | ||
Fig. 4 Comparison between our predicted and literature overall rate constants, ktot, of the oxazole + OH → products reaction at P = 100 torr as a function of temperature. Laboratory/computational data are collected from the works of Witte et al. (“Expt. (Witte 1986)”)2,3 and Shiroudi et al. (“Calc. (Shiroudi 2021)”),1 respectively. |
It is worth mentioning that the HIR correction plays a minor role in obtaining the reliable rate constants (e.g., it speeds up the calculated total rate by a factor of ∼1.5 in the temperature range of 299–468 K and P = 100 torr, cf. Table S4, ESI†). On the other hand, the tunneling treatment is found to play a more important role in kinetics and to be more relevant to the H-abstraction channels than OH-addition ones (e.g., the tunneling effect increases the individual rate constants by a factor of 30.0, 18.7, 16.9, 1.3, 1.3, 1.1, 3.6, 4.4 and 1.9 for channels viaTSx (x = 1–9), respectively, at T = 299 K, cf. Table S5, ESI†). Interestingly, a negative temperature dependence (NTD) manner, consistent with the experimental observation,2,3 of ktot(T) is also predicted for the title reaction in the considered conditions. Such behavior is similar to that of imidazole4/pyrrole5 + OH and mainly due to the faster re-dissociation of the reactant-complexes to the reactants at higher temperatures (cf. Fig. S3, ESI† for Keq(T) = kforward(T)/kreverse(T) plot as a function of temperature for the oxazole + OH → RC1 channel).
A comparison between our predicted ktot(T) and those calculated by Shiroudi and co-workers was carried out at P = 750 torr (∼1 bar), shown in Fig. 5. The calculated ktot(T) values at different temperatures and pressures (T = 299–468 K & P = 0.76–76000 torr) are tabulated in Table S2 (ESI†). Our values are found to be higher than those derived from the TST and RRKM models by Shiroudi and co-workers with an average factor of 3.1 and 125.3, respectively (cf.Fig. 5). Consequently, our calculated species branching ratios are significantly different from those suggested by Shiroudi et al. (cf.Fig. 6). Our model reveals that the production of adduct I3 (viaTS6) dominates the stabilization of I1 (viaTS4) and the bimolecular product formation from three H-abstraction channels P1/P2/P3 + H2O is insignificant in the considered conditions (T = 299–468 K & P = 100 torr). It is also observed that the yield of the main adduct I3 gradually decreases with the temperature increase (e.g., 94.9 and 92.4% at 299 and 468 K, respectively, at P = 100 torr).
With the PES explored at the same M06-2X/aug-cc-pVTZ level, the difference between our numbers and the values suggested by Shiroudi et al. (cf.Fig. 4) is mainly due to the use of two different kinetic models. Specifically, our models (i.e., models RC1&RC2 and RC1) explicitly include all species involved to mimic the dynamics of the title reaction without introducing the hypothesis of a fast pre-equilibrium between the reactants and the van der Waals pre-reactive complexes (e.g., ) as introduced in Shiroudi et al.'s model (as in the reaction scheme earlier proposed by Singleton et al.6). Such a fast equilibrium hypothesis is found to be invalid for this title system. For example, at T = 300 K, the [RC1]t/[OH]t × [oxazole]t ratio does not reach its equilibrium value until 10−8 seconds (cf.Fig. 7). Therefore, the model proposed by Shiroudi et al. is likely to be compromised, at least under the considered conditions. Note that the exclusion of RC1 certainly affects the detailed reaction mechanism; for example, the negative temperature dependence of the rate constants is not observed if the reaction complexes are excluded from our kinetic models.
The overall rate constants, ktot(T), based on two electronic structure methods (i.e., M06-2X/aug-cc-pVTZ and ROCBS-QB3//M06-2X/aug-cc-pVTZ) are shown in Fig. 8. Surprisingly, the calculated ktot(T) using the PES constructed at the M06-2X/aug-cc-pVTZ level is better than at the ROCBS-QB3//M06-2X/aug-cc-pVTZ method when compared to the experimental data with MADs of 9.6 and 83.7%, respectively. Since transition state TS6 (the main channel, cf.Fig. 1 and 3) is found at −1.6 and −2.2 kcal mol−1 (cf. Table S3, ESI†) at the M06-2X/aug-cc-pVTZ and ROCBS-QB3//M06-2X/aug-cc-pVTZ level, respectively, the calculated ktot(T) with the former method is slower than with the later one by an average factor of ∼2.0 in the range of T = 299–468 K & P = 100 torr. It is worth noting that the uncertainty of a kinetic model comes from different sources such as the PES and RRKM-ME model (including any hypothesis introduced); therefore, using a low level of theory might provide better kinetic data if the errors compensate the missing physics in such a model. Within this context, the M06-2X/aug-cc-pVTZ level is found to be suitable for investigating the kinetics of similar reactions (e.g., imidazole4/pyrrole5/aniline7 + OH).
In conclusion, the use of the rigorous rate model (i.e., the stochastic RRKM-based ME statistical model with the corrections of the HIR and Eckart tunneling treatments) on the M06-2X/aug-cc-pVTZ potential energy surface helps to resolve the large kinetic discrepancy between the theoretical1 and experimental2,3 studies for the title reaction without introducing any assumption as in the previous calculations. Such an achievement confidently reveals detailed mechanistic insights via the model to advance atmospheric applications related to the title reaction.
The hindered internal rotation (HIR) corrections were included in thermodynamic/kinetic analyses. The M06-2X/cc-pVDZ level was used to obtain the hindrance potentials of the rotation around the C–O and O–H “single” bonds (cf. Fig. S4, ESI†) via relaxed surface scans. The HIR parameters were determined with the use of the Multi-Species Multi-Channels (MSMC) graphical user interface.14,15 The detailed procedure of the HIR treatment was presented in our previous studies.16 The tunneling correction using the Eckart tunneling model17 was included in the rate constant calculations for all elementary reaction channels. Note that the spin–orbit splitting correction (∼139.7 cm−1)18 was also considered in calculating the partition function of OH radicals.
The microcanonical rate constants k(E) of the barrierless channels (e.g., the association of oxazole and OH radicals to form the reactant-complexes, RCs, cf.Fig. 1) were treated using the Inverse Laplace Transform (ILT) technique19 with the high-pressure-limit rate constant, k∞(T), of 4.0 × 10−10 cm3 per molecule per s, with the assumed temperature-independence based on the long-range transition state (LR-TST).20 It has been shown that this approach was successfully used to capture the rate coefficients for similar systems.4,5,7 The temperature-dependent expression of 〈ΔEdown〉 = 260.0 × (T/298)0.8 cm−1 was used for the bath gas of Ar4,5 regarding the energy-transfer manner. The Lennard-Jones parameters of ε/kB = 113.5 K and σ = 3.465 Å were used for Ar21 whereas ε/kB = 467.4 K and σ = 4.6 Å were adopted for the [oxazole–OH] complex and its adducts.1 The stochastic22,23 RRKM-based ME24 statistical rate model (with the trial number of 108), including the corrections of HIR and tunneling treatments, was used for the kinetic analyses. The MSMC code25 was used for all thermodynamic and kinetic calculations.
Footnote |
† Electronic supplementary information (ESI) available. See DOI: 10.1039/d1nj01020d |
This journal is © The Royal Society of Chemistry and the Centre National de la Recherche Scientifique 2021 |