Journal of Chemistry

Volume 2013 (2013), Article ID 659682, 11 pages

http://dx.doi.org/10.1155/2013/659682

Research Article

Prediction of Mechanism and Thermochemical Properties of O₃ + H₂S Atmospheric Reaction

Morteza Vahedpour,¹ Reza Baghary,² and Freshte Khalili²

¹Chemistry Department, University of Zanjan, P.O. Box 45195-313, Zanjan, Iran

²Chemistry Department, Payame Noor University, Abhar Branch, Abhar, Iran

Received 16 June 2012; Revised 19 September 2012; Accepted 9 October 2012

Academic Editor: Ai-Min Ren

Copyright © 2013 Morteza Vahedpour et al. This is an open access article distributed under the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

Abstract.

Ozone and hydrogen sulfide reaction mechanism including a complex was studied at the B3LYP/6-311++G(3df,3pd) and CCSD/6-311++G(3df,3pd)//B3LYP/6-311++G(3df,3pd) levels of computation. The interaction between sulfur atom of hydrogen sulfide and terminal oxygen atom of ozone produces a stable H₂S-O₃ complex with no barrier. With the decomposition of this complex, four possible product channels have been found. Intrinsic reaction coordinate, topological analyses of atom in molecule, and vibrational frequency calculation have been used to confirm the suggested mechanism. Thermodynamic data at T = 298.15 K and the atmospheric pressure have been calculated. The results show that the production of H₂O + SO₂ is the main reaction channel with ΔG° = −645.84 kJ/mol. Rate constants of H₂S + O₃ reaction show two product channels, SO₂ + H₂O and HSO + HOO, which compete with each other based on the temperature.

1. Introduction

Ozone and hydrogen sulfide are two important atmospheric pollutants. Ozone as a reactive species plays an important role in the upper atmosphere. In view of the relatively large ozone concentration in the middle atmosphere and the importance of ozone chemistry in the depletion of ozone layer, by a series of reactions with atmospheric species, its gas phase reactions are of interest. Recently, the reaction of ozone with atmospheric molecules and radicals has been studied extensively in the theory and experiment [1–6].

Sulfur compounds are released into the atmosphere in various forms by means of biogenic (biomass rotting and volcanic eruptions) and anthropogenic (chemical industry) sources [7]. Hydrogen sulfide is one of the sulfur compounds that play an important role in the chemistry of the atmosphere [8]. Also, it is a large source of sulfur in the atmosphere. Hydrogen sulfide is highly toxic and is used in a number of chemical industries [9]. Therefore, its atmospheric reaction to turn other compounds with lower toxicity is a topical subject in a number of studies.

The reaction of hydrogen sulfide with different molecules and radicals in the stratosphere was the subject of many experimental and theoretical studies [10–13]. In experimental studies, for many years, most of hydrogen sulfide is removed by the well-known Claus process [2, 12]. Also, H₂S reacts with ozone to produce SO₂ and H₂O, which potentially is one of the sources for acid rains [2, 3]. In theoretical work, kinetics and thermodynamics of hydrogen sulfur atmospheric reaction with some radical species such as CH₃, OH, and NO₂ are studied [1].

The reaction of ozone and H₂S can aid to remove the toxic pollutant which is important in atmospheric chemistry. Moreover a slow gas-phase reaction of ozone with hydrogen sulfide may be significant if the reaction proceeds to produce (HSOH) and (HSOOOH) that may be an important source of OH, HS, and HOO (hydroperoxyl) radicals which are well known as oxidizing agents and take part in the processes of the destruction of ozone. Therefore, in this work, we have computationally investigated the details of the mechanism and thermodynamics of the hydrogen sulfide and ozone reaction. Also, the hydrogen thioperoxide molecule is indicated as the other new product of the reactions.

Computational Methods. The geometries of the reactants, complexes (denoted as Cs), products, and transition states (denoted as TSs) involved in the reaction are optimized using density functional theory at the spin unrestricted B3LYP with the 6-311++G(3df,3pd) bases set. To obtain more reliable relative energy higher levels of electronic correlation method, CCSD(T)/6-311++G(3df,3pd)//B3LYP/6-311++G(3df,3pd), is employed in the single point energy calculation to improve the accuracy of energetic information on minimum energy path. The molecular energy at the CCSD(T) level and the zero-point energy (ZPE) which is obtained from the B3LYP level are reported. Harmonic vibrational frequencies were obtained at the B3LYP level for verification of the optimized geometries. The nature of the stationary points and transition states are determined according to the number of negative eigenvalue of the Hessian matrix. Additionally, intrinsic reaction coordinate (IRC) was performed to examine connections between all the species involved in the reaction. The counterpoise procedure is used to correct the interaction energy for basis set superposition error, and CCSD(T) method is used to calculate the diagnostic values for the structures.

At the level of this method, it is possible to generate a wave function in a suitable form to execute a topological analysis of atoms in molecules [14] using the AIM2000 series programs [15]. Moreover, thermodynamic data have been calculated using the statistical mechanics. The rate constants of the reaction channels for thermodynamic and kinetic products have been calculated using TST and RRKM theories due to the B3LYP/6-311++G(3df,3pd) method. All calculations reported in this work were carried out with the Gaussian 03 program package [16].

2. Result and Discussion

The optimized geometries of reactants, complexes (as denoted C and a number), transition states (as denoted TS and a number), and products with the corresponding available experimental data for some species are represented in Figure 1. The results of geometry optimization for all species that were obtained at the B3LYP/6-311++G(3df,3pd) level of theory in this work are in agreement with experimental values.

Figure 1: Optimized geometries of reactants, products, complexes, and transition states of studied reaction in the B3LYP/6-311++G (3df,3pd) level (bond lengths is in Å and bond angels and dihedral angles are in degree).

One reactant complex, H₂S-O₃ (denoted as C1), has been formed between the reactants on the singlet potential energy surface. In the complex, the terminal oxygen atoms of O₃ are associated with sulfur atom of H₂S and a four-member ring is formed. Molecular geometry indicates that the C1 complex has an asymmetric structure.

In the C1 complex, the bonds length of newly formed S–O1 and S–O2 are 3.371 and 3.235 Å at B3LYP level, respectively, and the bond of O–O is about 0.004 Å longer than that of the corresponding bond in the parent ozone molecule. To identify newly bonds formation and existence of ring structure, we have used the topological analysis of electronic charge density. The values of bond critical point parameters and electronic density gradient for newly formed bonds and ring critical point for ring structure are schematized in Figure 2. The Laplacians of electronic charge density value of S–O1 and S–O2 bonds are 0.0262 and 0.0265 au, respectively. Both low and positive values of Laplacian from AIM analysis for S–O bonds of C1 complex show that there is a weak van der Waals interaction between sulfur atom of H₂S and terminal oxygen atoms of ozone. The electronic charge density corresponds with the newly bonds are 0.0095 and 0.0685 au, respectively. The calculated electronic charge densities confirm the weak interaction of ozone and hydrogen sulfide in C1 complex.

Figure 2: The values of the AIM theory topological parameters for some selected bonds of species in reaction.

The binding energy of C1 is 5.78 kJ/mol below the original reactants (H₂S + O₃) at B3LYP level. The base set superposition error (BSSE) value for the complex C1 energy is 0.01 kJ/mol. No transition state has been found for its formation. So, the formation of C1 is a barrierless process. Subsequently, via variety of C1, three kinds of products are obtained. In these paths, we have found five transition states, TS1, TS2, TS3, TS4, and TS5, on the potential energy surface of H₂S + O₃ reaction and another one, TS6, which is described as the transformation of SOO to SO₂ and four other stationary points (C2, C3, HSOOOH, and SOO). As a matter of fact, the corresponding molecular graph (Figure 2) shows the existence of a BCP between bonded atoms. Also, ring structures are formed when ring critical points appear in all species. The total energies and relative energies of all species at the B3LYP/6-311++G(3df,3pd), under atmospheric pressure and 298.15 K are listed in Table 1 and CCSD(T) relative energies are shown in Figure 3. Higher levels of electronic correlation method CCSD(T)//B3LYP are employed in the single point energy calculation to improve the accuracy of energetic information on minimum energy path. Harmonic vibrational frequencies were obtained at the B3LYP level for verification of the optimized geometries. Also the results of the single point energy calculations at the CCSD/6-311++G(3df,3pd)//B3LYP/6-311++G(3df,3pd) level of computation are summarized in Table 1. All energies in Table 1 are corrected by BSSE.

Table 1: Total energies (Hartree), relative energies (kJ/mol) (in parentheses) including ZPE corrections at the B3LYP/6-311++G(3df,3pd) and CCSD(T)//B3LYP methods, ZPE (Hartree), and T1 diagnostic for various species.

species
B3LYP/6-311++G(3df,3pd)
CCSD(T)/6-311++G(3df,3pd)
ZPE
T1

H₂S + O₃
−624.9053 (0)
−624.0099 (0)
0.0224

TS1
−624.8647 (106.53)
−623.9442 (172.63)
0.0237
0.099

TS2
−624.8809 (63.65)
−623.9678 (110.69)
0.0250
0.047

TS3
−624.8880 (45.33)
−623.9768 (86.92)
0.0235
0.043

TS4
−624.9069 (−4.41)
−624.0015 (22.14)
0.0262
0.041

TS5
−624.8634 (110.05)
−623.9495 (158.65)
0.0246
0.047

C1
−624.9063 (−2.74)
−624.0121 (−5.78)
0.0239
0.023

C2
−624.8890 (42.76)
−623.9869 (60.40)
0.0274
0.042

C3
−624.9690 (−167.30)
−624.0811 (−186.98)
0.0293
0.029

H₂SO–O₂
−624.9468 (−109.12)
−624.0741 (−168.55)
0.0270
0.014

HSOOOH
−624.9551 (−130.76)
−624.0678 (−151.96)
0.0290
0.021

TS6 + H₂O
−624.8529 (137.47)
−623.9998 (26.46)
0.0248

TS6
−548.4097
−547.6707
0.0035
0.046

SOO
−548.5234
−547.7462
0.0056
0.027

SO₂
−548.7090
−547.9449
0.0070
0.022

H₂SO + ³O₂
−625.0074 (−268.10)
−624.1236 (−298.39)
0.0262

H₂SO + ¹O₂
−624.9461 (−107.25)
−624.0666 (−148.89)
0.0262

HSO + HOO
−624.9478 (−111.74)
−624.0812 (−187.24)
0.0245

H₂O + SO₂
−625.1521 (−648.19)
−624.2740 (−693.29)
0.0283

H₂O + SOO
−624.9666 (−160.91)
−624.0753 (−171.62)
0.0269

(a)

(b)

Figure 3: Profile of calculated potential energy surface of H₂S + O₃ reaction at CCSD(T)//B3LYP level of computation.

The vibrational frequencies of all stationary points along with the available experimental values are listed in Table 2. As can be seen, the calculated frequencies at B3LYP method are in acceptable agreement with available experimental data. All the reactants, products, and complexes have only real frequencies, while TS is identified with one and only one negative eigenvalue of the Hessian matrix.

Table 2: Vibrational frequencies of all species at B3LYP/6-311++G(3df,3pd) level of theory.

Species
Frequency (cm⁻¹)

H₂S
1210 (1183), 2695 (2615), 2708 (2626)

TS1
326i, 102, 198, 265, 355, 602, 694, 911, 1234, 1432, 1927, 2680

TS2
750i, 196, 344, 358, 424, 481, 639, 685, 763, 1281, 2016, 3774

TS3
1111i, 190, 331, 398, 476, 644, 775, 1031, 1113, 1162, 1484, 2701

TS4
1432i, 358, 411, 469, 501, 566, 695, 880, 1094, 1125, 1589, 3785

TS5
500i, 70, 167, 213, 272, 476, 535, 752, 1208, 1328, 2684, 3107

TS6
918i, 378, 1162

C1
45, 52, 68, 127, 169, 218, 752, 1201, 1211, 1231, 2696, 2709

C2
227, 248, 287, 393, 538, 575, 904, 1040, 1221, 2665, 3766

C3
58, 124, 142, 196, 236, 419, 502, 826, 1093, 1633, 3754, 3876

H₂SO–O₂
32, 1014, 49, 65, 2616, 99, 1619, 129, 1190, 755, 3791, 472

HSOOOH
71, 264, 392, 435, 519, 585, 770, 968, 1028, 1427, 2593, 3684

O₃
753 (716), 1213 (1089), 1264 (1135)

3O₂
1646 (1556)

H₂SO
487, 759, 1015, 1195, 2610, 3797

HSO
1013, 1091, 2426

HOO
1171, 1435, 3604

Our calculation led to the identification of the reaction paths on singlet potential energy surfaces which can be summarized in the four following paths:

P1(1): O₃ + H₂S→C1→TS1→C2→TS2→C3→SOO + H₂O→TS6 + H₂O→SO₂ + H₂OP1(2): O₃ + H₂S→C1→TS1→C2→TS2→C3→TS4→HSOOOH→HSO + HOOP2: O₃ + H₂S→C1→TS3→HSOOOH→HSO + HOOP3: O₃ + H₂S→C1→TS5→¹O₂ + HSOH→³O₂ + HSOH

where HSOOOH is a complex between the HOO and HSO molecules. The details of the reaction mechanism on the singlet potential energy surfaces will be discussed below.

2.1. Reaction Pathway Properties

We have only one collisional complex (C1) which was located on the potential energy surface for the H₂S + O₃ reaction which is formed when the sulfur atoms of hydrogen sulfide were attached to the terminal oxygen atoms of the ozone. For the H₂S + O₃ reaction, there are three possible products as SO₂ + H₂O, HSO + HOO, and HSOH + O₂.

2.1.1. SO₂ + H₂O Formation Channel

In path P1(1), C1 undergoes H–O1 and S–O2 bonds formation and S–O1, S–O3, and H–S bonds rupture process to form C2 via TS1 with the energy barriers of 178.41 kJ/mol. Then C2 transforms to C3 via TS2 with the energy barriers of 50.29 kJ/mol by O1–S and H–O3 bond formation and hydrogen transfers from sulfur atom to O1 atom.

In the four-member ring structure transition region (TS2) of C2→C3 process, the electronic charge density and its Laplacian are 0.0356 and 0.1718 au, respectively.

The energy height shows C1→C3 conversion is energetically feasible with two low energy transition states. The reaction path continues with the formation of water and SOO molecules from C3 without passing through any transition state. Finally, SOO undergoes O–S bond formation and O–O bond rupture process to the formation of SO₂ molecule via TS6 with the energy barriers of 26.46 kJ/mol relative to original reactants and 198.08 kJ/mol relative to SOO molecule. The height energy of the final product (H₂O + SO₂) is 693.29 kJ/mol lower than the original reactants (H₂S + O₃). Because the formation of SO₂ can be produced from the rearrangement of SOO and the number of their atoms is different from the others, we show this step of path P1(1) individually in Figure 3. The corresponding molecular graph of all species obtained from AIM calculations in path P1(1) in Figure 2 confirms that the formation and the rupture of the bonds are indicated by the appearance of BCP and the collapse of BCP, respectively. Also, the imaginary frequencies of TS1, TS2 and TS6, which denote the cleavage of some bonds and the formation of other bonds in the reaction coordinate, are 327i, 751i and 970i cm⁻¹, respectively. In all the steps of path P1(1), the formation of products from reactants is confirmed using related IRC curves.

2.1.2. HSO + HOO Formation Channels

For the HSO + HOO product, there are two pathways on the singlet potential energy surfaces denoted as P1(2) and P2 which are confirmed using vibrational frequencies of transition states, AIM analysis and IRC curves for transition states.

In path P1(2), the formation of C3, a five-member ring structure with a low value of electron charge density, 0.0081 au, and Laplacian of charge density, 0.0366 au, are similar to path P1(1). C3 undergoes H–O1 bond rupture process and disappear of ring structure to form HSOOOH via TS4 with the energy barriers of 209.12 kJ/mol. Finally HSO + HOO molecules are formed by the breaking of O2–O3 bond without undergoing any transition state. The formation of products from C3 is confirmed using AIM results and IRC curve.

In path P2, C1 transforms to the HSOOOH complex and its corresponding products (HSO and OOH) in the one step reaction via TS3 with the energy barriers of 92.70 kJ/mol. Its transformation, accompanied with H2–O3 bond formation and S–O3 bond rupture process to form a five-membered ring transition state, TS3, molecule. Then, ring structure of TS3 disappears and form HSOOOH specie. Ring existing, H2–O3 bond, formation and S–O3 bond cleavage are confirmed using AIM calculation which are schematized in Figure 2.

2.1.3. HSOH + O₂ Formation Channel

Another simple step of the reaction is path P3 which is an elementary reaction. In path P3, C1 undergoes to the products (HSOH + O₂) using oxygen abstraction of ozone by hydrogen sulfide and at the same time hydrogen shift from sulfur atom to newly oxygen atom via TS5 with the energy barriers of 164.43 kJ/mol and vibrational imaginary frequency of 500i . The oxygen abstraction and hydrogen shift lead to HSOH formation which are confirmed using intrinsic reaction coordinate calculation (IRC).

In summary, the H₂S + O₃ reaction is energetically feasible to the formation of all the three products throughout corresponding paths which all steps of reaction confirm using IRC and AIM calculations.

2.2. Topological Analysis of Electronic Density of the HSOH, HSO, and SOO Molecules

AIM topological analysis of electronic charge density has been elaborated using the AIM2000 program package. AIM topological analysis has been obtained using the electronic charge density integrated over atomic basins (up to 0.001 e/bohr³ level) as well as in terms of electron density, , density Laplacian, , and bond ellipticity, , at bond critical points (BCP) where and are the eigenvalues of the Hessian matrix of the BCP electron density.

According to the topological analysis of electronic charge density, in the theory of the atoms in molecules (AIM) [14], electron density () and Laplacian of the electron density, , are used to describe the strength and the characteristic of the bond, respectively. The Laplacian () is the sum of , , and , where is the th eigenvalue of Hessian matrix of the electronic density. If a critical point has two negative and one positive eigenvalues, it is called (3,−1) or the bond critical point (BCP). If a critical point has two positive and one negative eigenvalues it is called (3,+1) or the ring critical point (RCP), which indicates that a ring structure exists. According to the theory of AIM, the Laplacian of the electron density () describes the characteristic of the bond. In general, when < 0 the bond is covalent, but when > 0 the bond belongs to the electrostatic interaction.

The topological analysis of the electronic density were performed for all bonds in the HSOH, HSO, and SOO species as relatively stable products of our studied reaction. The topological characteristics in the BCP’s of S–O, O–H, and H–S bonds are tabulated in Table 3. Also, the values in Table 3 indicate that all the bonds in HSOH have covalent character ( < 0). Their values show the bond strength manner O–H > H–S > S–O is corresponding to the electronic charge densities 0.3781, 0.2225, and 0.1929 au, respectively. In HSO radical, the H–S bond has a covalent character while the S–O bond electrostatic interaction nature with electronic charge densities value 0.2174 and 0.2686 au, respectively. The electronic charge density of HS bond of HSOH in comparison with corresponding bond in HSO bond shows that the HSOH bond is stronger than the corresponding bond of HSO, while S–O bond in HSO has electron charge density higher than the corresponding bond in HSOH molecule which is similar to sulfur and oxygen double bond. In SOO species, the O–O bond has a covalent character while S–O bond has an electrostatic character with the electronic charge densities values of 0.4369 and 0.2048 au, respectively. The S–O bond strange in SOO molecule is weaker than the corresponding bond in HSOH molecule and stronger than the corresponding bond in HSO species. Its nature is the same with HSO molecule.

Table 3: Topological analysis of HSOH, HSO, and SOO molecules in atomic units.

HSOH
HSO
SOO

H–S bond
S–O bond
H–O bond
H–S bond
S–O bond
O–O bond
S–O bond

−0.47454
−0.28127
−1.9606
−0.46772
−0.45831
−1.1149
−0.24991

−0.39663
−0.22767
−1.9184
−0.42171
−0.4015
−1.7519
−0.24141

0.18834
0.390444
0.91642
0.225109
1.36699
1.9251
0.77949

0.22245
0.192901
0.37809
0.217379
0.26855
0.43688
0.20482

−0.68282
−0.1185
−2.9627
−0.66432
0.50718
−0.26501
0.28817

0.19643
0.235402
0.021986
0.10908
0.14150
0.036940
0.035514

The ellipticity is a measure of the ratio of the rate of density decrease in the two directions perpendicular to the bond path at the bond critical point, the general shape of the bond, and the degree of π-character. A value of zero indicates a symmetrical distribution of density about the bond path, such as that found in standard single and triple bonds, while large values indicate a preference for density buildup in a particular orientation. The values for the ellipticity of various bonds for HSOH molecule are listed in Table 3. The high ellipticity () for S–O bond in both HSOH and HSO molecules shows its high π bonding character, whereas low value belongs to the O–H and H–S bonds which are corresponding with the low π bond contributions. In SOO molecule, both bonds have almost the same ellipticity (0.0369 and 0.0352 au). This manner shows both bonds have almost the same symmetry in the distribution of bonding character and its distribution is almost symmetrical between the two bonds. The comparison of the two bonds shows that the O–O bond has a higher π-character compared with the S–O bond which is confirmed by the experimental predictions.

2.3. Molecular Orbital Analysis

Electron affinity, ionization potential, chemical potential, chemical hardness, and some other parameters for chemical reaction are predicted by the HOMO and the LUMO orbitals analysis. The HOMO (highest occupied molecular orbital) energy characterizes the ability of electron giving; the LUMO (lowest unoccupied molecular orbital) energy characterizes the ability of electron accepting [17]. We can determine the way the molecule interacts with other species; therefore, they are called the frontier orbitals. The energy gap between HOMO and LUMO characterizes that the molecules are chemically stable or not, and it is a critical parameter in determining molecular electrical transport properties because it is a measure of electron conductivity. In addition, the energy of the HOMO is directly related to the ionization potential, while the amount of electron affinity is related to the energy of the LUMO. The energy gaps are largely responsible for the chemical and spectroscopic properties of the molecules [18].

The features of the HOMO and LUMO are calculated by B3LYP/6-311++G(3df, 3pd) level in gaseous phase. The energy band gap () (transition from HOMO to LUMO), chemical potential (μ), chemical hardness (), ionization potential, and electron affinity of the HSOH, HSO, and SOO species can be seen in Table 4.

Table 4: Hardness (η), chemical potential (μ), HOMO energy (), and LUMO energy () (in atomic unit) of HSOH and HSO molecules.

μ
η

HSOH
−0.2561
−0.0509
0.1535
−0.1026
0.2053

HSO
−0.2651
−0.02933
0.1472
−0.1179
0.2358

SOO
−0.28904
−0.1953
0.2422
−0.04685
0.0937

The HOMO-LUMO energy gap in Table 4 predicts the HSO (0.2358 au) with a large energy gap () are more stable than the HSOH (0.2052 au) and SOO (0.0937 au) molecules.

The frontier molecular orbitals of HSO and HSOH are depicted in Figure 4. In HSO molecule, the LUMO is localized on the bond between oxygen-sulfur atoms. The HOMO is mainly composed on the p-type atomic orbitals of oxygen and sulfur atoms. In SOO molecule, the HOMO orbital is localized on the sulfur atom while the LUMO is composed in the π type orbitals of all atoms.

Figure 4: Frontier molecular orbitals of HSO, HSOH, and SOO molecules.

The chemical potential (μ) and chemical hardness () [18] have been determined using HOMO and LUMO energies which are calculated using the relations and , where is the ionization potential () and is the electron affinity () of the molecular system. According to the maximum hardness principal, the most stable anions should have maximum hardness value which is a minimum energy structure at a constant chemical potential [18]. The data in Table 4 show that the maximum absolute value of chemical hardness, −0.1179 au, belongs to HSO molecule. Also, the minimum chemical potential, 0.1472 au, of HSO molecule is in agreement with the predicted chemical potential.

Determination of molecular charge distribution is fundamental to the application of the molecular modeling method to a realistic chemical system. The variation in the molecular charge distribution with structure can be applied to explain to stability of species. Atomic charges found with our calculations inside atoms are shown in Figure 5. All atomic charges for HSOH, HSO, and SOO molecules are obtained by the Mulliken method at B3LYP/6-311++G(3df, 3pd) level. Results show that sulfur atoms have positive charge in three species and the oxygen atoms take negative charges in both HSOH and HSO molecules and, in SOO, the terminal oxygen atom has negative charge while the other atoms have positive charge as shown in Figure 5. But absolute value of charge distribution of oxygen atom in SOO is smaller than the HSO molecule and its oxygen negative charges smaller than the HSOH molecule that corresponds with chemical potential and HOMO-LUMO molecular orbital analysis.

Figure 5: HSO and HSOH calculated partial charge distribution with the Mulliken method.

Dipole moment shows the molecular charge distribution and is given as a vector in three dimensions. Therefore, it can be used as a descriptor to depict the charge movement across the molecule. The direction of the dipole moment vector in a molecule depends on the centers of positive and negative charges. In HSO species, dipole moment vector (2.298 D) is from sulfur atom to oxygen atom parallel with the S–O bond length but in HSOH molecule this vector (1.594 D) is to the outside of the molecule and maybe can cause the low stability in the molecule. For SOO molecule, this vector (1.525 D) is from terminal oxygen atom to sulfur atom.

Rotational constants of the HSOH molecule are 202.1, 15.2, and 14.7 GHz and for HSO are 298.3, 20.3, and 19.0 GHz and for SOO species are 95.2, 6.9, and 6.4 GHz. Because the experimental value of the rotational constants are not accessible, so our calculated results are reported in an absolute form.

2.4. The Thermodynamic Data in the O₃ and H₂S Reaction Process

The change of thermodynamic characteristics for each reaction channel is the difference between the corresponding thermodynamic properties of the products and the reactants. The thermodynamic data are calculated by G4 level of computation in gaseous phase and are corrected by ZPE for O₃ + H₂S reaction.

The calculated relative internal energies, enthalpies, Gibbs free energies, and entropies of all products of the reactions in gas-phase at atmospheric pressure and temperature of 298.15 K, have been summarized in Table 5.

Table 5: Relative energies , enthalpies , free energies and entropies (kJ mol⁻¹) for reaction of ozone and hydrogen sulfide with ZPE correction at B3LYP level of computation.

Reaction

H₂S + O₃ HSOH + ¹O₂
−117.09
−117.09
−1.25
115.85

H₂S + O₃ HSOH + ³O₂
−235.23
−235.23
1.47
−236.70

H₂S + O₃ HSO + HOO
−135.15
−135.16
7.98
−143.13

H₂S + O₃ SO₂ + H₂O
−655.63
−655.63 (−660.72)
−1.93 (2.29)
−653.70 (−718.4)

Table 5 data show that the Δ of all overall product channels are negative that means the reaction pathways are exothermic. This value for the SO₂ + H₂O path is greater than the others with high excess of heat. Also, Gibbs free energy values of all production channels are negative so, all product channels are spontaneous in gas phase at atmospheric pressure and 298.15 K. The SO₂ + H₂O path product channel thermodynamically is the most feasible in comparison with the other products by −653.70 kJ/mol in the Gibbs free energy. Our calculated data are compared with the experimental value for the SO₂ + H₂O reaction channel and the obtained results are acceptable in comparison with the experiment (9 percent relative error in Gibbs free energy).

2.5. Calculation of the Rate Constants

Rate constant is a key parameter in kinetic study of reaction that can be calculated from statistical mechanics and its equation is obtained from transition state theory (TST) based on statistical mechanics [17, 18]. For bimolecular reactions, one traditionally uses transition state theory, which requires knowledge of an accurate transition structure so the activation energy is known. The following is a familiar equation for calculating the rate constant for a bimolecular reaction, originally derived from Eyring [19]:

where is the tunneling correction; is the Avogadro constant; is the energy barrier between TS and reactants including zero point energy corrections; , , and are the partition functions of the transition structure and reactant species, respectively. The parameters for one species are the product of translational, rotational, vibrational, and electronic partition functions, is Boltzmann’s constant, and is Planck’s constant.

The translational partition function per unit volume of a molecule of mass in three dimensions is

The external rotational partition function for a nonlinear molecule has the following form:

here, the parameters are the moments of inertia and is a symmetry number. As two reactants approach each other, three (for nonlinear molecules) or two (for linear molecules) external rotations would transform to hindered rotations or harmonic vibrations as they approach the transition state.

Vibrations in the transition state are usually treated as harmonic oscillators. The vibrational partition function for one degree of freedom is calculated as

where, is vibrational frequency. According to classical mechanics, a particle can cross a potential energy barrier only if its energy is higher than the barrier height, while quantum mechanics predicts that there is a probability for the crossing at lower particle energies. This transmission probability is very sensitive to the particle energy and mass, the barrier height, and the shape of the barrier. The tunneling correction factor is defined as the quotient of the quantum mechanical rate to the classical rate. Although the tunneling process is a multidimensional phenomenon, for simplicity this is treated as a one-dimensional process. Also, one possible cause for a deviation from Arrhenius behavior is quantum mechanical tunneling of reactants through classical barrier which was assumed conservation of vibrational frequencies and moments of inertia of the fragments during the course of reaction. In this study, we used the simplest approach to assess the role of quantum tunneling which is the Wigner correction to the reaction rate [20]. It is related to the magnitude of the imaginary frequency along the reaction coordinates at the transition state. The rate constant is enhanced by a factor of

where is imaginary frequency with negative value that accounts for the vibrational motion along the reaction path, and is the speed of light.

The rate constant for unimolecular reaction channel is calculated by RRKM (Rice-Ramsperger-Kassel-Marcus) theory. In RRKM theory, the thermal rate coefficient is obtained by integrating over from to infinity of the following equation:

where , , , and are defined based on the Lindemann mechanism for the first-order chemical reactions. The various terms of the rate expression are now evaluated using statistical mechanics. In the high pressure limit, RRKM theory reduces to transition state theory and the results of the two theories coincide.

The rate constants of the reaction have been calculated on the singlet potential energy surfaces for two pathways which are defined as thermodynamic path (P1(1)) and kinetic path (P2) by the RRKM [20] (multichannel reaction) and TST [19] (one step reaction) theories implemented in the Ssumes [21] and Gpop [22] programs.

For multichannel reactions (path P1(1)), total rate constant are evaluated in the high pressure limit of the Lindemann mechanism. In the high pressure limit, RRKM theory reduces to transition state theory and the results of the two theories coincide. High pressure limit of total rate constant is defined as based on the Lindemann mechanism, where , , and are the rate expressions for the forward steps of the reaction that and and are the reverse steps of reaction follow unimolecular reactions. For the channels with one transition state, P2, total rate constant is calculated using the transition state theory for bimolecular reactions (Eyring equation (1)).

On the singlet potential energy surfaces, four reaction pathways, P1(1), P1(2), P2, and P3, are obtained. The rate constants have been calculated for two paths (P1(1) as thermodynamic path and P2 as kinetic path). The first pathways, P1(1), has the combination of two steps, so, it is a complex reaction. One is the formation of SOO and H₂O as products and, the next, SOO molecule rearranged to SO₂ molecule as a final product after passing TS6. In the rate constant calculation the final step (SOO→SO₂) was ignored because it is a unimolecular reaction. So, it is very faster than the formation of SOO steps. Also, the number of the atoms in the transition state is different from the number of the atoms in the other TSs which is a limitation for the software. The HSO and HOO adducts produce using two channels. The first channel for production of HSO and HOO molecules is passing TS3 and the second is passing TS1, TS2, and TS4. Between them the mechanism that includes TS3 is more reliable than the other three-step channel and the rate constant is calculated for the reliable path (P2).

For the multichannel reaction of H₂S and O₃, the calculated rate constant for two paths are denoted as and which correspond with P1(1) and P2 paths, respectively. Table 6 shows the individual rate constant of two paths at the temperature range of 200–2500 K at B3LYP/6-311++G(3df, 3pd) level. The rate constant values for the paths in Table 6 reveal that the path P2 is kinetically the main path in all the temperature ranges. Therefore the production of HSO and HOO molecules is the main product. Because of low level energy of transition state from TS3 channel, production of HSO + HOO is the most important adducts (path P2 and corresponding ) in comparison with the other paths. The value of rate constant shows in competition between P1(1) and P2 paths, path P2 is more reliable than the path P1(1). On the other word, the reaction of H₂S + O₃ kinetically proceeds HSO + HOO as the main products on the singlet potential energy surface through path P2 via the rate constants of . In thermodynamic aspect, the H₂O + SO₂ adducts are the main products on the singlet potential energy surface through path P1(1) via the rate constants of . Figure 6 shows the Arrhenius plots for the and rate constants.

Table 6: The overall rate constant and branching ratio of the reaction (in cm³ molecule⁻¹ s⁻¹), calculated at the temperature range of 200–2500 K at the B3LYP/6-311++G(3df,3pd).

T

/
/

200
9.99E − 38
9.17E − 24
9.17E − 24
1.09E − 12
1.00E + 02

300
4.66E − 29
3.63E − 21
3.63E − 21
1.28E − 06
1.00E + 02

400
9.33E − 25
1.79E − 19
1.79E − 19
5.21E − 04
1.00E + 02

500
3.49E − 22
2.40E − 18
2.40E − 18
1.45E − 02
1.00E + 02

600
1.83E − 20
1.56E − 17
1.56E − 17
1.18E − 01
9.99E + 01

700
3.15E − 19
6.47E − 17
6.51E − 17
4.84E − 01
9.95E + 01

800
2.71E − 18
2.02E − 16
2.04E − 16
1.33E + 00
9.87E + 01

900
1.48E − 17
5.14E − 16
5.28E − 16
2.80E + 00
9.72E + 01

1000
5.82E − 17
1.13E − 15
1.19E − 15
4.90E + 00
9.51E + 01

1400
2.22E − 15
1.08E − 14
1.30E − 14
1.71E + 01
8.29E + 01

1900
2.94E − 14
6.31E − 14
9.25E − 14
3.18E + 01
6.82E + 01

2500
1.95E − 13
2.59E − 13
4.54E − 13
4.29E + 01
5.71E + 01

Figure 6: Arrhenius plot for the (–∎–) and (−▲−) rate constants of H₂S + O₃ reaction.

Arrhenius parameters for the and rate constants are as follows:

The overall rate constant for of H₂S + O₃ reaction is denoted as which is the sum of rate constant of two paths ( + ). The branching ratios for the paths P1(1) and P2 are and , respectively. Table 6 shows the individual rate constant of paths, overall rate constants, and the branching ratios at the temperature range of 200–2500 K. The rate constant and branching ration values for several paths in Table 6 reveals that the path P2 is kinetically the main path in temperature lower than 1000 K. Therefore, the production of HSO and HOO molecules is the main adducts. In temperature higher than 1000 K, the production of SO₂ and H₂O molecules from TS1 and TS2 channel (path P1 (1) and corresponding and ) is the most important channels can kinetically compete with the HSO and HOO products. However, the production of SO₂ and H₂O molecules is thermodynamically important production. Figure 6 shows the Arrhenius plots for the and rate constants.

As shown in Table 6, the rate constant values of H₂S + O₃ under atmospheric pressure increase with increasing the temperature. Therefore, the main products of reaction depend on the temperatures and the reaction mechanism, so, this reaction can play different roles in the formation of acid rain.

3. Conclusion

The reaction of H₂S and O₃ has been investigated at the B3LYP/6-311++G(3df,3pd), CCSD/6-311++G(3df,3pd)//B3LYP/6-311++G(3df,3pd), and G4 levels of computation. By dissociation of H₂S-O₃ complex, four product channels have been found, two channels for production of HSO and HOO products, one channel for production of SO₂ and H₂O, and another one for the production of HSOH and O₂ molecules via corresponding transition states and complexes. AIM analyses were used to elucidate the reaction mechanism by comparing of BCP in complexes and transition states. Furthermore, AIM indicates ring structure for some transition structures and complexes. Additionally intrinsic reaction coordinate (IRC) calculation was performed to examine the correspondence of the calculated TS to the reactants and products. The energy profile of reaction mechanism is obtained. The thermodynamic results show that the SO₂ + H₂O product is the main product of H₂S + O₃ reaction with the most negative value of Δ. Also the rate constant and branching ratio are calculated and discussed. Results show two channels are competing with each other based on the temperature. The HSO + HOO adducts are the main kinetically products while the SO₂ + H₂O molecules are the main products thermodynamically.

References

A. Mansergas and J. M. Anglada, “The gas-phase hydrogen-bonded complex between ozone and hydroperoxyl radical. A theoretical study,” The Journal of Physical Chemistry A, vol. 111, no. 5, pp. 976–981, 2007.

Z. M. Xu and M. C. Lin, “Ab initio study on the kinetics and mechanisms for O₃ reactions with HO₂ and HNO,” Chemical Physics Letters, vol. 440, no. 1–3, pp. 12–18, 2007.

D. Bing, Y. Zhao, F. Hao et al., “Ab initio study on the reaction mechanism of ozone with bromine atom,” International Journal of Quantum Chemistry, vol. 107, no. 5, pp. 1085–1091, 2007.

W. T. Chan, C. Weng, and J. D. Godard, “Concerted and stepwise reaction mechanisms for the addition of ozone to acetylene: a computational study,” The Journal of Physical Chemistry A, vol. 111, no. 22, pp. 4792–4803, 2007.

A. Mansergas and J. M. Anglada, “Reaction mechanism between carbonyl oxide and hydroxyl radical: a theoretical study,” Journal of Physical Chemistry A, vol. 110, no. 11, pp. 4001–4011, 2006.

J. Yang, Q. S. Li, and S. Zhang, “Reaction-path dynamics and theoretical rate constants for the reaction CH₄ +  O₃→ HOOO + CH₃,” International Journal of Quantum Chemistry, vol. 107, no. 10, pp. 1999–2005, 2007.

S. M. Resende and F. R. Ornellas, “Atmospheric reaction between the HS radical and chlorine,” Chemical Physics Letters, vol. 318, no. 4-5, pp. 340–344, 2000.

S. H. Mousavipour, M. A. Namdar-Ghanbari, and L. Sadeghian, “A theoretical study on the kinetics of hydrogen abstraction reactions of methyl or hydroxyl radicals with hydrogen sulfide,” Journal of Physical Chemistry A, vol. 107, no. 19, pp. 3752–3758, 2003.

Y. A. Aleksandrov, A. V. Belyaev, A. Y. Sozin, and M. F. Churbanov, “Hydrogen sulfide oxidation by ozone on sorbents-catalysts,” Russian Journal of Applied Chemistry, vol. 83, no. 9, pp. 1548–1552, 2010.

R. J. A. Grim and J. M. Greenberg, “Photoprocessing of H₂S in interstellar grain mantles as an explanation for S₂ in comets,” Astronomy and Astrophysics, vol. 181, no. 1, pp. 155–168, 1987.

D. L. Singieton, G. Paraskevopoulos, and R. S. Irwin, “Mechanism of the atomic oxygen(3P) + hydrogen sulfide reaction. Abstraction or addition?” The Journal of Physical Chemistry, vol. 86, no. 14, pp. 2605–2609, 1982.

D.-K. Park, Y.-G. Cho, D.-W. Park, and H.-C. Woo, “Kinetic study on the selective oxidation of H₂S containing excess water and ammonia,” Reaction Kinetics and Catalysis Letters, vol. 74, no. 1, pp. 57–66, 2001.

C. Wilson and D. M. Hirst, “Reaction of H₂S with OH and a study of the HSO and SOH isomers. High-level ab initio calculations,” Journal of the Chemical Society, Faraday Transactions, vol. 90, no. 20, pp. 3051–3059, 1994.

R. F. W. Bader, Atoms in Molecules: A Quantum Theory, Oxford University Press, Oxford, UK, 1990.

F. Biegler-Konig and J. Schoenbohm, AIM2000, Buro fur Innovative Software, Bielefeld, Germany, 2nd edition, 2002.

M. J. Frisch, G. W. Trucks, H. B. Schlegel et al., Gaussian 03, Revision B.03. Gaussian, Inc., Pittsburgh, Pa, USA, 2003

P. W. Atkins, Physical Chemistry, Oxford University Press, Oxford, UK, 2001.

R. G. Pearson, Chemical Hardness, John Wiley & Sons, Weinheim, Germany, 1997.

H. Eyring, “The activated complex in chemical reactions,” The Journal of Chemical Physics, vol. 3, no. 2, p. 107, 1935.

W. Forst, Theory of Unimolecular Reactions, Academic Press, New York, NY, USA, 1973.

A. Miyoshi, Steady-State Unimolecular Master-Equation Solver (SSUMES), University of Tokyo, 2010

A. Miyoshi, Gaussian Post Processor (GPOP), University of Tokyo, 2010