A Nine-Dimensional Calculation of the Vibrational OH Stretching and

Apr 23, 2012 - The kinetic energy operator is a sum of the kinetic energy operators of the monomer units. We use the coupled cluster method with singl...
2 downloads 0 Views 795KB Size
Article pubs.acs.org/JPCA

A Nine-Dimensional Calculation of the Vibrational OH Stretching and HOH Bending Spectrum of the Water Trimer Teemu Salmi, Elina Sal̈ li, and Lauri Halonen* Laboratory of Physical Chemistry, Department of Chemistry, A.I. Virtasen aukio 1 (P.O. BOX 55), FI-00014, University of Helsinki, Finland S Supporting Information *

ABSTRACT: We have studied the vibrational high-frequency spectrum of the water trimer computationally. We expand an earlier study [J. Chem. Phys. A 2009, 113, 9124−9132] where we approximated the water trimer as three individually vibrating water monomer units. Some intramolecular potential energy coupling terms are now included in the previous model. The six OH bond lengths and the three HOH bending angles are used as the internal coordinates. The kinetic energy operator is a sum of the kinetic energy operators of the monomer units. We use the coupled cluster method with single, double, and perturbative triple excitations method [CCSD(T)] with augmented correlation consistent polarized valence triple-ζ (aug-cc-pVTZ) basis set to calculate the potential energy surface (PES). The counterpoise correction is included in the one-dimensional part of the PES. We calculate the vibrational energy eigenvalues using the variational method. The corresponding eigenfunctions are used to obtain the absorption intensities.

1. INTRODUCTION The water trimer consists of three water monomer units that are hydrogen bonded to each others (see Figure 1). In the

has considered the vibrations of the deuterated water trimer species. The other earlier experimental methods to study the water trimer include helium droplet experiments11−13 and gasphase methods such as far-infrared laser absorption, terahertz laser spectroscopy, infrared cavity ring-down and laser difference frequency spectroscopy in supersonic expansions, and infrared molecular beam depletion and fragment spectroscopy.14−21 Most of the experimental studies have concentrated on the intermolecular vibrations. The vibrational properties of the water trimer have also been studied computationally.22−34 These works concentrate mostly on the tunneling motion between the potential energy minima or the fundamental vibrational transitions. There are only a few computational studies on the high-frequency overtones.24,32 This work is based on our earlier study on the water trimer.32 We treated the water trimer as three individually vibrating monomer units. The model worked mostly fine but failed to reproduce the correct band structure for the hydrogen-bonded OH stretching modes.9 In this work, we expand the earlier model by including some intermolecular potential energy coupling terms into the Hamiltonian as was also suggested in ref 9.

Figure 1. Water trimer.

cyclic minimum energy structure, the oxygen atoms form an almost equilateral triangle. All three monomer units act as both hydrogen donor and acceptor in a hydrogen bond. To minimize the tension, the free hydrogen atoms are oriented to maximize the distance from each other. The bonded hydrogen atoms lie close to the plane containing the oxygen atoms. The matrix isolation spectroscopic method is an efficient way to study the water trimer experimentally. Even the overtones can be observed. There has been experimental work since the 1970s,1−8 but the topic is still actual as there are two recent studies on high-frequency modes9,10 with which we shall compare our results. In the first one,9 the fundamentals and some overtone 10 18 transitions of (H16 2 O)3 and (H2 O)3 are given. The second one © 2012 American Chemical Society

2. VIBRATIONAL MODEL The starting point for this work is our previous article concerning the water trimer.32 In the present contribution, as well as in the previous work, we use the CCSD(T)/AVTZ + CP32 Received: February 22, 2012 Revised: April 23, 2012 Published: April 23, 2012 5368

dx.doi.org/10.1021/jp3017584 | J. Phys. Chem. A 2012, 116, 5368−5374

The Journal of Physical Chemistry A

Article

(counterpoise corrected,35,36 coupled cluster method with single, double, and perturbative triple excitations with the aug-cc-pVTZ37,38 basis set) potential energy surface (PES) for the one-dimensional slice of the PES and the CCSD(T)/AVTZ surface for the two-dimensional part. The definitions and the parameters can be found in ref 32. (There has been some reported cases where the CP correction has not resulted in a significant improvement when using the CCSD(T) method.39 However, we were able to improve the hydrogen-bonded OH stretching energies by 12−14 cm−1 by using the CP correction in the one-dimensional part of the PES).32 Earlier, we used the exact kinetic energy operator for a triatomic molecule (see refs 40, 41, and references therein). The kinetic energy operator includes a term, the pseudopotential term, which does not have any wave function derivatives. It depends on the Jacobian of the coordinate transformation from the Cartesian to the internal coordinates. Hence, the triatomic Jacobian differs from the nine atomic one which we use in this work. Therefore, we ignore the pseudopotential term in this work. We have run test calculations for a single monomer unit both with and without the pseudopotential term. Its effect on the fundamental and overtone transitions turned out to be insignificant when compared with the experimental work.9 Now, we have a Hamiltonian operator

quantum number for the HOH bending modes is 14 and for the sum of the OH stretching modes is 16. 2. Second, we take 30 lowest energy eigenstates from each monomer calculation. We use a product of the monomer eigenstates as the basis set for the trimer calculation. The maximum of the sum of the monomer quantum numbers is 30. After performing these steps, we are able to reduce the number of basis functions to be small enough for a feasible computational time. This procedure produces a good convergence at least up to the first OH stretching overtone region. The intermolecular potential energy coupling surface (eq 2) is calculated using the CCSD(T)/AVTZ method around the corresponding minimum energy geometry.32 All ab initio calculations are performed using Molpro 2010.1.42 The grid points are obtained displacing two rbi distances at a time independently from rebi − 0.20 Å to rebi + 0.35 Å with a 0.05 Å step size and a maximum total displacement of 0.40 Å. The coefficients of the potential energy function in eq 2 are optimized with the least-squares method. These intermolecular coupling coefficients are presented in Table 1. We tried to improve our model by adding the potential energy terms that couple rb with rf and θ. However, this did not produce a significant change in the results. Therefore, we decided to ignore these couplings to keep our model as simple as possible. We calculated absorption intensities as the dimensionless oscillator strengths using the formula43

Ĥ 0 = H1̂ (rb1 , rf1 , θ1) + Ĥ 2(rb2 , rf2 , θ2) + H3̂ (rb3 , rf3 , θ3) (1)

where subscripts 1, 2, and 3 refer to the labeling of the monomer units, rb and rf are the hydrogen bonded and the free OH distances, respectively, and θ is the HOH bending angle. All three parts contribute only to a single water monomer unit. The eigenvalues of these three Hamiltonians can be solved independently as we did in our previous work.32 In this work, we have improved the model by adding potential energy coupling terms between the OHb stretching modes of the monomer units as follows ̂ (rb , rb , rb ) = V123 1 2 3

∑ i ,j,k

fijk i! j! k!

feg = 4.70·10−7[cm D−2 ]νeg̃ |⟨e|μ ⃗ |g|⟩|2

where ν̃eg is the transition wavenumber from the ground state, g, to an excited state, e, and μ⃗ is the dipole moment vector operator in Debye (1 D = 3.33564 × 10−30 Cm). We use our earlier dipole moment vector from ref 32, where the definitions of the dipole moment axes are given, for the intramolecular degrees of freedom. The intermolecular dipole moment coupling vector, which is calculated in this work, is added to the intramolecular part. The two-dimensional grid is calculated using the CCSD(T)/AVDZ method. We obtain the dipole moment grid by displacing two rbi distances at a time between ±0.20 Å with a 0.05 Å step size. The intermolecular dipole moment function coefficients are presented in Table 2.

(Δrb1)i (Δrb2) j (Δrb3)k (2)

where the quantity f ijk is a potential energy parameter which describes the strength of interaction, Δr = r − re, and re is the equilibrium OH bond length. The improved Hamiltonian is ̂ Ĥ = H1̂ + Ĥ 2 + H3̂ + V123

(4)

(3)

3. RESULTS We may consider the planar structure as a reference when we evaluate the vibrational transitions of the water trimer. The planar reference structure belongs to the symmetry point group C3h. The fundamental bonded OH stretches (OHb) span A′ + E′ irreducable species. Hence, the IR transition from the ground state to the symmetric A′ OHb stretching state is group theoretically forbidden. The two asymmetric E′ modes are degenerate. Even though the minimum energy structure of the water trimer belongs to the C1 point group, it is close to C3h. A similar energy level structure can be derived also for the free OH stretches (OHf) and the HOH bends. However, the coupling of these modes is weak, and the symmetric/ asymmetric difference is small. The variationally calculated vibrational band centers of the water trimer are presented in Table 3. We use a notation where δ is a HOH bend and νs/aOHb is a hydrogen-bonded OH stretch, where the subscripts s and a refer to symmetric and asymmetric modes, respectively (even though there is no

There is no direct intermolecular coupling in the kinetic energy operator. The coupling would exist indirectly via the lowfrequency modes. As these modes are excluded from our model, no intermolecular kinetic energy coupling is taken into account. We solve the Schrödinger equation using the variational method. When the intermolecular coupling term is added to Ĥ 0, the monomer units cannot be treated independently anymore. There are nine degrees of freedom in the Hamiltonian Ĥ . Therefore, the variational problem consists of nine variables. If we use a product of the nine one-dimensional basis functions as the basis set for the trimer, we would need to diagonalize a very large matrix. That is an unnecessarily demanding task. Therefore, we need to choose the most important basis functions and discard the idle ones. We reduce the basis set in two steps. 1. First, we solve the uncoulped Schrödinger equation for each monomer unit independently. The maximum 5369

dx.doi.org/10.1021/jp3017584 | J. Phys. Chem. A 2012, 116, 5368−5374

The Journal of Physical Chemistry A

Article

Table 1. Intermolecular Potential Energy Coupling Coefficients (CCSD(T)/AVTZ) for the Water Trimera f rb1rb2 (aJ Å−2)

−0.08176

f rb1rb2rb2 (aJ Å−3)

−0.11990

μrb1rb2 (D Å−2)

−0.15606

μrb1rb2rb2 (D Å−3)

f rb1rb2rb2rb2 (aJ Å−4) −5

f rb1rb2rb2rb2rb2 (aJ Å )

x

−3

μrb1rb1rb2 (D Å )

−0.18200

−4

f rb1rb1rb2 (aJ Å )

−0.19656

μrb1rb1rb2rb2 (D Å )

f rb1rb1rb2rb2 (aJ Å−4)

−0.24601

μrb1rb3 (D Å−2)

−3

−3

y

z

0.05518

−0.67153

−0.09520

−0.44519

−1.08614

−0.08917

0.29146

−1.81770

−0.24336

−0.51010

−2.95151

−0.23479

0.57137

0.41606

0.04801

1.25569

0.12400 −0.02302

f rb1rb1rb2rb2rb2 (aJ Å )

−0.68249

μrb1rb3rb3 (D Å )

1.48916

f rb1rb1rb1rb2 (aJ Å−4)

−0.35087

μrb1rb1rb3 (D Å−3)

1.17848

0.22912

f rb1rb1rb1rb2rb2 (aJ Å−5)

−0.52961

μrb1rb1rb3rb3 (D Å−4)

2.87120

1.22620

−0.01468

0.26591

−0.00840 0.08838

−5

−2

μrb2rb3 (D Å )

−0.61237

f rb1rb3 (aJ Å−2)

−0.08720

μrb2rb3rb3 (D Å−3)

−0.78313

0.82815

f rb1rb3rb3 (aJ Å−3)

−0.21401

μrb2rb2rb3 (D Å−3)

−1.76886

0.57290

0.03970

−0.37192

μrb2rb2rb3rb3 (D Å−4)

−2.51570

1.66140

0.21311

−5

f rb1rb1rb1rb1rb2 (aJ Å )

−4

f rb1rb3rb3rb3 (aJ Å ) f rb1rb3rb3rb3rb3 (aJ Å−5) f rb1rb1rb3 (aJ Å−3)

0.24360

0.48301 −0.12928

−4

f rb1rb1rb3rb3 (aJ Å )

−0.28198

f rb1rb1rb3rb3rb3 (aJ Å−5)

−0.60708

f rb1rb1rb1rb3 (aJ Å−4)

−0.17279

−5

f rb1rb1rb1rb3rb3 (aJ Å )

−0.58409

f rb1rb1rb1rb1rb3 (aJ Å−5)

−0.68210

f rb2rb3 (aJ Å−2)

−0.08024

−3

f rb2rb3rb3rb3 (aJ Å−4)

−0.17571

f rb2rb3rb3rb3rb3 (aJ Å−5)

−0.24918

f rb2rb2rb3 (aJ Å−3)

−0.18528

−4

f rb2rb2rb3rb3 (aJ Å )

−0.25656

f rb2rb2rb3rb3rb3 (aJ Å−5)

−0.65961

f rb2rb2rb2rb3 (aJ Å−4)

−0.32058

−5

f rb2rb2rb2rb3rb3 (aJ Å )

−0.62431

f rb2rb2rb2rb2rb3 (aJ Å−5)

−0.20568

1 Å =10

−10

Table 3. Calculated Vibrational Band Centers and Intensities for the Water Trimerc

−0.12664

f rb2rb3rb3 (aJ Å )

a

Table 2. Intermolecular Dipole Moment Coupling Coefficients (CCSD(T)/AVDZ) for the Water Trimer

m.

symmetry between the water monomer units, we shall denote these states as symmetric and asymmetric ones in the text and the tables; the main contribution to the symmetric wave function consists of parts of the same sign; in the asymmetric parts, not all the signs are the same), and νOHf is a free OH stretch. We have calculated spectra for both (H16 2 O)3 and (H18 2 O)3. The band centers are calculated as the weighted average of the transitions corresponding to a given band. The calculated intensity is used as the weighing function. The presented intensities are the sums of the individual transitions. In Table 3, we have scaled the intensities so that the strongest band νaOHb gets the value 100. In this way, computed intensities can be compared with the scaled observed intensities.9 The energies and the intensities of the individual transitions are presented in Table 4. In sections 3.1−3.6, we shall discuss the transitions of (H16 2 O)3. In sections 3.7 and 3.8, we shall present vibrational transitions for (H18 2 O)3 and the deuterated species, respectively. 3.1. Fundamentals of the Bonded OHb Stretching Modes. There are three fundamental OHb stretching transitions in the water trimer. In our earlier work, we calculated the spectrum without the intermolecular coupling terms.32 We obtained three energies: 3445, 3453, and 3467 cm−1 with

a

Overlap with dimer bands makes intensity measurement inaccurate. Bands which consist of combinations of the OHb stretching fundamentals: vOHbi + vOHbj. sh = shoulder. cThe calculated intensities are scaled so that the total intensity of the most intense transition, νaOHb, gets the value 100. b

oscillator strengths 6.6 × 10−5, 6.3 × 10−5, and 5.7 × 10−5, respectively. The total intensity is 18.6 × 10−5. When we include the intramolecular coupling terms into our model, we obtain a band structure which is closely associated with the planar C3h transition state: a weak transition to a symmetric state at 3386 cm−1 and two strong transitions to almost degenerate asymmetric states at 3451 and 3463 cm−1. Even though there is a significant difference between the coupled and the uncoupled models, the total intensity remains the same, 18.6 × 10−5. 5370

dx.doi.org/10.1021/jp3017584 | J. Phys. Chem. A 2012, 116, 5368−5374

The Journal of Physical Chemistry A

Article

3472.2, 3517.5, and 3529.5 cm−1.9 Our improved method reproduces the observed splittings and transition intensities well. However, as expected, our model still overestimates the redshift in the OHb stretching modes. This problem is discussed in more detail in refs 32, 40, and 44. Even though it may seem that the difference between the observed and the calculated energy levels is larger in the new than in the old model, this is actually not the case. What should be compared instead of the plain numbers is the band centers. The weighted average calculated using the new model OHb stretching fundamentals is 3455.7 cm−1, whereas with the old model it is 3454.5 cm−1. Therefore, the inclusion of the intermolecular coupling has little effect on the band centers. 3.2. Bonded OHb Stretching Overtones. As usually, the intensity of the first hydrogen-bonded OHb stretching overtone is smaller than the one of the free OHf stretch.45 There are six states at the OHb stretching overtone region. Again, according to our calculations, there are a weak symmetric (6638 cm−1) and two stronger asymmetric transitions (6666 and 6700 cm−1) accompanied with three combination transitions. Only one peak (at 6793 cm−1, assigned as 2νsOHb) was found in the matrix isolation experiment.9 Our calculated intensity for 2νsOHb at 6638 cm−1 is an order of magnitude smaller than the observed one at 6793 cm−1. A similar behavior was also observed when we compared our water dimer calculations40 with the matrix isolation observation46 (0.05 vs 0.36). Therefore, it seems that our model overestimates the intensity decrease in the OHb band in the overtone region when compared with the matrix isolation experiments. However, the total intensity of 2νsOHb and the two 2νaOHb transitions is 0.08 which is close to the observed one. Therefore, we cannot rule out the possibility that the observed feature at 6793 cm−1 could also be 2νaOHb instead of 2νsOHb. The intensities are calculated using eq 4. Small intensities are difficult to calculate accurately as small changes in the eigenfunctions can have large effects on them. The three combination states are centered at 6878 cm−1 with a relative intensity of 0.054. When the wavenumber error within our model, which is about 70 cm−1, is taken into account, we suggest that it could be observed at (6878 + 70 + 70) cm−1 ≈ 7018 cm−1. However, the water dimer band located at 7018 cm−1 46 could make the observation difficult with the matrix isolation technique. 3.3. Free OHf Stretching Modes. The three computed OHf stretching fundamentals lie close to each other centered at 3709 cm−1, which is in good agreement with the experiment. The difference is 15 cm−1. The computed first overtones are found between 7266 and 7276 cm−1 (Table 4). The difference between the calculated and the experimental values is 41 cm−1. On the basis of our experience with the water dimer,40 we assume that these results could be improved with a more advanced ab initio method. In ref 40, the difference in the OHf stretching fundamental decreased from 24.1 cm−1 to 6.6 cm−1 when the basis set was increased from AVTZ to AVQZ. The intensities of the fundamentals and the overtones are in good agreement with the experiment. 3.4. νOHf + νOHb Combination Bands. Of the nine νOHf + νOHb states, we compute three groups consisting of three states. The first group consists of energy levels at 7082, 7083, and 7089 cm−1. Their wave functions are symmetric in a sense discussed earlier in section 3.1. The remaining six states form two groups of asymmetric levels; the first one consists of states which are formed approximately as 2|1⟩fi|1⟩bi −|1⟩fi|1⟩bj − |1⟩fi|1⟩bk,

Table 4. Calculated Vibrational Energy Levels and Oscillator Strengths for the Water Trimer Using the Coupled and Uncoupled Models uncoupled model40

this work ν̃/cm−1 δ



νOHb

νOHf νOHb + δ

νOHf + δ

νOHb + 2δ

2νOHb

νOHf + 2δ

νOHf + νOHb

2νOHf

1610 1610 1615 3186 3187 3198 3386 3451 3463 3705 3709 3710 4983 4986 4994 5050 5050 5061 5066 5070 5072 5297 5300 5305 6533 6542 6564 6612 6615 6627 6645 6646 6647 6638 6666 6700 6831 6871 6884 6855 6858 6869 7082 7083 7089 7141 7145 7148 7161 7165 7171 7266 7271 7276

osc. strength 8.3 1.0 9.2 8.1 1.5 1.2 2.9 9.5 8.8 6.2 1.3 1.5 9.4 7.0 4.5 6.3 5.5 6.7 2.4 5.4 1.5 1.2 1.2 1.3 3.9 3.8 2.4 8.6 4.1 5.9 1.5 6.7 2.0 1.4 7.3 6.0 3.9 3.1 6.5 2.2 2.3 7.0 5.7 1.5 6.7 1.2 1.2 2.0 4.7 1.2 3.5 1.7 1.7 1.6

× × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × × ×

−6

10 10−5 10−6 10−7 10−6 10−6 10−6 10−5 10−5 10−6 10−5 10−5 10−7 10−7 10−7 10−7 10−7 10−7 10−7 10−8 10−8 10−6 10−6 10−6 10−9 10−9 10−9 10−9 10−9 10−9 10−9 10−10 10−9 10−8 10−8 10−8 10−9 10−8 10−8 10−8 10−8 10−8 a 10−8 10−7 10−8 10−7 10−7 10−7 10−9 10−8 10−9 10−7 10−7 10−7

ν̃/cm−1

osc. strength

1611 1611 1615 3188 3189 3199 3445 3453 3467 3707 3710 3711 5032 5041 5062

9.3 9.0 9.1 1.3 1.2 1.2 6.6 6.3 5.7 1.1 1.1 1.1 1.5 1.3 1.3

× × × × × × × × × × × × × × ×

10−6 10−6 10−6 10−6 10−6 10−6 10−5 10−5 10−5 10−5 10−5 10−5 10−6 10−6 10−6

5299 5302 5307 6577 6589 6618

1.2 1.2 1.3 5.8 4.8 6.2

× × × × × ×

10−6 10−6 10−6 10−9 10−9 10−9

6704 6723 6756

6.4 × 10−8 7.5 × 10−8 6.0 × 10−8

6858 6860 6872 7127 7129 7144

2.0 2.2 2.2 2.3 2.4 2.6

7270 7276 7281

1.7 × 10−7 1.7 × 10−7 1.6 × 10−7

× × × × × ×

10−8 10−8 10−8 10−7 10−7 10−7

a

The state is strongly coupled with the 2vaOH(comb) state. Hence, the transition intensity is unexpectedly large.

In a recent Ne matrix isolation spectroscopic experiment, the OHb stretching fundamental transitions were observed at 5371

dx.doi.org/10.1021/jp3017584 | J. Phys. Chem. A 2012, 116, 5368−5374

The Journal of Physical Chemistry A

Article

and the second one is |1⟩fi|1⟩bj − |1⟩fi|1⟩bk, where the subscripts i, j, and k refer to the three monomer units in the water trimer. The three computed band centers are at 7084, 7145, and 7165 cm−1. The splittings are almost identical with those of the fundamental OHb stretching region. However, the computed intensity of the symmetric transition is not weak compared to the asymmetric ones. The first group of asymmetric transitions (centered at 7145 cm−1) is an order of magnitude stronger than the second one (7165 cm−1). Therefore, we expect the second group to appear as a shoulder to the first one. Only two peaks were observed in the matrix isolation experiment. A stronger one was at 7182 cm−1, and a shoulder was at 7197 cm−1. The one at the lower energy was assigned as the symmetric and the higher one as the asymmetric stretch. The differences between the observed and the calculated wavenumbers would be 98 cm−1 and 52 cm−1, respectively. Also, the splitting between the symmetric and asymmetric states should be larger on the basis of the OH stretching fundamentals. This does not seem to be correct. Assigning the two observations as the asymmetric stretches would sound tempting. However, the observed minus calculated values assuming that the observations are the two asymmetric modes would be 37 cm−1 and 32 cm−1. These separations are too small compared with the errors in the fundamental OH stretching region. Therefore, we suggest that the transitions at 7182 cm−1 and 7197 cm−1 are from the ground to the symmetric states. The two calculated peaks at 7082 cm−1 and 7083 cm−1 form a single peak, and the one at 7089 cm−1 is the shoulder. Doing this, we obtain a splitting of 6 cm−1, which is slightly smaller than the observed one, 15 cm−1. 3.5. HOH Bending Modes. As the HOH bending modes are not strongly coupled, there are only small changes between our current and previous work.32 The calculated energies and the transition intensities are in good agreement with the experiment in the case of the fundamental and overtone transitions.9 Similarly, with the stretching modes, combination bands appear in the overtone region. However, these are weak and are likely not to be observed. 3.6. OH Stretching + HOH Bending Combination Bands. As there are three stretching and three bending modes, there are in total nine νOH + δ states. Three of these consist of νsOH + δ, and the remaining six consist of νaOH + δ. The center of the νsOHb + δ band is located at 4986 cm−1. The observed νsOHb + δ band was found at 5078.5 cm−1, and therefore, the difference is 92.5 cm−1. Our model overestimates the redshift of νsOHb by 86.2 cm−1, and hence, the agreement with the experiment is good. The asymmetric combination band, νaOHb + δ, consists of six transitions. These are presented in Table 4. Two bands are observed in the matrix isolation experiment. Their splitting is 23.5 cm−1. Combining our calculated transitions leads into two peaks with a slightly smaller splitting. The νOHb + 2δ as well as the νOHf + 2δ region consists of 18 vibrational states. In Table 4, we report the most intense ones. 3.7. Vibrations of (H18 2 O)3. We have also calculated the vibrational transitions for (H18 2 O)3. The results as well as the matrix isolation observations from ref 9 are presented in Table 3. Replacing 16O with the heavier isotope 18O results in a slight lowering of all calculated vibrational transitions. Our model produces shifts that are in a good agreement with the observations. Our calculated intensities remain basically unchanged. The only observable difference is found in νOHf + 2δ for which the intensity is halved. This is due to coupling of one νOHf + 2δ state with a 2νOHb state. Hence, the intensity of the 2νOHb

band increases slightly. In the matrix isolation experiment, the intensities of the 2νsOHb and νOHf + νOHb bands decrease by more than a half. We do not see this in our calculations. 3.8. Deuterated Species. There are six unequivalent hydrogen atoms in the minimum energy structure of the water trimer. If we replace 1−6 of these with deuterium atoms, we obtain 63 different isotopologues: 6 × (H2O)2HDO, 12 × H2O(HDO)2, 3 × (H2O)2D2O, 9 × (HDO)3, 12 × H2OHDOD2O, 12 × (HDO)2D2O, 3 × H2O(D2O)2, 6 × HDO(D2O)2, and (D2O)3. The isotopologues were studied experimentally in ref 10. In the matrix isolation conditions (T < 10 K), the partially deuterated monomer units are usually oriented so that the deuterium atom forms the hydrogen bond.10 However, the energy difference between the two isomers (for instance, (H2O)2DbHfO and (H2O)2HbDfO) is less than 100 cm−1. Therefore, when the temperature is increased, for instance, to 50 K, there is a significant amount of both isomers present. We have calculated the vibrational spectra of all 63 isotopologues using the method described earlier. The transitions that do not contain bonded OH stretching character remain almost unchanged from the uncoupled model.32 Therefore, we do not present these kinds of transitions here. The full listing of the calculated transition wavenumbers and intensities of the transititons that include OHb or ODb stretching character is presented in the Supporting Information. To analyze the results, we have divided the isotopologues into four groups according to the masses of the bonded H/D atoms. The isotopologues belonging to the first group have three bonded and 0−3 free H atoms. The three OHb stretches are strongly coupled, and they form a band structure that was described earlier in the case of (H2O)3: a weak symmetric and two strong asymmetric fundamental modes. When free D atoms are introduced to the system, the OHb transition energies increase slightly. There is also a small increase in most of the transition intensities. There is one bonded D atom and two bonded H atoms in the group two isotopologues. The coupling between the OHb and ODb stretches is weaker than between OHb and OHb. Therefore, the stretching modes can be divided into one ODb and two OHb stretches. Despite the weakness of the coupling, the ODb stretches are red-shifted from the uncoupled model. For instance, the fundamental ODb stretching energies for (H2O)2D2O are 2524, 2528, and 2538 cm−1 in this work and 2539, 2543, and 2552 cm−1 in ref 32. In the overtones, the redshift increases harmonically (15 cm−1 → 30 cm−1). Similar shifts are found in the modes where the vibration is isolated in one monomer unit, for example, νOD + δDOD. The two OH stretching modes are a symmetric and an asymmetric combination of the OH oscillators. In the fundamental OH stretching region, their splitting is ≈50 cm−1. The transition to the asymmetric state is about 3 times stronger than to the symmetric one. However, the νsOHb + δ transitions are stronger than the νaOHb + δ ones. The third group (2 × Db + 1 × Hb) would resemble closely the previous group if the labels H and D were reversed. There are two OD stretches and one OH stretch. There are three bonded deuterium atoms in the molecules belonging to the last group. The results are qualitatively similar with the first group. Finally, we compare our results with a spectroscopic experiment. Recently, Ceponkus et al. have measured the bonded OH/D stretching fundamentals of various isotopologues of the water trimer in the Ne matrix.10 For each peak in the observed spectrum, we calculate multiple energy levels. We have calculated the weighted averages for each mode reported in ref 10. We assume that all partially deuterated species are oriented so that a deuterium bond is formed. The comparison 5372

dx.doi.org/10.1021/jp3017584 | J. Phys. Chem. A 2012, 116, 5368−5374

The Journal of Physical Chemistry A

Article

arguments, we suggest that the unobserved symmetric ODb state of (DHO)3 should be found at 2565 cm−1. 3.9. Comparison with Other Computational Work. A few computational studies of the high-frequency fundamentals of the water trimer have been published recently. Both Wang and Bowman34 and Watanabe et al.31 calculate three HOH bending fundamentals: two transitions close to 1600 cm−1 and a third one at approximately 1625 cm−1. In our case, the splitting is smaller (1610 and 1615 cm−1). We tried to improve our model by adding the PES terms that couple the bending modes. However, it did not make a difference in the results and, therefore, we decided to ignore these terms. Only a single band at 1608.3 cm−1 was observed in a Ne matrix experiment.9 However, the gas-phase observations show two bands with similar intensities (2:1) located at 1609 and 1638 cm−1.21 Wang does not report the transition intensities. The calculations by Watanabe predict that the intensity ratio between the low- and high-frequency modes is 40:1. Therefore, the interpretation of the higher energy band is not fully clear. As mentioned earlier, our calculated asymmetric OHb stretching fundamentals are approximately 65 cm−1 lower than the experimental ones.9 Watanabe et al. are able to calculate these modes within a few wavenumbers from the Ne matrix experimental values, and Wang’s numbers are approximately 15 cm−1 lower. We are unable to reproduce this, but the splitting between our calculated symmetric and asymmetric modes is in a good agreement with the experiment. The high-frequency spectrum of the water dimer has been calculated earlier using a model which is similar with the one used in this work excluding the intermolecular coupling.40 The bonded OH stretching motion was found to be lower than the experimental one by 31 cm−1, whereas the other fundamentals were calculated within 7 cm−1 from the experimental values.46 In the second overtone region, the difference was 105 cm−1. This difference was reduced to 25 cm−1 by fitting some of the one-dimensional PES parameters using observed energies from the fundamental and the first overtone region. Therefore, we assume that the error lies mostly in the description of the OHb stretching coordinate and could be fixed by coupling the high- and lowfrequency modes.47 However, this is not feasible in our current model. Moreover, the aim of this study was not to solve the problem involving the overestimated redshift in the bonded OH stretching mode but to improve the shape of the spectral band. Hardly any split has been either observed9 or calculated31,34 in the free OH stretching modes. As discussed earlier, our current calculations predict these modes in a slightly lower energy, but we expect to improve our values by using more advanced ab initio calculations.

is presented in Table 5. Again, our model overestimates the redshift in the OHb/Db stretching energies. However, our Table 5. Comparison of the Calculated Band Centers and Experimental OHb(Db) Stretching Energies (in cm−1)a this work

expt10

expt−comp

2597

39

2561 2591 2595

54 38 33

2555 2587 2594

50 38 35

2568b 2585

65 38

(DbHfO)3 νsODb νaODb νaODb

2510 2558 2565

νsODb νaODb νaODb

2507 2553 2562

νsODb νaODb νaODb

2505 2549 2559

νsODb νaODb νaODb

2503 2547 2553

νsODb νaODb νOHb

2524 2560 3440

(DbHfO)2D2O

DbHfO(D2O)2

(D2O)3

H2O(DbHfO)2

νsODb νaODb νOHb νsODb νaODb νOHb νODb νsOHb νaOHb νODb νsOHb νaOHb

2577 2599 3508 H2ODbHfOD2O 2519 2567b 2554 2592b 3440 H2O(D2O)2 2516 2561 2548 2585 3440 3502 (H2O)2DbHfO 2539 2589 3413 3490 3461 3522 (H2O)2D2O 2529 2574 3413 3488 3461 3521

53 39 68 48 38

45 37 62 50 77 61 45 75 60

a

The band centers are calculated as the weighted averages of all the transition energies of the isomers belonging to the given group. The transition intensities are used as the weighing coefficients. bThe average of two observed bands.

calculated results are consistent with the observations. The difference between the observed and the calculated remains almost constant. The best agreement is obtained in the case of asymmetric OD stretches. In the seven studied cases, the difference varies between 32 and 39 cm−1. Our model seems to slightly overestimate the splitting between the symmetric and asymmetric energy levels. In most cases of νsODb, our calculated energy is about 50 cm−1 smaller than the observed one. One exception is (D2O)3, where the difference is 65 cm−1. On the basis of our calculations, the energy should decrease when the free hydrogen atoms are replaced with the deuterium atoms from (DHO)3 to (D2O)3. Hence, we suggest that the symmetric ODb stretching state of (D2O)3 should be found at 2550 cm−1. However, the intensity of the transition is small which makes the observation difficult. On the basis of similar

4. CONCLUSIONS AND SUMMARY We have calculated the vibrational high-frequency spectrum of the water trimer using a nine-dimensional model. In general, our earlier model succeeded well in reproducing most of the observed features of the water trimer IR spectrum. The aim of this work was to improve the model so that it would reproduce the observed band structure for the hydrogen-bonded OH stretching modes. We introduced the intramolecular coupling to our earlier model. In general, the band positions remain basically unchanged when the coupling terms are added to the potential energy operator. However, the point of this work was not to improve the description of the band centers but to obtain the experimentally detected splittings in the hydrogenbonded OH stretching bands. The new model succeeds in this 5373

dx.doi.org/10.1021/jp3017584 | J. Phys. Chem. A 2012, 116, 5368−5374

The Journal of Physical Chemistry A

Article

(20) Brown, M. G.; Viant, M. R.; McLaughlin, R. P.; Keoshian, C. J.; Michael, E.; Cruzan, J. D.; Saykally, R. J.; van der Avoird, A. J. Chem. Phys. 1999, 111, 7789−7800. (21) Paul, J. B.; Provencal, R. A.; Chapo, C.; Petterson, A.; Saykally, R. J. J. Phys. Chem. A 1999, 103, 2972−2974. (22) Schütz, M.; Bürgl, T.; Leutwyler, S.; Bürgl, H. B. J. Chem. Phys. 1993, 99, 5228−5238. (23) Estrin, D. A.; Paglieri, L.; Corongiu, G.; Clementi, E. J. Phys. Chem. 1996, 100, 8701−8711. (24) Low, G. R.; Kjaergaard, H. G. J. Chem. Phys. 1999, 110, 9104− 9115. (25) Milet, A.; Moszynski, R.; Wormer, P. E. S.; van der Avoird, A. J. Phys. Chem. A 2003, 103, 6811−6819. (26) Taketsugu, T.; Wales, D. J. Mol. Phys. 2002, 100, 2793−2806. (27) Anderson, J. A.; Crager, K.; Fedoroff, L.; Tschumper, G. S. J. Chem. Phys. 2004, 121, 11023−11029. (28) Dunn, M. E.; Evans, T. M.; Kirschner, A.; Shields, G. C. J. Phys. Chem. A 2006, 110, 303−309. (29) van der Avoird, A.; Szalewicz, K. J. Chem. Phys. 2008, 128, 014302. (30) Wang, Y.; Carter, S.; Braams, B. J.; Bowman, J. M. J. Chem. Phys. 2008, 128, 071101. (31) Watanabe, Y.; Maeda, S.; Ohno, K. J. Chem. Phys. 2008, 129, 074315. (32) Salmi, T.; Kjaergaard, H. G.; Halonen, L. J. Phys. Chem. A 2009, 113, 9124−9132. (33) Richardson, J. O.; Althorpe, S. C.; Wales, D. J. J. Chem. Phys. 2011, 135, 124109. (34) Wang, Y.; Bowman, J. M. J. Chem. Phys. 2011, 134, 154510. (35) Boys, S. F.; Bernardi, F. Mol. Phys. 1970, 19, 553−566. (36) Valiron, P.; Mayer, I. Chem. Phys. Lett. 1997, 275, 46−55. (37) Dunning, T. H., Jr. J. Chem. Phys. 1989, 90, 1007−1023. (38) Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J. J. Chem. Phys. 1992, 96, 6796−6806. (39) Lane, J. R.; Kjaergaard, H. G. J. Chem. Phys. 2009, 131, 034307. (40) Salmi, T.; Hänninen, V.; Garden, A. E.; Kjaergaard, H. G.; Tennyson, J.; Halonen, L. J. Phys. Chem. A 2008, 112, 6305−6312; see erratum ref 41. (41) Salmi, T.; Hänninen, V.; Garden, A. E.; Kjaergaard, H. G.; Tennyson, J.; Halonen, L. J. Phys. Chem. A 2012, 116, 796−797. (42) Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M.; Celani, P.; Korona, T.; Lindh, R.; Mitrushenkov, A.; Rauhut, G.; et al. Molpro, version 2010.1, a package of ab initio programs, 2010; see http://www.molpro.net. (43) Atkins, P. V. Molecular Quantum Mechanics, 2nd ed.; Oxford University Press: Oxford, 1983. (44) Sälli, E.; Salmi, T.; Halonen, L. J. Phys. Chem. A 2011, 115, 11594−11605. (45) Kjaergaard, H. G.; Low, G. R.; Robinson, T. W.; Howard, D. L. J. Phys. Chem. A 2002, 106, 8955−8962. (46) Bouteiller, Y.; Tremblay, B.; Perchard, J. P. Chem. Phys. 2011, 386, 29−40. (47) Kjaergaard, H. G.; Garden, A. L.; Chaban, G. M.; Gerber, R. B.; Matthews, D. A.; Stanton, J. F. J. Phys. Chem. A 2008, 112, 4324−4335.

as the bonded OH stretching band profile gets a correct shape after the intermolecular coupling is added. Our model is also able to describe and even to predict some of the band positions for the deuterated water trimer species. We believe that we have made a significant improvement to our earlier work by including a potential energy coupling between the water monomer units in the water trimer. In spite of the substantial increase of the size of the Hamiltonian matrix, we have been able to obtain fully converged eigenvalues for the high-frequency mode fundamentals and first overtones.



ASSOCIATED CONTENT

S Supporting Information *

The calculated vibrational energy levels for the fully and partially deuterated isotopologues of the water trimer. This material is available free of charge via the Internet at http:// pubs.acs.org.



AUTHOR INFORMATION

Corresponding Author

*E-mail: lauri.halonen@helsinki.fi. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS We thank the Academy of Finland and the University of Helsinki for financial support. We also thank CSC Scientific Computing Ltd. for providing computing time.



REFERENCES

(1) Ayers, P. G.; Pullin, A. D. E. Spectrochim. Acta 1975, 32A, 1629− 1639. (2) Ayers, P. G.; Pullin, A. D. E. Spectrochim. Acta 1975, 32A, 1641− 1650. (3) Ayers, P. G.; Pullin, A. D. E. Spectrochim. Acta 1975, 32A, 1689− 1693. (4) Ayers, P. G.; Pullin, A. D. E. Spectrochim. Acta 1975, 32A, 1695− 1704. (5) Perchard, J. P. Chem. Phys. 2001, 266, 109−124. (6) Perchard, J. P. Chem. Phys. 2001, 273, 217−233. (7) Engdahl, A.; Nelander, B. J. Mol. Struct. 1989, 193, 101−109. (8) Ceponkus, J.; Uvdal, P.; Nelander, B. J. Phys. Chem. A 2008, 112, 3921−3926. (9) Tremblay, B.; Madebène, B.; Alikhani, M. E.; Perchard, J. P. Chem. Phys. 2010, 378, 27−36. (10) Ceponkus, J.; Uvdal, P.; Nelander, B. J. Chem. Phys. 2011, 134, 064309. (11) Fröchtenicht, R.; Kaloudis, M.; Koch, M.; Huisken, F. J. Chem. Phys. 1996, 105, 6128−6140. (12) Burnham, C. J.; Xantheas, S. S.; Miller, M. A.; Applegate, B. E.; Miller, R. E. J. Chem. Phys. 2002, 117, 1109−1122. (13) Slipchenko, M. N.; Kuyanov, K. E.; Sartakov, B. G.; Vilesov, A. F. J. Chem. Phys. 2006, 124, 241101. (14) Liu, K.; Brown, M. G.; Viant, M. R.; Cruzan, J. D.; Saykally, R. J. Mol. Phys. 1996, 89, 1373−1396. (15) Huisken, F.; Kaloudis, M.; Kulcke, A. J. Chem. Phys. 1996, 104, 17−25. (16) Paul, J. B.; Collier, C. P.; Saykally, R. J.; Scherer, J. J.; O’Keefe, A. J. Phys. Chem. A 1997, 101, 5211−5214. (17) Viant, M. R.; Cruzan, J. D.; Lucas, D. D.; Brown, M. G.; Liu, K.; Saykally, R. J. J. Phys. Chem. A 1997, 101, 9032−9041. (18) Paul, J. B.; Provencal, R. A.; Petterson, A.; Saykally, R. J. J. Chem. Phys. 1998, 109, 10201−10206. (19) Viant, M. R.; Brown, M. G.; Cruzan, J. D.; Saykally, R. J.; Geleijns, M.; van der Avoird, A. J. Chem. Phys. 1999, 110, 4369−4380. 5374

dx.doi.org/10.1021/jp3017584 | J. Phys. Chem. A 2012, 116, 5368−5374