ARTICLE pubs.acs.org/JPCA
CCSD(T) Study of Dimethyl-Ether Infrared and Raman Spectra M. Villa† Departamento de Química, UAM-I Purísima y Michoacan, s/n, CP 09340, Mexico D.F.
M. L. Senent* Departamento de Química y Física Teoricas, Instituto de Estructura de la Materia, IEM-CSIC, Serrano 121, Madrid 28006, Spain
R. Dominguez-Gomez‡ Departamento de Ingeniería Civil, Catedra de Química, E.UIT Obras Publicas, Universidad Politecnica de Madrid, Spain
O. Alvarez-Bajo§ and M. Carvajal|| Departamento de Física Aplicada, Facultad de Ciencias Experimentales, Universidad de Huelva, 21071 Huelva, Spain ABSTRACT: CCSD(T) state-of-the-art ab initio calculations are used to determine a vibrationally corrected three-dimensional potential energy surface of dimethyl-ether depending on the two methyl torsions and the COC bending angle. The surface is employed to obtain variationally the lowest vibrational energies that can be populated at very low temperatures. The interactions between the bending and the torsional coordinates are responsible for the displacements of the torsional overtone bands and several combination bands. The effect of these interactions on the potential parameters is analyzed. Second order perturbation theory is used as a help for the understanding of many spectroscopic parameters and to obtain anharmonic fundamentals for the 3N 9 neglected modes as well as the rotational parameters. To evaluate the surface accuracy and to verify previous assignments, the calculated vibrational levels are compared with experimental data corresponding to the most abundant isotopologue. The surface has been empirically adjusted for understanding the origin of small divergences between ab initio calculations and experimental data. Our calculations confirm previous assignments and show the importance of including the COC bending degree of freedom for computing with a higher accuracy the excited torsional term values through the Fermi interaction. Besides, this work shows a possible lack of accuracy of some available experimental transition frequencies and proposes a new assignment for a transition line. As an example, the transition 100 f 120 has been computed at 445.93 cm1, which is consistent with the observed transition frequency in the Raman spectrum at 450.5 cm1.
’ INTRODUCTION Astrophysical surveys contain many spectral lines whose identification requires complete molecular catalogs.1 These lines derive from sources where very well-known molecules coexist with new uncharacterized species and with previously detected molecules not fully described. This is the case of dimethyl-ether (DME), for which only the vibrational ground state has been well explored.27 DME is an astrophysical relevant molecule first detected in emission from the Orion Nebula by Snyder et al. in 1974.8 Later on, it was found to be abundant in star-forming regions where highly excited rotational lines have been observed.9 DME represents the most abundant ether. It is generally accepted that ethers and alcohols coexist in gas phase sources where they are involved in common chemical reactions, coming from gasgrain processes which play an important role in their formation.10,11 Many facts concerning their abundances in hot core regions are expected to be clarified with the new ALMA observatory that will r 2011 American Chemical Society
be a key instrument for their studies.12 For this reason, various laboratories have demonstrated a great interest in DME spectroscopic studies. They try to describe well the millimeter and submillimeter regions for isotopic varieties containing deuterium and 13C, focusing on the monodeuterated species. In addition, the rotational spectrum of the most abundant isotopologue needs to be recorded at the first excited vibrational states that can be populated at hot core temperatures. The expectation of these new experiments has motivated the present study performed with state-of-the art ab initio calculations. From the theoretical point of view, DME is a nonrigid molecule whose potential energy surface contains nine minima.1315 The two methyl group torsions that intertransform and connect Received: July 1, 2011 Revised: October 14, 2011 Published: October 15, 2011 13573
dx.doi.org/10.1021/jp2062223 | J. Phys. Chem. A 2011, 115, 13573–13580
The Journal of Physical Chemistry A the multiple minima by tunneling effects show very low energy levels. The two corresponding fundamentals lie below 250 cm1. A third large amplitude mode, the COC-bending, has been observed in the Raman spectrum at 412 cm1. In 1977, Groner and Durig13 recorded and assigned the far-infrared region spectra of various isotopic DME varieties. For their assignments, they used an effective two-dimensional Hamiltonian adequate for the tunneling splitting description. One year before, in 1976, the microwave spectrum was recorded,2 and more recently, the submillimeter spectrum corresponding to the vibrational ground state has been measured.57 In 1995, some of us performed two ab initio studies of DME using M€ollerPlesset theory and a double-ζ basis set.14,15 One of these studies14 was performed with a two-dimensional model (2D) first developed for the study of acetone.16 The two independent variables were the two methyl group torsions, whereas the remaining 3N 8 vibrational coordinates were optimized in all the conformations selected for the potential energy surface determination. Generally, this procedure describes well the interactions among the torsions and the other vibrational degrees of freedom if they are relatively small. However, in DME, the coupling between the COC bending and the torsions was found to be too large to be accurately described by the relaxation, obliging us to solve, in a second paper, a three-dimensional Hamiltonian (3D) where the COC bending coordinate was added explicitly.15 We concluded that, whereas the fundamental frequencies calculated with the 2D are enough accurate, overtones and several combination bands involving excited torsional levels are displaced by non-negligible Fermi interactions. The actual computational resources that allow us to perform more accurate calculations than those of 1995 and the fresh interest coming from the development of Herschel and the ALMA observatory are our motivations for recovering the theoretical DME study determining a new three-dimensional potential energy surface (3D-PES) using highly correlated ab initio methods. The present paper focuses on the new surface calculation that we test, using it for the study of the most abundant DME isotopologue, 12CH316O12CH3, for which available experimental data can be used to discuss the accuracy. In the future, once new experimental data will be available, we will use this surface to predict many spectroscopic properties of various isotopologues to help assignments. For this purpose, we will combine recent techniques previously employed for acetic acid and methyl-formate isotopologues.17,18 The new potential energy surface has been calculated with CCSD(T)/aug-cc-pVTZ, optimizing the geometry of 126 molecular conformations. Methodological improvements are included to improve old results (i.e., we use a more proper definition of the methyl group torsional coordinates and add the zero point vibrational energy correction to contemplate isotopic effects). Vibrational levels are determined variationally, although we also use perturbation theory as a help to understand the molecular spectroscopic properties.
’ COMPUTATIONAL DETAILS All the ab initio calculations have been performed with the Gaussian 09 package.19 The geometries and the harmonic force field have been calculated with coupled cluster theory using both single and double substitutions.20 To obtain very accurate energies, single point calculations were performed by adding triple excitations noniteratively (CCSD(T)) on the CCSD
ARTICLE
Figure 1. DME at the equilibrium geometry. Independent coordinates α, θ1, and θ2.
Table 1. Structural Parameters (in Å, degrees), Dipole Moment (in Debyes), and Rotational Constants (in MHz) Corresponding to the Minimum Energy Geometry of DME Calculated with CCSD/aug-cc-pVTZ O1C2dO1C3
1.4090
H4C2dH7C3
1.0879
H5C2dH6C2dH8C3dH9C3
1.0962
C3O1C2
111.4
H4C2O1dH7C3O1
107.5
H5C2O1dH6C2O1dH8C3O1dH9C3O1 H4C2O1C3dH7C3O1C2
111.2 180.0
H5C2O1H4dH9C3O1H7
119.4
H6C2O1H4dH8C3O1H7
119.4
μ (MP2/aug-cc-pVTZ)
1.4882
Ae
Be
Ce
38973.4024
10145.2871
8961.8951
geometries.21 To describe well long-range effects, such as the intramolecular nonbonding interactions, the augmented aug-ccpVTZ set was employed in all the computations.22 The vibrational levels corresponding to the three vibrations (the two torsions and the COC bending) have been calculated variationally from a CCSD(T)/CCSD three-dimensional, vibrationally correcte potential energy surface and the ENEDIM code.23 To determine spectroscopic parameters for the 3N 9 remaining coordinates, the full-dimensional anharmonic spectroscopic analysis was performed with second order perturbation theory implemented in FIT-ESPEC.23 The starting point was a CCSD quadratic force field and cubic and quartic force fields determined with second order M€ollerPlesset theory (MP2).24 The vibrational corrections of the 3D-PES are also determined with MP2.
’ RESULTS AND DISCUSSION DME Molecular Structure. At the ground electronic state, the most abundant isotopologue of DME shows a C2v geometry. One hydrogen atom of each methyl group is anti to the other carbon atom (H4C2O1C3dH7C3O1C2 = 180°; see Figure 1). In Figure 1, our three independent coordinates, θ1, θ2, and α, corresponding to the two methyl group torsions and the COCbending modes, are also shown. These coordinates are defined below. At the equilibrium geometry used as reference, θ1 = 0°, θ2 = 0°, and α = 0°. The structural parameters and the rotational 13574
dx.doi.org/10.1021/jp2062223 |J. Phys. Chem. A 2011, 115, 13573–13580
The Journal of Physical Chemistry A
ARTICLE
Table 2. Torsional Energy Barrier V060 and 2D Relative Maximum V6060 (in cm1) Dependence with the Bending Coordinate (in degrees) MP4/6-31G(d,p) α
V060
V6060
6
1260
3012
3
1124
2509
0 3
1006 905
6 9
CCSD(T)/aug-cc-pVTZ
V6060 2V060
V060
V6060
492
1190
2898
261
1046
2390
2082 1734
70 76
921 816
819
1454
184
748
1238
258
n
n
n
n
n
þ
n
n
n
n
518
1210
2919
499
298
1065
2410
280
1966 1618
124 14
939 832
1984 1634
105 30
729
1343
115
742
1353
132
656
1127
185
668
1133
204
1
∑ ∑ ∑ ∑ fijkl qi qjqkql i ¼ 1 j ¼ 1 k l ¼ 1 24 1
V6060 2V060
V6060
∑ ∑ fij qi qj þ t∑¼ 1 j∑¼ 1 k∑¼ 1 6 fijk qi qj qk t ¼1 j¼1 2 1
V6060 2V060
V060
constants calculated with CCSD/aug-cc-pVTZ and the MP2/ aug-cc-pVTZ dipole moment are shown in Table 1. At the equilibrium geometry, the COC bending angle was evaluated to be 111.4° and the rotational constants to be Ae = 38973.4024 MHz, Be = 10145.2871 MHz, and Ce = 8961.8951 MHz. In principle, they can be compared with the experimental values for the vibrational ground state (A0 = 38788.18 MHz, B0 = 10056.48 MHz, and C0 = 8886.83 MHz).5,7 The MP2 dipole moment, μ = 1.4882 D, differs from the experimental one of μ = 1.302 D given in ref 7. DME is a nonrigid molecule whose most abundant isotopologue can be classified in the G36 Molecular Symmetry Group (ref 16 and references therein). The internal rotation of the two methyl groups intertransforms nine minima through feasible energy barriers. Table 2 shows the one-dimensional V060 and the two-dimensional V6060 barriers and the V6060 2V060 difference as a function of the COC bending coordinate (V060 = E(α,0,60) E(α,0,0); V6060 = E(α,60,60) E(α,0,0)), and E(0,0,0) is the energy of the reference equilibrium geometry. V6060 2V060 represents a measure of the two methyl group interactions. When the COC bending angle is augmented, the nonbonding interactions among different methyl group hydrogen atoms decrease. For small values of α, V6060 2V060 > 0, whereas, for large values, V6060 2V060 < 0. All the interaction terms of the torsional Hamiltonian depend strongly on α. Full Dimensional Anharmonic Analysis. A preliminary anharmonic analysis using second order perturbation theory (PT2) has been performed for evaluating the number of internal coordinates that could be considered as independent variables in the variational calculation of the CH3 torsional levels. Although PT2 represents a deficient theoretical model for nonrigid molecules, it is taken as a first approximation for the determination of the torsional frequencies as well as for the analysis of possible Fermi interactions among our three independent modes and the remaining vibrations. At the same time, this theory allows us to obtain the fundamental frequencies for the remaining vibrational modes. For this analysis, the quadratic, cubic, and quartic force field terms expressed in V ðq1 , q2 , :::, qn Þ ¼
CCSD(T)/aug-cc-pVTZ+ZPVE
ð1Þ
have been evaluated using CCSD and MP2 calculations and the subroutines implemented in Gaussian.19 We use the FIT-ESPEC code23 to transform the force field in normal coordinates
Table 3. Harmonic (ω) and anharmonic (ν) fundamental frequencies of DME-h6 (in cm‑1) MP2/aug-cc-pVTZ C2v A1
A2
B1
B2
a
mode
ω
ν
CCSD/aug-cc-pVTZ
ν + Δν
a
ω
ν
υ1
3179
3043
3050
3147
3014
υ2
3024
2888
2827
3011
2755
υ3
1533
1489
1490
1537
1492
υ4
1498
1461
1461
1511
1474
υ5
1275
1244
1244
1291
1260
υ6
958
931
931
971
946
υ7 υ8
417 3090
410 2961
405 2978
421 3059
413 2958
υ9 υ10
1507 1173
1461 1149
1461 1151
1511 1180
1465 1156
υ11
210
202
202
202
189
υ12
3082
2954
2950
3055
2939
υ13
1518
1471
1470
1521
1474
υ14
1201
1173
1173
1212
1183
υ15
260
246
246
252
235
υ16
3178
3043
3045
3145
3013
υ17 υ18
3019 1519
2903 1471
2898 1471
3002 1522
2956 1475
υ19
1467
1434
1434
1478
1445
υ20
1211
1176
1176
1231
1198
υ21
1131
1105
1105
1139
1113
Δν = Fermi displacements (emphasized in bold if Δν > 10 cm1).
produced by Gaussian, in a force field defined in Cartesian and internal coordinates. All the 3N 6 frequencies are localized and identified in Table 3. Bands displaced by strong Fermi interactions are emphasized in bold type. The CCSD calculations are quite time demanding; consequently, the quadratic force field is considered at the CCSD level, whereas the cubic and quartic fields are evaluated with MP2. The MP2 anharmonic effects included in the CCSD case are a good approximation due to the small value of the anharmonic correction. Then, the anharmonic frequencies of Table 3 were determined from MP2 theory (4th column) and combining CCSD and MP2 theories (7th column). Harmonic frequencies allow us to compare MP2 and CCSD theories. In general, for the hydrogen stretching modes and the torsional modes, CCSD frequencies are lower than the MP2 ones whereas for the remaining modes CCSD produces higher values than MP2. Anharmonic fundamentals follow the same tendency 13575
dx.doi.org/10.1021/jp2062223 |J. Phys. Chem. A 2011, 115, 13573–13580
The Journal of Physical Chemistry A
ARTICLE
Table 4. Predicted Spectroscopic Properties Using Second Order Perturbation Theory MP2/aug-cc-pVTZ CCSD/aug-cc-pVTZ
obs7
Fermi Displacements of Torsional Bands (in cm‑1) 2ν11 (torsion)
401 f 398
370 f 370
2ν15 (torsion) ν7 (COC bending)
491 f 498 410 f 405
467 f 479 413 f 401
Equilibrium Rotational Constants (in MHz) Ae
38549.07
38973.40
Be
10192.86
10145.29
Ce
8974.20
8961.90
following our experience in this type of studies, the calculated variation of the rotational constant with the vibrational energy (α-rovibrational parameters) is trustworthy. 3D-Variational Calculations. The energy levels corresponding to the two large amplitude vibrations and the COC bending mode of DME have been determined by solving variationally the following three-dimensional Hamiltonian:26 ! ! 3 3 ∂ ∂ ^ 3D ¼ ∑ ∑ H Bij þ V ðα, θ1 , θ2 Þ ∂qi ∂qj i j þ V 0 ðα, θ1 , θ2 Þ þ V ZPVE ðα, θ1 , θ2 Þ
i, j ¼ α, θ1 , θ2
ð2Þ
Rotational Constants (in MHz) A0
38405.04
38966.84
38788.18
B0
10038.45
10000.01
10056.48
C0 Aν7
8851.91 39030.36
8846.44 39517.12
8886.83
Bν7
10014.68
10030.27
Cν7
8817.27
8810.80
Aν11
38417.92
38205.88
Bν11
10009.19
9987.19
Cν11
8828.64
8801.94
Aν15
38419.36
38539.24
Bν15 Cν15
9967.78 8824.86
9961.91 8806.22
The first term, which depends on the Bij parameters (the G matrix elements in cm1), is a 3D-kinetic energy operator T3D; V(α,θ1,θ2), V0 (α,θ1,θ2), and VZPVE(α,θ1,θ2) represent the potential energy surface, the “Podolsky pseudopotential”, and the zero point vibrational energy correction, respectively (ref 26 and references therein). As the 3D-PES transforms as the totally symmetric representation of the G36 Molecular Symmetry Group,1416 we use the following analytical expression for the potential energy surface: V ðα, θ1 , θ2 Þ ¼
4
þ
Centrifugal Distortion Constants (in MHz) ΔJ
0.0095
0.0091
0.0091
ΔK
0.3346
0.3328
0.3420
ΔJK
0.0313
0.0292
0.0269
δJ
0.0019
0.0018
0.0018
δK
0.0087
0.0117
0.0138
rule, with the exception of ν7 (COC bending) and ν17 (H-asymmetric stretching), strongly displaced by Fermi interactions. With MP2, ω7 and ν7 are calculated to be 417 and 405 cm1, whereas, with CCSD, they are calculated to be 421 and 401 cm1 (see Tables 3 and 4). The tendency rule infringement, which produces an overestimation of the Fermi displacements, may be explained from the mix of levels of theory used for anharmonic CCSD calculations. To our knowledge, DME fundamentals of the medium and high frequency modes have been measured in the solid phase,25 except for the ν6, ν17, and ν21 modes, which were measured in the gas phase.6 We cannot compare our results with available experimental data, especially when possible errors are too small. The COC bending (ν7) interacts strongly with the 2ν15 torsional overtone (see Table 4) and increases it an amount of 712 cm1. Perturbation theory confirms the need for a 3-dimensional (3D) model to obtain torsional excited CH3 vibrational states above the fundamentals. Fortunately, interactions with the remaining 3N 9 modes are negligible. A 3D model is sufficient. In Table 4, some rotational constants and centrifugal distortion constants, derived from second order perturbation theory, are shown. For the vibrational ground state, they are compared with experimental data.7 As far as we know, there are not published experimental data for excited vibrational levels. There is a good agreement for the centrifugal distortion constant, and
2
2
∑ ∑ ∑ ½AMNL αM cosð3Nθ1 Þ cosð3Lθ2 Þ M¼0 N ¼0 L¼0 BM11 αM sinð3θ1 Þ sinð3θ2 Þ ∑ M¼0
ð3Þ
In this equation AMNL = AMLN. Formally identical expressions can be employed for V0 , VZPVE, and Bθ1θ2. Although the kinetic term T3D is totally symmetric (A1) (ref 16 and references therein), some kinetic parameters present low symmetry properties in several conformations. The Hamiltonian parameters have been determined from the energies and optimized geometries of a set of 126 selected conformations. They have been chosen for different values of the COC bending angle (104.676° f 119.676°, Δα = 3°), and the H4C2O1C3 and H7C3O1C2 dihedral angles (0, 30, 90, 150, 180, 30, 90, and 150°). The remaining 3N 9 parameters were allowed to be relaxed. Whereas the geometries were calculated with CCSD/aug-cc-pVTZ, the energies were obtained by performing single point CCSD(T)/aug-cc-pVTZ calculations on the partially optimized geometries. For each molecular structure, the G matrix elements and V0 were computed with the MATRIZG subroutine of the ENEDIM code from the internal coordinates,23 and the ZPVE correction was determined from the MP2/aug-cc-pVTZ harmonic frequencies ωi: EZPVE ¼
i ¼ 3N 6
∑ i¼n þ 1
ωi 2
In this equation, the sum neglects two torsional modes and the COC bending modes (n = 3). V is isotopically invariable, but VZPVE, V0 , and Bij are different for the different isotopologues. For this reason, the effective potential VEF (VEF = V + VZPVE + V0 ) and the effective barriers depend on the nuclear masses. The vibrational correction effect on the DME-h6 potential parameters is shown in Table 2. For α = 0, V060 varies from 921 to 939 cm1 and V6060 from 1966 to 1984 cm1 when ZPVE is considered. The correction augments the barrier height, allowing prediction of a displacement of the 13576
dx.doi.org/10.1021/jp2062223 |J. Phys. Chem. A 2011, 115, 13573–13580
The Journal of Physical Chemistry A
ARTICLE
spectra to high frequencies and a reduction of the torsional splittings. The geometry optimization is a partial way to consider vibrational interactions. Thus, during the motions of the independent variables (torsions and COC bending), the remaining coordinates are forced to lie at their minima. Later on, once the ZPVE correction is added, they are fixed at their zero point vibrational energy. This procedure leads to accurate results but introduces numerical errors that need to be minimized. For example, the tree dihedral angles H4C2O1C3 (D11), H5C2O1C3 (D12), and H6C2O1C3 (D13) of each methyl group (or H7C2O1C3 (D21), H8C2O1C3 (D22), and H9C2O1C3 (D23)) are not identically treated. One is fixed (D11 or D21), and the other two are relaxed, producing an artificial distortion and symmetry breaking. To obtain a proper definition of torsional coordinates, we used linear combinations of internal coordinates transforming as the C3v symmetry group representation for defining the vibrational coordinate:27 θ1 ¼ ðD11 þ D12 þ D13 2πÞ θ2 ¼ ðD21 þ D22 þ D23 2πÞ
Table 5. Expansion Coefficients for the CCSD(T)/aug-ccpVTZ Potential Energy Surface, V(α,θ1,θ2), the Vibrationally Corrected VEF(α,θ1,θ2) (VEF = V + V0 + VZPVE), and the Refined VREF Surfaces (in cm1/degreesM)a coeff
a
Not totally symmetric linear combinations are used to define four deformation coordinates: 1 β1 ¼ ðD11 D12 Þ 2 1 β3 ¼ ðD21 D22 Þ 2 1 1 β2 ¼ ðD11 D12 Þ þ D13 4 2 1 1 β4 ¼ ðD21 D22 Þ þ D23 4 2
ð5Þ
To obtain the 3D-PES, the total electronic energies of the 126 structures are fitted to a 7D-PES, given by the function of eq 3 plus some additional terms depending on the deformation coordinates: Eα, θ1 , θ2 , β1 , β2 , β3 , β4 ¼
1
4
½CMi αM sinð3θ1 Þ sinð3θ2 Þβi ∑ ∑ M¼0 i¼1 þ
1
4
4
∑ ∑ ∑ ½CijαM sinð3θ1Þ sinð3θ2Þβiβj M¼0 i¼1 j¼i þ 1 ð6Þ
The Cij coefficients corresponding to the additional terms depending on βi are very small, and they are neglected.17 It has to be pointed out that errors derived from the combination of different levels of theory (CCSD(T), CCSD, and MP2 used for energy, geometry, and vibrational corrections. respectively) are really small (
