Computational Studies on the Thermodynamic and Kinetic Parameters of Oxidation of 2-Methoxyethanol Biofuel via H-Atom Abstraction by Methyl Radical

In this work, a theoretical investigation of thermochemistry and kinetics of the oxidation of bifunctional 2-Methoxyethanol (2ME) biofuel using methyl radical was introduced. Potential-energy surface for various channels for the oxidation of 2ME was studied at density function theory (M06-2X) and ab initio CBS-QB3 levels of theory. H-atom abstraction reactions, which are essential processes occurring in the initial stages of the combustion or oxidation of organic compounds, from different sites of 2ME were examined. A similar study was conducted for the isoelectronic n-butanol to highlight the consequences of replacing the ϒ CH2 group by an oxygen atom on the thermodynamic and kinetic parameters of the oxidation processes. Rate coefficients were calculated from the transition state theory. Our calculations show that energy barriers for n-butanol oxidation increase in the order of α ‹ O ‹ ϒ ‹ β ‹ ξ, which are consistent with previous data. However, for 2ME the energy barriers increase in the order α ‹ β ‹ ξ ‹ O. At elevated temperatures, a slightly high total abstraction rate is observed for the bifunctional 2ME (4 abstraction positions) over n-butanol (5 abstraction positions).


computational Details
All calculations were carried out using Gaussian-16W program 44 . Geometry optimizations of reactants, transition states, and products have been performed using Density function theory (DFT) Minnesota M06-2X hybrid meta functional with 54% HF exchange 45 with the 6-31 + G(d, p) basis set. M06-2X function was designed by Zhao et al. 45 to give accurate thermokinetic calculations and tested in many advanced publications [46][47][48][49][50] . The high ab initio CBS-QB3 51-53 method was also used for providing accurate energies at low computational cost. The CBS-QB3 methodology includes geometry optimization and frequency calculation at the B3LYP/6-311 G(d, p) level followed by CCSD(T)/6-31 + G (d), MP4SDQ/6-31 + G (d, p), and MP2/6-311 + G (2df, 2p) single point energy calculations with CBS extrapolation. Transition states for H-atom abstraction pathways were located using the eigenvector-following (EF) optimization technique which is implemented in the Gaussian suit of programs. Linear Synchronous Transit (LST) method is used to search for the saddle point (maximum energy) on the linear path between reactants and products. Synchronous Transit-Guided Quasi-Newton (STQN) 54,55 method is a techniques of Synchronous transit (ST) methods that use the quadratic synchronous transit (QST) approach which chooses the intermediate point in a perpendicular direction to the LST trajectory to get closer to the quadratic region of the transition state then uses a quasi-Newton or eigenvector-following algorithm to complete the optimization process. The STQN methods can be obtained by invoked keywords QST2 which requires input file contain two molecule specifications of the reactant and product, and QST3 that requires three molecule specifications of the reactant, the product, and an initial structure for the transition state. For accurate transition state location, the STQN methods 54,55 (QST2 and QST3 keywords) were used to check the efficiency of that obtained by Berny algorithm (OPT = TS). Vibrational frequency calculations were conducted at the same level of theory to characterize the nature of those points as minima (real frequencies) or transition state (only one imaginary frequency) and to correct energies for zero-point (ZPE) and thermal contributions at 298 K. Vibrational modes of different structures were visualized using the ChemCraft program 56 . For step by step verifying of transition states existence, the reactants connected with desired products minimum energy paths (MEP) were computed through intrinsic reaction coordinates(IRC) 57,58 at M06-2X/6-31 + G(d, p).
The highly accurate kinetic program Kisthelp 59 was used for rate constants (k) calculations of chemical reactions over 200-2000 K using Classical transition state theory 60 and Wigner correction 61 as follows: where k B , h, R, and T are Boltzmann, Planck, ideal gas constant, and the system's temperature in Kelvin, respectively; σ is the reaction path degeneracy; P° is the standard pressure = 1 atm; Δ ± G°(T) is the standard Gibbs free energy of activation for reaction, while Δn can takes integer values zero for unimolecular decomposition and one for bimolecular oxidation.
To account for tunneling effects, the transmission coefficient χ (T) along the reaction coordinate and, thus, the rate constant k TST/W (T) including tunneling correction is given by The transmission coefficient χ (T) was calculated using the Wigner correction 61 which is the simplest form and assumes a parabolic potential for the nuclear motion near the transition state. The Wigner transmission coefficient is given by where v is the imaginary frequency of the transition state.
Tukey honestly significant difference 62 (Tukey HSD) is a statistical method which is used to find whether the relation between two groups is significantly different or not. Tukey criterion (T) can be obtained from the formula where q α(c, n−c) is studentized range distribution based on the degree of freedoms (df) of c (number of columns), n is total sample size, MSE is the mean square error which is obtained from the analysis of the variance "ANOVA" table, and n i is the sample size in a specific column with the smallest number of observation.

Results and Discussion
Both of n-butanol and 2ME have three main dihedral angles which, in case of n-butanol and 2ME (in parenthe- Each dihedral angle can be exist in trans (T, t), gauche (G, g) or anti-gauche (G-, g-) form. n-Butanol conformations were being a subject for many past discussions at different levels of theory 41,42,63,64 . Results confirmed the existence of fourteen conformers of n-butanol with a relatively high stability of tGt conformer with narrowed energy range ≤2 kcal/mol at the G3 41 and CBS-QB3 63 methods. Our results of conformational analysis of 2ME are tabulated in Table 1, where two computational methods are used to calculate the relative conformer energies. The most stable conformer is tGg-with an intramolecular hydrogen bond (2.329 Å, see Table S1 in the SI file), which lies below the least stable conformer, gGt, by 4.38 kcal/mol. The next most stable conformer is gGg-which is less than 2 kcal/mol (Table 1) above the global minimum structure. Our results at CBS-QB3 and G3 are in good agreement with the reported data of 2ME conformers stability 65 . The higher stability of tGg-conformers enhanced its high yield at different combustion temperature in contrast to n-butanol. Optimization and energies of different 2ME conformers are presented in the Supporting Information (SI). Our study will be conducted on the most stable structures tGt, tGg-for n-butanol and 2ME, respectively at CBS-QB3. Figure 1 shows bond dissociation energies (BDEs) of n-butanol (tGt), 2ME (tGg-) at 298 K and their optimized structures at CBS-QB3. The results obtained for tGt conformer n-butanol are in good agreement with literatures 9,12 and our previous results 66 at the same level. From Fig. 1, it is obvious that C β -C ϒ (n-butanol) and C β -O β (2ME) have nearly equal bond energies of 90.3 and 90.0 kcal/mol, respectively. The O-H bond has the highest bond energy in the two compounds with the values of 105.3 and 108.2 kcal/mol for n-butanol and 2ME, respectively. The C α -H bonds have the weakest bond energies among different abstraction positions with 95.6 kcal/mol for n-butanol and 96.2 kcal/mol for 2ME. Replacing C ϒ with O β lowers bonds energies for C β -H and C ξ -H with 4.0 and 4.4 kcal/mol, respectively. For bifunctional 2ME, H-atom abstraction from α and β positions is easier and faster compared to the ξ position. This trend could be explained by the effect of the adjacent two active groups 67 .
Based on similar studies of Katsikadakos et al. 42,43 for n-butanol oxidation and by comparing all transition states for different H-abstractions from the same site (see SI Tables S3 and S4). The results indicated that all hydrogen atoms linked to the same carbon atom can be considered equivalent; the notations of carbon atoms are named relative to the alcoholic oxygen. n-Butanol has five abstraction sites α, β, ϒ, ξ, and alcoholic hydrogen. These formed radicals are CH 3   www.nature.com/scientificreports www.nature.com/scientificreports/ Optimized structures of transition states for H-atom abstraction from 2ME (tGg-) and n-butanol (tGt), by the • CH 3 radical are displayed in Fig. 2. For n-butanol, α transition state proceeds through stretching the broken C α -H bond by 17.1% and elongation of the formed C methyl -H bond by 33.3% relative to reactants and products, respectively. The bond breaking/bond forming of the other transition states TS β-nbuOH , TS ϒ-nbuOH , TS ξ-nbuOH , and TS O-nbuOH are 20/28.4%, 20.2/27.7%, 21.2/27%, and 24.6/19.9%, respectively. For 2ME, bond breaking/bond forming are 17.2 / 32.7%, 18/ 32.3%, 18.7/31.1%, and 25.8/18.9% for TS α-2ME , TS β-2ME , TS ξ-2ME , and TS O-2ME , respectively. The calculated data are consistent with the Hammond postulate 68 which stated that for exothermic reactions, the transition state structure and energy are close to reactants rather than products, while for endothermic reactions (alcoholic H-atom abstraction), the structure and energy of the transition state are close to products rather than reactants. Tables 2, 3 summarize energy barriers and reaction energies for the bimolecular oxidation of n-butanol, and 2ME with methyl radical calculated at different levels of theory. We have recalculated the energy barriers for the oxidation of n-butanol with methyl radical at CBS-QB3 to compare 2ME at the same level of theory. Comparing the obtained transition states for H-atom abstraction from α and β sites for 2ME and n-butanol by Berny algorism with that of STQN methods indicated that the obtained transition states from the two methods are similar (barrier height difference is less than 0.1 kcal/mol). So the calculations based on Berny algorism (OPT = TS) can give accurate barrier heights.
For n-butanol, the estimated barrier energies at CBS-QB3 agree well with the reported data at ROCBS-QB3 42 with very small energy differences being maximum for ϒ channel of 0.6 kcal/mol. Our results at CBS-QB3 illustrate less agreement with CCSD(T)/CBS 43 having the highest energy difference of 1.9 kcal/mol for H-atom abstraction from oxygen. For 2ME, the β and ξ products are much stable relative to those of n-butanol due to the ability to form a delocalized π bond with a lone pair of etheric oxygen. The results of M06-2X are in good agreement with that obtained using the acceptable CBS-QB3 level. From Tables 3 and S7, the computed reaction energies and enthalpies for n-butanol and 2ME at CBS-QB3 and M06-2 × /6-31 + G(d, p) suggest that all hydrogen atom abstraction channels from the carbon atoms are exothermic processes and only the abstraction reaction from the oxygen atom is endothermic. The most exothermic process is the α hydrogen atom abstraction and the least exothermic one is the H-atom abstraction from ξ position for both compounds at the two computational methods. The potential energy diagram of 2ME oxidation by • CH 3 radical at CBS-QB3 appears in Fig. 3.
The calculated barrier heights order of 2ME oxidation is consistent with the corresponding product radical stability. α H-atom abstraction represents the most preferable decomposition pathway both thermodynamically and kinetically with barrier energy of 10.1 kcal/mol and reaction energies of −9.7 and −9.3 kcal/mol for n-butanol and 2ME, respectively. All channels are exothermic, except that of alcoholic H-abstraction with reaction energies of 0.1 and 2.6 kcal/mol for n-butanol and 2ME, respectively, at CBS-QB3. Alcoholic H-abstraction pathway is more kinetically preferable for n-butanol oxidation compared to that for 2ME. This could be easily understood on the basis of intramolecular hydrogen bond formation in 2ME. Table 4 presents enthalpies of formation of radicals derived from n-butanol and 2ME oxidation and their stabilities relative to the least stable www.nature.com/scientificreports www.nature.com/scientificreports/ n-butoxy radical CH 3 CH 2 CH 2 CH 2 O • and methoxyethoxy radical, CH 3 OCH 2 CH 2 O • , respectively. According to Table 4, the calculated enthalpies of formation for n-butanol and its radicals are in an excellent agreement with theoretical values obtained by Black et al. 9 and with the available experimental results.
Estimation of ionization energies (IE) and electron affinities (EA) for chemical compounds are a crucial step for the determination of many chemical properties such as softness, hardness, electronegativity and the chemical potential 69,70 . IEs and EAs of n-butanol, 2ME, and their radicals calculated theoretically using the adiabatic and vertical approaches at CBS-QB3 level and the available experimental data are presented in Table 5. For further confirmation of the accuracy of our theoretical calculations at CBS-QB3 level, the Enthalpies of formation (AE) and adiabatic ionization energies (AIEs) for a group of oxygenated compounds which previously experimentally detected are collected in Table 6, while Table 7 tests the theoretically estimated barrier heights against the experimental activation energy (E a ) for their reaction with methyl radical.  Tables 4-7 show a good agreement between the theoretical values of AE, IE, EA, and E a and the available experimental data indicating the suitability of the employed level of theory.
Tukey test is used to examine the difference between theoretical and experimental results presented in Table 6. Our results (Table S8) reveal that the absolute difference between the averages of the experimental and theoretical data│X exp . − X theo .│is 0.53 which is less than T of 129 indicating a high constancy between the theoretical and experimental results.
Branching ratio analysis of n-butanol (Table S11) indicates domination of α abstraction with around 16% contribution of the alcoholic H-atom abstraction at temperature up to 500 K. At T > 500 K, the contribution of the ϒ channel increases with the decline of contribution from α and alcoholic channels. At high temperature (T  G(d, p)  , kcal/mol) for H-atom abstraction from n-butanol (tGt) and 2ME (tGg-) by the • CH 3 radical at different levels of theory. a Ref. 43 , b ref. 42 , and c current study.  G(d, p) Table 3. Reaction energies (E 0 , kcal/mol) for H-atom abstraction from n-butanol (tGt) and 2ME (tGg-) by the • CH 3 radical at different levels of theory. a Ref. 43 , b ref. 42 , and c current study.   Table 4. Enthalpies of formation using atomization energy approach (AE, ∆H f,298 ), and relative radicals' stabilities (∆E) derived from n-butanol and 2ME (kcal/mol) at CBS-QB3. a Ref. 71 , b ref. 72 , c ref. 73 , d ref. 74 .
≥ 1100 K), the ϒ channel becomes the main preferable pathway of n-butanol oxidation which agrees with previously reported data 43 . On the other hand, branching ratios of 2ME (Table S12) indicate domination of β H-atom abstraction at the applied temperature range with a remarked competition of ξ H-atom abstraction, especially at higher temperatures. A significant contribution of H-atom abstraction from the O atom of n-butanol compared to 2ME at the applied temperatures was observed. Figures 4, 5 depict rate constants (cm 3 /mol/s) of all abstraction sites of n-butanol and 2ME, respectively, versus temperature at CBS-QB3. Figure 6 displays a comparison of the contribution of similar H-atom abstraction sites from n-butanol and 2ME at 200-2000 K. Figure 6a shows a decrease of α H-atom abstraction contribution for the two molecules with rising of temperature, a noticeable decrease of α abstraction contribution for n-butanol (from 86% at 200 K to 21% at 2000 K) compared to a moderate decrease for 2ME (from 54.5% at 200 K to 21% at 2000 K). Figure 6b illustrates β channels for the two selected biofuels. For n-butanol, the graph indicates a gradual increase in channel contribution with rising of temperature from 1% at 200 K to 21% at 2000 K. However, 2ME has a stable branching contribution with small variation. For the ξ position shown in graph 6c, an increase of contribution of H-atom abstraction is observed for the two molecules with rising of temperature. This channel contributes almost zero at 200 K and rises to 19% and 35% at 2000 K for n-butanol and 2ME, respectively. Figure 6d reveals a noticeable high contribution from the alcoholic H-atom abstraction for n-butanol compared to that for 2ME. The branching ratio for n-butanol shows a maximum at 400 K with 17% which decreases to 9% at 2000 K, while for 2ME, a regular branching ratio increases from zero at 200 K to 5% at 2000 K.
Based on the relation between ln(k) and 1000/T (Figs S1, S2), kinetics of bimolecular oxidation reactions usually adopt non-Arrhenius behavior that fitted to a modified three-parameter Arrhenius relation of A, n, and E a .  Table 7. Theoretical barrier heights and the experimental activation energy (E a ) for H-abstraction by the • CH 3 radical from some oxygenated compounds (kcal/mol) at CBS-QB3. a Ref. 94 , b ref. 95 , c ref. 96 , d ref. 97 . Table 8 collects rate constants of individual positions and total abstraction rate for n-butanol, and 2ME molecules. The results indicate that the calculated rate constant of individual sites and total abstraction of n-butanol at CBS-QB3 are in good agreement with those obtained by Katsikadakos et al. 43 . Figure 7 shows a comparison of total abstraction rate constants for 2ME and n-butanol at CBS-QB3 which indicates a slightly higher total H-atom abstraction rate for 2ME (four abstraction sites) than that for n-butanol (five abstraction sites) as a result of presence of etheric oxygen atom.
Based on Fig. 7, at T ≥ 1300 K, the total rate for H-abstraction by methyl radical from 2ME was preceded that of n-butanol. This can be explained based on the TST Eq. (1) where the rate k TST value is inversely proportional to the Δ ǂ G°(T) value.
Similar to the free energy change equation, the Gibbs energy of the activated complex (transition state) can be obtained from: where Δ ǂ G°, Δ ǂ H°, and Δ ǂ S° are the transition state standard Gibbs free energy, the standard enthalpy, and the standard entropy, respectively. T is the absolute temperature.

Summary and conclusions
The current study describes the main thermodynamic and kinetic features of H-atom abstraction from 2ME by • CH 3 radical based on a comparison with n-butanol at the accurate ab initio CBS-QB3 methodology. The results can be summarized as follows:  www.nature.com/scientificreports www.nature.com/scientificreports/ 1. There are some agreements between n-butanol and 2ME regarding: • All investigated channels are exothermic except the abstraction of the alcoholic hydrogen atom.   Table 8. Modified three-parameter Arrhenius expression (cm 3 /mol/s) for individual sites and total abstraction of 2ME and n-butanol at CBS-QB3 over 200-2000 K.
2. The results of barrier heights and reaction energies of 2ME oxidation at CBS-QB3 are in excellent agreement with the corresponding M06-2X values with energy discrepancy less than 1 kcal/mol. 3. Energy barriers and reaction energies of 2ME oxidation increase in the order α < β < ϒ < O. 4. Replacing ϒ CH 2 group in n-butanol with etheric oxygen lowers C-H bond dissociation energies for β and ξ hydrogen atoms which enhances oxidation of the bifunctional biofuel compared to the uni-functional one. 5. Branching ratio of n-butanol indicates domination of the α channel up to 1100 K. Above 1000 K, the ϒ channel becomes the main abstraction route. Our study of 2ME oxidation confirms domination of the β channel at the applied temperatures.

Data availability
All data generated through this study are collected in this manuscript and the Supporting Information file.