4810
J. Phys. Chem. 1996, 100, 4810-4814
Structure, Energy, Vibrational Frequencies, and Potential Energy Curve of 2,3,7,8-Tetrachlorinated Dibenzo-p-dioxin: Ab Initio MO Studies Toshihiro Fujii,* Kenji Tanaka, Hiroaki Tokiwa,† and Yuko Soma National Institute for EnVironmental Studies, 16-2 Onogawa, Tsukuba, Ibaraki 305, Japan ReceiVed: September 15, 1995; In Final Form: January 3, 1996X
Optimized geometries and total energies for the 2,3,7,8-tetrachlorinated dibenzo-p-dioxin (TCDD) molecules were calculated using the Gaussian 92 system of programs at two ab initio MO levels: RHF/6-31G and RHF/6-31G*. In addition to the above levels, the corresponding basis sets, including diffuse function, were also used. The fully-optimized geometry is given for the TCDD molecule. The results were compared and agree closely with the observed structures, which were obtained by X-ray spectrometry. Harmonic vibrational frequencies are calculated at the above-mentioned levels of theory on the basis of the optimized geometries. Comparison with experimental IR results is discussed. It is predicted from vibrational analyses that the butterfly flapping motion of two benzo-planes is a motion having a very low fundamental frequency (3.6 cm-1). For those modes which are IR inactive, the predicted Raman frequencies are reported. Some ambiguity concerning the dynamic behavior has also been investigated in terms of the potential energy curve (PEC) as a function of the folding angle (θ), to specify the butterfly flapping motion of the two benzo-planes. It is shown that the potential has little curvature near the minimum. The small curvature of the PEC appears to be consistent with the far-infrared 3.6 cm-1 frequency of the flapping motion.
Introduction Highly toxic polychlorinated dibenzo-p-dioxins (PCDDs) have been detected as trace contaminants in incinerator fly ash,1 herbicide formulations, and the synthesis of several important commercial chlorophenolic compounds.2 As a result, a number of poisoning fears have occurred in recent years due to release of these substances into the environment.3 Among 75 possible PCDDs having from one to eight chlorine atoms, 2,3,7,8-tetrachlorinated dibenzo-p-dioxin (2,3,7,8-TCDD) shows the highest degrees of toxicity.4 This has triggered a large number of scientific studies5 aimed at elucidating the chemical and biological properties of TCDD-related compounds. Research into the behavior of PCDDs requires a data base of molecular properties. Unfortunately, due to the great number of chlorinated dioxins and the high toxicity of certain isomers, experimental information about these molecular parameters is hard to obtain. In the absence of experimental data, molecular orbital calculations may be used to provide vital information about these molecules. Earlier theoretical studies on PCDD were unreliable, and some of the results may actually be misleading, due to the technical limitations of the day. MMP2, a computational force field method, was employed to derive molecular geometries of some PCDDs. A preliminary report of structural features and a correlation chart of dihedral angles against IR data are presented.6 A semiempirical investigation of PCDDs, which was associated with the study of the photodegradation mechanism, has been reported by Makino et al.7 But a geometrical consideration has not been made. Koester and Hites8 used the MNDO method to calculate physical properties of many PCDDs such as the heats of formation, ionization potentials, dipole moments, and LUMOs. Also the geometries and energies have not been reported. * To whom correspondence should be directed. † Guest Research Scientist from Rikkyo University, Toshimaku, Tokyo 171, Japan. X Abstract published in AdVance ACS Abstracts, March 1, 1996.
0022-3654/96/20100-4810$12.00/0
An ab initio study of the electronic structures of a series of 25 PCDDs has been made by Cheney and Tolly.9 However, they did not report the theoretical geometry. Recently, Schaefer and Sebastian10 reported ab initio calculations but only for the nonchlorinated dibenzo-p-dioxin (DD). Furthermore, they reported only the 4-31G MO geometry and did not report energies at all. Hence, a systematic reexamination of this system is opportune. The ab initio MO approach to larger molecules is made possible by recent advances in supercomputer technology and methodology. The mechanism of the toxicity attributed to the PCDDs, presumably, depends on structurally related properties.11 From such a point of view, molecular geometry and conformational and stereoelectronic features such as planarity, dipole moment, or polarizability are essential subjects to investigate. The potential energy for the ring folding angle, θ, which gives the dihedral angle between the benzo planes, is also an interesting subject to study.12 The geometry of TCDD is based upon X-ray crystallographic data that exhibit near-planar structures in the crystal lattice. Few structural investigations have been conducted in solution or in the gas phase. The valance force equation approximations of C-O-C bond angles from infrared data are also near-planar conformations12 for TCDD geometries in the gas phase. Bond angles for nonchlorinated-dibenzodioxin (DD) and related compounds, determined in solution by a variety of techniques, have generated structures ranging from folded to near-planar. The dipole moment measurement of DD13 and its 13C nuclear magnetic resonance14 suggest a near-planar folded geometry. Dielectric absorption measurements on the DD molecule in a polystyrene matrix15 and in benzene solution16 are consistent with two shallow potential minima (θ ) 180°) or with one minimum characterized by a very flat potential. The birefringence studies17 imply a planar molecule in the solution. Both extended Huckel and ab initio calculations find the planar molecule to be most stable (θ ) 180°). This is in contradiction to experiments. Such an ambiguity concerning its conformational behavior must be solved. © 1996 American Chemical Society
2,3,7,8-Tetrachlorinated Dibenzo-p-dioxin We have performed ab initio calculations for 2,3,7,8-TCDD at a relatively high level of theory because only semiempirical studies were previously done and no complete ab initio MO results were reported. First, the equilibrium geometries were fully optimized on the RHF/6-31G* potential energy surface with diffuse functions partially included (+O,Cl). This permits direct geometrical and energetic comparison between experimental and theoretical results. On the basis of the optimized geometries, we have calculated harmonic vibrational frequencies of IR- and Raman-active modes at the RHF/6-31G* level. The latter are compared with Patterson’s IR experimental results. Finally, geometry-optimized ab initio MO computations at various folding angles (dihedral angles between the benzene planes) are reported, for which some equilibrium geometry controversies exist.
J. Phys. Chem., Vol. 100, No. 12, 1996 4811
Figure 1. Structure of the TCDD molecule and its notation with bonding scheme. Refer to Table 1, for definition of structural parameters for the TCDD molecule.
TABLE 1: Optimized Geometriesa of the TCDD Molecule from Various ab Initio RHF Levels
Computational Methods All computations were performed with an NEC SX-3 supercomputer system and the Gaussian 9218 program. The computational resources used in this study include Kubota, Inc., a Vistra workstation, and a CONVEX C-3J superminicomputer. Full geometry optimization procedures allowed for variations in bond angles and lengths. The constraints were, however, that the molecule have D2h symmetry for θ ) 180° and C2 symmetry at other values of θ. Each benzene moiety was not required to be planar. Computations were made at the following levels: RHF/6-31G and RHF/6-31G*. Calculations including diffuse functions were also performed. This was proposed by Clark et al.19 for use with the 6-31G* basis sets. This addition is denoted as 6-31G*(+O,Cl). It includes a single additional diffuse sp shell on oxygen and chlorine atoms only. The total energies (Et) as well as the first ionization potential (IP) are obtained from a complete geometry optimization. Estimates for the effect of electron correlation, which require an AO/MO integral transformation process, were not obtained because of the file (read/write) size limitation of the Gaussian 92 program package. Harmonic frequency analysis using analytical second derivatives was carried out at the RHF/6-31G*(+O,Cl) level at a stationary point to confirm the equilibrium structure and to provide zero-point vibrational energy (ZPE) corrections. The harmonic frequency analysis also provided, with calculated intensities, a prediction of the IR and Raman spectra. A calculation was done for the dependence of Et of the TCDD molecule on the angle θ (potential energy curve for the ring puckering). The calculations at the above-mentioned levels were done for every 1° increment of θ between 150 and 180. Again, no constraint was given with the exception of the C2 symmetry at values of θ. Results and Discussion 1. Geometry. The geometry and its notation for the TCDD structure are shown in Figure 1. The optimized structural data for dioxin at different basis set levels of ab initio MO calculations are presented in Table 1. Listed in Table 2 are the total energies (Et) and ionization potentials (IPs). The geometrical parameters of TCDD, determined experimentally from an X-ray diffraction (crystallography) study by Boer et al.,20 are also listed in Table 1. Examination of Table 1 indicates that there is generally good agreement between the calculated and experimentally determined bond distances. The symmetry-related C-C bonds denoted as b, g, l, and q have a computed length of 1.375-1.376 Å, which is about 0.01 Å less than all other C-C bond distances of 1.380-1.386 Å. It also indicates some degree of conjugation (localized), which is
a (f, k, p) b (g, l, q) c (h, m, r) d (i, n, s) e (j, o, t) u (w) v (x) ab () fg, kl, pq) bc () gh, lm, qr) cd () hi, mn, rs) be () gj, lo, qt) au () fw, kw, pu) bu () gw, jw, qu) cv () hx, mx, rv) ak () fp) aub (fwg, kwl, puq) bec (gjh, lom, qtr) ecd (ihj, omn, srt) cdv (rsv, mnx, hix) dvs (nxi)
6-31G
6-31G*
6-31G*(+O,Cl)
exptlb
1.385 1.375 1.384 1.791 1.069 1.383 1.380 118.7 119.4 117.9 119.5 120.7 120.4 120.2 118.2 180.0 180.0 180.0 180.0 180.0
1.360 1.375 1.386 1.731 1.072 1.384 1.385 118.4 120.2 118.2 119.2 121.6 120.0 119.8 116.8
1.361 1.376 1.386 1.732 1.072 1.384 1.386 118.4 120.2 118.2 119.3 121.6 120.1 119.8 116.9
1.379 1.377 1.384 1.728 1.01 1.387 1.385 117.6 121 118.9 119 122.2 115.7
a Refer to Figure 1 for definition of parameters; e.g., a ) bond length (A), ab ) bond angle (deg), aub ) dihedral angle (deg). b From ref 20, averaged bond distances and bond angles from the two independent TCDD molecules.
TABLE 2: Total Energiesa and IP of TCDD Calculated at Various ab Initio Levels E0 RHF/6-31G -2444.034 64 RHF/6-31G* -2444.371 45 RHF/6-31G*(+O,Cl) -2444.376 88
E(ZPE)b 88.6 87.7 87.6
Et
IP (eV)
-2443.893 45 7.9 (9.1)c -2444.231 69 7.3 (8.5) -2444.237 28 7.4 (8.6)
a Total energy, Et ) E0 + E(ZPE), in atomic units. b Unscaled zeropoint energy, in kcal/mol. c Koopmans’ IP in parentheses.
also found in the 4-31G geometry of nonchlorinated dibenzop-dioxin.10 This result seems reliable, compared with experimental data. The experimentally determined C-C bond length in b, g, l, and q is 1.377 Å, while the other C-C bond distances, such as c, h, m, and r, are 1.384 Å. To assess the importance of basis set size on geometry and to highlight trends in the calculations, the results from different levels of calculations may be compared, since experimental data were available. Throughout the calculations, we have noticed the effect of ab initio levels on the bond length. Compared with the experimental result, the RHF/6-31G calculations severely overestimate the C-Cl bond (d, i, n, s) length by as much as 0.065 Å. Calculations with the larger 6-31G* basis set do much to improve this discrepancy. It is noticeable that some of the RHF/6-31G structural data are fortuitously in close agreement with the observed structure. The optimized C-C bond length, noted as c, h, m, and r, is
4812 J. Phys. Chem., Vol. 100, No. 12, 1996
Fujii et al.
TABLE 3: Analysis of Harmonic Vibrational Frequencies (in Units of cm-1) of TCDD experimentala assignment COC (sym stretch) ring breathing CCC (ring bending) CH (in-plane deform) COC (asym stretch) skeletal vibrations
RHF/6-31G*(+O,Cl) freq, int 895, vw 916, vw 932, vw 1103, m 1114, m 1173, w 1248, vw 1255, vw 1306, m 1321, m (1391, w) 1465, vs 1489, s 1566, m
freqb (sym,c IR intd)
assignment
875 (B2u, 210) 915 (B3u, 103) 923 (B1u, 28) 1126 (B2u, 247)
CH in plane + COC str CH out of plane skeletal vib + CCl str + COC str CH in plane + skeletal vib
1173 (B1u, 92) 1243 (B1u, 27)
COC str + CCC bend CH in plane
1332 (B2u, 444)
COC str
1391 (B1u, 51) 1499 (B2u, 1377)
skeletal vib (ring str) skeletal vib
1587 (B2u, 138)
skeletal vib
a
From the IR spectra; ref 12. Abbreviations used in the reference are as follows: s, strong; m, medium; w, weak; v, very; sym, symmetric; asym, asymmetric. b Frequencies (cm-1) are scaled by 0.9. c Symmetry of the vibrational mode (refer to Figure 2). d Infrared intensities are calculated by the atomic polar method as implemented in GAUSSIAN 92, in km/mol.
1.384 Å (RHF/6-31G) vs 1.384 Å (observed). However, it is lengthened to 1.386 Å (RHF/6-31G*) when a higher level method with polarized basis sets is used. Worse results can be obtained for the C-O bond length when extended basis sets incorporating d orbital polarization are included in the calculation. The calculated C-O bond distances at the 6-31G*(+O,Cl) level are 0.017 Å shorter than the experimental value. This is a frequently encountered problem with ab initio SCF calculations regardless of the basis set used.21 Perhaps, significant improvement is realized with the post-HF method. The structural parameters calculated at the RHF/6-31G* and RHF/6-31G*(+O,Cl) levels are practically the same, although, as expected, there is a little lengthening in all the bond distances at the 6-31G*(+O,Cl) level of theory. The difference for all the bond lengths is less than (0.001 Å. This trend holds for all the bond angles. The major discrepancy between the calculated (6-31G) and experimentally determined bond angles is observed in the COC angle (ak, fp). The COC bond angle is computed as 118.2° while the experimental value is 115.7°. As expected, this difference can be corrected when higher levels are used with polarized basis set and d functions. The other calculated bond angles at the 6-31G level are also considerably different from the measured values with no specific pattern. However, some improvement is realized with the 6-31G* or 6-31G*(+O,Cl) calculations. All the dihedral angles are calculated to be 180.0; suggesting no deformations of benzo planes. The agreement between calculated and measured values is not considered because of a lack of precise experimental data. 2. Energy. The total energies and IPs at various levels are listed in Table 2. Among the examined levels the lowest total energy of -2444.376 88 au for D2h TCDD is obtained at the RHF/6-31G*(+O,Cl) level. The difference between the RHF/ 6-31G* and RHF/6-31G*(+O,Cl) levels is 0.005 57 au (3.5 kcal/ mol). According to Koopmans’ theorem, the negative of the oneelectron energies of filled molecular orbitals conventionally gives vertical IPs. A value for the first IP of the molecule is derived from the energy of the highest occupied molecular orbital (HOMO). The first IP is also computed as the energy difference between the cation and the neutral molecule since it is defined as the energy required to remove an electron. Here we denote this as the ∆SCF’s IP. The results of both IPs (listed in the last column of Table 2) show the fluctuation of the ionization energy among the calculation levels when the calculation is changed from RHF/6-31G to RHF/6-31G*-
(+O,Cl). The difference between the ∆SCF’s and Koopmans’ IPs remains unchanged at 1.2 eV. In the case of RHF/6-31G*(+O,Cl), comparing the total energies of the positive ion and the neutral (∆SCF’s IP) leads to IP ) 7.4 eV, whereas the orbital energy of the neutral (Koopmans’ theorem) leads to IP ) 8.6 eV. This is typical. The ∆SCF procedure leads to too small an IP. Koopmans’ theorem assumes that the orbitals do not relax upon ionization, leading to too large a value. The average value of 8.0 eV may be a good prediction22 as a TCDD IP. No comparison with experimental results can be made because of nonavailability. Calculations indicate that the HOMO is due to ejection of an electron from an antibonding π-orbital mainly centered on the C-O bond.9 Approximately 1.08 eV above the first HOMO (5B3u), the second HOMO (4B2g) is a delocalized π-orbital around the O-C-C bond, whose C-C bond is not in the oxine ring but in the benzene ring. 3. Vibrational Frequency Analysis. Harmonic vibrational frequencies are calculated for TCDD at the RHF/6-31G, RHF/ 6-31G*, and RHF/6-31G*(+O,Cl) levels of theory on the basis of optimized geometries. However, they are barely affected by the calculation levels, although there is some improvement. For example, the frequency for the skeletal vibration mode at the RHF/6-31G*(+O,Cl) level is smaller, in closer agreement with the experimental value of 1566 cm-1. Only the results at the RHF/6-31G*(+O,Cl) level, together with computed IR intensities at more than 25 km/mol and their detailed assignment, are listed in Table 3. The animated vibrational modes listed in Table 3 are shown in Figure 2, because the description, based on animation of the normal modes for the equilibrium geometries, is approximate. The most intense bands are due to the skeletal vibrations at 1499 cm-1. Grainger et al.23 have reported the infrared spectrum (4000600 cm-1) of the 2,3,7,8-TCDD molecule in the gas phase. It was obtained by GC/FTIR measurements on a Nicolet 170 SX Fourier transform infrared spectrometer and a Hewlett-Packard Model 4882A gas chromatographer. Only the peaks appearing clearly in the FTIR spectrum were listed in Table 3 for comparison. Comparison with experimental results is discussed below. Overall, there is good agreement between the observed and calculated frequencies for the TCDD molecules. The agreement is slightly less satisfactory for some frequencies. For instance, there is a difference of more than 20 wavenumbers between the observed and calculated values for the skeletal vibration. When considering assignments of mode and the intensities of
2,3,7,8-Tetrachlorinated Dibenzo-p-dioxin
J. Phys. Chem., Vol. 100, No. 12, 1996 4813
Figure 3. Raman-active normal vibrational modes of 2,3,7,8-TCDD illustrated together with symmetry, Raman scattering activities (Å4/ amu), and frequency (cm-1). Figure 2. Normal vibrational modes of 2,3,7,8-TCDD having comparatively strong IR intensities (refer to Table 3).
absorptions, we note that the ab initio results are in excellent agreement with the experiments. The largest difference is the broad doublet peaks (1103, 1114; 1248, 1255; 1306, 1321; 1465, 1489 cm-1) seen in the observed spectrum. The calculated results predict that these four doublet bands have a corresponding singlet peak. A possible explanation might be (1) that these observed doublets are due to the isomer of TCDD because of poor gas chromatographic separation, (2) that rotation-vibration interaction causes splitting of peaks into doublets, or (3) that when symmetry is weakened by molecular distortions from D2h to D2 symmetry, the IRRaman mutual exclusion rule is relaxed and some Raman-active bands can be observed in the IR spectrum. A good agreement of all other modes gives confidence that the singular peak is real. Figure 3 gives the animated normal vibrational modes of 2,3,7,8-TCDD having comparatively strong Raman scattering activities with an intensity of more than 10 Å4/amu. These modes should be prominent in the Raman spectrum of TCDD. We hope these results will promote upcoming experimental investigations of the Raman spectrum of this interesting molecule. 4. PEC and Puckering. Figure 4 shows a plot of the relative energy difference, located on the potential energy surface [calculated with full geometry optimization subject to the constraints of D2h for the planar ring and C2s symmetry for bended conformations at the RHF/6-31G, RHF/6-31G*, and RHF/6-31G*(+O,Cl) levels], as a single structural parameter of the folding angle, θ, which is the dihedral angle specifying the puckering of the two benzo-planes. Computations were made at every 1° increment of θ.
Figure 4. Ab initio MO potential, E(θ), displayed against the folding angle θ for TCDD at various calculation levels of 6-31G (bold solid line), 6-31G* (dotted line), and 6-31G*(+O,Cl) (solid line). θ ) 180° corresponds to a planar molecule, in which the total energy is minimum. In the graph, the normal vibrational mode with a frequency of 3.6 cm-1 is also illustrated.
All the basis sets imply the planar form of TCDD as most stable. They do not indicate two shallow potentials (a double minimum potential). The barrier to ring puckering motion is, as expected, lower when higher ab initio levels are used with extended basis sets and diffuse functions. A difference could arise from a less than adequate description of oxygen and chlorine atoms. Since the RHF/6-31G model pictures oxygen and chlorine lone pairs as being too tightly bound, this has an overestimated effect on the calculated potential curve.
4814 J. Phys. Chem., Vol. 100, No. 12, 1996 Compared to thermal energies, the potential has low curvature near 180°, especially at the RHF/6-31G*(+O,Cl) level. The six energies at the RHF/6-31G*(+O,Cl) level are given as follows: At 180° the total energy is -2444.237 29 au and the relative energies at other angles are, in kcal/mol (deg), 0.006 (175), 0.030 (170), 0.089 (165), 0.210 (160), 0.438 (155), and 0.834 (150). At ambient temperatures, nonplanar conformations will be significantly populated for this potential. This prediction accounts for the experimental results that the electric dipole moment of the DD molecules is reported to be 0.64 D.24 The structural changes due to the bending conformers show some interesting features. The bond length changes during the ring puckering are very small (
