Infrared Spectra of Protonated Water Clusters, H ... - ACS Publications

Mar 9, 2017 - The CPU time for the oc-VSCF calculation was 10.4 min for ..... For a full list see Figure S4 in the Supporting Information. ... Such ch...
0 downloads 0 Views 4MB Size
Subscriber access provided by HACETTEPE UNIVERSITESI KUTUPHANESI

Article +

2

4

Infrared Spectra of Protonated Water Clusters, H(HO), in Eigen and Zundel Forms Studied by Vibrational Quasi-Degenerate Perturbation Theory Kiyoshi Yagi, and Bo Thomsen J. Phys. Chem. A, Just Accepted Manuscript • DOI: 10.1021/acs.jpca.6b11189 • Publication Date (Web): 09 Mar 2017 Downloaded from http://pubs.acs.org on March 16, 2017

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a free service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are accessible to all readers and citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

The Journal of Physical Chemistry A is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Infrared Spectra of Protonated Water Clusters, H+(H2O)4, in Eigen and Zundel Forms Studied by Vibrational Quasi-Degenerate Perturbation Theory Kiyoshi Yagi,*,†,¶ and Bo Thomsen† †

Theoretical Molecular Science Laboratory and ¶iTHES, RIKEN, 2-1 Hirosawa, Wako, Saitama

351-0198, Japan.

ACS Paragon Plus Environment

1

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 2 of 54

ABSTRACT

The infrared spectrum of H+(H2O)4 recently observed in a wide spectral range has shown a series of bands in a range of 1700-2500 cm-1, which can not be understood by the standard harmonic normal mode analysis. Here, we theoretically investigate the origin of these bands with a focus on (1) the possibility of co-existence of multiple isomers in the Eigen [H3O+(H2O)3] and Zundel [H5O2+(H2O)2] forms and (2) the effect of anharmonic coupling that gives rise to non-zero intensities for overtones and combination bands. Anharmonic vibrational calculations are carried out for the Eigen and Zundel clusters by the second-order vibrational quasi-degenerate perturbation theory (VQDPT2) based on optimized coordinates. The anharmonic potential energy surface and the dipole moment surfaces are generated by a multiresolution approach combining one-dimensional (1D) grid potential functions derived from CCSD(T)-F12, 2D and 3D grid potential functions derived from B3LYP for important coupling terms, and a quartic force field derived from B3LYP for less important terms. The spectrum calculated for the Eigen cluster is in excellent agreement with the experiment, assigning the bands in the range of 17002500 cm-1 to overtones and combination bands of a H3O+ moiety in line with recent reports [J. Phys. Chem. A 2015, 119, 9425; Science 2016, 354, 1131]. On the other hand, characteristic OH stretching bands of the Zundel cluster is found to be absent in the experimental spectrum. We therefore conclude that the experimental spectrum originates solely from the Eigen cluster. Nonetheless, the present calculation for the Eigen cluster poorly reproduces a band observed at 1765 cm-1. A possible nature of this band is discussed.

ACS Paragon Plus Environment

2

Page 3 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

1. INTRODUCTION Interplay between vibrational spectroscopy and computational chemistry has greatly advanced our understanding on the structure and dynamics of hydrated proton. Although the proton in solution has been classically considered to be in two forms, i.e., the Eigen form1-2 (a hydrated hydronium cation, H3O+(H2O)3) and the Zundel form3 (a proton shared between two water molecules, H2O···H+···OH2 or H5O2+), ab initio molecular dynamics (MD) simulations4-6 have suggested that the Eigen and Zundel forms are both meta-stable in the presence of surrounding water molecules, and that the dynamics of water molecules in the second shell is of key importance for an interconversion of Eigen/Zundel motifs that facilitates the transport of protons through liquid water in the Grotthuss mechanism. Vibrational spectroscopic studies7-22 have revealed such a remarkable structural complexity in protonated water clusters, H+(H2O)n. Headrick et al.9 have found that the Eigen and Zundel forms interconvert in a size range of 2 ≤ n ≤ 11 by monitoring an infrared (IR) signal of a proton that dramatically changes from 1000 to 3000 cm-1. Miyazaki et al.16 have shown that the excess charge greatly affects the solvation structure. The water molecules surrounding the proton grow in a two-dimensional, sheet-like structure up to n=20, in contrast to neutral water clusters which start to form a three-dimensional structure already at a size as small as hexamer. The fundamental, spectroscopic information on H+(H2O)n has provided a firm basis for elucidating the dynamics of a Zundel complex in liquid water.23 These studies have utilized theoretical calculations to interpret the experimental spectrum. Density function theory (DFT) offers a cost-effective way to explore wide geometric space and to find stable isomers. Subsequently, the energetic rank of the isomers is refined using accurate electron correlation methods, and the vibrational spectrum is calculated for the lowest energy isomers. Finally, the calculated spectrum is compared with the observed one to identify

ACS Paragon Plus Environment

3

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 4 of 54

Figure 1. (a) Eigen and (b) trans-Zundel isomers of H+(H2O)4. the structure of the cluster. However, such a strategy is now reaching a limit due to the accuracy of vibrational analysis methods used to calculate the spectrum. A showcase example is recent discussion on the structure of H+(H2O)4. Historically, the IR spectrum of H+(H2O)4 was investigated in a high frequency, OH stretching region.15, 24-25 The spectrum in this region matched well with that of an Eigen cluster [H3O+(H2O)3] (Figure 1a) obtained by the harmonic normal mode analysis, and thus the structure of n=4 has been considered to be the Eigen form. The advances in experimental techniques have made it feasible to measure the spectrum in a wider window (600 – 3800 cm-1) as well as with higher resolution.11-13 Surprisingly, the spectrum in the whole range has found not only a signal of fundamental transitions of the Eigen cluster that are characterized by the harmonic analysis, but also a series of congested bands in a range of 1700-2500 cm-1. This feature has brought about an intense debate. Kulig and Agmon26 have suggested that the spectrum may be an indication of another isomer in a Zundel form, H5O2+(H2O)2 (Figure 1b), based on ab initio MD calculations using DFT that yielded an intense band around 1750 cm-1. However, the MD calculation was unable to reproduce other bands observed in a range of 18002500 cm-1. On the other hand, McCoy et al.12 have carried out three-dimensional quantum mechanical calculations for H3O+Ar3 using two bending modes and one frustrated rotation of H3O+, and have shown that a combination band acquires sizable IR intensity through anharmonic

ACS Paragon Plus Environment

4

Page 5 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

coupling. From this result, a strong band observed at 2250 cm-1 has been suggested to be such a combination band of H3O+(H2O)3 (the so-called α band). Recently, an anharmonic calculation has been carried out for H+(H2O)414,

22

and its deuterated species22 using the second-order

vibrational perturbation theory (VPT2).27 The VPT2 result for the Eigen cluster has not only reinforced the above assignment but also featured other combination bands around 1850 cm-1. Together with the recent IR-IR double resonance measurement, which has revealed that the spectrum originates from a single isomer,14 it has been suggested that the observed isomer is solely the Eigen cluster. [We note that n=5 is also a single isomer, whereas n=6 and 7 have multiple isomers.17-21] Nevertheless, the agreement between the VPT2 and experimental spectrum was not sufficient, leaving some of the bands unassigned in the range of 1700-2500 cm1

. Furthermore, VPT2 calculation for the Zundel cluster (Fig. S1 of Ref. 14 in the supporting

information) is too different from the harmonic spectrum, which could be due to a divergence of the perturbative expansion. Therefore, a thorough vibrational analysis of the Eigen and Zundel clusters using a higher level of theory is still necessary. Quantum vibrational calculations on H+(H2O)n are still few due to a challenge in solving the vibrational Schrödinger equation (VSE). Park et al.28 have calculated the IR spectrum of Zundel (n = 2 and 6) and Eigen (n = 4) clusters using the correlation corrected vibrational self-consistent field (cc-VSCF) method.29-31 Torrent-Sucarrat and Anglada32 have studied the Eigen clusters (n = 3, 4, and 21) using VPT2. Although VPT2 and cc-VSCF are cost-efficient, these methods may diverge in the presence of strong vibrational resonance, because they are both based on the perturbation theory. Vibrational configuration interaction (VCI)33-35 based on VSCF36 is more stable and accurate, and has been applied to the Zundel cluster (n=2).10 Meyer and coworkers37-39 also studied the Zundel cluster (n=2) using the multiconfiguration time-dependent Hartree

ACS Paragon Plus Environment

5

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 6 of 54

(MCTDH) method.40 VCI and MCTDH have successfully assigned the IR signal of a shared proton in the Zundel cluster. However, applications to larger clusters are so far limited due to unfavorable scaling of the number of configurations with respect to system size. Additionally, MCTDH requires a derivation of the kinetic energy term tailored for each molecular system. We have recently developed efficient and accurate vibrational analysis methods applicable to general polyatomic molecules based on a vibrational Hamiltonian in terms of rectilinear coordinates. Optimized coordinate VSCF (oc-VSCF)41 provides a “good” set of coordinates to describe the system, which are hereafter denoted optimized coordinates. We have demonstrated that optimized coordinates of water clusters are localized to each water molecules, and that the locality, which reduces the anharmonic coupling, dramatically accelerates the convergence of vibrational calculations.41-42 The second-order vibrational quasi-degenerate perturbation theory (VQDPT2)43-44 is an efficient solver of the VSE, which treats the strong resonance between quasi-degenerate states by configuration interaction (CI) and the other weak interaction by perturbation theory (PT). The combination of CI and PT alleviates the divergence of perturbation expansion, while retaining its cost-efficiency. Besides the difficulties in solving the VSE, an intensive cost is generally required for generating the anharmonic potential energy surface (PES). For this purpose, we have developed a multiresolution method,45-46 which generates strongly coupled, important terms of the PES using accurate electronic structure theory and potential functions, and other weakly coupled terms by less accurate, but cost-effective methods. The use of multiple resolutions as well as a massively parallel implementation has made an efficient generation of the PES feasible. These methods have been recently applied to a pentapeptide to analyze the spectrum in the OH and NH stretching region.47

ACS Paragon Plus Environment

6

Page 7 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

In this study, we have calculated the IR spectrum of H+(H2O)4 in both Eigen and Zundel forms using VQDPT2 with optimized coordinates and a multiresolution PES. The methods are described in Section 2. Section 3 presents the spectrum calculated for the Eigen and Zundel clusters. The spectra are analyzed by investigating the component of the vibrational wavefunction and by introducing simplified, reduced dimensional models. The analysis has revealed that the bands observed in 1700-2500 cm-1 originates from overtones and combination bands of a H3O+ moiety of the Eigen cluster. The possibility of co-existence of the Eigen/Zundel clusters is discussed in Section 4.

ACS Paragon Plus Environment

7

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 8 of 54

2. METHOD 2.1. Optimized coordinate VSCF. The oc-VSCF method41 solves the following equation using a Taylor expansion PES around the equilibrium geometry, (1)

(2)

where f is the vibrational degrees of freedom, and

and

denote the s-th coordinate

and a one-mode function, respectively. These quantities are represented as, (3) M

(s) (s) φn(s) = ∑ dmn χm ,

(4)

s

m=0

where Qi and χ m(s) are the normal coordinates and the harmonic oscillator (HO) basis functions, respectively. The coefficients, U and d (s) , are simultaneously optimized by minimizing the VSCF energy using a Jacobi sweep algorithm.41 The coordinates,

, thus obtained are the

optimized coordinates. The behavior of optimized coordinates is known in two extreme cases. The coordinates become localized for non-interacting anharmonic systems, whereas they are (delocalized) normal coordinates for a harmonic system, i.e., when the PES in Eq. (2) is truncated at the second order. Therefore, optimized coordinates are a natural extension of normal coordinates for anharmonic systems, and are generally expected to have a “good” locality. The locality provides a compact

ACS Paragon Plus Environment

8

Page 9 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

representation of the PES and a superior basis for vibrational structure calculations compared to normal coordinates. In this study, we employed a quartic force field (QFF) for the oc-VSCF calculation. The equilibrium geometry and the QFF coefficients in terms of normal coordinates were computed using the B3LYP48-49 functionals with aug-cc-pVTZ50 basis sets. The Cartesian coordinates are given in Table S1 of the supporting information (SI). The second-order terms were obtained from the analytical Hessian, while the third- and fourth-order terms were obtained by numerical differentiations of the Hessian.51-52 The stepsize was set to

with δ = 0.5 for each

coordinate, where ωi is the harmonic frequency. The numerical differentiation required the computation of Hessian at 2,179 points. The 9 and 7 lowest frequency modes were neglected for the Eigen and Zundel clusters, respectively, to avoid a divergence of the VSCF cycle. 6 and 4 of them, respectively, corresponded to wagging and twist motions of terminal water molecules (shown in Fig. S1). The number of modes is different, because the Eigen cluster has three such water molecules, whereas the Zundel cluster has only two. The remaining three modes were collective, slow motions with the frequency less than 100 cm-1 (also shown in Fig. S1). The effect of these modes on the resulting IR spectrum is unknown; yet we note that the IR spectrum of neutral water clusters has been calculated reasonably well with similar settings, for example, for a trimer53 and a hexamer.54 Note also that the present model includes O···O (or water···water) stretching modes that modulate the hydrogen bond distance and are known to affect high frequency motions significantly.55 The number of HO functions was set to 11 for all modes, i.e. M=10. The target VSCF state was the vibrational ground state, i.e., ns = 0 for all s. The Jacobi sweep was carried out using a parameter p = 2 for the order of the Fourier expansion, a pair selection44 with a threshold value of 500 cm-1, and a convergence criterion of 10-5 cm-1 for

ACS Paragon Plus Environment

9

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 10 of 54

the VSCF energy. The CPU time for the oc-VSCF calculation was 10.4 min for the Eigen cluster and 32.2 min for the Zundel cluster using a single core of Intel Xeon5650 2.66 GHz. Further details are given in the SI. 2.2. Generation of a multiresolution PES. Although QFF already incorporates anharmonicity up to the fourth order, the higher order terms are needed to achieve spectroscopic accuracy. Here, we use a multiresolution method45 to construct the PES in terms of optimized coordinates. In this method, the PES is expanded in terms of mode coupling up to the n-th order (n-mode coupling representation: nMR),56 (5) and the resolution of coupling terms is varied according to their strength. Important mode coupling terms are generated by an accurate, grid method, while the other, less important terms are represented by QFF or even neglected. Although the cost of grid PES generation is expensive, the method takes advantage of the fact that the mode coupling strength (MCS) decays rapidly so that the number of important coupling terms remains moderate. For example, the two-mode potential, Vst , may be negligibly weak when

and

have little spatial overlap and/or when

their associated frequencies are vastly different. Based on such concept, we have previously proposed an estimate of MCS, which can be readily calculated given the harmonic frequency and the QFF coefficient.46 Note that the scheme assumed that cst = 0 (s ≠ t) , which is true for normal coordinates but not necessarily for optimized coordinates. We have therefore taken into account an additional term, (6)

ACS Paragon Plus Environment

10

Page 11 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

where

ω s is the harmonic frequency associated with the s-th coordinate, i.e,

. Eq. (6)

is added to two-mode MCS, Eq. (13) of Ref. [46]. The important terms are identified from the magnitude of MCS and subsequently generated by the grid method. The dipole moment surface (DMS), which is needed for computing the IR intensity, is constructed in a similar way, (7) Note that electronic structure methods such as DFT compute the dipole moment together with the energy, so that the DMS is generated with no additional cost. In the present study, the multiresolution PES was constructed combining the grid PES and QFF. The grid PES used 11 grid points along each coordinate, which were determined by the discrete variable representation (DVR) method.57 The 1MR grid PES was generated at the CCSD(T)-F1258-60 level with aug-cc-pVDZ basis sets for all active modes, i.e., 24 and 26 modes for Eigen and Zundel clusters, respectively. The grid PES was generated for 2MR and 3MR terms with MCS larger than 1.0 cm-1 at the B3LYP/aug-cc-pVTZ level of theory. Other 2MR and 3MR terms (i.e., with MCS smaller than 1.0 cm-1) and all the 4MR terms were represented by QFF generated at the level of B3LYP/aug-cc-pVTZ. The 5MR and higher order terms were neglected. The correction for higher electronic correlation effects using CCSD(T)-F12, though limited to 1MR-PES, substantially refined the quality of the PES. For example, the multiresolution PES yielded the O···O, O···H and O-H distances of the Eigen cluster as 2.5535, 1.5504, and 1.0056 Å, respectively, which were shifted from those of B3LYP (2.5625, 1.5502, and 1.0143 Å, respectively) to agree with those of CCSD(T)-F12 (2.5544, 1.5507, and 1.0058 Å, respectively) within 0.001 Å. The harmonic frequencies of the Eigen cluster obtained from the

ACS Paragon Plus Environment

11

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 12 of 54

multiresolution PES agreed with those of CCSD(T)-F12 with a mean absolute deviation (MAD) of 4.6 cm-1 and a maximum deviation of 16.6 cm-1. The result contrasts to those of B3LYP, which deviated from CCSD(T)-F12 by a MAD of 35.2 cm-1 and a maximum deviation of 132.1 cm-1. A similar trend was observed for the Zundel cluster as well. The Cartesian coordinates and harmonic frequencies are listed in Tables S1 and S2, respectively. The DMS was constructed using the same grid for PES and computed at the B3LYP/aug-cc-pVTZ level of theory. The electronic structure calculations were performed using Gaussian0961 and Molpro2012.1.62-63 The total number of grid points for generating the grid PES was 143,041 and 296,761, for the Eigen and Zundel clusters, respectively. The energy and the dipole moment at these grid points were calculated in parallel using supercomputers of RIKEN and Nagoya University. Note that optimized coordinates yield a more compact representation of the PES than normal coordinates. The number of grid points increases to 287,041 and 479,761 for the Eigen and Zundel clusters in terms of normal coordinates. See Table S3 in the SI for more details. The saving in this step more than compensates the cost of the oc-VSCF calculation. 2.3. VQDPT2 calculation. Finally, we solve the vibrational Schrödinger equation using the multiresolution PES, (8) Note that the vibrational angular momentum terms are not considered in this method. The largest rotational constant is 0.09 and 0.34 cm-1 for the Eigen and Zundel clusters, respectively, so that their effect is expected to be on the order of 1 cm-1. The VSCF calculation is first performed for the vibrational ground state to obtain the one-mode functions for ns = 0 and the virtual functions for ns > 0. These functions are used to construct a VSCF configuration function,

ACS Paragon Plus Environment

12

Page 13 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(9) In VQDPT2,43-44 the effective Hamiltonian reads,

( Hˆ ) (2) eff

  p Hˆ q q Hˆ p′ 1 1  = p Hˆ p′ + ∑  (0) + , (0) (0)  pp′ 2 Ep(0) ′ − Eq  q∈Q  Ep − Eq

(10)

where p is a component of P space, in which VSCF configurations are close in energy, and Q is a complementary space containing non-degenerate VSCF configurations. Ep(0) is the zeroth order energy, f

E

(0) p

= ∑ε (s) ps .

(11)

s=1

Note that the method reduces to cc-VSCF29 or the second-order vibrational Møller-Plesset perturbation method (VMP2),64-65 if the P space contains only one VSCF configuration, i.e., if there are no configurations close in energy. The diagonalization of the effective Hamiltonian matrix, Eq. (10), yields the energy and the wavefunction for a state n, (2) Hˆ eff C = E(2)C,

(12)

Ψ (2) n = ∑ Cnp p .

(13)

p∈P

Then, the IR spectrum is obtained as,

I(ν ) = ∑ I n Γ(ν − ν n ),

(14)

n

where νn and In are the transition frequency and the IR intensity, respectively,

In =

,

(15)

2 2πν n nMR Ψ (2) 0 , n d 2 3c

(16)

ACS Paragon Plus Environment

13

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 14 of 54

with the speed of light c. Γ is a Lorentzian line shape function, Γ(ν − ν n ) =

 2 γ  , π  4(ν − ν n )2 + γ 2 

(17)

where γ is a full width at half maximum (FWHM) of the vibrational band. In the present study, the P and Q space was constructed using an algorithm in Ref.[44] with parameters set to k = 4 and Ngen = 3. The target configurations for constructing the P space were set to the fundamentals and combination bands for the Eigen cluster, and to the fundamentals for the Zundel cluster. The FWHM was set to 10 cm-1 when constructing the spectrum, though the appearance of the spectrum is insensitive to this parameter as far as it is chosen in a range of 5 – 20 cm-1 (see in Fig. S2 of SI).

ACS Paragon Plus Environment

14

Page 15 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

3. RESULTS 3.1. Eigen cluster. Optimized coordinates obtained from the oc-VSCF calculation are depicted in Figure 2 (a). The 24 active coordinates can be classified into three different types of motion: (1) nine internal vibrations of water molecules, (2) six internal vibrations and three frustrated rotations of H3O+, and (3) six intermolecular motions. The internal vibrations of H2O are localized to each water molecule, in striking contrast to normal coordinates, which are delocalized over the three water molecules according to the C3 symmetry of the cluster. This result is consistent with the previous work on water clusters,41-42 where oc-VSCF yielded OH stretching and bending modes localized to each water molecule. The optimized coordinates of H3O+ look similar to normal coordinates except that OH stretching modes are localized to each OH bond. As for the intermolecular motions, the rocking vibrations are localized to each water molecule, while H3O+ ··· H2O stretching vibrations remain more or less the same as normal coordinates. We emphasize that the locality of these coordinates is automatically derived by ocVSCF without requiring any input (or insight) on how local the coordinates should be. Further comparison of normal and optimized coordinates is given in the SI. VQDPT2 calculations have been performed using the 24 optimized coordinates (denoted 24D). In addition, we have carried out two types of reduced dimensional calculations to interpret the spectrum. In the first model, denoted 18D, VQDPT2 was carried out in 18 dimensions neglecting six intermolecular modes. In the second model, denoted 9D/3D×3, the 18-dimensional system was further divided into three subsystems; nine vibrations of H3O+ and three internal vibrations of each water molecule. Note that the coupling between H3O+ – H2O and H2O – H2O is turned off in this model, and that each molecule is treated independently. VQDPT2 calculations were performed for each molecule and the resulting IR spectra were summed to obtain the total

ACS Paragon Plus Environment

15

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 16 of 54

Figure 2. (a) Optimized coordinates of the Eigen cluster. For a full list see Figure S4 in SI. (b) IR spectrum of the Eigen cluster obtained by VQDPT2 calculations based on 24D, 18D, and 9D/3D×3 models, and the harmonic approximation at the B3LYP/aug-cc-pVTZ level together with the experiment.13 The fundamental transition is labeled with a number. The labels, m×n, m+m’, and pg, denote an (n – 1)-th overtone of mode m, a combination band of mode m and m’, and progression bands, respectively.

spectrum. The calculated spectra are shown in Figure 2 (b) together with the experiment.13 The spectrum obtained from the harmonic approximation at the B3LYP/aug-cc-pVTZ level exhibits characteristic two bands (1 and 2) that originate from OH stretching modes of water molecules, and another two bands (3 and 4) from OH stretching modes of H3O+. Note that 4

ACS Paragon Plus Environment

16

Page 17 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

carries a large IR intensity, so that its intensity is scaled to half to make other bands visible. In a lower frequency region, the bending motion of H2O (6), and the umbrella motion and frustrated rotation of H3O+ (7 and 8, respectively) are obtained. The bending motion of H3O+ and the other component of the frustrated rotation (5 and 9) have negligibly small intensity, but their positions are shown in Fig. 2 (b) for later convenience. The harmonic spectrum is obviously different from the experiment in a region of 1700-2800 cm-1, completely missing the observed bands. Therefore, a high-level vibrational structure theory is needed to elucidate the experimental spectrum in the whole range. In the 9D/3D×3 model, despite its simplicity, the spectrum is considerably improved over the harmonic one, not only correcting the band positions to the right place but also showing six new bands in a range of 1700 – 3000 cm-1. The wavefunction analysis (cf. Table S8) has revealed that ×2 and 8+7, are the these bands are highly excited states of H3O+. Two bands around 2000 cm-1, 8× first overtone of the frustrated rotation and a combination band of the frustrated rotation and the umbrella motion, respectively. Three bands in a range of 2300-2750 cm-1, 9+5, 8+5, and 7+5, are combination bands involving the bending motion, and a band at 2795 cm-1, 8× ×3, is the second overtone of the frustrated rotation. These vibrational states acquire intensity due to an intensity borrowing from the fundamental transitions of three OH stretching motion of H3O+ (denoted , and

,

). The two most intense bands, 4 and 7+5, share the component of 71 , 81 , and 91

by 0.79 and 0.53, respectively, in weight, where m1 denotes a VSCF configuration function for the fundamental transition of the m-th mode. The weight is much smaller in the wavefunction of other states: 8× ×2, 8+7, 9+5, 8+5, and 8× ×3 have a weight of 0.04, 0.14, 0.18, 0.19, and 0.09, respectively. Nevertheless, the component dominantly contributes to their absorption, since 71 ,

ACS Paragon Plus Environment

17

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 18 of 54

81 , and 91 carry a large IR intensity. To confirm this mechanism, we carried out a calculation removing one-mode components of

,

, and

from the DMS to find that these bands were

all diminished. Although the 9D/3D×3 model brings in much insight, the drawback is that each molecule is decoupled and thus the vibrational wavefunctions are restricted to be local. The 18D model, in which the coupling is explicitly considered, recovers delocalized nature of the wavefunctions. Such character is prominent for the bending motion of H3O+ (

,

) and H2O (

,

,

).

The wavefunctions of two components of 5 are obtained as

Ψ15 = 0.74 111 − 0.50 121 − 0.42 131

(18)

Ψ 52 = 0.82 101 − 0.39 141

(19)

and a similar mixing is seen for 6. Thus, the bending motions of H3O+ and H2O are heavily mixed. Although the mixing affects the fundamental transitions (5 and 6) only slightly, the effect on the higher excited states is significant. An enlarged view shown in Figure 3 indicates that the combination bands, n+5, in the 9D/3D×3 model split into two in the 18D model for n = 7, 8, and 9. By analyzing the wavefunction (cf. Table S9), these bands are assigned to n+5 and n+6 (n = 7, 8, and 9). The incorporation of intermolecular motions in the 24D model brings about a heavy mixing of VSCF configurations, because the low-frequency modes constitute a dense density of states that manifests vibrational resonance. Therefore, the intensity is distributed over many vibrational states. 56 states are found to carry the IR intensity more than 10 km mol-1 in a range of 2250 – 2750 cm-1. For example, the following states,

Ψ14 = 0.23 91 − 0.25 232 211121 − 0.27 211181101 ,

(20)

ACS Paragon Plus Environment

18

Page 19 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

2655

(a) Exp. 2240 α 2300 2340

2485

2990

4 7+6 7+5, 8×3, 3 8+6 8+5 pg2

(b) VQDPT2 (24D) 9+6 9+5

(c) VQDPT2 (18D)

4 7+6

9+6 9+5

8+6 8+5

7+5 8×3 3

(d) VQDPT2 (9D/3D×3) 4 9+5

8+5

7+5 8×33

wavenumber / cm-1 Figure 3. Enlarged view of Figure 2 (b) in a range of 2200 – 3000 cm-1. (a) The experimental spectrum, and the VQDPT2 spectra obtained from (b) 24D, (c) 18D, and (d) 9D/3D×3 models.

Ψ 42 = −0.25 151121 − 0.28 232 211121 − 0.24 232191111 + 0.24 232 201131 ,

(21)

Ψ 34 = 0.24 71 − 0.23 151141 ,

(22)

Ψ 44 = 0.28 221231211121 + 0.23 221231191101 − 0.28 221231191141 ,

(23)

are obtained in a range of 2680-2710 cm-1 with an IR intensity larger than 200 km mol-1 at 2695.6, 2695.6, 2704.6, and 2705.2 cm-1, respectively. The wavefunction of other states are listed in Table S10. The strong mixing hampers a clear-cut assignment of the calculated bands. Nevertheless, as shown in Figure 3, the band shape of the 24D model resembles to that of the 18D model, and thus the prominent bands can be assigned to 9+6/9+5, 8+6/8+5, 7+6 and 4. The

ACS Paragon Plus Environment

19

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 20 of 54

assignment is consistent with the previous studies12, 14 in which the band at 2240 cm-1 (the α band) has been assigned to a combination band of a frustrated rotation and bending motions of H3O+. The three bands, 3, 7+5, and 8× ×3, are hard to identify in the 24D spectrum, but are likely to be around 2740 cm-1. A unique feature in the 24D spectrum is three new bands, designated as pg1, pg2, and pg3. Hydrogen bonds, D-H···A, are known to exhibit a progression of slow, intermolecular D···A stretching motion upon an excitation of fast D-H stretching vibration in the IR spectrum in analogy to the vibrational progression in the electronic excitation spectrum. We refer to Ref. [55], in which the mechanism has been thoroughly discussed. The three bands are manifestation of such a progression, because they emerge upon the inclusion of intermolecular motions. Since the frequency of the intermolecular motions are about 300 cm-1, it is plausible to assign pg1 and pg2 to a progression upon exciting 4 and 8+6/8+5, respectively. The origin of pg3 is less clear, but seems to arise from 8+7. The calculated vibrational frequencies and IR intensities are listed in Table 1 together with the experimental frequencies obtained by using D2,13 Ar,9, 11, 15 and Ne18 for a tagging species. The ×2 are in excellent agreement with experiment with OH stretching frequencies of H2O, 9+5, and 8× a deviation of ~5 cm-1. On the other hand, the main band, 4, and the combination bands, 9+6 and 8+7, show larger deviations of 35 – 50 cm-1. The bending motion of H2O is lower than the experiment by 25 – 30 cm-1. The accuracy is ~20 cm-1 on average for the sharp bands. As for the broad bands, pg1, 8+6/8+5, and 7/8 are calculated well within the envelope of the bands observed around 2990, 2485, and 1045 cm-1 (see Fig. 2 (b)), though the peak positions do not necessarily match well.

ACS Paragon Plus Environment

20

Page 21 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

It is also notable in Table 1 that the harmonic frequencies of 3 and 4 are obtained as 2978 and 2850 cm-1 using B3LYP/aug-cc-pVTZ, which are ~130 cm-1 lower than those by CCSD(T)F12/aug-cc-pVDZ (3079 and 2982 cm-1, respectively). In contrast, the multiresolution PES yields 3096 and 2977 cm-1, respectively, in better agreement with the CCSD(T)-F12/aug-ccpVDZ result. The result indicates the importance and effectiveness of the 1MR correction using CCSD(T)-F12/aug-cc-pVDZ for achieving high accuracy in the final anharmonic frequencies. In view of the quantitative agreement, we are tempted to say that the observed spectrum originates solely from the Eigen cluster. However, there remains one problem regarding a band observed at 1765 cm-1. The band was assigned to the bending motion of H3O+ in the previous work.12 However, the present calculation does not support this assignment, because 5 deviates by 126 cm-1, which is beyond the expected accuracy, and also because the intensity is too weak. ×2 based on a VPT2 calculation. In the Recently, Fournier et al.14 have assigned this band to 8× present study, 8× ×2 better matches the frequency with a band at 1840 cm-1, which is very weak in the measurement using D2 tagging, but clearly discernible in a spectrum obtained by using Ar tagging.12 On the other hand, Kulig and Agmon26 have suggested that this band is a signal of the Zundel cluster. Thus, we next investigate the spectrum of the Zundel cluster.

ACS Paragon Plus Environment

21

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

Page 22 of 54

Table 1. Vibrational frequencies (ν in cm-1) and IR intensities (I in km mol-1) of the Eigen cluster obtained by the harmonic approximation and VQDPT2 calculations based on 9D/3D×3, 18D, and 24D models, and the experimental value. banda

assignment

Harmonicb

18D

9D/3D×3

24D

Exp.

ν

I

ν

I

ν

I

ν

I

ν

1

asym. OH str. of H2O

3870 (3914)

482

3723

475

3720

475

3724

460

3723/3733c 3730d, 3732e, 3732f

2

sym. OH str. of H2O

3780 (3828)

149

3662

138

3657

138

3652

131

3636/3648c 3644d, 3646e, 3647f

pg1

pg g for 4

-

-

-

-

-

-

3096

46

2990c, 3050d, 3053e

pg2

pg g for 8+6 and/or 8+5

-

-

-

-

-

-

2877

12

-

3

sym. OH str. of H3O+

133

2823

130

2814

125

~2740

-

-

8×3

2nd overtone of 8

0

2795

271

2784

520

~2740

-

-

7+5

combination band of 7 and 5

0

2753

1634

2743

1049

~2740

-

-

4

asym. OH str. of H3O+

2978 (3096) 2921 (2894) 2858 (2894) 2850 (2977) 2789 (2825) 2669 (2686) 2600 (2617) 2426 (2429)

6009

2711

2353

2719

1899

2699

1613

2655c, 2665d, 2670e, 2657f

0

-

-

2684

770

2667

581

-

0

2551

554

2550

267

2560

183

~2485c, ~2500d

0

-

-

2516

207

2507

189

~2485c, ~2500d

0

2319

438

2321

213

2340

59

2340c

-

-

-

-

-

-

2302

35

2300c, 2307d

2357 (2360)

0

-

-

2292

226

2276

255

2240c, 2245d

7+6 8+5 8+6 9+5

9+6

combination band of 7 and 6 combination band of 8 and 5 combination band of 8 and 6 combination band of 9 and 5

combination band of 9 and 6

ACS Paragon Plus Environment

22

Page 23 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

The Journal of Physical Chemistry

pg3

pg g for 8+7

-

-

-

-

-

-

2133

93

8+7

combination band of 8 and 7

2136 (2137)

0

2032

427

2012

417

1970

94

-

-

-

-

-

-

1952

149

1905c, 1904d

1947 (1930)

0

1902

94

1882

85

1845

87

1840c, 1847d

-

-

-

-

-

-

-

-

1765c, 1750d

22

1624

10

1635

27

1639

38

-

120

1604

88

1583

112

1590

105

1615c,d

258

1120

247

1112

242

1095

222

1045c,d

113

954

107

943

106

940

103

1045c,d

0

722

0

717

0

703

0

-

8×2

1st overtone of 8

?h 5

bending of H3O+

6

bending of H2O

7

umbrella of H3O+

8

frustrated rot. of H3O+

9

frustrated rot. of H3O+

1695 (1721) 1627 (1652) 1163 (1173) 974 (965) 730 (708)

a

The label of vibrational bands in Figure 2 (b).

b

Frequencies from B3LYP/auc-cc-pVTZ (without parentheses) and the multiresolution PES (with parentheses).

c

D2 tag, Fournier et al. (2014).

d

Ar tag, Headrick et al. (2005), Olesen et al. (2011).

e

Ar tag, Douberly et al. (2010).

f

Ne tag, Mizuse and Fujii (2012).

g

Progression of intermolecular vibrations.

h

Unassigned band found in the experiment.

-

ACS Paragon Plus Environment

23

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 24 of 54

3.2. Zundel cluster. Optimized coordinates of the Zundel cluster are depicted in Figure 4 (a). The coordinates can be classified into three types: (1) six internal vibrations of outer water molecules (W1 and W2), (2) 16 vibrations of the Zundel core (six internal vibrations of WZ1 and WZ2, three motions of a proton, and seven motions of WZ1 and WZ2) and (3) four intermolecular motions of W1 and W2. The internal vibrations of W1 and W2 are localized to a single water molecule as in the Eigen cluster. On the other hand, the vibrations of the Zundel core have various types of locality. The OH stretching motions are localized to a single OH bond, whereas the bending motions and the proton motions are delocalized over the Zundel core. The wagging motion is mainly a motion of a single OH bond, though neighboring atoms also contribute to this motion. The rocking motion involves one water molecule and a proton, resembling the frustrated rotation of H3O+. Further analysis of optimized coordinates is given in the SI. VQDPT2 calculations have been performed using the 26 vibrational coordinates (denoted 26D). In addition, we have carried out two types of reduced dimensional calculations. The first model, denoted 16D/3D×2, treated the Zundel core, W1, and W2 independently, in 16, 3, and 3 dimensions, respectively, neglecting the four intermolecular vibrations. Note that the 16D model incorporates 15 internal vibrations and one overall rotation of the Zundel core. In the second model, denoted 12D/3D×2, the Zundel core was further reduced to 12 dimensions neglecting the rocking motions and two wagging motions of a dangling OH bond of WZ1 or WZ2. The calculated spectra are shown in Figure 4 (b) together with the harmonic spectrum. The calculated band position and intensity are listed in Table 2. The harmonic spectrum obtained by the B3LYP/aug-cc-pVTZ method gives three bands around 3800 cm-1 (νa, νz, and νs) and two intense bands around 3000 cm-1 (1 and 2). The former originates from OH bonds free from hydrogen bonds. The highest and lowest frequency bands,

ACS Paragon Plus Environment

24

Page 25 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

νa and νs, are asymmetric and symmetric stretching modes of W1 and W2, respectively, whereas the one in the middle, νz, arises from dangling OH bonds of WZ1 and WZ2. The latter arises from OH bonds of WZ1 and WZ2 that donate a hydrogen bond to W1 and W2. An asymmetric stretching combination of these OH bonds carries a large intensity of 2582 km mol-1. Four bands in a range of 1400-1800 cm-1 are bending modes (3 and 4) and motions of proton perpendicular to the O –

(a)

(b)

H2O (3D×2)

OH stretch×2

11 910

bending

H5O2+ (16D)

7

9+8, 6, 5

4 3

8, 11+9

OH stretch×4

bending×2

pg

3 11×2 6 9+8

11 910

7

11+10×2 54

7×2

stretch

proton transfer

9+8

wagging×4

3

rocking×2

1

νs

2

νz νa

Harmonic

7 10 9

stretch×2

VQDPT2 (12D/3D×2)

8

Intermolecular motions (4D)

rocking×2

654

7

9

VQDPT2 2 1, 10×2+3 (16D/3D×2) 10×2+9+8 νs νz νa 8+7×2

2

8 proton×2

VQDPT2 (26D) νa νz νs

2, 1

11×2

8, 11+9

654

3

1

νs

νz νa

wavenumber / cm-1

Figure 4. (a) Optimized coordinates of the Zundel cluster. For a full list see Figure S6 in SI. (b) IR spectrum of the Zundel cluster obtained by VQDPT2 calculations based on 26D, 16D/3D×2, and 12D/3D×2 models, and the harmonic approximation at the B3LYP/aug-ccpVTZ level. The labels are defined in the same way as in Figure 2 (b).

ACS Paragon Plus Environment

25

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

Page 26 of 54

Table 2. Vibrational frequencies (ν in cm-1) and IR intensities (I in km mol-1) of the Zundel cluster obtained by the harmonic approximation and VQDPT2 calculations based on 12D/3D×2, 16D/3D×2, and 26D models. banda

assignment

Harmonicb ν

νa

asym. OH str. of W1/W2

νz

dangling OH str. of WZ1/WZ2

νs

sym. OH str. of W1/W2

1

sym OH str. of WZ1/WZ2

2

asym OH str. of WZ1/WZ2

pg

progression of intermolecular motions for 3

7×2

1st overtone of 7

11+10×2

combination band of 11 and 1st overtone of 10

3

bending of WZ1/WZ2

4

bending of W1/W2

9+8

combination band of 9 and 8

5

proton motion

6

proton motion

3877 (3920) 3825 (3867) 3786 (3834) 3042 (3184) 3016 (3159) 2010 (1969) 1537 (1530) 1737 (1770) 1638 (1660) 1471 (1422) 1539 (1575) 1458 (1490)

12D/3D×2

26D

16D/3D×2

I

ν

I

ν

I

ν

I

315

3727

308

3727

308

3730

308

292

3691

292

3692

293

3695

224

101

3667

87

3667

87

3656

92

657

2846

674

2878

728

~2845

889

2582

2830

2588

2852

968

~2815

872

-

-

-

-

-

1906

250

0

-

-

1838

77

-

-

0

-

-

1749

253

-

-

793

1753

730

1710

882

1742

144

160

1599

184

1599

184

1606

140

0

1464

93

1335-1410

361

1250-1400

59

1493

60

1488

30

1250-1400

241

1418

215

1435

206

1250-1400

30

ACS Paragon Plus Environment

26

Page 27 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48

The Journal of Physical Chemistry

11×2 11+9 7

1st overtone of 11 combination band of 11 and 9 wagging of H-bonded OH bonds of WZ1/WZ2

8

proton transfer motion

9

WZ1-WZ2 str.

10

rocking of WZ1/WZ2

11

wagging of dangling OH bonds of WZ1/WZ2

617 (605) 934 (945) 1005 (985) 846 (780) 626 (643) 615 (613) 308 (303)

0

-

-

1070

169

1107

142

0

-

-

1044

834

1087

158

1325

976

155

958

151

950

139

1848

982

2571

994

912

1053

434

4

581

2

571

2

557

3

192

-

-

617

12

614

15

143

-

-

460

162

450

72

a

The label of vibrational bands in Figure 4 (b).

b

Frequencies from B3LYP/auc-cc-pVTZ (without parentheses) and the multiresolution PES (with parentheses).

ACS Paragon Plus Environment

27

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 28 of 54

(a) 7

1005 cm-1 1325 km mol-1

8

846 cm-1 1848 km mol-1

(b) 7

976 cm-1 155 km mol-1

8

982 cm-1 2571 km mol-1

Figure 5. (a) Normal and (b) optimized coordinates corresponding to the wagging motion (7) and the PT motion (8). The associated harmonic frequencies and IR intensities are indicated.

O axis of WZ1 and WZ2 (5 and 6). Below 1000 cm-1, two prominent bands are obtained, which are assigned to wagging motions of a hydrogen bonded OH bond of WZ1 or WZ2 (7) and a proton transfer (PT) motion (8). We note in passing that the harmonic spectrum suffers from the accuracy of B3LYP/aug-cc-pVTZ from a quantitative point of view. The harmonic frequencies of 1 and 2 deviate more than 140 cm-1 compared to those of CCSD(T)-F12/aug-cc-pVDZ. In striking contrast, the multiresolution PES yields the harmonic frequency in agreement with CCSD(T)-F12/aug-cc-pVDZ with a deviation of < 1 cm-1 for 1 and 2, and with a MAD of 6.7 cm-1 (see Table S2). As in the Eigen cluster, the result substantiates the effectiveness of 1MR correction for quantitative purposes. VQDPT2 based on the 12D/3D×2 model gives a single intense band below 1000 cm-1, which contrasts to the two bands (7 and 8) in the harmonic spectrum. The spectral change is caused for

ACS Paragon Plus Environment

28

Page 29 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Figure 6. Contour plots of (a) harmonic and (b) anharmonic PES of the Zundel cluster as a function of

(PT motion) and

(WZ1-WZ2 stretching motion). The contours are drawn

every 0.001 Hartree.

two reasons. The first is due to a change in the character of vibrational coordinates of 7 and 8. The wagging and PT motions are somehow mixed in the harmonic approximation as shown in Figure 5 (a). Thus, a strong IR intensity associated with the PT motion is distributed to 7 and 8, rendering an intensity of 1325 and 1848 km mol-1, respectively. On the other hand, the wagging and PT motions are separated in optimized coordinates (see Figure 5 (b)). Consequently, the intensity of 8 increases to 2571 km mol-1, while that of 7 decreases to 155 km mol-1. The other reason is that 8 makes a blue-shift from 846 to 982 cm-1, so that 8 overlaps with 7. The blue shift implies a steep rise of the PES in a direction of the PT motion ( calculation using

gives an even larger blue shift to 1313 cm-1. The frequency relaxes back to

982 cm-1 mainly due to a WZ1-WZ2 stretching motion ( without

). In fact, a one-dimensional

); an 11 dimensional calculation

gives 1187 cm-1. Figure 6 (a) and (b) compares a section of harmonic and

anharmonic PESs, respectively, as a function of only a stiff curve along

and

. The anharmonic PES shows not

but also strong deformation in a coupling region. In a positive region

ACS Paragon Plus Environment

29

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

of

Page 30 of 54

, where WZ1 and WZ2 are stretched out, the proton is destabilized in the middle of two

water molecules and prefers to attach to either WZ1 or WZ2. Such interaction is obviously missing in the harmonic PES. In a region above 1000 cm-1, the 12D/3D×2 spectrum appears similar to the harmonic one. However, the band positions shift in a complex way. νa, νz, and νs are red-shifted by 100 – 150 cm-1 (3 – 4% of the frequency), while the shift in 1 and 2 is more significant, as many as 200 cm1

(6 % of the frequency). The shift in 3 – 6 is relatively smaller (1 – 3 % of the frequency); yet 3

makes a blue-shift. Therefore, the common practice to multiply the harmonic frequency by a constant factor does not work in this case. Including the rocking and wagging vibrations of WZ1 and WZ2 significantly affects the overall appearance of the spectrum. In the 16D/3D×2 model, 2 and 8 are in strong resonance and distribute their intensity to background overtones and combination bands. Consequently, other bands appear relatively larger than those in the 12D/3D×2 model. The wavefunction analysis (see Table S11) reveals that the PT motion is in strong resonance with a combination band of the wagging vibrations and the WZ1-WZ2 stretching vibration, giving rise to a sharp doublet at 994 and 1044 cm-1 (8 and 11+9, respectively). This characteristic signature was also reported in the previous studies on a Zundel monomer (H5O2+) by MCTDH (calculated at 918 and 1033 cm-1)37, 39

and by experiment (observed at 928 and 1047 cm-1).10 Furthermore, the experiment observed a

strong band at 1763 cm-1 and two relatively small bands at 1470 and 1878 cm-1 for the Zundel monomer. The 16D/3D×2 model yields a strong band at 1710 cm-1 (3) and a series of weaker bands around 1410 cm-1 (9+8). The assignment of these bands is consistent with the MCTDH result. However, the present calculation did not yield any band around 1878 cm-1, which may be

ACS Paragon Plus Environment

30

Page 31 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

a difference between H5O2+ and H5O2+(H2O)2. The result attests to the quality of the present theoretical treatment of the Zundel core. The 26D model shows a further change in the spectrum, in particular, in a region below 2000 cm-1. An enlarged view in Figure 7 (a) shows a series of weak bands in a range of 1100-1300 cm1

and a strong band around 1900 cm-1. It is also notable that 3 significantly loses its intensity.

Although the wavefunction analysis was attempted for these states, the character was unclear due to extremely strong mixing of VSCF configuration functions (cf. Table S12). Instead, we have performed additional reduced dimensional calculations to reveal the motion responsible for the

(a)

26D 18D/3D×2

8, 11+9 11×2 12

7

9+8, 6, 5

11

4

pg

3

9 10



(b)

16D/3D×2 22D

8, 11+9

3

11 9 10



7

6 11×2 9+8

11+10×2 5

4

7×2

wavenumber / cm-1

Figure 7. The IR spectrum of the Zundel cluster obtained by VQDPT2 based on (a) 26D and 18D/3D×2 models, and (b) 16D/3D×2 and 22D models. The 18D model adds stretching motions of W1 and W2 to the 16D model. The 22D model couples the motions of a Zundel core, W1, and W2.

ACS Paragon Plus Environment

31

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 32 of 54

spectral change. As shown in Fig. 7 (a), the 26D spectrum is well reproduced by an 18D/3D×2 model, in which W1-WZ1 and W2-WZ2 stretching motions are added to the 16D model. Also shown in Fig. 7 (b) is the result of a 22D model, which couples the motion of the Zundel core (16D) and the water molecules (3D×2). The spectrum shows no essential change from the 16D/3D×2 model, indicating that the coupling between the Zundel core and the internal vibrations of W1 and W2 is negligibly small. From these results, it is plausible to assign the strong band around 1900 cm-1 to a progression of W1-WZ1 and W2-WZ2 stretching motions upon exciting 3. The origin of weak bands in a range of 1100-1300 is less clear; yet it may be a progression for the PT motion (e.g., 8, 11+9), a red-shift of 9+8, 6 and 5, or both. Although many combination bands are seen below 2000 cm-1, the spectrum remains blank in a range of 2000-2700 cm-1. This feature contrasts to the spectrum of the Eigen cluster, where many combination bands appear in this range. In a high frequency region, the OH stretching bands remain essentially the same as in the 16D/3D×2 model.

ACS Paragon Plus Environment

32

Page 33 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

4. DISCUSSION Having established the IR spectrum of Eigen and Zundel clusters, we now discuss on the coexistence of two clusters. Figure 8 compares the experimental spectrum with those calculated for Eigen and Zundel clusters. In the middle panel, the experiment shows a strong band at 2655 cm-1, which agrees with 4 of the Eigen cluster (at 2699 cm-1), but deviates from 1 and 2 of the Zundel cluster by about 200 cm-1. Furthermore, the high frequency region, enlarged in the right panel, shows that the experimental spectrum gives two bands without νZ of the Zundel cluster. Note that a slight splitting seen in the experiment is caused by a tagging molecule, D2, which breaks the C3 symmetry of the Eigen cluster.25 The characteristic OH stretching bands indicate that the Eigen cluster is a major product in the experiment. However, the Eigen cluster agrees poorly with the experiment in a low frequency range, as shown in the left panel of Figure 8, missing a strong band at 1765 cm-1 and a broad feature around 1045 cm-1. It is noteworthy that the Zundel cluster gives 3 and pg (at 1742 and 1906 cm-1, respectively), and a strong, broad band of a proton transfer motion (8, 11+9, etc.), which agree surprisingly well with the experiment. In fact, the previous MD work conjectured the existence of the Zundel cluster on this basis.26 To test this hypothesis, we have constructed a spectrum by mixing the Eigen and Zundel spectra by 9:1 and 1:1, which is shown in Fig 8 (b) and (c), respectively. The mixing ratio of 1:1 reproduces the experiment in the low frequency range fairly well. However, the mixed spectrum gives notable peaks of 1/2 and νZ of the Zundel cluster in the OH stretching range, which are absent in the experiment. On the other hand, the ratio of 9:1 reduces these OH stretching bands, but also makes the bands at 1765 and 1045 cm-1 too weak. Therefore, the present calculation does not support the mixing of Eigen and Zundel clusters. The

ACS Paragon Plus Environment

33

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1045

1840 1765 1905 1615

(a) Exp.

Page 34 of 54

2655

(b) E:Z=9:1

(c) E:Z=1:1

8

pg 4

7 8

3

8⋅ 2 8+7 6 5

(d) Zundel

2, 1

8

(e) Eigen

νZ

νZ

4

8 7

wavenumber / cm-1

Figure 8. Comparison of (a) the experimental spectrum with theoretical ones for (d) Zundel and (e) Eigen clusters, and a mixture of Eigen and Zundel by (b) 9:1 and (c) 1:1. The left and right panels show enlarged view in a range of 800-2400 cm-1 and 3600-3900 cm-1, respectively.

result is consistent with a recent double resonance experiment,14 which has indicated that the spectrum originates from a single isomer. If the mixing of Eigen and Zundel is not the case, the discrepancy in the low frequency range should be attributed to shortcomings of the present theoretical treatment: the neglect of nine low frequency modes, the accuracy of the PES and DMS (both in terms of the functional form and the electronic structure level), the truncation of a higher level of correlation in the vibrational wavefunction, the lack of finite temperature effects, the lack of a tagging species, etc. The band at 1045 cm-1 seems to originate from 7 and/or 8, though the mechanism of broadening needs to

ACS Paragon Plus Environment

34

Page 35 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

be clarified. The origin of the band at 1765 cm-1 remains an open question. We note that a protonated pentamer (n=5), which takes an Eigen form, also gives a strong band around 1900 cm-1 and a shoulder around 1800 cm-1. Interestingly, the Zundel tetramer (n=4) and the Eigen pentamer (n=5) have some similarities. The structure of the pentamer is related to the tetramer by adding a water molecule to a dangling OH bond of WZ1. The additional water molecule breaks the symmetry and shifts the proton (P) to WZ1 to form an Eigen core; yet it is somehow distorted due to a strong hydrogen bond between P and WZ2 compared to the other two. Consequently, the OH stretching bands of the Eigen core are observed at two distinctly different frequencies: 1889 cm-1 for the WZ1-P bond, and 2837/2882 cm-1 for the other two.21 Furthermore, the pentamer yields hardly any band in a range of 1900-2800 cm-1. It is intriguing that the pattern resembles the spectrum calculated for the Zundel tetramer. This observation implies that these bands may arise from a distorted Eigen form. Although no such minima is found on the PES, dynamical fluctuations that are induced upon inclusion of low frequency modes and finite temperature effects may enable to sample such distorted structures. A thorough theoretical investigation will be the scope of future works. Finally, we comment on future perspectives toward larger clusters. As for Eigen clusters, an extension to the 9D/3D×3 model is a promising direction, because it only requires to solve a 9dimensional problem for H3O+ and 3-dimensional problems for each water molecules, and thus is readily extendable to large clusters. The local coordinates for H3O+ and H2O can be obtained using localization techniques,66-70 bypassing the oc-VSCF step. Note that VCI based on such local coordinates has been recently carried out to predict the IR spectrum of H3O+(H2O)3Cl-.71 Similarly, Zundel clusters may be treated by extending the models, 16D/3D×2 or 12D/3D×2. However, the treatment of the H5O2+ moiety needs a careful attention, since a regular

ACS Paragon Plus Environment

35

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 36 of 54

perturbation theory severely fails whereas VCI is costly. In the present study, the decoupling of wagging and PT motions (see Fig. 5) seems to be of great importance for the perturbative expansion in VQDPT2 to work reasonably well, though a stringent test is still needed in comparison with an experimental spectrum of Zundel clusters (e.g., n=2 and 6). The vibrational calculation based on local coordinates combined with a fragment approach for the electronic structure calculation72-73 will make feasible to compute the IR spectrum of large clusters with high accuracy and bring about much insight into the motion of proton in cluster and solution.

ACS Paragon Plus Environment

36

Page 37 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

ASSOCIATED CONTENT Supporting Information. The Supporting Information is available free of charge on the ACS Publications website at DOI: xxx. Cartesian coordinates (Table S1) and harmonic frequencies (Table S2) of the Eigen and Zundel clusters obtained by B3LYP/aug-cc-pVTZ, CCSD(T)-F12/aug-cc-pVDZ, and the multiresolution PES. Normal coordinates not considered in this work (Figure S1). Plots of the IR spectrum of the Eigen cluster using FWHM of 5, 10, and 20 cm-1 for the Lorentz functions to construct the spectrum (Figure S2). Computational details on oc-VSCF calculations. The number of important coupling terms in normal and optimized coordinates (Table S3). Comparison of normal and optimized coordinates: Normal and optimized coordinates (Figure S3 and S4) and harmonic frequency (Table S4) of the Eigen cluster, Weight of H3O+ and water molecules (Table S5), Normal and optimized coordinates (Figure S5 and S6) and harmonic frequency (Table S6) of the Zundel cluster, Weight of H5O2+ and water molecules (Table S7). Vibrational wavefunction analysis: 9D/3D×3, 18D, and 24D (Table S8, S9, and S10) for the Eigen cluster, and 12D/3D×2, 16D/3D×2, and 26D (Table S11, S12, and S13) for the Zundel cluster.

AUTHOR INFORMATION Corresponding Author *[email protected]

ACS Paragon Plus Environment

37

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 38 of 54

ACKNOWLEDGMENT We are grateful to Mark Johnson for providing us with the experimental data for Figure 2 and 8. We are also grateful to Mark Johnson and Ken Jordan for their comments on this work. BT is supported by the Special Postdoctoral Researchers Program at RIKEN. This research is partially supported by the Center of innovation Program from Japan Science and Technology Agency, JST, and JSPS KAKENHI Grant No. JP16H00857 (to K. Y.). We used computational resources provided by the HPCI System Research Project (Project ID: hp140105, hp150023), and the RIKEN Integrated Cluster of Clusters (RICC).

ACS Paragon Plus Environment

38

Page 39 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

REFERENCES 1. Eigen, M.; de Maeyer, L., Self-Dissociation and Protonic Charge Transport in Water and Ice. Proc. R. Soc. London 1958, 247, 505-533. 2. Eigen, M., Proton Transfer, Acid-Base Catalysis, and Enzymatic Hydrolysis. Angew. Chem. Int. Ed. 1964, 3, 1-19. 3. Zundel, G., Hydrogen Bonds with Large Proton Polarizability and Proton Transfer Processes in Electrochemistry and Biology. Adv. Chem. Phys. 2000, 111, 1-217. 4. Marx, D.; Tuckerman, M. E.; Hutter, J.; Parrinello, M., The Nature of the Hydrated Excess Proton in Water. Nature 1999, 397, 601-604. 5. Marx, D., Proton Transfer 200 Years after von Grotthuss: Insights from Ab Initio Simulations. ChemPhysChem 2006, 7, 1848-1870. 6. Berkelbach, T. C.; Lee, H.-S.; Tuckerman, M. E., Concerted Hydrogen-Bond Dynamics in the Transport Mechanism of the Hydrated Proton: A First-Principles Molecular Dynamics Study. Phys. Rev. Lett. 2009, 103, 238302. 7. Asmis, K. R.; Pivonka, N. L.; Santambrogio, G.; Brümmer, M.; Kaposta, C.; Neumark, D. M.; Wöste, L., Gas-Phase Infrared Spectrum of the Protonated Water Dimer. Science 2003, 299, 1375-1377. 8. Shin, J.-W.; Hammer, N. I.; Diken, E. G.; Johnson, M. A.; Walters, R. S.; Jaeger, T. D.; Duncan, M. A.; Christie, R. A.; Jordan, K. D., Infrared Signature of Structures Associated with the H+(H2O)n (n=6 to 27) Clusters. Science 2004, 304, 1137-1140. 9. Headrick, J. M.; Diken, E. G.; Walters, R. S.; Hammer, N. I.; Christie, R. A.; Cui, J.; Myshakin, E. M.; Duncan, M. A.; Johnson, M. A.; Jordan, K. D., Spectral Signatures of Hydrated Proton Vibrations in Water Clusters. Science 2005, 308, 1765-1769. 10. Hammer, N. I.; Diken, E. G.; Roscioli, J. R.; Johnson, M. A.; Myshakin, E. M.; Jordan, K. D.; McCoy, A. B.; Huang, X.; Bowman, J. M.; Carter, S., The Vibrational Predissociation Spectra of the H5O2+ RGn (RG=Ar, Ne) Clusters: Correlation of the Solvent Perturbations in the Free OH and Shared Proton Transitions of the Zundel Ion. J. Chem. Phys. 2005, 122, 244301. 11. Olesen, S. G.; Guasco, T. L.; Roscioli, J. R.; Johnson, M. A., Tuning the Intermolecular Proton Bond in the H5O2+ 'Zundel Ion' Scaffold. Chem. Phys. Lett. 2011, 509, 89-95.

ACS Paragon Plus Environment

39

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 40 of 54

12. McCoy, A. B.; Guasco, T. L.; Leavitt, C. M.; Olesen, S. G.; Johnson, M. A., Vibrational Manifestations of Strong Non-Condon Effects in the H3O+·X3 (X = Ar, N2, CH4, H2O) Complexes: A Possible Explanation for the Intensity in the “Association Band” in the Vibrational Spectrum of Water. Phys. Chem. Chem. Phys. 2012, 14, 7205-7214. 13. Fournier, J. A.; Johnson, C. J.; Wolke, C. T.; Weddle, G. H.; Wolk, A. B.; Johnson, M. A., Vibrational Spectral Signature of the Proton Defect in the Three-Dimensional H+(H2O)21 Cluster. Science 2014, 344, 1009-1012. 14. Fournier, J. A.; Wolke, C. T.; Johnson, M. A.; Odbadrakh, T. T.; Jordan, K. D.; Kathmann, S. M.; Xantheas, S. S., Snapshots of Proton Accommodation at a Microscopic Water Surface: Understanding the Vibrational Spectral Signatures of the Charge Defect in Cryogenically Cooled H+(H2O)n=2–28 Clusters. J. Phys. Chem. A 2015, 119, 9425-9440. 15. Douberly, G. E.; Walters, R. S.; Cui, J.; Jordan, K. D.; Duncan, M. A., Infrared Spectroscopy of Small Protonated Water Clusters, H+(H2O)n (n = 2-5): Isomers, Argon Tagging, and Deuteration. J. Phys. Chem. A 2010, 114, 4570-4579. 16. Miyazaki, M.; Fujii, A.; Ebata, T.; Mikami, N., Infrared Spectroscopic Evidence for Protonated Water Clusters Forming Nanoscale Cages. Science 2004, 304, 1134-1137. 17. Mizuse, K.; Fujii, A., Infrared Photodissociation Spectroscopy of H+(H2O)6·Mm (M = Ne, Ar, Kr, Xe, H2, N2, and CH4): Messenger-Dependent Balance between H3O+ and H5O2+ Core Isomers. Phys. Chem. Chem. Phys. 2011, 13, 7129-7135. 18. Mizuse, K.; Fujii, A., Tuning of the Internal Energy and Isomer Distribution in Small Protonated Water Clusters H+(H2O)4-8: An Application of the Inert Gas Messenger Technique. J. Phys. Chem. A 2012, 116, 4868-4877. 19. Heine, N.; Fagiani, M. R.; Rossi, M.; Wende, T.; Berden, G.; Blum, V.; Asmis, K. R., Isomer-Selective Detection of Hydrogen-Bond Vibrations in the Protonated Water Hexamer. J. Am. Chem. Soc. 2013, 135, 8266-8273. 20. Heine, N.; Fagiani, M. R.; Asmis, K. R., Disentangling the Contribution of Multiple Isomers to the Infrared Spectrum of the Protonated Water Heptamer. J. Phys. Chem. Lett. 2015, 6, 2298-2304. 21. Fagiani, M. R.; Knorke, H.; Esser, T.; Heine, N.; Wolke, C. T.; Gewinner, S.; Schöllkopf, W.; Gaigeot, M.-P.; Spezia, R.; Johnson, M. A.; Asmis, K. R., Gas Phase Vibrational

ACS Paragon Plus Environment

40

Page 41 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

Spectroscopy of the Protonated Water Pentamer: The Role of Isomers and Nuclear Quantum Effects. Phys. Chem. Chem. Phys. 2016, 18, 26743-26754. 22. Wolke, C. T.; Fournier, J. A.; Dzugan, L. C.; Fagiani, M. R.; Odbadrakh, T. T.; Knorke, H.; Jordan, K. D.; McCoy, A. B.; Asmis, K. R.; Johnson, M. A., Spectroscopic Snapshots of the Proton-Transfer Mechanism in Water. Science 2016, 354, 1131-1135. 23. Thämer, M.; Marco, L. D.; Ramasesha, K.; Mandal, A.; Tokmakoff, A., Ultrafast 2D IR Spectroscopy of the Excess Proton in Liquid Water. Science 2015, 350, 78-82. 24. Schwarz, H. A., Gas Phase Infrared Spectra of Oxonium Hydrate Ions from 2 to 5 μ. J. Chem. Phys. 1977, 67, 5525-5534. 25. Okumura, M.; Yeh, L. I.; Myers, J. D.; Lee, Y. T., Infrared Spectra of the Cluster Ions H7O3+⋅H2 and H9O4+⋅H2. J. Chem. Phys. 1986, 85, 2328-2329. 26. Kulig, W.; Agmon, N., Both Zundel and Eigen Isomers Contribute to the IR Spectrum of the Gas-Phase H9O4+ Cluster. J. Phys. Chem. B 2014, 118, 278-286. 27. Barone, V., Anharmonic Vibrational Properties by a Fully Automated Second-Order Perturbative Approach. J. Chem. Phys. 2005, 122, 014108. 28. Park, M.; Shin, I.; Singh, N. J.; Kim, K. S., Eigen and Zundel Forms of Small Protonated Water Clusters: Structures and Infrared Spectra. J. Phys. Chem. A 2007, 111, 10692-10702. 29. Norris, L. S.; Ratner, M. A.; Roitberg, A. E.; Gerber, R. B., Møller-Plesset Perturbation Theory Applied to Vibrational Problems. J. Chem. Phys. 1996, 105, 11261-11267. 30. Chaban, G. M.; Jung, J. O.; Gerber, R. B., Ab Initio Calculation of Anharmonic Vibrational States of Polyatomic Systems: Electronic Structure Combined with Vibrational SelfConsistent Field. J. Chem. Phys. 1999, 111, 1823-1829. 31. Chaban, G. M.; Jung, J. O.; Gerber, R. B., Anharmonic Vibrational Spectroscopy of Hydrogen-Bonded Systems Directly Computed from Ab Initio Potential Surfaces: (H2O)n, n=2, 3; Cl-(H2O)n, n=1, 2; H+(H2O)n, n=1, 2; H2O-CH3OH. J. Phys. Chem. A 2000, 104, 2772-2779. 32. Torrent-Sucarrat, M.; Anglada, J. M., Anharmonicity and the Eigen-Zundel Dilemma in the IR Spectrum of the Protonated 21 Water Cluster. J. Chem. Theory Comput. 2011, 7, 467472.

ACS Paragon Plus Environment

41

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 42 of 54

33. Bowman, J. M.; Carrington, T.; Meyer, H.-D., Variational Quantum Approaches for Computing Vibrational Energies of Polyatomic Molecules. Mol. Phys. 2008, 106, 21452182. 34. Bowman, J. M., The Self-Consistent-Field Approach to Polyatomic Vibrations. Acc. Chem. Res. 1986, 19, 202-208. 35. Christoffel, K. M.; Bowman, J. M., Investigations of Self-Consistent Field, SCF CI and Virtual State Configuration-Interaction Vibrational Energies for a Model 3-mode System. Chem. Phys. Lett. 1982, 85, 220-224. 36. Bowman, J. M., Self-Consistent Field Energies and Wavefunctions for Coupled Oscillators. J. Chem. Phys. 1978, 68, 608-610. 37. Vendrell, O.; Gatti, F.; Meyer, H.-D., Dynamics and Infrared Spectroscopy of the Protonated Water Dimer. Angew. Chem., Int. Ed. 2007, 46, 6918-6921. 38. Vendrell, O.; Gatti, F.; Lauvergnat, D.; Meyer, H.-D., Full-Dimensional (15-Dimensional) Quantum-Dynamical Simulation of the Protonated Water Dimer. I. Hamiltonian Setup and Analysis of the Ground Vibrational State. J. Chem. Phys. 2007, 127, 184302. 39. Vendrell, O.; Gatti, F.; Meyer, H.-D., Full-Dimensional (15-Dimensional) QuantumDynamical Simulation of the Protonated Water Dimer. II. Infrared Spectrum and Vibrational Dynamics. J. Chem. Phys. 2007, 127, 184303. 40. Beck, M. H.; Jäckle, A.; Worth, G. A.; Meyer, H.-D., The Multiconfiguration TimeDependent Hartree (MCTDH) Method: a Highly Efficient Algorithm for Propagating Wavepackets. Phys. Rep. 2000, 324, 1-105. 41. Yagi, K.; Keçeli, M.; Hirata, S., Optimized Coordinates for Anharmonic Vibrational Structure Theories. J. Chem. Phys. 2012, 137, 204118. 42. Thomsen, B.; Yagi, K.; Christiansen, O., Optimized Coordinates in Vibrational Coupled Cluster Calculations. J. Chem. Phys. 2014, 140, 154102. 43. Yagi, K.; Hirata, S.; Hirao, K., Vibrational Quasi-Degenerate Perturbation Theory: Applications to Fermi Resonance in CO2, H2CO, and C6H6. Phys. Chem. Chem. Phys. 2008, 10, 1781-1788. 44. Yagi, K.; Otaki, H., Vibrational Quasi-Degenerate Perturbation Theory with Optimized Coordinates: Applications to Ethylene and Trans-1,3-Butadiene. J. Chem. Phys. 2014, 140, 084113.

ACS Paragon Plus Environment

42

Page 43 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

45. Yagi, K.; Hirata, S.; Hirao, K., Multiresolution Potential Energy Surfaces for Vibrational State Calculations. Theor. Chem. Acc. 2007, 118, 681-691. 46. Seidler, P.; Kaga, T.; Yagi, K.; Christiansen, O.; Hirao, K., On the Coupling Strength in Potential Energy Surfaces for Vibrational Calculations. Chem. Phys. Lett. 2009, 483, 138142. 47. Otaki, H.; Yagi, K.; Ishiuchi, S.; Fujii, M.; Sugita, Y., Anharmonic Vibrational Analyses of Pentapeptide Conformations Explored with Enhanced Sampling Simulations. J. Phys. Chem. B 2016, 120, 10199-10213. 48. Becke, A. D., Density-Functional Thermochemistry. III. The Role of Exact Exchange. J. Chem. Phys. 1993, 98, 5648-5652. 49. Lee, C.; Yang, W.; Parr, R. G., Development of the Colle-Salvetti Correlation-Energy Formula into a Functional of the Electron Density. Phys. Rev. B 1988, 37, 785-789. 50. Kendall, R. A.; Dunning, T. H., Jr.; Harrison, R. J., Electron Affinities of the First-Row Atoms Revisited. Systematic Basis Sets and Wave Functions. J. Chem. Phys. 1992, 96, 6796-6806. 51. Yagi, K.; Hirao, K.; Taketsugu, T.; Schmidt, M. W.; Gordon, M. S., Ab Initio Vibrational State Calculations with a Quartic Force Field: Applications to H2CO, C2H4, CH3OH, CH3CCH, and C6H6. J. Chem. Phys. 2004, 121, 1383-1389. 52. Yagi, K.; Taketsugu, T.; Hirao, K.; Gordon, M. S., Direct Vibrational Self-Consistent Field Method: Applications to H2O and H2CO. J. Chem. Phys. 2000, 113, 1005-1017. 53. Salmi, T.; Salli, E.; Halonen, L., A Nine-Dimensional Calculation of the Vibrational OH Stretching and HOH Bending Spectrum of the Water Trimer. J. Phys. Chem. A 2012, 116, 5368-5374. 54. 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. 55. Nibbering, E. T. J.; Elsaesser, T., Ultrafast Vibrational Dynamics of Hydrogen Bonds in the Condensed Phase. Chem. Rev. 2004, 104, 1887-1914. 56. Carter, S.; Culik, S. J.; Bowman, J. M., Vibrational Self-Consistent Field Method for ManyMode Systems: A New Approach and Application to the Vibrations of CO Adsorbed on Cu (100). J. Chem. Phys. 1997, 107, 10458-10469.

ACS Paragon Plus Environment

43

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 44 of 54

57. Light, J. C.; Hamilton, I. P.; Lill, J. V., Generailzed Discrete Variable Approximation in Quantum Mechanics. J. Chem. Phys. 1985, 82, 1400-1409. 58. Adler, T. B.; Knizia, G.; Werner, H.-J., A Simple and Efficient CCSD(T)-F12 Approximation. J. Chem. Phys. 2007, 127, 221106. 59. Knizia, G.; Adler, T. B.; Werner, H.-J., Simplified CCSD(T)-F12 Methods: Theory and Benchmarks. J. Chem. Phys. 2009, 130, 054104. 60. Werner, H.-J.; Knizia, G.; Manby, F. R., Explicitly Correlated Coupled Cluster Methods with Pair-Specific Geminals. Mol. Phys. 2011, 109, 407-417. 61. Frisch, M. J.; Trucks, G. W.; Schlegel, H. B.; Scuseria, G. E.; Robb, M. A.; Cheeseman, J. R.; Scalmani, G.; Barone, V.; Mennucci, B.; Petersson, G. A., et al., Gaussian 09 Revision C.01. Gaussian, Inc. Wallingford CT, 2009. 62. Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M., Molpro: A GeneralPurpose Quantum Chemistry Program Package. WIREs Comput. Mol. Sci. 2012, 2, 242-253. 63. Werner, H.-J.; Knowles, P. J.; Knizia, G.; Manby, F. R.; Schütz, M.; others, MOLPRO, version 2012.1, a package of ab initio programs, see http://www.molpro.net. (Date of access: 3-Mar-2017). 64. Christiansen, O., Møller-Plesset Perturbation Theory for Vibrational Wave Functions. J. Chem. Phys. 2003, 119, 5773-5781. 65. Yagi, K.; Hirata, S.; Hirao, K., Efficient Configuration Selection Scheme for Vibrational Second-Order Perturbation Theory. J. Chem. Phys. 2007, 127, 034111. 66. Wang, Y.; Bowman, J. M., Towards an Ab Initio Flexible Potential for Water, and PostHarmonic Quantum Vibrational Analysis of Water Clusters. Chem. Phys. Lett. 2010, 491, 110. 67. Jacob, C. R.; Reiher, M., Localizing Normal Modes in Large Molecules. J. Chem. Phys. 2009, 130, 084106. 68. Panek, P. T.; Jacob, C. R., On the Benefits of Localized Modes in Anharmonic Vibrational Calculations for Small Molecules. J. Chem. Phys. 2016, 144, 164111. 69. Cheng, X.; Steele, R. P., Efficient Anharmonic Vibrational Spectroscopy for Large Molecules Using Local-Mode Coordinates. J. Chem. Phys. 2014, 141, 104105. 70. Cheng, X.; Talbot, J. J.; Steele, R. P., Tuning Vibrational Mode Localization with Frequency Windowing. J. Chem. Phys. 2016, 145, 124112.

ACS Paragon Plus Environment

44

Page 45 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

71. Mancini, J. S.; Bowman, J. M., Isolating the Spectral Signature of H3O+ in the Smallest Droplet of Dissociated HCl Acid. Phys. Chem. Chem. Phys. 2015, 17, 6222-6226. 72. Hirata, S., Fast Electron-Correlation Methods for Molecular Crystals: An Application to the α, β1, and β2 Modifications of Solid Formic Acid. J. Chem. Phys. 2008, 129, 204104. 73. Gordon, M. S.; Fedrov, D. G.; Pruitt, S. R.; Slipchenko, L. V., Fragmentation Methods: A Route to Accurate Calculations on Large Systems. Chem. Rev. 2012, 112, 632-672.

ACS Paragon Plus Environment

45

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Page 46 of 54

Table of Contents Graphic

ACS Paragon Plus Environment

46

Page 47 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(b)

(a)

P

W3 W2

W1

WZ1

WZ2 W2

W1

Figure 1. (a) Eigen and (b) trans-Zundel isomers of H+(H2O)4. (84×35 mm)

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(a)

(b)

Exp.

H2O (3D×3)

1840 1765 1905 1615

1045

VQDPT2 (24D) OH stretch

bending

H3O+ (9D)

8

9

7

VQDPT2 (18D) OH stretch×3

bending×2

umbrella

8

9

7

2340 2300 2240 2485

Intermolecular motions (6D)

9

7

2990

4 7+6 9+6 8+6 7+5, 8×3, 3 8+7 9+5 8+5 pg2 pg1 6 5 8×2 pg3 4 7+6 7+5 9+6 8+6 8×3 9+5 8+5 8+7 6 5 8×2 3 4

6 5 8×2

8+7 9+5 8+5

Harmonic

9 rocking×3

8

2655

α

VQDPT2 (9D/3D×3)

frustrated rotation×3

Page 48 of 54

8

3723/3733 3636/3648

1 2

2

1

7+5 8×3 2

3

1

4 ×1/2

7

3

6 5

stretch×3

2

1

wavenumber / cm-1

Figure 2. (a) Optimized coordinates of the Eigen cluster. For a full list see Figure S4 in SI. (b) IR spectrum of the Eigen cluster obtained by VQDPT2 calculations based on 24D, 18D, and 9D/3D×3 models, and the harmonic approximation together with the experiment.13 The fundamental transition is labeled with a number. The labels, m×n, m+m’, and pg, denote an (n – 1)-th overtone of mode m, a combination tone of mode m and m’, and progression bands, respectively. (177×126 mm)

ACS Paragon Plus Environment

Page 49 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

2655

(a) Exp. 2240 α 2300 2340

2485

(b) VQDPT2 (24D) 9+6 9+5

2990

4 7+6

8+6 8+5 4

(c) VQDPT2 (18D) 9+6 9+5

8+6 8+5

7+6

(d) VQDPT2 (9D/3D×3) 9+5

7+5, 8×3, 3 pg2 7+5 8×3 3

4

8+5

7+5 8×33

wavenumber / cm-1

Figure 3. Enlarged view of Figure 2 (b) in a range of 2200 – 3000 cm-1. (a) The experimental spectrum, and the VQDPT2 spectra obtained from (b) 24D, (c) 18D, and (d) 9D/3D×3 models. (84×89 mm)

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

(a)

Page 50 of 54

(b) H2O (3D×2)

OH stretch×2

bending

H5O2+ (16D)

11 910

7

9+8, 6, 5

4 3

8, 11+9

bending×2

OH stretch×4

11 910

7

pg

3 11×2 6 9+8

11+10×2 54

7×2

stretch

proton transfer

9+8

wagging×4

3

VQDPT2 (12D/3D×2) 1 νs 2

rocking×2

8

Intermolecular motions (4D) 10 9 rocking×2

654

7

9

VQDPT2 2 1, 10×2+3 (16D/3D×2) 10×2+9+8 νs νz νa 8+7×2 2

8 proton×2

VQDPT2 (26D) νa νz νs

2, 1

11×2

8, 11+9

7 654

stretch×2

3

νz νa

Harmonic 1

νs

νz νa

wavenumber / cm-1

Figure 4. (a) Optimized coordinates of the Zundel cluster. For a full list see Figure S6 in SI. (b) IR spectrum of the Zundel cluster obtained by VQDPT2 calculations based on 26D, 16D/3D×2, and 12D/3D×2 models, and the harmonic approximation. The labels are defined in the same way as in Figure 2. (177×121 mm)

ACS Paragon Plus Environment

Page 51 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(a)

7

8

1005 cm-1 1325 km mol-1

846 cm-1 1848 km mol-1

(b) 7

8

976 cm-1 155 km mol-1

982 cm-1 2571 km mol-1

Figure 5. (a) Normal and (b) optimized coordinates corresponding to the wagging motion (7) and the PT motion (8). The associated harmonic frequencies and IR intensities are indicated. (84×81 mm)

ACS Paragon Plus Environment

The Journal of Physical Chemistry

40

Q!17 / bohr emu1/2

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

20

40

(a)

20

0

0

-20

-20

-40 -60

-40

-20

0

20

40

60

Page 52 of 54

(b)

-40 -60

-40

Q!18 / bohr emu1/2

-20

0

20

40

60

Q!18 / bohr emu1/2

Figure 6. Contour plots of (a) harmonic and (b) anharmonic PES of the Zundel cluster as a function of Q!17 (PT motion) and Q!18 (WZ1-WZ2 stretching motion). The contours are drawn every 0.001 Hartree.

(84×34 mm)

ACS Paragon Plus Environment

Page 53 of 54

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

The Journal of Physical Chemistry

(a)

26D 18D/3D×2

8, 11+9 11×2 12

7

11

9+8, 6, 5

4

3

pg

9 10

(b)

16D/3D×2 22D

8, 11+9

3 11×2 7

11 910

6

11+10×2

9+8 5

4

7×2

wavenumber / cm-1

Figure 7. The IR spectrum of the Zundel cluster obtained by VQDPT2 based on (a) 26D and 18D/3D×2 models, and (b) 16D/3D×2 and 22D models. The 18D model adds stretching motions of W1 and W2 to the 16D model. The 22D model couples the motions of a Zundel core, W1, and W2. (84×96 mm)

ACS Paragon Plus Environment

The Journal of Physical Chemistry

1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1045

1840 1765 1905 1615

(a) Exp.

Page 54 of 54

2655

(b) E:Z=9:1

(c) E:Z=1:1

8 4 7 8

3

pg

8" 2 8+7 6 5

(d) Zundel

2, 1

8

!Z

!Z

4

(e) Eigen

8 7

wavenumber / cm-1

Figure 8. Comparison of (a) the experimental spectrum with theoretical ones for (d) Zundel and (e) Eigen clusters, and a mixture of Eigen and Zundel by (b) 9:1 and (c) 1:1. The left and right panels show enlarged view in a range of 800-2400 cm-1 and 3600-3900 cm-1, respectively. (177×114 mm)

ACS Paragon Plus Environment