Local-Monomer Calculations of the Intramolecular IR Spectra of the

Jul 10, 2014 - ABSTRACT: Dilute mixtures of HOD in pure H2O and D2O ices and liquid have been used by experimentalists to focus on the spectrum and ...
0 downloads 0 Views 310KB Size
Download PDF
Article pubs.acs.org/JPCB

Local-Monomer Calculations of the Intramolecular IR Spectra of the Cage and Prism Isomers of HOD(D2O)5 and HOD and D2O Ice Ih Hanchao Liu, Yimin Wang, and Joel M. Bowman* Cherry L. Emerson Center for Scientific Computation and Department of Chemistry, Emory University, Atlanta, Georgia 30322, United States ABSTRACT: Dilute mixtures of HOD in pure H2O and D2O ices and liquid have been used by experimentalists to focus on the spectrum and vibrational dynamics of the local OH and OD stretches and bend of HOD in these complex and highly heterogeneous environments. The hexamer version of the mixture is HOD(D2O)5. The cage isomer of this cluster was recently studied and analyzed theoretically using local-mode calculations of the IR spectrum by Skinner and co-workers. This and the further possibility of experimental investigation of this cluster have stimulated us to study HOD(D2O)5 using the three-mode, localmonomer model, with the ab initio WHBB dipole moment and potential energy surfaces. Both the cage and prism isomers of this cluster are considered. In addition to providing additional insight into the HOD portion of the spectrum, the spectral signatures of each D2O are also presented in the range of 1000−4000 cm−1. The OH stretch bands of both the prism and cage isotopomers exhibit rich structures in the range of 3100−3700 cm−1 that are indicative of the position of the HOD in these hexamers. A preliminary investigation of the site preference of the HOD is also reported for both cage and prism HOD(D2O)5 using an harmonic zero-point energy analysis of the entire cluster. This indicates that the energies of free-OH sites are lower than the ones of H-bonded OH sites. Finally, following our earlier work on the IR spectra of H2O ice models, we present IR spectra of pure D2O and HOD.

1. INTRODUCTION The water hexamer has attracted the interest of chemists for more than 2 decades because it is the “smallest droplet” of water.1 The several low-lying isomers with only small differences in energy and entropy have stimulated a series of state-of-the-art spectroscopic and theoretical studies. The seminal work using vibration−rotation tunneling spectra has identified the cage as the most stable structure at 6 K.1 Rare gas tagging vibrational action spectra measured at higher temperature show that the book isomer is dominant.2,3 In helium nanodroplets, the isomer form seen is the higher-energy cyclic ring, evidently because the rapid cooling of monomers in the helium surrounding prohibits the cyclic structure from isomerizing to the low-lying forms.4,5 Very recently, the prism isomer of (H2O)6 has been observed for the first time from broad-band rotational spectroscopy;6 however, the cage isomer was reported as the dominant isomer, in agreement with previous microwave experiments.1 Theoretically, coupled-cluster calculations with single, double, and perturbative triple excitations, CCSD(T), using complete-basis-set extrapolation, shows the that prism isomer has the lowest electronic energy.7−10 The energy of the cage isomer is approximately 0.25 kcal/mol higher than the prism. A recent study using the parametric two-electron reduced density matrix (2-RDM) method shows that the cage isomer has a 0.07 kcal/mol higher electronic energy than the prism isomer.11 The © XXXX American Chemical Society

small electronic energy difference between isomers can be compensated for with the zero-point energy (ZPE) and at nonzero temperature entropy, resulting in a different freeenergy ordering at different temperatures. At 0 K, including ZPE, full-dimensional diffusion Monte Carlo (DMC) calculations using the ab initio WHBB water potential show that the prism and cage are nearly isoenergetic.12 Replica-exchange path-integral molecular dynamics (RE-PIMD) shows that the fraction of cage versus prism isomers increases as the temperature rises.12 Up to 60 K, the book isomer has a negligible population but becomes dominant when the temperature is above 150 K due to the its floppier structure and thus higher entropy.12 A recent PIMD simulation by Babin and Paesani,13 using WHBB and their new related ab initio PES, denoted HBB2-pol, shows that for (D2O)6, the prism is dominantly populated at temperatures below at least 30 K. This important finding, which is understandable based on the diminished effect of ZPE, will hopefully stimulate new experiments on the isotopologues. Special Issue: Spectroscopy of Nano- and Biomaterials Symposium Received: June 19, 2014 Revised: July 8, 2014

A

dx.doi.org/10.1021/jp5061182 | J. Phys. Chem. B XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry B

Article

There have been a number of theoretical simulations of the IR spectra of hexamer clusters, most recently by Tainter and Skinner14 and Wang and Bowman.15 The latter considered the overtone of the monomer bend and showed its importance in the “OH stretch” portion of the spectrum. With inclusion of these bands, agreement between theory and experiment for the cyclic ring and book isomers is very good. Dilute mixtures of HOD in D2O and H2O have been used to decouple and thereby isolate the spectroscopy and vibrational relaxation dynamics of the OD and OH stretches in HOD; see, for example, refs 16−18 and references therein for liquid and ref 19 and references therein for ice. Simulations of the vibrational spectra of dilute HOD in water systems have been reported.20−23 We recently reported quantum calculations of the vibrational relaxation pathways of the OH and OD stretches for dilute HOD in D2O and H2O ice Ih.24 In recent work by Skinner and co-workers, “dilute” HOD in the water hexamer cage was investigated.23 Specifically, the IR spectrum of the cage isotopologue of HOD(D2O)5 in the range of 3000−3800 cm−1 was calculated using that group’s combined quantum and classical approach, with each OH stretch treated as a local mode. The OH stretch frequencies of the 12 unique sites for HOD in this cluster were shown to be correlated with the hydrogen-bonding class, based on both the donor and acceptor H-bonding environment. By contrast, for HOD in liquid and ice D2O, the calculated and experimental spectra show two seemingly simple broad peaks centered on the OD stretch and the OH stretch. Stimulated by this work and also by the recent prediction of Babin and Paesani13 that the prism form of (D2O)6 is lower in free energy at low temperatures than the cage isomer, we consider both the prism and cage isomers of HOD(D2O)5 as it is likely that both are significantly populated at low temperatures. We also extend the spectral range down to 1000 cm−1 in order to encompass all of the monomer intramolecular spectral features. In addition, we examine the energetic site dependence of the HOD in these clusters to shed light on whether there is a possible significant nonuniform distribution among the 12 sites. The answer to this question may be relevant to the site distribution of dilute HOD in the liquid and ice environments mentioned above. Finally, following our earlier calculations on the IR spectrum of pure H2O ice Ih,25 we present the analogous spectra for pure D2O and pure HOD ice Ih.

In this work, as noted already, we consider the prism and cage isomer and use the structures from WHBB, which were reported in ref 26. We systematically hydrogenate the 12 unique positions of the cage and prism (D2O)6 and therefore consider 12 HOD(D2O)5 isotopomers for both the cage and prism. The Local-Monomer Model30 (LMon) is used to calculate the IR spectra. In this approach, a normal-mode analysis is done for each monomer, with all other monomers fixed at a reference configuration. In the present application, this is done at the minimum of the cage and prism configuration. Therefore, for monomer m in the cluster, the 9 × 9 mass-scaled Hessian is determined in Cartesian coordinates of that monomer. The potential is the full WHBB potential; however, only the coordinates of the mth monomer are varied. The diagonalization of the Hessian is done to get the “local” normal modes and harmonic frequencies of each monomer. These modes are used to obtain the coupled anharmonic vibrational eigenvalues and eigenfunctions of the Watson Hamiltonian for each monomer. Specifically, the Schrödinger equation for monomer m [Tm̂ + Um(Q m) − Ei(m)]ϕi(m)(Q m) = 0

where Qm are the local normal modes, T̂ m is the full kinetic energy operator, and Um is the potential, is solved using MULTIMODE.31−33 In the present calculations, we dynamically couple the three high-frequency local modes, that is, the monomer bend and two stretches. The eigenfunctions are then used to obtain the dipole transition matrix and IR intensities for each monomer (see below for more details). This “minimalist” three-mode approach has been used successfully to obtain the IR spectra of various isomers of (H2O)6 (see refs 15 and 30) as well as numerous other water clusters, up to the 20-mer,30 and for 192-mer models of ice.25 Inclusion of intermolecular modes and their effect on the high-frequency ones can easily be included in these calculations. This was done in a recent application where six-mode calculations were done to model the vibrational relaxation dynamics of an HOD monomer in ice Ih.24 (It is perhaps worth noting for the isolated HOD molecule that the LMon treatment is exact. The OH and OD stretch frequencies and excellent agreement with experiment were shown in ref 24.) For the present application to HOD(D2O)5, the LMon calculation is straightforwardly performed for each monomer for the cage and prism isomers of HOD(D2O)5. For each isomer, the LMon calculations are done for all 12 isotopomers, with respect to the different H position. The details of the vibrational self-consistent field (VSCF)/virtual-state configuration interaction (VCI) calculations are, in brief, the following. For each mode, 6 numerical basis functions, contracted from 9, 16, and 16 harmonic oscillator functions for the bend, OD, and OH, respectively, are used. These contracted functions are eigenfunctions of 1d Hamiltonians in each mode. From these, optimized “HEG” quadrature grids are obtained; we used 8, 10, and 10 optimized quadrature points for the bend, OD stretch, and OH stretch modes, respectively. Numerical quadrature of matrix elements is done with these quadrature points. The VSCF equations are solved for the ground vibrational state, and from each effective single-mode Hamiltonian, a virtual excitedstate basis is obtained and used in the VCI calculations. The excitation space for each mode goes up to a maximum of 5 quanta. Then, a strategy is applied to restrict the excitation space, such that the maximum sum of excitation in the one-

2. THEORY AND COMPUTATIONAL DETAILS The ab initio WHBB potential energy surface (PES) and dipole moment surface (DMS) are used to calculate the IR spectrum of the water hexamers. The PES and DMS are described in detail in ref 26. To briefly summarize here, the PES is a manybody expansion of one-body, two-body, three-body, and fourand higher-body terms. The one-body potential uses a spectroscopically accurate PES. The two-body potential is a permutationally invariant fit of roughly 30 000 CCSD(T)/augcc-pVTZ energy points. The three-body potential is fitted by roughly 40 000 MP2/aug-cc-pVTZ points. The four- and higher-body terms include the long-range polarization interactions. The PES has been tested and applied broadly. This includes the successful prediction of De from the water dimer to the 22-mers,26 the accurate benchmark calculated D0 for the dimer,27 the trimer,28,29 and several low-lying hexamers.12 The DMS is represented by one-body and two-body terms, fitted from about 30 000 MP2/aug-cc-pVTZ dipole moment calculations.26 B

dx.doi.org/10.1021/jp5061182 | J. Phys. Chem. B XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry B

Article

mode, two-mode, and three-mode representations of the potential31 is 5, 8, and 5, respectively. The final size of the VCI Hamiltonian matrix is small, 92, but as verified by larger calculations, this gives frequencies of the bend fundamental and overtone and the OD and OH stretch fundamentals that are converged to within a couple wavenumbers or less. The procedure of calculating the infrared intensity is standard. First, the dipole matrix elements between vibrational eigenstates ϕi(m) (Qm) are evaluated. Then, the intensity is given by I(vi → j) = vi → j

∑

|R αij|2

α=x ,y,z

where the sum is over the three components of the transition dipole vector. Rαij is the αth vector component of the dipole matrix element from state i to state j. These matrix elements are evaluated numerically using the vibrational eigenfunctions and the WHBB DMS. In the present calculations, the initial state is the ground vibrational state, and the resulting IR spectra from this procedure are a series of “sticks”. Thus, this spectrum corresponds to 0 K and also ignores any relaxation of these eigenstates due to coupling to other modes in the cluster. In recognition of the relaxation, we apply a Lorentzian line shape function to each stick. The full width at half-maximum (fwhm) is chosen to be 5 cm−1, which corresponds to a reasonable guesstimate of the lifetime on the order of a picosecond, which is based on lifetimes reported for HOD in ice; see ref 24 and references cited therein. The resulting spectra should look a bit more like experimental ones at low temperature, where rotational broadening of the band is quite narrow. Stick spectra are also shown in separate figures just in the OH stretch region to assign each feature to a site and the nature of the H-bonding. Finally, to estimate the H atom site energetics, harmonic ZPE analysis is done for each of the 12 isotopomers. This was done by performing a full-dimensional, that is, 6 × 9 = 54, normal-mode analysis for each isotopomer from which the harmonic ZPE is then trivial to obtain for each isotopomer. The full normal-mode analysis is efficient using the PES, taking roughly 5 s per isotopomer on a single CPU of a 5 year old workstation. It is perhaps worth noting that it is not necessary to diagonalize the Hessian to get the total harmonic ZPE, as reported previously by Higgs.34 We did not use this interesting approach, however.

Figure 1. Structure and atom numbering for the cage (top) and prism (bottom) isomers of the water hexamer.

3. RESULTS AND DISCUSSION 3.1. IR Spectra of Prism and Cage HOD(D2O)5. The structures and atom number assignments of the cage and prism hexamers are shown in Figure 1. As noted, there are 12 unique hydrogen sites for both cage and prism isomers; therefore, there are 12 isotopomers of HOD(D2O)5. These are numbered in the figure in ascending order by the OH stretch frequency, based on our calculations. The complete IR spectrum for all 12 isotopomers of the cage HOD(D2O)5 in the range of 1000− 4000 cm−1 from LMon calculations is shown in Figure 2. The contributions from each “dopant” HOD and the five “solvent” D2O are also shown in the bottom two panels. From the bottom up, the spectrum of the HOD monomers have five bands: the HOD bends clustered at roughly 1500 cm−1, the OD stretches in the range of 2300−2750 cm−1, the bend overtone at roughly 3000 cm−1, the OH stretches in the range of 3150− 3750 cm−1, including weak features at around 3700 cm−1 due to the HOD bend plus the H-bonded OD stretch. For D2O, the

Figure 2. IR spectra of cage HOD(D2O)5.

bending modes are at roughly 1200 cm−1, and the bands in the range of 2300−2750 cm−1 include the OD stretches and small contributions from the overtone of the D2O bend in the range of 2350−2410 cm−1. The total spectrum of HOD(D2O)5, obtained by combining the bands of both the HOD dopant and D2O solvent, is shown in the top panel of Figure 2. Notice that the intensities for D2O bend and OD stretch bands are larger C

dx.doi.org/10.1021/jp5061182 | J. Phys. Chem. B XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry B

Article

hydrogen sites, spanning a range from 3150 to 3700 cm−1. (Note the HOD bend overtone is below 3000 cm−1 and not shown in Figure 4, unlike the case of pure H2O where the bend overtone is embedded in the OH stretch region of the spectrum.25) As seen, each stick is color-coded according to the well-known donor (D)/acceptor (A) classification.2 We discuss this in detail below, but first, we make a more “coarse-grained” observation about these stick spectra. It is based on identification of two types of monomers. One is with one free-OH(D) and one H-bonded OH(D). There are four such monomers in the cage and three in the prism. The second monomer type is with no free-OH(D) stretches and thus with two H-bonded OH(D) stretches. There are two of these monomers in the cage and three in the prism. We now note that for the first type of monomer, that is, with one freeOH(D), the two OH stretches’ spectral sticks are at the extremes of the spectral range. The free-OH stretches, shown in blue, are the highest-frequency bands in the spectrum, and the H-bonded OH stretches, shown in red and orange, are the lowest-frequency bands of the spectrum. The OH bands for the second type of monomer, which has two H-bonded OH stretches, shown in green, are between the free-OH bands and the lowest-frequency H-bonded OH bands, which are assigned to the first type of monomer. (Note that the same correlation of OH stretch frequency and monomer type is seen in the corresponding pure (H2O)6.) To go beyond this coarse-grained analysis, the D/A classification is used. Recall that a monomer with two donor H atoms and two acceptor H atoms is labeled as DDAA. We use this notation to indicate only the H-bonded OH stretches and continue to use the term “free” to indicate free-OHs. This labeling is used in Table 1, where each labeled spectral line is

than those of the corresponding bands of HOD monomers; this is mainly due to the ratio of dopant and solvent, which is 1:5. Note that to obtain the overall spectrum, an equal population of the 12 isotopomers of HOD(D2O)5 is assumed. In section 3.3, we discuss this in detail. The corresponding spectra of the prism HOD(D2O)5 are shown in Figure 3. Compared with the spectra of the cage in

Figure 3. IR spectra of prism HOD(D2O)5.

Figure 2, the prism spectra show the same bands in similar frequency ranges. However, a noticeable difference from the cage spectrum is seen in the structure of the OH stretch portion of the spectrum. Detailed discussion of these OH stretch bands of both isomers will be presented below. The same assumption that populations of the 12 isotopomers are equal is applied to the prism as well. Next, we focus on the OH stretch region of the IR spectra for the cage and prism HOD(D2O)5. These bands are of most interest because the OD stretch bands in HOD are largely masked by the D2O solvent stretch bands. To analyze this portion of the spectra, we present the stick spectra of the cage and the prism from LMon calculations in Figure 4. There are 12 sticks for each isomer, which correspond to the 12 unique

Table 1. Monomer Classification of the 12 HOD(D2O)5 Isotopomers of the Cage and Prism Isomersa hydrogen site

cage

prism

1 2 3 4 5 6 7 8 9 10 11 12

DAA DAA DA DA DDA DDA DDA DDA free free free free

DAA DAA DAA DDA DDA DDA DDA DDA DDA free free free

a

Hydrogen sites follow the atomic numbering in Figure 1. For H-bond type, D denotes one donor, and A denotes one acceptor.

identified by the site numbering shown in Figure 1. For the cage isomer, it is straightforward to see that the H-bonded OH stretches on the first type of monomer are either DAA or DA. The H-bonded OH stretches on the second type of monomer are DDA. We see that the DAA OH stretches, shown as red sticks, are more red-shifted than the DA ones, shown by orange sticks. For the prism, all H-bonded OH stretches on the first type of monomer are DAA, and the ones on the second type are DDA. To summarize, the OH stretch frequencies are in the order of DAA < DA < DDA < free. This order, seen in the cage and prism isomers, is consistent with the results of the book

Figure 4. OH stretch band of IR spectra of the cage and prism HOD(D2O)5. D

dx.doi.org/10.1021/jp5061182 | J. Phys. Chem. B XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry B

Article

isomer in pure (H2O)6.2 Finally, it is worthwhile to note that Ohno et al.35 and later Tainter et al.23 analyzed the OH frequencies by considering the H-bond configurations not only on the donor monomer but also on the acceptor monomer, using an “M index”. This provides additional insights into the pattern of spectral features and also in making connections with the condensed phase and surfaces, and we refer the reader to those papers for this insightful and detailed analysis. Next, we compare our numerical results for the cage isomer with previous local-mode calculations of the OH stretch fundamental of Skinner and co-workers23 in Table 2. For

fundamental to the OH stretch and bend overtone region.25 These ices were modeled by a cluster of 192 monomers with a focus on a core subcluster of 105 monomers, with no free-OH stretches.25 The initial configuration of the proton-disordered ice Ih 192mer was taken from the paper by Hayward and Reimers (H−R).36 These configurations were then optimized using the WHBB potential.25 (We discuss the “optimization” further below.) Analogous calculations were done for the IR spectra of HOD and D2O ice Ih, and the results are shown in Figure 5a and b, along with experimental results for dilute

Table 2. OH Stretch Frequencies of the 12 Sites of the Cage and Prism HOD(D2O)5 along with the OH Stretch Frequencies of the Cage HOD(D2O)5 with the Local-Mode Calculations from Reference 23a cage

prism

peak

this work

ref 23

this work

1 2 3 4 5 6 7 8 9 10 11 12

3163 3279 3338 3439 3466 3536 3537 3542 3702 3703 3719 3726

3148 3274 3330 3392 3421 3517 3547 3576 3709 3714 3727 3738

3161 3251 3386 3446 3514 3539 3558 3582 3617 3694 3699 3738

The frequencies from ref 23 were determined by us from digitized versions (Figure 1 in that paper).

Figure 5. IR spectra of HOD and D2O ice Ih, compared with the experiment. Panel (a) is the optimized HOD ice Ih. Panel (b) is the optimized D2O ice Ih. Panel (c) is the experimental spectrum at 210 K, obtained from ref 37. Panel (d) is the H−R HOD ice Ih. Panel (e) is the H−R D2O ice Ih. See the text for details.

completeness, we also give these frequencies for the Prism isotopolog. For consistency, we show the current OH stretch fundamentals; however, we remind the reader that in this energy range, there are also combination states of the HOD bend plus H-bonded OD stretch. As seen, the two calculations are consistent with each other, but there are some differences too. First, peaks 4 and 5 from our calculations are blue-shifted by 40 cm−1 compared with the local-mode ones. Second, in our calculations, sticks 6−8 are very close in frequencies, whereas the local-mode peaks are more dispersed. Overall, it is gratifying to see good agreement between these two calculations given that there are significant differences in their details. Perhaps the biggest difference is the treatment of the electronic structure, which in the local-mode calculations is DFT-B3LYP/6-311+ +G** for the OH stretch and a semiempirical three-body interaction in the LMon ones, using WHHB, which is based on CCSD(T)/aug-cc-pVTZ for the intrinsic two-body and MP2/ aug-cc-pVTZ for the intrinsic three-body interactions. A significant point in common is the explicit consideration of the important three-body interactions. Finally, we compare the spectra of the cage and prism isomers. The most significant difference between the two isomers is perhaps in the region above 3500 cm−1. In this region, the spectrum of the prism is more dispersed than that of the cage. This suggests that IR spectroscopy should be able to clearly distinguish the cage and prism isomers of HOD(D2O)5. 3.2. IR Spectra of HOD and D2O Ice Ih. Previously, we reported LMon-WHBB calculations of the IR spectra of H2O ice Ih and an amorphous ice in the region of the bend

HOD in D2O at 210 K, shown in Figure 5c.37 We also present spectra from the H−R configurations for HOD and D2O in Figure 5d and e. All of the calculated spectra were obtained by smoothing the stick spectra in the same way as reported in ref 25. Consider first the calculated spectra at the previously adopted optimized structures in panels (a) and (b). For HOD ice, in Figure 5a, we see three bands, the OD stretch, the (weak) HOD bend overtone, and the OH stretch. Given as pairs of numbers, the frequency at the peak centers and the fwhm of the three bands are (2480 cm−1, 200 cm−1) for the OD stretch, (2872 cm−1, 196 cm−1) for the HOD bend overtone, and (3353 cm−1, 259 cm−1) for the OH stretch. For D2O ice, the spectrum shown in Figure 5b has one broad band in the range of 2150−2750 cm−1, with the main peak at 2476 cm−1 and a shoulder at 2389 cm−1 and fwhm of 303 cm−1. This band is mainly due to the OD stretch and also contains a small portion of intensities from the D2O bend overtone. The experimental spectrum at 210 K of dilute HOD in D2O ice Ih37 is shown in Figure 5c. This spectrum was digitized by us from Figure 1 in ref 37. In this spectrum, we see two major bands. One has a peak at approximately 2440 cm−1 with a shoulder at 2370 cm−1. The fwhm of this band is roughly 180 cm−1. The other band has a peak at about 3300 cm−1. Note that this band has significant intensity in the “wings”, and therefore, the usual fwhm estimate of the bandwidth is not an accurate characterization of the band. It is evident that the extent of this band, including the wings, would be comparable to the fwhm of

a

E

dx.doi.org/10.1021/jp5061182 | J. Phys. Chem. B XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry B

Article

Table 3. Harmonic ZPEs and ΔEH of the 12 HOD(D2O)5 Isotopomers of the Cage Hexamer, Relative to the Cage (D2O)6, along with the Harmonic Frequency for Each OHa

the calculated widths but somewhat closer to the width of the H−R HOD ice. To compare theory with experiment thus far, we first note that the peak positions of the two major bands agree well, that is, 2480 and 3353 cm−1 from theory and 2440 and 3300 cm−1 from experiment. However, the fwhm from theory is larger than that from experiment for both bands. To gain insight into this apparent discrepancy, now consider the calculated spectra using the H−R configuration, shown in Figure 5d and e. From the former panel, the peak and fwhm of the three bands are (2464 cm−1, 125 cm−1) for the OD stretch, (2885 cm−1, 109 cm−1) for (the weak) HOD bend overtone, and (3303 cm−1, 157 cm−1) for the OH stretch. For the latter, the D2O ice, the center of the peak is at 2468 cm−1, with a shoulder at 2357 cm−1. A rough estimate of the entire bandwidth is 150 cm−1. Thus, the major difference with the optimized spectra is in the widths of the bands, with the widths of the H−R configuration narrower and in better agreement with experiment. This apparent better agreement with experiment may be due to the moderate-size cluster used in the optimized simulations. The optimization of this cluster may introduce more structural disorder, especially for the region away from the center of the cluster, than would be found in a much larger cluster. Calculations with a larger cluster are considerably more computationally intensive and therefore have not been done. There are several additional caveats worth noting about the comparisons with experiment. First, as noted already, our calculations represent the IR spectra at 0 K, whereas experiment is done at 210 K. There is evidence that the change in temperature does not introduce major effects in the peak center and width. For example, based on the experimental IR spectra of pure H2O ice from 16 to 100 K, the peak center has a roughly 35 cm−1 blue shift, and the width changes from roughly 130 to 160 cm−1. (The spectra at these two temperatures are shown in ref 25.) Second, we notice that in HOD ice, the OH band is more intense than the OD band, whereas in experiment, the OD band is more intense than the OH band. This is because in such dilute HOD in D2O ice, the OD band is dominated by the D2O solvent monomers rather than the HOD dopant monomers. For the comparisons of dilute HOD in D2O solvent, it is important to note that the calculated spectra assume equal populations of all HOD sites in the cluster. If this is not the case, then the calculated spectral lines shapes could be affected significantly. A detailed investigation of this possible site dependence is difficult (and temperature-dependent) and will hopefully be done in the future. Finally, it is worthwhile to briefly compare the spectrum of HOD ice with HOD(D 2 O) 5 . Overall, there is good correspondence between the two. The hexamer’s spectral bandwidth is larger and obviously sparser than the bandwidth of the ice spectrum. This can be explained, roughly at least, by the more heterogeneous H-bonding environments in the two hexamer isomers, that is, DA, DAA, DDA, and free, compared to ice Ih, where nominally only DDAA monomers are present. 3.3. Harmonic ZPEs of Isotopomers of Prism and Cage HOD(D2O)5. It is of interest to investigate the site dependence of the ZPEs of the cage and prism isotopomers of HOD(D2O)5. These were obtained in the harmonic approximation, as described above, and the results are given in Tables 3 and 4, respectively. As seen, these range from 71.30 to 71.50 kcal/mol for the cage and 71.56 to 71.73 kcal/mol for the prism. Note that these energies are relative to the potential

hydrogen site

frequency (cm−1)

ZPE (kcal/mol)

ΔZPE (kcal/mol)

1 2 3 4 5 6 7 8 9 10 11 12

3382 3503 3549 3635 3671 3727 3732 3737 3893 3899 3906 3918

71.50 71.48 71.48 71.49 71.50 71.46 71.46 71.47 71.33 71.34 71.35 71.30

2.18 2.16 2.16 2.17 2.18 2.14 2.14 2.15 2.00 2.01 2.03 1.97

a

The site labels are given in Figure 1.

Table 4. Harmonic ZPEs and ΔEH of the 12 HOD(D2O)5 Isotopomers of the Cage Hexamer, Relative to the Prism (D2O)6, along with the Harmonic Frequency for Each OHa hydrogen site

frequency (cm−1)

ZPE (kcal/mol)

ΔZPE (kcal/mol)

1 2 3 4 5 6 7 8 9 10 11 12

3385 3473 3596 3646 3715 3730 3748 3776 3805 3891 3904 3921

71.73 71.72 71.71 71.73 71.69 71.69 71.70 71.67 71.66 71.56 71.58 71.59

2.18 2.16 2.15 2.18 2.13 2.13 2.14 2.12 2.11 2.01 2.02 2.03

a

The site labels are given in Figure 1.

minimum of the corresponding isomer. It is well-established that the prism minimum is below the cage one, as discussed in more detail below. From these ZPEs, we define the “hydrogenation energy” as the ZPE difference between the HOD(D2O)5 and the reference energies of cage or prism (D2O)6 ΔE H = E[HOD(D2 O)5 ] − E[(D2 O)6 ]

We see that this energy is on the order of 2 kcal/mol for any of the 12 sites in the cage or prism isomers. However, at the freeOD sites, ΔEH is about 0.15−0.2 kcal/mol lower than ΔEH at hydrogen-bonded positions, for both the cage and prism. This suggests that at 0 K to about roughly 100 K, hydrogenation is energetically favored at the free-OD positions. We must add an important caveat to this conclusion though. As is well-known, especially by us, harmonic ZPEs are not accurate to within 0.15−0.2 kcal/mol for a cluster this large; however, the estimate of the difference in harmonic ZPEs is probably at least of the right sign, based on cancellation of errors. Therefore, it should be a reasonable guide to the site preference of H atom substitution. Also note that this suggests that the relative intensities of IR spectral bands might be altered due to the nonequal populations of the 12 sites. A similar analysis could be done for dilute HOD in the condensed phase; however, that would be a fairly major F

dx.doi.org/10.1021/jp5061182 | J. Phys. Chem. B XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry B

Article

(7) Dahlke, E. E.; Olson, R. M.; Leverentz, H. R.; Truhlar, D. G. Assessment of the Accuracy of Density Functionals for Prediction of Relative Energies and Geometries of Low-Lying Isomers of Water Hexamers. J. Phys. Chem. A 2008, 112, 3976−3984. (8) Santra, B.; Michaelides, A.; Fuchs, M.; Tkatchenko, A.; Filippi, C.; Scheffler, M. On the Accuracy of Density-Functional Theory Exchange−Correlation Functionals for H Bonds in Small Water Clusters. II. The Water Hexamer and Van Der Waals Interactions. J. Chem. Phys. 2008, 129, 194111. (9) Bates, D. M.; Tschumper, G. S. CCSD(T) Complete Basis Set Limit Relative Energies for Low-Lying Water Hexamer Structures. J. Chem. Phys. A 2009, 113, 3555−3559. (10) Góra, U.; Podeszwa, R.; Cencek, W.; Szalewicz, K. Interaction Energies of Large Clusters from Many-Body Expansion. J. Chem. Phys. 2011, 35, 224102. (11) Foley, J. J.; Mazziotti, D. A. Cage versus Prism: Electronic Energies of the Water Hexamer. J. Phys. Chem. A 2013, 117, 6712− 6716. (12) Wang, Y.; Babin, V.; Bowman, J. M.; Paesani, F. The Water Hexamer: Cage, Prism, or Both. Full Dimensional Quantum Simulations Say Both. J. Am. Chem. Soc. 2012, 134, 11116−11119. (13) Babin, V.; Paesani, F. The Curious Case of the Water Hexamer: Cage vs. Prism. Chem. Phys. Lett. 2013, 580, 1−8. (14) Tainter, C. J.; Skinner, J. L. The Water Hexamer: Three-Body Interactions, Structures, Energetics, and OH-Stretch Spectroscopy at Finite Temperature. J. Chem. Phys. 2012, 137, 104304. (15) Wang, Y.; Bowman, J. M. IR Spectra of the Water Hexamer: Theory, With Inclusion of the Monomer Bend Overtone, and Experiment Are in Agreement. J. Phys. Chem. Lett. 2013, 4, 1104− 1108. (16) Nibbering, E. T. J.; Elsaesser, T. Ultrafast Vibrational Dynamics of Hydrogen Bonds in the Condensed Phase. Chem. Rev. 2004, 104, 1887−1914. (17) Rey, R.; Møller, K. B.; Hynes, J. T. Ultrafast Vibrational Population Dynamics of Water and Related Systems: A Theoretical Perspective. Chem. Rev. 2004, 104, 1915−1928. (18) Bakker, H. J.; Skinner, J. L. Vibrational Spectroscopy as a Probe of Structure and Dynamics in Liquid Water. Chem. Rev. 2010, 110, 1498−1517. (19) Perakis, F.; Borek, J. A.; Hamm, P. Three-Dimensional Infrared Spectroscopy of Isotope-Diluted Ice Ih. J. Chem. Phys. 2013, 139, 014501. (20) Bergren, M. S.; Schuh, D.; Sceats, M. G.; Rice, S. A. The OH Stretching Region Infrared Spectra of Low Density Amorphous Solid Water and Polycrystalline Ice Ih. J. Chem. Phys. 1978, 69, 3477. (21) Wojcik, M. J.; Buch, V.; Devlin, J. P. Spectra of Isotopic Ice Mixtures. J. Chem. Phys. 1993, 99, 2332. (22) Li, F.; Skinner, J. L. Infrared and Raman Line Shapes for Ice Ih. I. Dilute HOD in H2O and D2O. J. Chem. Phys. 2010, 132, 204505. (23) Tainter, C. J.; Ni, Y.; Shi, L.; Skinner, J. L. Hydrogen Bonding and OH-Stretch Spectroscopy in Water: Hexamer (Cage), Liquid Surface, Liquid, and Ice. J. Phys. Chem. Lett. 2013, 4, 12−17. (24) Liu, H.; Wang, Y.; Bowman, J. M. Ab Initio Deconstruction of the Vibrational Relaxation Pathways of Dilute HOD in Ice Ih. J. Am. Chem. Soc. 2014, 136, 5888−5891. (25) Liu, H.; Wang, Y.; Bowman, J. M. Quantum Calculations of Intramolecular IR Spectra of Ice Models Using Ab Initio Potential and Dipole Moment Surfaces. J. Phys. Chem. Lett. 2012, 3, 3671−3676. (26) Wang, Y.; Huang, X.; Shepler, B. C.; Braams, B. J.; Bowman, J. M. Flexible, Ab Initio Potential, and Dipole Moment Surfaces for Water. I. Tests and Applications for Clusters up to the 22-Mer. J. Chem. Phys. 2011, 134, 094509. (27) Shank, A.; Wang, Y.; Kaledin, A.; Braams, B. J.; Bowman, J. M. Accurate Ab Initio and “Hybrid” Potential Energy Surfaces, Intramolecular Vibrational Energies, and Classical IR Spectrum of the Water Dimer. J. Chem. Phys. 2009, 130, 144314. (28) Wang, Y.; Bowman, J. M. Communication: Rigorous Calculation of Dissociation Energies (D0) of the Water Trimer, (H2O)3 and (D2O)3. J. Chem. Phys. 2011, 135, 131101.

computational effort for a future study. However, on the basis of the very limited study for the hexamer, which shows the largest difference for the free-OH, and the fact that there are no free-OH stretches in the core regions of the ice or liquid, it is likely that the differences in ZPEs are quite small from site to site. Finally, we attempt to answer the question of which isomer of HOD(D2O)5 is more stable at 0 K, cage or prism? The CCSD(T)/CBS calculations show that the electronic energy of the prism is lower than the cage by 0.25 kcal/mol.9 For the harmonic ZPE, all 12 isotopomers of the prism are higher than the highest-energy isotopomer of the cage. The relative harmonic ZPE of the prism over cage ranges from 0.06 to 0.43 kcal/mol. Therefore, this suggests that both the cage and the prism could be the most stable isomer, depending on the hydrogen sites. These preliminary predictions may be tested by new experiments on the isotopologues of water hexamers.

4. SUMMARY AND CONCLUSIONS We presented local-monomer IR spectra of the cage and prism isomers of HOD(D2O)5 in the range of 1000−4000 cm−1. The appearance of the bend overtone in the region of 2700−3000 was noted. The OH stretch bands of both isomers show rich structures due to the abundant H-bonding environments. Also, the spectra of HOD and D2O ice Ih were presented, and these are in good quantitative agreement with experiment. Finally, we reported preliminary site energetics of the 12 isotopomers of the cage and prism HOD(D2O)5 based on harmonic ZPE analysis.

■

AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS We thank the National Science Foundation (CHE-1145227) for financial support. J.M.B. also thanks Brooks Pate for discussions.

■

REFERENCES

(1) Liu, K.; Brown, M. G.; Carter, C.; Saykally, R. J.; Gregory, J. K.; Clary, D. C. Characterization of a Cage Form of the Water Hexamer. Nature 1996, 381, 501−503. (2) Steinbach, C.; Andersson, P.; Melzer, M.; Kazimirski, J. K.; Buck, U.; Buch, V. Detection of the Book Isomer from the OH-Stretch Spectroscopy of Size Selected Water Hexamers. Phys. Chem. Chem. Phys. 2004, 6, 3320−3324. (3) Diken, E. G.; Robertson, W. H.; Johnson, M. A. The Vibrational Spectrum of the Neutral (H2O)6 Precursor to the “Magic” (H2O)6− Cluster Anion by Argon-Mediated, Population-Modulated Electron Attachment Spectroscopy. J. Phys. Chem. A 2003, 108, 64−68. (4) Nauta, K. Formation of Cyclic Water Hexamer in Liquid Helium: The Smallest Piece of Ice. Science 2000, 287, 293−295. (5) Burnham, C. J.; Xantheas, S. S.; Miller, M. A.; Applegate, B. E.; Miller, R. E. The Formation of Cyclic Water Complexes by Sequential Ring Insertion: Experiment and Theory. J. Chem. Phys. 2002, 117, 1109. (6) Pérez, C.; Muckle, M. T.; Zaleski, D. P.; Seifert, N. A.; Temelso, B.; Shields, G. C.; Kisiel, Z.; Pate, B. H. Structures of Cage, Prism, And Book Isomers of Water Hexamer from Broadband Rotational Spectroscopy. Science 2012, 336, 897−901. G

dx.doi.org/10.1021/jp5061182 | J. Phys. Chem. B XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry B

Article

(29) Ch’ng, L. C.; Samanta, A. K.; Wang, Y.; Bowman, J. M.; Reisler, H. Experimental and Theoretical Investigations of the Dissociation Energy (D0) and Dynamics of the Water Trimer, (H2O)3. J. Phys. Chem. A 2013, 117, 7207−7216. (30) Wang, Y.; Bowman, J. M. Ab Initio Potential and Dipole Moment Surfaces for Water. II. Local-Monomer Calculations of the Infrared Spectra of Water Clusters. J. Chem. Phys. 2011, 134, 154510. (31) Carter, S.; Culik, S. J.; Bowman, J. M. Vibrational SelfConsistent Field Method for Many-Mode Systems: A New Approach and Application to the Vibrations of CO Adsorbed on Cu(100). J. Chem. Phys. 1997, 107, 10458. (32) Carter, S.; Bowman, J. M.; Handy, N. C. Extensions and Tests of “Multimode”: A Code to Obtain Accurate Vibration/Rotation Energies of Many-Mode Molecules. Theor. Chem. Acc. 1998, 100, 191−198. (33) Bowman, J. M.; Carter, S.; Huang, X. MULTIMODE: a Code to Calculate Rovibrational Energies of Polyatomic Molecules. Int. Rev. Phys. Chem. 2003, 22, 533−549. (34) Higgs, P. W. A Method for Computing Zero-Point Energies. J. Chem. Phys. 1953, 21, 1300. (35) Ohno, K.; Okimura, M.; Akai, N.; Katsumoto, Y. The Effect of Cooperative Hydrogen Bonding on the OH Stretching-Band Shift for Water Clusters Studied by Matrix-Isolation Infrared Spectroscopy and Density Functional Theory. Phys. Chem. Chem. Phys. 2005, 7, 3005− 3014. (36) Hayward, J. A.; Reimers, J. R. Unit Cells for the Simulation of Hexagonal Ice. J. Chem. Phys. 1997, 106, 1518. (37) Iglev, H.; Schmeisser, M.; Simeonidis, K.; Thaller, A.; Laubereau, A. Ultrafast Superheating and Melting of Bulk Ice. Nature 2006, 439, 183−186.

H

dx.doi.org/10.1021/jp5061182 | J. Phys. Chem. B XXXX, XXX, XXX−XXX