Quantum Chemical Study of Potential Energy Surface in the Formation of Atmospheric Sulfuric Acid^{†}
Ⅰ. INTRODUCTION
Atmospheric reactions of sulfur compounds, particularly for the formation of sulfuric acid, have been a topic of considerable study in recent years [1, 2]. Sulfuric acid, a major component of acid rain, is also known to influence atmospheric nucleation processes for the formation of atmospheric aerosols [36]. Because virtually all sulfur in the atmosphere is either emitted as or eventually converted into sulfur dioxide (SO$_2$), many studies focus on the oxidation of SO$_2$ as an initial step for sulfuric acid formation.The general reaction process of sulfur dioxide oxidation is fairly wellestablished, beginning with the reaction of SO$_2$ with hydroxyl radical (OH) [1, 2, 7]. Sulfur compounds emitted to the atmosphere are oxidized to form SO$_2$, which reacts first with OH to form the HOSO$_2$ radical. The HOSO$_2$ radical then forms SO$_3$ and HO$_2$ through collisions with O$_2$. Finally SO$_3$ goes on to react with water via a water/acid catalyzed reaction to form sulfuric acid:
(1) SO$_2$ + OH + M $\rightarrow$ HOSO$_2$ + M
(2) HOSO$_2$ + O$_2$ $\rightarrow$ SO$_3$ + HO$_2$
(3) SO$_3$ + H$_2$O $\rightarrow$ H$_2$SO$_4$
The first step of this mechanism has been more or less confirmed by theoretical and experimental evidence; its product, the HOSO$_2$ radical, has been detected directly in the atmosphere [1]. Step (2) has been proposed as a dominant channel for the reaction of the HOSO$_2$ radical in the atmosphere [7]. While both the thermochemical estimation [8] and theoretical calculation [911] showed that the reaction was endothermic with an energy change of 8 kcal/mol or larger, a recent study revealed it to be slightly exothermic with an enthalpy of 2.3 kcal/mol [12]. Some studies propose an alternate step involving the formation of a HOSO$_4$ radical instead of SO$_3$ and HO$_2$ [13]. However, this step, while a likely atmospheric process, is of little importance in the overall reaction mechanism because the HOSO$_4$ radical may form SO$_3$ and HO$_2$ on its own through a fast hydrogen transfer step [2].
Reaction (3) also has been a subject of several theoretical studies [14, 15]. Sulfuric acid is unlikely to form directly from an elementary reaction of SO$_3$ and H$_2$O due to the high energy barrier involved (about 30 kcal/mol) [16]. Although various studies [15, 17] show that reaction of SO$_3$ in a cluster of two, three, or four waters significantly lowers the barrier, making this step appear more favorable, clusters of multiple waters are in fact highly unstable in the gas phase under normal atmospheric conditions. Therefore, it is still uncertain that such a process actually occurs in the atmosphere, and an alternate reaction pathway should be considered.
This study attempts to find an alternate pathway to better describe the sulfuric acid formation mechanism. One possible pathway for consideration is the direct reaction of a water molecule with the intermediate product SO$_3$$\cdot$HO$_2$ from reaction (2), rather than releasing the separate intermediate products SO$_3$ and HO$_2$. The reaction steps can be revised as follows.
(4) HOSO$_2$ + O$_2$ $\rightarrow$ SO$_3$$\cdot$HO$_2$
(5) SO$_3$$\cdot$HO$_2$ + H$_2$O $\rightarrow$ H$_2$SO$_4$ + HO$_2$
The SO$_3$$\cdot$HO$_2$ complex is known to be highly stable. A binding energy of 12.5 kcal/mol was found for SO$_3$$\cdot$HO$_2$ with respect to the separate SO$_3$ and HO$_2$ species [11]. Water, abundant under normal atmospheric conditions, is widely available to facilitate the hydrolysis of SO$_3$$\cdot$HO$_2$ to form sulfuric acid and HO$_2$. The final outcome of the new mechanism is identical to that of the original mechanism in the literature, except that the production and complete release of SO$_3$ is avoided in the new mechanisms. Other possible pathways could also be considered. For example, a radical complex HOSO$_2$$\cdot$H$_2$O could be formed prior to the reaction with O$_2$. The radical complex HOSO$_2$$\cdot$H$_2$O is also very stable with a binding energy of 10.6 kcal/mol [18].
In this study, the atmospheric reaction mechanism for the production of sulfuric acid from SO$_2$ is investigated with quantum chemical calculations using density functional theory and highlevel ab initio molecular orbital theory. The new reaction pathway outlined above is the main focus of this study. All reactant complexes, transition states, reaction intermediates, product complexes, and separate products are characterized for their molecular and thermodynamic properties. A new reaction mechanism will be established based on the results, and the atmospheric implications will be discussed.
Ⅱ. THEORETICAL METHODS
Equilibrium geometries and energies of the reactant complexes, transition states, and reaction intermediates involved in reactions (1) to (5) were calculated using density functional theory (DFT) and ab initio molecular orbital theory. These geometries and energies represent the stationary points on a new potential energy surface in the atmospheric formation of H$_2$SO$_4$ from OH+SO$_2$. The DFT method in this study employed Becke's threeparameter functional with nonlocal correlation provided by Lee, Yang, and Parr (B3LYP) [19, 20], along with three different basis sets, 631+G$^*$ [21], augccpVDZ [22], and augccpVTZ [22]. Two ab initio methods for different levels of electron correlation involving valence electrons only were employed: secondorder MøllerPlesset perturbation approximation (MP2) [23] and coupledcluster method with single, double, and noniterative triple excitations (CCSD(T)) [2427], both using the augccpVDZ basis set. The CCSD(T) calculations were carried out using singlepoint energy calculations on the MP2/augccpVDZ geometries. Frozencore approximation was used in both MP2 and CCSD(T) calculations. Level of consistency in the results from the diverse theoretical methods and basis sets may indicate the accuracy and reliability in the results of calculation. Transition state geometries were found using the synchronous transitguided quasiNewton (STQN) method [28, 29]. Transition states were verified by harmonic frequency calculations and by Intrinsic Reaction Coordinate (IRC) calculations [30, 31].
Relative electronic energy ($\Delta E$), enthalpy ($\Delta H$), and Gibbs free energy ($\Delta G$) were determined as the difference between the energy of the given molecular complex (or transition state) and the sum of the energies of the respective reactant monomers HOSO$_2$, H$_2$O, and O$_2$ (in the triplet ground state), that is,
$
\Delta E=E(E_{\textrm{HOSO}_2}+E_{\textrm{H}_2\textrm{O}}+E_{\textrm{O}_2})
$

(1) 
$
\Delta H=H(H_{\textrm{HOSO}_2}+H_{\textrm{H}_2\textrm{O}}+H_{\textrm{O}_2})
$

(2) 
$
\Delta G=G(G_{\textrm{HOSO}_2}+G_{\textrm{H}_2\textrm{O}}+G_{\textrm{O}_2})
$

(3) 
where $E$, $H$, and $G$ are the total electronic energy, enthalpy, and Gibbs free energy of the molecular complex (or transition state), and $E_m$, $H_m$, and $G_m$ are the corresponding energies of the reactant monomer (HOSO$_2$, H$_2$O, or O$_2$) which the molecular complex is composed of. The relative enthalpy ($\Delta H$) and Gibbs free energy ($\Delta$$G$), reported at 298.15 K and 1 atm, were calculated using harmonic vibrational frequencies. All calculations in this study were carried out on the Gaussian 09 computational chemistry program [32] running on a UNIX computer cluster.
Included in our calculations are all species involved in reactions (1), (2), and (3), particularly the reactant complex HOSO$_2$$\cdot$O$_2$, product complex SO$_3$$\cdot$HO$_2$, and the corresponding transition state TS0. Instead of the dissociation channel of SO$_3$$\cdot$HO$_2$ into SO$_3$ and the HO$_2$ radical, our focus is on the continuing reaction of SO$_3$$\cdot$HO$_2$ with water which is abundant in the atmosphere. As a new reaction pathway on the potential energy surface, a water molecule is introduced to the complex SO$_3$$\cdot$HO$_2$, resulting in the trimolecular complex SO$_3$$\cdot$HO$_2$$\cdot$H$_2$O as involved in reaction (5). It should be noted that several conformations are expected for the complex SO$_3$$\cdot$HO$_2$$\cdot$H$_2$O and the present study primarily focuses on the most stable conformation as well as the one which would directly facilitate the hydrolysis reaction to form sulfuric acid in the following step.
Ⅲ. RESULTS AND DISCUSSION
FIG. 1 shows the equilibrium geometries of all complexes and transition states at the stationary points along the new reaction pathway on the potential energy surface. It should be noted that selected bond lengths in the figure were obtained from MP2/augccpVDZ calculations. The equilibrium geometries from B3LYP calculations with the three basis sets, 631+G$^*$, augccpVDZ, and augccpVQZ respectively, are similar to those from MP2/augccpVDZ calculations. The reaction pathway is composed of three reaction steps: (ⅰ) the radical complex HOSO$_2$ reacts with O$_2$ to form the radical complex SO$_3$$\cdot$HO$_2$, (ⅱ) SO$_3$$\cdot$HO$_2$ combines H$_2$O to form the radical complex SO$_3$$\cdot$HO$_2$$\cdot$H$_2$O which is rearranged to SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$, and (ⅲ) SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$ reacts to form the product complex H$_2$SO$_4$$\cdot$HO$_2$. Each of the reaction steps involves a transition state. For reaction step (ⅰ), the reactant complex, transition state, and the product complex are denoted HOSO$_2$$\cdot$O$_2$, TS0, and SO$_3$$\cdot$HO$_2$, respectively. Similarly, the corresponding species for step (ⅱ) are SO$_3$$\cdot$HO$_2$$\cdot$H$_2$O, TS1, and SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$, and those for step (ⅲ) are SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$, TS2, and H$_2$SO$_4$$\cdot$HO$_2$.
The total electronic energies of the species of all species at the stationary points along the reaction pathways on the potential energy surface, including those of separate reactant and product molecules, are available from the authors directly. These energies are from B3LYP/631+G$^*$, B3LYP/augccpVDZ, B3LYP/ augccpVTZ, MP2/augccpVDZ, and CCSD(T)/augccpVDZ calculations, respectively. Table Ⅰ gives the relative electronic energies ($\Delta E$), relative enthalpies ($\Delta$$H$), and relative Gibbs free energies ($\Delta G$) with respect to the separate reactants, HOSO$_2$, O$_2$, and H$_2$O. Note that only electronic energies ($\Delta E$) are available from CCSD(T)/augccpVDZ calculations. FIG. 2 shows a schematic of the corresponding potential energy profile for the reactions of HOSO$_2$, O$_2$, and H$_2$O to form H$_2$SO$_4$ and HO$_2$.
Table Ⅰ Relative electronic energies ($\Delta E$), enthalpies ($\Delta H$), and Gibbs free energies ($\Delta G$) (all in kcal/mol) from B3LYP, MP2, and CCSD(T) calculations with different basis sets.
As shown in Table Ⅰ and FIG. 2, the reaction HOSO$_2$+O$_2$$\rightarrow$SO$_3$+HO$_2$ has a very small energy barrier from the reactant complex HOSO$_2$$\cdot$O$_2$: about 4.0, 0.1, and 1.4 kcal/mol at the B3LYP, MP2, and CCSD(T) levels, respectively. The corresponding activation energies are either slightly negative or near zero, which would suggest the reaction is fast and kinetically favorable. However, the separate products SO$_3$ and HO$_2$ are about 10 kcal/mol higher in energy than the reactants HOSO$_2$ and O$_2$. The result is consistent with the thermal estimation [8] and the reported theoretical calculations [911]. It indicates that the reaction HOSO$_2$+O$_2$$\rightarrow$SO$_3$+HO$_2$ is not thermodynamically favorable. The similar, positive values of $\Delta G$ at the various levels confirm the thermodynamic unfavorability of the reaction. The strong intermolecular energy between SO$_3$ and HO$_2$, as in the product complex SO$_3$$\cdot$HO$_2$, is mainly responsible for this result. It can thus be concluded that the reaction HOSO$_2$+O$_2$ may take place, but it only produces the SO$_3$$\cdot$HO$_2$ complex.
In the new reaction pathway, a water molecule is added to the SO$_3$$\cdot$HO$_2$ complex formed from the reaction between HOSO$_2$ and O$_2$. At least two conformations were found for the SO$_3$$\cdot$HO$_2$$\cdot$H$_2$O complex. In the most stable conformation, H$_2$O accepts a hydrogen bond from HO$_2$ while donating a hydrogen bond to an oxygen atom of SO$_3$. In the second conformation, SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$, the water is inserted between the O atom of HO$_2$ and the S atom of SO$_3$. In other words, HO$_2$ and H$_2$O exchange positions on SO$_3$. The SO$_3$$\cdot$H$_2$O unit is an electron donoracceptor complex involved as a prereaction complex for reaction (3). There is a transition state, TS1, for the conversion between the two complexes, SO$_3$$\cdot$HO$_2$$\cdot$H$_2$O and SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$. The fairly low value of the relative energy at TS1, 1.6 kcal/mol at the MP2 level or 2.7 kcal/mol at CCSD(T), indicates that the conversion should be kinetically favorable. The further next step of the reaction is the hydrolysis of SO$_3$ in the complex SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$. The transition state, TS2, involves the transfer of an H atom from H$_2$O to HO$_2$ as the latter returns its own H atom to SO$_3$. The relative energy at TS2, 6.8 kcal/mol at the MP2 level or 4.7 kcal/mol at CCSD(T), is again very low, and so the reaction is also expected to be fast. The product complex H$_2$SO$_4$$\cdot$HO$_2$ is highly stable with a low relative energy of 13.4 kcal/mol at the MP2 level or 17.3 kcal/mol at CCSD(T). Finally, the H$_2$SO$_4$$\cdot$HO$_2$ complex is dissociated into the separate products H$_2$SO$_4$ and HO$_2$. The overall reaction, HOSO$_2$+O$_2$+H$_2$O$\rightarrow$SO$_3$+HO$_2$, has a value of $\Delta E$=5.0 kcal/mol at the CCSD(T) level. The reaction is expected to be slightly endothermic or exothermic as shown by small values of $\Delta H$, 8.1 and 1.3 kcal/mol from B3LYP/augccpVTZ and MP2/augccpVDZ calculations, respectively.
The hydrolysis of SO$_3$, with the transition state TS2, is a critical step in the reaction pathway. The relative energy of the transition state is surprisingly low. The corresponding energy barrier, about 6 kcal/mol from both MP2 and CCSD(T) methods, is much smaller than the energy barriers for the hydrolysis of SO$_3$ in small water clusters SO$_3$$\cdot$(H$_2$O)$_n$, ca. 21 and 14 kcal/mol with one and two waters, respectively [14, 15]. The barrier height is comparable to the one found for the cluster of three water molecules [15]. Clearly, HO$_2$ plays a more efficient role in transferring the proton/hydrogen atom in the hydrolysis reaction than a single H$_2$O molecule or a water dimer (H$_2$O)$_2$. At the transition state TS2, the H$\cdots$OO$\cdots$H chain appears to be equivalent to a stable molecule of hydrogen peroxide (H$_2$O$_2$), which might be responsible for the stability of TS2 and the low energy barrier for the hydrolysis of SO$_3$ in the complex SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$. Note that the progressively large unit of the hydrated proton was responsible for the stepwise decrease in the energy barrier in the clusters SO$_3$$\cdot$(H$_2$O)$_n$ [15].
The results should be reliable as shown by generally consistent values of the relative energies ($\Delta E$, $ \Delta H$, and $\Delta G$) from different methods and basis sets in Table Ⅰ. There are some noticeable fluctuations in these values at the B3LYP level among the three basis sets, 631+G$^*$, augccpVDZ, and augccpVTZ. It appears that the proposed new reactions are more favored by the use of the largest basis set (augccpVTZ). On the other hand, MP2 and B3LYP calculations using the same basis set (augccpVDZ) give surprisingly consistent results despite of the two very different theoretical methods. Furthermore, CCSD(T) calculations using the same basis set have nearly reproduced the results from MP2 calculations. It is expected that the energetic values for the proposed reactions would become more favorable if CCSD(T) calculations were performed using a larger basis set beyond augccpVDZ. Such calculations were attempted but were not complete because of limitations in our computer resources.
Finally, it is important to point out that the spincontamination problem is very mild throughout the calculations and it does not severely affect the reliability of the results in the present study. Typically unrestricted DFT methods such as B3LYP are subject to less spin contamination than molecular orbital methods such as UHF and UMP2. Coupledcluster theory is even more effective in reducing the spincontamination problem and CCSD(T) nearly completely eliminates the two major contaminants (S+1 and S+2). This is another reason why three different methods, B3LYP, MP2, and CCSD(T), were performed in our calculations. Almost all molecular species (reactant complexes, transition states, and product complexes) in the study are openshell with the exact expectation eigenvalue of 0.75 for the spinsquared expectation operator $\langle S^2\rangle$. The unrestricted B3LYP, MP2, and CCSD(T) procedures were used for the calculations of these openshell species. The B3LYP calculations with the different basis sets give the $\langle S^2\rangle$ values between 0.7550 and 0.7563 before annihilation (0.7500 after annihilation). The MP2 calculations give the values between 0.7639 and 0.7816 before annihilation (0.75010.7507 after annihilation). As expected, no spincontamination problems were found in the CCSD(T) calculations. It is clear that the results in the present study are not affected by the spincontamination problem.
Ⅳ. CONCLUSION AND ATMOSPHERIC IMPLICATIONS
Density functional theory and highlevel ab initio molecular orbital theory have been used to explore the new potential energy surface for the atmospheric formation of sulfuric acid from OH+SO$_2$ by simultaneous inclusion of O$_2$ and H$_2$O molecules. Instead of producing the separate products of SO$_3$ and HO$_2$, as for the reaction HOSO$_2$+O$_2$ discussed in the literature, the intermediate complex SO$_3$$\cdot$HO$_2$ directly reacts with a water molecule and produces the final products H$_2$SO$_4$ and HO$_2$. The entire reaction pathway is composed of three consecutive elementary steps: (1) HOSO$_2$+O$_2$$\rightarrow$SO$_3$$\cdot$HO$_2$, (2) SO$_3$$\cdot$HO$_2$+H$_2$O$\rightarrow$SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$, (3) SO$_3$$\cdot$H$_2$O$\cdot$HO$_2$$\rightarrow$ H$_2$SO$_4$+HO$_2$. Molecular geometries at the stationary points of PES were calculated using B3LYP and MP2 methods with several basis sets, corresponding to the reactant complex, transition state, and product complex involved in each reaction step. The results consistently show that all three elementary steps have very low activation energies and are exothermic or nearly exothermic. It implies that the reactions are favorable both kinetically and thermodynamically. The results have been further confirmed by singlepoint energy calculations with highlevel CCSD(T) method using MP2 geometries.
The new reaction pathway for the formation of sulfuric acid from SO$_2$ is more direct and favorable than existing mechanisms in the literature which assume the release of SO$_3$ from the initial steps of reaction. No additional species are required in the new pathway to initiate or complete the reaction other than OH, O$_2$, and H$_2$O, which are abundant in the atmosphere and known to be the principal ingredients in the conversion of SO$_2$ into sulfuric acid. There are several major advantages of the new pathway over existing mechanisms in the literature. First, the hydrolysis of SO$_3$ takes place before SO$_3$ is released from the SO$_3$$\cdot$HO$_2$ complex. In other words, no release of SO$_3$ is necessary in the new reaction pathway. There is very strong intermolecular binding in the SO$_3$$\cdot$HO$_2$ complex, about 10 kcal/mol, which is thermodynamically unfavorable to dissociate into SO$_3$ and HO$_2$. Second, the radical HO$_2$ in the SO$_3$$\cdot$HO$_2$ complex assists the hydrolysis by exchanging its hydrogen between H$_2$O and an O atom of SO$_3$ and lowering the energy barrier of reaction; without HO$_2$, additional water molecules, at least two, must be added to serve the same purpose of lowering the energy barrier. Third, the hydrogen of HO$_2$ generated in the new mechanism originates from the water molecule was directly used for the hydrolysis of SO$_3$, instead of coming from HOSO$_2$ (originally the OH radical). Note that the production of HO$_2$ is important to the overall chemistry of SO$_2$ oxidation [1]. In the presence of NO, the HO$_2$ reacts to regenerate OH by the reaction HO$_2$+NO$\rightarrow$OH+NO$_2$. As a result, the oxidation of SO$_2$ is a chain reaction without a net loss of OH. The production of HO$_2$ is predicted in the new mechanism as well as in existing mechanisms in the literature.
Clearly, there are several significant atmospheric implications with the new mechanism. First, no SO$_3$ is released or can be observed in the conversion of SO$_2$ into sulfuric acid. As an intermediate, SO$_3$ is predicted in the existing mechanism. To our best knowledge, however, no direct evidence has been reported in support of the generation of SO$_3$ in the oxidation of SO$_2$ by the OH radical in the atmosphere. Second, the rate of conversion in the new mechanism is expected to weakly depend on the humidity of air because only one water molecule is required for each SO$_2$. In contrast, the rate of conversion in the existing mechanism would be strongly dependent on the humidity because multiple water molecules are required for the conversion process of each SO$_2$.
Ⅴ. ACKNOWLEDGMENTS
This work was partially funded by National Science Foundation of the United States (No.1012994) and by California State University, Fullerton.