Environ. Sci. Technol. 2009, 43, 6220–6227
Development of a Group Contribution Method To Predict Aqueous Phase Hydroxyl Radical (HO•) Reaction Rate Constants DAISUKE MINAKATA,† KE LI,‡ PAUL WESTERHOFF,§ AND J O H N C R I T T E N D E N * ,† School of Civil and Environmental Engineering, Georgia Institute of Technology, 800 West Peachtree Street, Suite 400 F-H, Atlanta, Georgia 30332, Faculty of Engineering, University of Georgia, Athens, Georgia 30602, and Department of Civil, Environmental and Sustainable Engineering, Arizona State University, P.O. Box 5306, Tempe, Arizona 85287
Received April 1, 2009. Revised manuscript received July 1, 2009. Accepted July 6, 2009.
The hydroxyl radical (HO•) is a strong oxidant that reacts with electron-rich sites of organic compounds and initiates complex chain mechanisms. In order to help understand the reaction mechanisms, a rule-based model was previously developed to predict the reaction pathways. For a kinetic model, there is a need to develop a rate constant estimator that predicts the rate constants for a variety of organic compounds. In this study, a group contribution method (GCM) is developed to predict the aqueous phase HO• rate constants for the following reaction mechanisms: (1) H-atom abstraction, (2) HO• addition to alkenes, (3) HO• addition to aromatic compounds, and (4) HO• interaction with sulfur (S)-, nitrogen (N)-, or phosphorus (P)-atom-containing compounds. The GCM hypothesizes that an observed experimental rate constant for a given organic compound is the combined rate of all elementary reactions involving HO•, which can be estimated using the Arrhenius activation energy, Ea, and temperature. Each Ea for those elementary reactions can be comprised of two parts: (1) a base part that includes a reactive bond in each reaction mechanism and (2) contributions from its neighboring functional groups. The GCM includes 66 group rate constants and 80 group contribution factors, which characterize each HO• reaction mechanism with steric effects of the chemical structure groups and impacts of the neighboring functional groups, respectively. Literature-reported experimental HO• rate constants for 310 and 124 compounds were used for calibration and prediction, respectively. The genetic algorithms were used to determine the group rate constants and group contribution factors. The group contribution factors for H-atom abstraction and HO• addition to the aromatic compounds were found to linearly correlate with the Taft constants, σ*, and electrophilic substituent parameters, σ+, respectively. The best calibrations for 83% (257 rate constants) and predictions for 62% (77 rate constants) of the rate constants were within 0.5-2 times * Corresponding author phone: (404)894-5676, fax: (404)894-7896, e-mail:
[email protected]. † Georgia Institute of Technology. ‡ University of Georgia. § Arizona State University. 6220
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 43, NO. 16, 2009
the experimental values. This accuracy may be acceptable for model predictions of the advanced oxidation processes (AOPs) performance, depending on how sensitive the model is to the rate constants.
Introduction The hydroxyl radical (HO•) is a reactive electrophile that reacts rapidly and nonselectively with most electron-rich sites of organic contaminants. It is the active species that potentially leads to complete mineralization of emerging contaminants in advanced oxidation processes (AOPs) (e.g., O3/H2O2, UV/H2O2, UV/TiO2, etc.) and natural water (1-3). Because of concerns of emerging contaminants and the large number of chemicals that are in use or being produced (4), there is a need for kinetic models that can quickly assess their removal by AOPs. The three critical components for building a kinetic model to predict AOPs performance are (1) numerical methods that solve ordinary differential equations (ODEs), (2) algorithms that can predict reaction pathways, and (3) algorithms that can predict reaction rate constants. Many kinetic models have been built for known reaction pathways (5-8). Recently, a model that can generate reaction pathways for the aqueous phase AOPs has been developed (9). Yet the capability of generating rate constants on-the-fly for the aqueous phase radical reactions is still limited because of the complexity of radical chemistry. A number of studies have been conducted to predict the HO• reaction rate constants in the gaseous phase, including a group contribution method (GCM) (10-12), bond dissociation energies (BDEs) (13), computational molecular orbital OH radical (MOOH) methods (14, 15), extrapolation of transition state theory (16), and neural networks (NNs) (17). When the experimental rate constants for compounds with different functional groups were available, the GCM was proven to be robust for the prediction of gas phase rate constants for compounds with a wide range of functional groups. Atkinson’s GCM is one of the most widely accepted methods and is implemented in the U.S. Envirnomental Protection Agency software, AOPWIN (18). In their method, the HO• rate constant was determined by the reaction mechanism and effect of the neighboring functional group of the compound of interest. For each reaction mechanism such as H-atom abstraction, HO• addition, and interactions of HO• with S-, N-, or P-atom-containing compounds, there were “group rate constants” that represented the reaction mechanism and “substitutent factors” that represented the effects of neighboring and next-nearest neighboring functional groups. Using this GCM, the gaseous-phase rate constants at 298 K of ∼90% of approximately 485 organic compounds were predicted within a difference of a factor of 2 from the experimental values (12). Although the GCM is successful in the gaseous phase, applying it to the aqueous phase requires carefully discerning the mechanistic differences between aqueous and gaseous phase reactions. The hydrogen bond and polarity of water molecules drastically change the reaction mechanisms for radicals. For example, the 1,2-H shift of the oxyl radical in the aqueous phase can proceed with a rate constant larger than 106 s-1, yet it has not been observed in the gaseous phase (19). Another example is the solvent cage effect (20). In addition, the solvation process affects the free energy of activation of reactions (21) and, thus, changes the rate constant. The impact varies significantly for different radicals and functional groups. There are only a few studies concerned with the prediction of aqueous phase reaction rate constants in AOPs. Monod 10.1021/es900956c CCC: $40.75
2009 American Chemical Society
Published on Web 07/17/2009
et al. (22) applied the GCM to estimate the aqueous phase HO• rate constants for the oxygenated aliphatic compounds. They reported that 84% of their 128 calculated data were within a difference of a factor of 5 from the experimental values. However, the parameters that represented the property of functional groups in their GCM did not follow trends that would be expected on the basis of the reaction mechanisms. Herrmann (23) used the correlation between the bond dissociation energy (BDE) (uncertainty of (8.4 kJ/ mol) of the weakest C-H bond broken and Arrhenius activation energy, Ea, and obtained the Evans and Polanyi relationship for 16 oxygenated compounds with a 0.75 correlation coefficient. Dutot et al. (24) used an artificial neural network and multilayer perception (MLP) to estimate HO• rate constants, and 87% of the MLP predictions were within a factor of 2 from the experimental data. In this study, a GCM is developed to predict the aqueous phase HO• rate constants that integrate the reaction mechanisms in the aqueous phase and the essential features of Atkinson’s GCM. The potential errors and limitations of the GCM will be discussed herein. Development of the Group Contribution Method (GCM). The GCM is based on Benson’s thermochemical group additivity (25). It hypothesizes that an observed experimental rate constant for a given organic compound is the combined rate of all elementary reactions involving HO•, which can be estimated using Ea. For each reaction mechanism, there is a base activation energy, E a0, and a functional group contribution of activation energy, E Ra i, due to the neighboring (i.e., R-position) and/or the next-nearest neighboring (i.e., β-position) functional group (i.e., Ri). These contributions to the rate constant can be parametrized and determined empirically when sufficient rate constant data are available. The GCM considers four reaction mechanisms that HO• initiates in the aqueous phase, which include (1) H-atom abstraction, (2) HO• addition to alkenes, (3) HO• addition to aromatic compounds, and (4) HO• interaction with sulfur (S)-, nitrogen (N)-, or phosphorus (P)-atom-containing compounds. Accordingly, an overall reaction rate constant, koverall, may be given in eq 1 koverall ) kabs + kadd-alkene + kadd-aromatic + kint
(1)
where, kabs, kadd-alkene, kadd-aromatic, and kint are the rate constants for the aforementioned reaction mechanisms, 1-4, respectively. The manner in which these rate constants are determined is discussed in the following section. Hydrogen-Atom Abstraction. For H-atom abstraction, the active bond is a C-H bond. According to the functional groups on the C atom, there are primary, secondary, and tertiary C-H bond(s), except in the special case of methane. Therefore, the fragments of a molecule are CH3R1, CH2R1R2, and CHR1R2R3, where Ri is a functional group (i ) 1-3). Each of the fragments corresponds to a partial rate constant kCH3R1, kCH2R1R2, and kCHR1R2R3, respectively. The Ea for the reaction of HO• is affected by the C-H bond itself and adjacent functional group(s). The contribution that results from the C-H bond to the Ea is defined as base activation energy, Ea0, while the contribution of the functional groups is defined as Ri , due to the functional a group contribution parameter, E a,abs group Ri for H-atom abstraction. For example, the base activation energy for H-atom abstraction from one of the i indicates the electronprimary C-H bonds is E0a, prim. The E Ra,abs donating and -withdrawing ability of the functional group. An electron-donating functional group decreases the Ea and, hence, increases the overall reaction rate constant and vice versa. Accordingly, the partial rate constant for the fragmented part, CH3R1, can be written as below R1
kCH3R1 ) 3Aprime-(Ea,prim+Ea,abs)/RT 0
(2)
where 3 is the amount of primary C-H bonds, Aprim denotes the Arrhenius frequency factor for the reaction of HO• with CH3R1, R is the universal gas constant, and T denotes absolute temperature. Similarly, the partial rate constants for other fragmented parts CH2R2 and CHR1R2R3 are expressed in eqs 3 and 4, using the corresponding frequency factors, Asec and Atert, and group contribution parameters, E Ra,2 abs and E Ra,3 abs, respectively. R1
R2
kCH2R1R2 ) 2Asece-(Ea,sec+Ea,abs+Ea,abs)/RT 0
R1
R2
(3)
R3
kCHR1R2R3 ) Aterte-(Ea,tert+Ea,abs+Ea,abs+Ea,abs)/RT 0
(4)
However, for eqs 2-4, the functional group contribution is ignored for cases where the neighboring functional groups -H have no effect on the H-atom abstraction (i.e., Ea, abs is zero, where a valence bond of the H atom is expressed as a line before H). The group rate constants that represent H-atom abstraction from the primary, secondary, and tertiary C-H bond 0 0 0 , ksec , and ktert and are expressed in eqs are defined as kprim 5-7. 0 0 kprim ) Aprim e-(Ea,prim)/RT
(5)
0 0 ksec ) Asec e-(Ea,sec)/RT
(6)
0 0 ktert ) Atert e-(Ea,tert)/RT
(7)
0
0
0
In addition, the group rate constant kR4 is defined for the HO• interaction with the functional group R4 (e.g., -OH and -COOH). The detailed discussions for R4 will be given in the following sections. The group contribution factor, XRi, that represents the influence of functional group Ri is defined in eq 8 Ri
XRi ) e-(Ea,abs)/RT
(8)
Because each reaction is independent from one another, the rate constant for H-atom abstraction, kabs, may be written as the sum of the partial rate constants I
kabs ) 3
∑k 0
J
0 primXR1
+2
∑k 0
0 secXR1XR2
+
K
∑k 0
0 tertXR1XR2XR3
+ kR4 (9)
where, I, J, and K denote the number of the fragments CH3R1, CH2R2, and CHR1R2R3, respectively. Equation 10 shows an example of the rate constant calculation for 1,2-dichloro3-bromopropane (CH2Cl-CHCl-CH2Br). The detailed definitions of the functional groups will be given in the Results and Discussion. 0 0 koverall ) 2ksec X-ClX-CHCl- + ktert X-ClX-CH2ClX-CH2Br + 0 2ksec X-BrX-CHCl- (10)
Rate constant additivity can be an important concern. To investigate this, the experimental rate constants for the linear and longer chain compounds with some functional groups were compared with the number of -CH2- alkyl functional groups as shown in Figure SI1-1 of Supporting Information 1. The linear relationship was observed, and it verified the rate constant additivity. For nonlinear compounds, nonadditive effects arise from different contributions to entropy of activation (26) and hence affect A. Generally, A can be VOL. 43, NO. 16, 2009 / ENVIRONMENTAL SCIENCE & TECHNOLOGY
9
6221
assumed to be constant for the same reaction mechanism because intrinsic entropy due to translation and rotation can be assumed to be constant. Strictly speaking, however, moments of inertia of primary, secondary, and tertiary C-H bonds would be different and thus affect entropy (27). Although this difference is insignificant when compared to other factors such as the impacts of neighboring functional groups (26), the GCM considers different Aprim, Asec, and Atert. Therefore, this should be sufficient to treat nonlinear compounds. HO• Addition to Alkenes. The detailed mechanisms of HO• addition to alkene in the aqueous phase have not been elucidated yet. It is expected that molecular solvation affects the Ea, and water molecules interact with the HO• approaching the carbon double bonds. Nevertheless, it can be a reasonable first effort assumption that there is little difference in the major reaction mechanisms between the aqueous and gaseous phases (28-30). Considering a basic CdC double bond structure (i.e., one σ-bond and one π-bond, >CdCCdCCdCCdC< H(trans), and H2CdCC)C< Y-ClY-ClY-ClY-Cl
(14)
HO• Addition to Aromatic Compounds. Although the HO• preferentially adds to the aromatic ring at the close to diffusion-controlled rates, the HO• addition to the aromatic ring can be expected to be highly regioselective due to the electrophilic HO•. Therefore, the electron-donating and -withdrawing functional groups on the aromatic ring can significantly affect the rate constants and the ratio of ortho-, meta-, para-, and ipso-positions. For the HO• addition to aromatic compounds, the following points are considered: (1) Probability for the symmetrical HO• addition to the benzene ring is identical. 6222
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 43, NO. 16, 2009
(2) Addition to the ipso-position is negligible due to the significant steric effect (31-33). Although some studies report the possibility for the addition to the ipso-position (34), it is quite negligible for the aromatic compounds with single functional groups, which are used for the calibration (e.g., C< were assumed to be identical because of the following reasons: (1) limited data availability for the >C< functional group and (2) identical BDEs affected by the corresponding functional group (44). When the alkyl functional groups had the ether functional groups at both sides of the R-positions, another group contribution factor X-O-(second) was considered. Group contribution factors for carboxylic and ester functional groups, -COOH and -COOR, were assumed to be identical because of their similar electron-withdrawing ability. Because of their strong electron-withdrawing ability, the influence of the β-position that resulted from the halogen functional groups was considered. This represented X-CHn(halogen)m where (n,m) ) (1,1), (1,2), or (2,1) for Cl or Br atom, (i.e., X-CH2Cl ≈ X-CHCl2 ≈ X-CHCl- ≈ X-CH2Br ≈ X-CHBr2 ≈ X-CHBr- ), respectively, which was assumed to be identical for the purpose of reducing the number of group contribution factors. The values of X-F, X-I, and X-CF2- were not available because no experimental data were available for these functional groups. Ring structures were expected to increase Ea because of the extra energy to form the ring structure; consequently, group contribution factors for the 3-ring XRS3 and the 5-ring XRS5 for saturated cyclic compounds were considered. 6224
9
ENVIRONMENTAL SCIENCE & TECHNOLOGY / VOL. 43, NO. 16, 2009
The group contribution factors for the H-atom abstraction linearly correlate with the Taft constant, σ* (45), as shown in Figure 2. Because the alkyl functional groups contribute to weakening the C-H bond with release of the steric compression as the alkyl functional group moves apart to form a planar radical, they increase in the HO• reactivity in the H-atom abstraction reactions. Therefore, X-CH3 and X-CH2≈ X>CH- ≈ X>C< values are greater than 1.0, which correspond to negative values of the Taft constant. All of the group contribution factors for the oxygenated and halogenated functional groups show negative values, which indicate the electron-withdrawing ability of the functional groups. HO• Addition to Alkenes. Twenty-eight alkenes were categorized based on the base CdC double-bond structures and position(s) of hydrogen atom(s) adjacent to those base structures. Table SI1-4 of Supporting Information 1 summarized the structural configurations for eq 13 on the basis of the compounds that were available in the literature. There were few rate constants reported for the conjugated and unconjugated dienes. Therefore, these compounds were excluded from calibration. As shown in Table SI1-1 of Supporting Information 1, the total DF were 14. A total of 79% of the calibrated data was within the EG. The SD was 0.4. Group Rate Constants and Group Contribution Factors for HO• Addition to Alkenes. Group rate constants and group contribution factors are summarized in Table SI1-5 of Supporting Information 1. It was found that the group contribution factors did not linearly correlate with the Taft constant (data not shown). Two reasons can be considered. First, the functional group contribution to the Ea does not follow the general inductive effect (i.e., Taft constant). For example, the group contribution factor for alkane does not indicate the electron-donating value (i.e., Y-alkane ) 0.171). Because of limited data availability and gaseous phase reaction mechanisms (46, 47), the rate constant expression as shown in eq 13 did not consider the different effects of functional groups to two unsaturated carbons (one being attacked and the neighboring one). For example, the effects of three chlorine functional groups for 1,1,2-trichloroethylene were treated identical in association with the two unsaturated carbons in the >CdCH- base structure. Although Peeters et al. (48) successfully applied the GCM to the gaseous phase HO• rate constants with nonconjugated and conjugated
(poly) alkenes on the basis of the number of functional groups on the neighboring carbon, their approach was limited to nonpolar functional groups, and it was not clear if it was applicable to the heteroatom functional groups and aqueous phase reactions. Second, the experimental rate constants do not seem to follow the inductive effect (e.g., vinyl chloride > ethylene > vinyl alcohol) because of experimental errors or the existence of unknown reaction mechanisms. We suggest more experimental studies be conducted for alkenes to confirm this. Despite the observation of the nonlinear correlation between the group contribution factors and the Taft constant, 79% of the calibrated rate constants were within the EG, and this might be acceptable for a rate constant estimator. It should be addressed that more quantitative investigations are required to examine the effect of the functional groups in the aqueous phase. HO• Addition to Aromatic Compounds. Table SI1-6 of Supporting Information 1 summarized the structural configurations for eq 17 on the basis of the compounds that were available in the literature. As shown in Table SI1-1 of Supporting Information 1 for the overall HO• addition to aromatic compounds, the total DF were 35. 88% of the rate constants from the calibration was within the EG. The SD was 0.73. Once the group rate constants and group contribution factors were calibrated, the rate constants for 46 compounds were predicted as shown in Table SI1-2 of Supporting Information 1. A total of 64% of the rate constants from the prediction was within the EG. The SD was 0.53. Group Rate Constants for HO• Addition to Aromatic Compounds. Table SI1-7 of Supporting Information 1 summarized the group rate constants for aromatic compounds. The approximate magnitude of the group rate constants can be compared with theoretical studies because there is no experimental data for Ea. According to Ashton et al. (49), Ea between 13 and 21 kJ/mol of the net reaction was observed for various aromatic compounds. From their conclusion that the electron-withdrawing functional groups increased the Ea by approximately 2.1-8.4 kJ/mol, the magnitude of group rate constants fell in the range of 109 orders. This is consistent with the group rate constants from the calibration. Because the group rate constants for pyridine, furan, imidazole, and triazine compounds include the hetero atoms (e.g., N, NH, and O) that affect the HO• addition to aromatic ring, the calibrated group rate constants vary in magnitude of order from 106 to 109 M-1 s-1. For example, the group rate constant, k0(2-pyr)-3,6, that represents one functional group on 2-position and the addition of HO• to either the 3- or 6-position is 9.9 × 108 M-1 s-1. This value and other group rate constants for pyridine are smaller than those of benzene compounds because of the lower reactivity and lower electrophilic nature of the pyridine nucleus toward the HO• (50). For triazine compounds, the group rate constant k0(1,3,5-triaz)-2,4,6 (i.e., functional groups located on the 1,3,5positions and HO• interacts with each N-atom located on the 2,4,6-positions) represents the reactivity of HO• with the triazine base structure other than the functional groups. However, this single group rate constant does not seem to represent the reaction mechanisms well. Additional group rate constants may be required and can be calibrated when more experimental data are available. Group Contribution Factors for HO• Addition to Aromatic Compounds. Table SI1-7 of Supporting Information 1 summarized the group contribution factors for the HO• addition to aromatic compounds. Figure 2 plots those values against electrophilic substituent constants, σ+ (40), for benzene (r ) 0.89), pyridine (r ) 0.93), and furan (r ) 0.65) compounds. Because only one σ+ for imidazole functional groups is available, Figure 2 does not include the plot for imidazole. These observations validate that the group contribution
factors that are empirically derived from the experimental rate constants linearly correlate with the general electrondonating and -withdrawing property. For furan compounds, the weak correlation may be due to insufficient number of data. HO• Interactions with Sulfur-, Nitrogen-, or PhosphorusContaining Compounds. As shown in Table SI1-1 of Supporting Information 1, the DF were 54. 74% of the calibrated data, which was within the EG. The SD was 1.5. Once the group rate constants and group contribution factors were determined, 18 of the rate constants were predicted and compared with the experimental rate constants as shown in Table SI1-2 of Supporting Information. A Total of 44% of the rate constants was within the EG. The SD was 2.2. Group Rate Constants for HO• Interaction. The group rate constants and group contribution factors were summarized in Table SI1-8 of Supporting Information 1. Because no direct interactions were experimentally observed (51, 52), the group rate constants, k-N-NO and k-N-NO2 were not considered. The group rate constants k-CN and k-NH2 can be compared with the rate constants for compounds that react with HO• via only interaction such as cyanogen and thiourea, respectively. The rate constant for thiourea is approximately twice k-NH2 because the electron positive -CS- functional group does not significantly affect the electron density of the N atom. The reaction of HO• with urea may be different because two amine functional groups in urea bond to the electron-negative functional group, -CO-. As a result, another group rate constant k-N-CO-N- was considered for methylurea, tetramethyl urea, and 1,3-dimethylurea. The magnitude of most group rate constants for the S-containing compounds was in the same order as for the amine-related compounds but approximately 1 order of magnitude larger than those for the amide-related compounds. This is probably because of the electronegative -CO- functional group that is a part of amide functional groups. S-, N-, or P-Atom-Containing Group Contribution Factors. The S-, N-, or P-atom-containing group contribution factors apparently play the same role as the functional groups for R H-atom abstraction, i.e., XRi ) e-(Ea,iabs)/RT. However, it is anticipated that S-, N-, or P-atom-containing functional groups may have different effects on H-atom abstraction. The group contribution factors for -S, -S-S-, and -SH, and -NH2, -NH-, and -N