Subscriber access provided by ONDOKUZ MAYIS UNIVERSITESI
Article
Predicting Gaseous Reaction Rates of Short Chain Chlorinated Paraffins with #OH: Overcoming the Difficulty in Experimental Determination Chao Li, Hong-bin Xie, Jingwen Chen, Xianhai Yang, Yifei Zhang, and Xianliang Qiao Environ. Sci. Technol., Just Accepted Manuscript • DOI: 10.1021/es504339r • Publication Date (Web): 05 Nov 2014 Downloaded from http://pubs.acs.org on November 10, 2014
Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.
Environmental Science & Technology is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.
Page 1 of 26
Environmental Science & Technology
1
Predicting Gaseous Reaction Rates of Short Chain Chlorinated Paraffins with ⋅OH:
2
Overcoming the Difficulty in Experimental Determination
3
Chao Li†, Hong-Bin Xie†, Jingwen Chen†∗, Xianhai Yang†, Yifei Zhang‡, Xianliang Qiao†
4
†
5
Environmental Science and Technology, Dalian University of Technology, Dalian 116024,
6
China.
7
‡ State Key Laboratory of Fine Chemicals, Dalian University of Technology, Dalian 116024,
8
China.
Key Laboratory of Industrial Ecology and Environmental Engineering (MOE), School of
9
Table of Contents
10 11
ABSTRACT
12
Short chain chlorinated paraffins (SCCPs) are under evaluation for inclusion in the Stockholm
13
Convention on persistent organic pollutants. However, information on their reaction rate
14
constants with gaseous ⋅OH (kOH) is unavailable, limiting the evaluation of their persistence in
15
the atmosphere. Experimental determination of kOH is confined by the unavailability of
16
authentic chemical standards for some SCCP congeners. In this study, we evaluated and
17
selected density functional theory (DFT) methods to predict kOH of SCCPs, by comparing the
18
experimental kOH values of 6 polychlorinated alkanes (PCAs) with those calculated by the
19
different theoretical methods. We found that the M06-2X/6-311+G(3df,2pd)//B3LYP/6-311 ∗
Corresponding author phone/fax: +86-411-84706269; e-mail:
[email protected]. 1
ACS Paragon Plus Environment
Environmental Science & Technology
Page 2 of 26
20
+G(d,p) method is time-effective and can be used to predict kOH of PCAs. Moreover, based on
21
the calculated kOH of 9 SCCPs and available experimental kOH values of 22 PCAs with low
22
carbon chain, a quantitative structure-activity relationship (QSAR) model was developed. The
23
molecular structural characteristics determining the ⋅OH reaction rate were discussed. logkOH
24
was found to negatively correlate with the percentage of chlorine substitutions (Cl%). The
25
DFT calculation method and the QSAR model are important alternatives to the conventional
26
experimental determination of kOH for SCCPs, and are prospective in predicting their
27
persistence in the atmosphere.
28 29
INTRODUCTION
30
Chlorinated paraffins (CPs) are commercially produced polychlorinated alkanes (PCAs)
31
having carbon chain lengths from C10 to C30 and chlorine content from 30% to 70% by
32
mass.1-3 Among them, short chain CPs (SCCPs) with carbon chain lengths ranging from C10 to
33
C13 have received much attention due to their widespread occurrence and persistence in the
34
environment,4-15 their potential
35
bioaccumulate8 and high toxicity.16,17 SCCPs show long-term toxicity to algae, aquatic
36
invertebrates and fish at concentrations as low as 19.6, 8.9 and 3.1 µg/L, respectively.18 There
37
is also sufficient evidence for carcinogenicity to humans from SCCPs of C12 carbon chain
38
length group with an average chlorine content of 60%.17 The vapor pressures of SCCPs are in
39
the range of other chlorinated organics with the same molecular weight range such as
40
polychlorinated biphenyls and toxaphenes, suggesting that they have the potential to undergo
41
long-range atmospheric transport through the grasshopper effect.4 SCCPs have been detected
42
in environmental media including air,10,11,13 water19 and sediment;9,12 in humans20 and
43
organisms like fish,14 birds,15 marine mammals;21 and even in remote areas like Canadian
44
Arctic lakes, the Hazen lake and the European Arctic.5,8,15,21 As SCCPs represent a potential
for
long-distance
transport,12,15
2
ACS Paragon Plus Environment
their
tendency
to
Page 3 of 26
Environmental Science & Technology
45
“new” category of persistent organic pollutants (POPs), several countries (e.g., the United
46
States and Canada) and the European Union have imposed regulations on the use of SCCPs.22
47
Additionally, SCCPs are now under evaluation for inclusion in the Stockholm Convention on
48
POPs.
49
Despite the efforts to characterize the extent of SCCP contamination in the environment,
50
little is known about their atmospheric degradation kinetics that governs the fate and
51
persistence of SCCPs. Among the various atmospheric reactions, the reactions initiated by
52
hydroxyl radicals (⋅OH) produced by atmospheric photochemical reactions, are of paramount
53
importance. Thus, the ⋅OH reaction rate constant (kOH) is a key parameter for the fate
54
assessment of pollutants in the troposphere. However, there are no kOH values available for
55
any SCCPs.
56
To date, several direct and indirect experimental methods have been developed to measure
57
kOH values.23-26 However, the experimental determination of kOH for SCCPs is difficult as the
58
quantities and kinds of authentic chemical standards for some SCCPs, which are necessary for
59
the determination, are very limited. Moreover, there are a large number of isomers for
60
SCCPs,22 which could further increase the workloads and difficulties encountered in the
61
experimental determination. Another solution is to predict kOH based on quantitative
62
structure-activity relationship (QSAR) models, on which the Organization for Economic
63
Co-operation and Development (OECD) has issued guidelines on model development and
64
validation.27 Nevertheless, the utility of QSAR models is constrained by the applicability
65
domain of the models.28,29 The currently available QSAR models on kOH do not cover SCCPs.
66
For example, two latest QSAR models on kOH prediction cover only a limited number of
67
PCAs with chain lengths ranging from C1 to C6.30,31 Thus, seeking a new solution for
68
predicting kOH of SCCPs is necessary.
69
In recent years, the increasing computational power has enabled quantum chemical 3
ACS Paragon Plus Environment
Environmental Science & Technology
70
calculations for molecules or systems comprising up to a few hundred atoms.32 Some studies
71
in recent years showed that the density functional theory (DFT) can accurately predict kinetics,
72
pathways and products for the reactions of chemicals with ⋅OH.33--35 For example, Zhou et al.
73
employed DFT to calculate the reaction of 4,4'-dibromodiphenyl ether with ⋅OH, and found
74
that the predicted products and the overall rate constants (kOH) at 298 K agreed well with the
75
experimental results.33 Ci et al. employed DFT to compute the hydrogen abstraction reactions
76
of CF3CH2CHO with ⋅OH, also found that the calculated kOH values were consistent with the
77
experimental data.35 Furthermore, the theoretical approaches have the advantage that they can
78
provide details involving reaction sites and favorable reaction pathways.33-36 As far as we
79
know, there are no previous quantum chemical calculations on reactions of SCCPs with ⋅OH.
80
It was the purpose of this study to find a solution for kOH prediction of SCCPs. By
81
employing 6 molecules (CH3CCl3, CH2ClCHCl2, CH2ClCHClCH3, CH3(CHCl)2CH3,
82
CH2Cl(CH2)3CH3 and CH2Cl(CH2)4CH3) with available experimental kOH values as test cases,
83
we compared the prediction accuracy of several DFT methods. The expensive ab initio
84
second-order Moeller-Plesset (MP2)37 method, which consistently describes all types of
85
correlation energy including intermolecular correlation and intramolecular correlation terms,
86
and therefore is reliable for geometry optimization,38-40 was also adopted as a reference. Based
87
on the established method, kOH values for 9 SCCPs with carbon chain ranging from C10 to C13
88
(Figure 1), were calculated. A QSAR model for predicting kOH of more SCCPs was developed
89
based on the DFT calculated and available experimental kOH values for PCAs. To the best of
90
our knowledge, this is the first study concerning kOH of SCCPs. The DFT calculation method
91
and the QSAR model can be employed to predict kOH of other SCCPs.
4
ACS Paragon Plus Environment
Page 4 of 26
Page 5 of 26
Environmental Science & Technology
92 93
Figure 1. Atom numbers of nine SCCPs (Dark gray: C; Light gray: H; Green: Cl).
94 95
COMPUTATIONAL DETAILS
96
Global Minimum Search. As the molecular structures of PCAs consist of sp3-hybridized
97
C-X (X = H, Cl) bonds, they have a range of conformations. Thus, the global minimum search
98
is needed to identify the optimal structures of the reactants under study, and the optimal
99
structures were then included in further structure calculations followed by dynamic
100
calculations of kOH.
101
The calculations were carried out using TURBOMOLE.41 Ab initio molecular dynamics
102
(AIMD) was used to generate reasonable gas-phase geometries by employing the BLYP
103
functional42 along with a triple-ξ valence polarized basis set (TZVP) within the 5
ACS Paragon Plus Environment
Environmental Science & Technology
104
resolution-of-the-identity (RI) approximation.41 In order to sample enough conformations in a
105
short time and increase to some extent the possibility for getting the real global minimum, the
106
temperature was set to 700 K with the premise that the molecules would not dissociate in the
107
AIMD calculation (detailed in the SI). A time step of 0.96752 fs was used, and the total length
108
of each AIMD simulation was 20000 time steps. We selected the conformations from the
109
AIMD run as the starting point for geometry optimization at the B3LYP/6-311+G(d,p) level
110
using the Gaussian 09 program suite.43 The minimum with the lowest energy was identified as
111
the global minimum used to investigate its reaction with ⋅OH.
112
Electronic Structure Calculations. A gradual scheme was used to screen reliable
113
theoretical methods for the geometry optimization and frequency analysis of the studied
114
systems. For all the systems, the single point energy calculation was performed at the
115
M06-2X/6-311+G(3df,2pd) level, since the M06-2X functional gives the best performance
116
(mean signed errors of -0.51 kcal/mol) without increasing the computational time for the
117
hydrogen-transfer barrier height calculation.44 As CH3CCl3 has a C3v symmetry that can
118
reduce the computational time, its reaction with ⋅OH was taken as an initial test to screen the
119
different computational methods (the functional BHandHLYP,42,45 MPW1K,46 B3LYP,42,47
120
TPSSh,48 B97-149 and ab initio MP237 in conjunction with the 6-311+G(d,p) basis set) by
121
comparing the calculated zero-point corrected reaction barriers (Ea). Subsequently, the
122
selected time-effective and reliable methods (B3LYP, TPSSh and B97-1 functionals) were
123
evaluated with CH2ClCHCl2 and CH2ClCHClCH3 as test cases, followed by using
124
CH3(CHCl)2CH3, CH2Cl(CH2)3CH3 and CH2Cl(CH2)4CH3 for further method validation.
125
Finally, an optimum method (B3LYP/6-311+G(d,p)) was selected for the SCCPs. All the
126
stationary points including reactants, products, pre-complexes and product-like complexes
127
have real frequencies, whereas the transition states (TSs) have only one imaginary frequency.
128
The minimum energy paths (MEP) were obtained by the intrinsic reaction coordinate (IRC) 6
ACS Paragon Plus Environment
Page 6 of 26
Page 7 of 26
129
Environmental Science & Technology
theory for each reaction channel.50
130
Dual-level Direct Dynamic Calculation. By means of the POLYRATE 2010-A
131
program,51 the rate constant for each reaction channel (ki) was calculated using the improved
132
canonical variational transition-state theory (ICVT) with small-curvature tunneling (SCT)
133
correction51 over the temperature range 200-500 K. The frequencies of the selected points
134
along the MEP were calculated at the B3LYP/6-311+G(d,p) level to obtain the Cartesian
135
coordinates, gradient and Hessian matrix; and the single-point energy of these points was
136
calculated by the M06-2X/6-311+G(3df,2pd) method. kOH was calculated by summing up the
137
ki values of different channels.
138
QSAR Modeling. The DFT calculated kOH values (T = 298 K) of 9 SCCPs, together with
139
those of 22 low carbon chain PCAs (C1~C6) with available experimental kOH values,31 were
140
employed for QSAR modeling. According to previous studies,30,31,52-55 the following quantum
141
chemical descriptors were selected to construct the model, including the atomic Mulliken
142
charges [the most negative net C-atom charge (QC-), the most positive net H-atom charge and
143
the most negative net Cl-atom charge], polarizability (α), the energy of the highest occupied
144
molecular orbital (EHOMO), the energy of the lowest unoccupied molecular orbital (ELUMO),
145
chemical potential (µ = (ELUMO + EHOMO)/2)56, and absolute hardness (η = (ELUMO -
146
EHOMO)/2)56. These descriptors were calculated from the geometries optimized by the
147
B3LYP/6-311+G(d,p) method using the Gaussian 09 program suite. Additionally, two
148
descriptors that describe the chlorine substitution characteristics [the number of chlorine
149
atoms (nCl) and percentage of chlorine substitutions (Cl%)], 3D-MoRSE (3D-molecular
150
representation of structure based on electron diffraction) and geometrical descriptors were
151
also considered, for which the DRAGON software was employed for the calculation.57
152
For the modeling, the data were randomly split into a training set and a validation set
153
with a ratio of 4:1. Stepwise multiple linear regression (MLR) analysis was employed to 7
ACS Paragon Plus Environment
Environmental Science & Technology
154
construct the QSAR model. The applicability domain of the QSAR model was assessed by the
155
Euclidean distance based method, and was conducted using the AMBIT Discovery (version
156
0.04).58
157 158
RESULTS AND DISCUSSION
159
Evaluation of the Different Calculation Methods. The reaction of PCAs with ·OH
160
proceeds via H-atom abstraction, as the abstraction of Cl atoms has much higher Ea and is
161
endothermic.59 For the reaction CH3CCl3+·OH, we considered only one H-abstraction process
162
from -CH3, as the three H atoms in -CH3 are equivalent (Figure S1). The Ea values calculated
163
from the different methods, together with the experimental values are listed in Table S1. We
164
found that the Ea values from the geometries optimized by the TPSSh, B3LYP and B97-1
165
functionals are close to those calculated by the more-reliable MP2/6-311+G(d,p) method, and
166
are closer to the experimental values25,60-63 than the other functionals under investigation.
167
Thus, the TPSSh, B3LYP and B97-1 functionals were employed for geometry optimization
168
and frequency analysis of the other two tested PCAs.
169
For CH2ClCHCl2 and CH2ClCHClCH3, the possible reaction pathways for the ·OH
170
reaction are described in Figure S1. The potential energy surface profiles for the 3 tested
171
PCAs above are shown in Figure S2. As H-abstraction by ·OH can occur from -CHCl2 and
172
CH2Cl-, CH2ClCHCl2 has three TSs (TS1b, TS2b and TS3b). Similarly, CH2ClCHClCH3 has six
173
TSs (TS1c TS2c, TS3c, TS4c, TS5c and TS6c).
174
The calculated and experimental25,26, 60-67 kOH values for the 3 PCAs are shown in Figures
175
2-4. Several experimental studies reported the kOH values of CH3CCl3,25,60-66 which makes it
176
possible to estimate the deviation of the experimental determinations. At 298 K and at 95%
177
confidence level, the deviation of experimental kOH values of CH3CCl3 is within a factor of
178
0.4-2.8 (detailed in the SI). As for CH2ClCHCl2, only two studies reported the kOH values,25,26 8
ACS Paragon Plus Environment
Page 8 of 26
Page 9 of 26
179
Environmental Science & Technology
and the biggest ratio between them is 1.7 (T = 295 K).
180 181
Figure 2. Calculated and experimental kOH values for CH3CCl3 over different temperature
182
ranges (a), and the enlarged view of (a) at T = 200-330 K is shown in (b).
183 184
Figure 3. Calculated and experimental kOH values for CH2ClCHCl2 over different temperature
185
ranges (a), and the enlarged view of (a) at T = 200-330 K is shown in (b).
186 187
Figure 4. Calculated and experimental kOH values for CH2ClCHClCH3 over different
188
temperature ranges (a), and the enlarged view of (a) at T = 200-320 K is shown in (b).
189 190
By comparing the calculated with the available experimental values (Figures 2-4), it can 9
ACS Paragon Plus Environment
Environmental Science & Technology
Page 10 of 26
191
be concluded that the TPSSh functional overestimates the kOH values of CH3CCl3 in the range
192
T = 200-720 K except for the experimental kOH value of Chang et al.,61 and underestimates kOH
193
values of CH2ClCHClCH3 in the range T = 200-372 K. Thus, the TPSSh functional was
194
excluded for its inability for accurately predicting kOH of the PCAs.
195
There are no significant differences between the kOH values calculated from the B97-1 and
196
B3LYP functionals in the low temperature range (200-330 K). However, in the range of T >
197
340 K, the kOH values calculated from the B3LYP functional are generally more consistent
198
with the experimental values than those calculated with the B97-1 functional. The kOH data in
199
the high temperature range are important in understanding the incineration process of
200
PCAs.25,65 Thus, the B3LYP method is more universal than the B97-1 functional for
201
predicting kOH in both the low and high temperature ranges. Furthermore, Wang et al.
202
calculated the kOH of CH3CH2CH2Cl and CH3CHClCH3, and also concluded that the
203
B3LYP/6-311G(d,p) method is reliable for geometries optimization.68 As the studied system
204
involves radicals and Cl atoms, inclusion of a diffuse function to the basis set should increase
205
the
206
M06-2X/6-311+G(3df,2pd)//B3LYP/6-311+G(d,p) method for predicting kOH values of
207
SCCPs.
accuracy
of
kOH
prediction.
Thus,
we
selected
the
208
The selected method was further verified with the 3 PCAs with C4-C6 carbon chain
209
lengths, and the results are listed in Table S2-S5. The deviations between the calculated kOH
210
values and the corresponding experimental values are within a factor of 0.6-1.7, 0.4-0.7,
211
0.6-0.7,
212
CH3(CHCl)2CH3, CH2Cl(CH2)3CH3 and CH2Cl(CH2)4CH3, respectively, which are generally
213
acceptable.69,70 Thus, the deviation of the DFT calculated kOH values of the PCAs is estimated
214
to be within a factor of 0.4-1.7.
215
0.5,
0.4-0.6
and
0.5-0.6
for
CH3CCl3,
CH2ClCHCl2,
CH2ClCHClCH3,
Reactions of the SCCPs with ⋅OH. For convenience, all the atoms of the 9 SCCPs are 10
ACS Paragon Plus Environment
Page 11 of 26
Environmental Science & Technology
216
numbered (Figure 1). Due to the C1 symmetry of the SCCPs, all the channels of H-abstraction
217
for the SCCPs were considered. Similar to the 3 PCAs with low carbon chains, all the
218
H-abstraction channels proceed via pre-complexes, TSs, product-like complexes and products
219
(Figure S3). The calculated activation free energy (∆G) and Ea for each channel are listed in
220
Tables S6-S14.
221
By comparing the Ea values of the different H-abstraction channels, we found the channels
222
with six-membered ring TSs formed by H-bonds either between H and O or between H and Cl
223
atoms, tend to have low Ea and ∆G values. Generally, the conformations with cis-H-C-C-Cl
224
sub-structural units facilitate the formation of the six-membered ring TSs. Thus, chlorine
225
substitutional positions and steric structures of SCCPs have great effects on the reaction rates
226
with ⋅OH. Beyond that, there seem no regularities for the H-abstraction channels. Thus, it is
227
necessary to take every H atom into account when calculating the kOH values.
228
The DFT calculated kOH values for the 9 SCCPs over the temperature range 200-500 K are
229
shown in Figure 5 and Table S15. In general, the kOH values display a negative temperature
230
dependence in the low temperature range. The reason can be ascribed to the negative Ea of the
231
reactions (Tables S6-S14). In principle, for reaction pathways with negative Ea, the reaction
232
rates at low temperatures are mainly controlled by the kinetics, and at high temperatures are
233
mainly affected by the thermodynamics.36 The similar temperature dependence of kOH was
234
also reported in the reaction of other compounds with ⋅OH.33,35,36,71
235 236
Figure 5. Plot of the DFT calculated kOH values for the 9 SCCPs versus temperature 11
ACS Paragon Plus Environment
Environmental Science & Technology
237 238
Driven by curiosity, we compared the DFT calculated kOH values at 298 K with those
239
predicted by QSAR models, although the training sets of the QSAR models did not cover
240
SCCPs. The QSAR predicted kOH values are all lower than the DFT calculated values. The
241
deviations between the DFT calculated kOH values and those predicted by the QSAR models
242
of Li et al.,31 Wang et al.,30 Gramatic et al.53 and Roy et al.,55 are within a factor of 3.0-6.8,
243
2.8-22.5, 1.2-15.9 and 40.5-397.4, respectively. Thus, it is necessary to develop a new QSAR
244
model that covers SCCPs in the applicability domain.
245 246 247
QSAR Modeling. The following optimum QSAR model for logkOH prediction (T = 298 K) was obtained: logkOH = 5.057SPH + 1.167QC- -13.991η - 12.186
248
The compounds and the descriptor values involved in the QSAR model are listed in Table S16.
249
For this model, there are 25 PCAs (7 SCCPs) included in the training set. The variable inflation
250
factors (VIF) for the predictor variables are all