Nuclear Quantum Effects and Thermodynamic Properties for Small

3 hours ago - This agreement for the solvation energies motivated us to extend the current study to larger (N = 6–21) clusters to explore their ther...
0 downloads 7 Views 2MB Size
Article Cite This: J. Phys. Chem. A XXXX, XXX, XXX−XXX

pubs.acs.org/JPCA

Nuclear Quantum Effects and Thermodynamic Properties for Small (H2O)1−21X− Clusters (X− = F−, Cl−, Br−, I−) Joel D. Mallory and Vladimir A. Mandelshtam* Department of Chemistry, University of California, Irvine, California 92697, United States S Supporting Information *

ABSTRACT: We carried out accurate diffusion Monte Carlo (DMC) studies for small (H2O)NX− clusters (N = 1−5; X− = F−, Cl−, Br−, I−) and their D2O isotopologues. We found remarkably good agreement (i.e., ∼0.1 kcal/mol or better) with no exceptions between the DMC solvation energies and the corresponding harmonic approximation (HA) estimates, due, apparently, to massive error cancellations. This is surprising, in particular, because HA does not account for a substantial (i.e., ∼ 3%) increase of the mean O−O distances, caused by the anharmonicity in conjunction with the nuclear quantum effects, although the other distances in the system are affected to a much lesser extent. This agreement for the solvation energies motivated us to extend the current study to larger (N = 6−21) clusters to explore their thermodynamic properties using the harmonic superposition method (HSM). The HSM results for the solvation free energies in turn reveal that at finite temperatures the nuclear quantum effects (including the isotope effects) in these systems are miniscule.



H2OX− dimers.56,61−63 However, not until recently have a few studies emerged that explicitly probe the quantum effects observed in positional isotope exchange reactions for the singly deuterated isotopomers of H2OX− systems as individual dimers64 and in the condensed phase65 using finite temperature path-integral molecular dynamics (PIMD) simulations. Notably, a similar PIMD analysis was carried out for protonated HX(H2O)N clusters,66 and a new joint experimental and theoretical study has also centered on the analysis of IR spectra for deuterated isotopologues of (H2O)2I− clusters and their various isotopomers using two-color, IR-IR photodissociation spectroscopy.67 Moreover, the underlying potential energy functions (PEFs) previously developed for (H2O)NX− clusters were unable to provide a fully comprehensive representation of the halide ion− water and water−water interactions. Some studies focused on analyzing the structural and spectroscopic properties of small halide ion−water clusters from ab initio electronic structure calculations at the HF, DFT, MP2, and MP4 levels,68−73 while further studies sought to develop semiempirical74−76 or ab initio-based77 PEFs. Nevertheless, PEFs for the halide ion− water interactions were often fit to MP4 electronic structure data, and as such, they lack a rigorous characterization of shortrange quantum effects. Still another ab initio-based PEF for the H2OCl− cluster was fit to CCSD(T) energies using Morse functions.55,56 Ultimately, even more complex PEFs were developed that use the WHBB potential78 for the water−water

INTRODUCTION Halide ion solvation in aqueous media plays a central role in many crucial phenomena with physical, chemical, biological, and atmospheric implications.1−5 The details of the solvation process have been explored for clusters of finite size and for condensed phase systems from a plethora of experimental2,6−17 and theoretical6,8−10,12−15,18−23 angles. Experimental studies have focused primarily on the structural analysis of small halide ion−water clusters using IR spectroscopy,24−33while theoretical studies have sought to gain molecular-level insight into the structural, thermodynamic, dynamical, and spectroscopic aspects of the solvation process34−45 with particular attention devoted to preferential ion solvation at the air−water interface.46−49 To this end, classical molecular dynamics (MD) simulations in conjunction with thermodynamic integration methodologies have revealed that there is a propensity for the fluoride ion to be repelled from the air− water interface,50 but that all of the larger, more polarizable halide ions tend to reside at the surface of the air−water interface.50−53 Yet, the influence of the nuclear quantum (and isotope) effects on the structure of the underlying hydrogen bond network has proven elusive and challenging to establish definitively even for small (H2O)NX− clusters (X− = F−, Cl−, Br−, I−). Several theoretical studies have attempted to elucidate the quantum effects by producing vibrational frequencies and/ or IR spectra from the local monomer model (LMM), VSCF/ VCI, or even full vibrational CI calculations that were subsequently compared to corresponding experimental results.54−60 In addition, experimental and theoretical IR spectra were generated to examine the quantum isotope effects for the deuterated isotopologues and the various isotopomers of © XXXX American Chemical Society

Received: January 28, 2018 Revised: March 7, 2018

A

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

The Journal of Physical Chemistry A

Article



COMPUTATIONAL DETAILS We use the DMC method originally proposed by Anderson,114,115 i.e., without importance sampling, to examine the ground state properties of the smallest halide ion−water clusters considered in this work. This methodology has already been applied on numerous occasions to a diverse range of weakly interacting systems in the past by us94,104,116,117 and by others.118−126 DMC is based on the isomorphism of the Schrödinger equation in “imaginary time” (or “projection time”) and the diffusion equation with an additional source-sink term. In order to obtain the ground state energy numerically, the instantaneous wave function Ψ(r;τ) at projection time τ is represented as an ensemble of NW = NW(τ) random walkers that samples the configuration space of the system:

interactions and exploit the many-body expansion such that permutationally invariant polynomials for the 2-body terms were fit to high-level CCSD(T) data, while the 3-body terms (if any) were fit to MP2-level data. This family of WHBB-based PEFs have been employed to generate equilibrium structural properties57,79,80 and IR spectra59,60,79 in good agreement with experimental data for small halide ion−water clusters and for various clusters of pure water and ice.23,81−87 However, these PEFs were developed only for specific ions, i.e., F−, Cl−, and Na+, and were not transferable to the other halide (or alkali metal) ions in the series. Therefore, we deploy the MB-nrg PEF to model the interactions between the water molecules and the halide ions X−.88 The full suite of MB-nrg PEFs (which also includes alkali metal ion−water systems)89,90 is the most sophisticated ab initio-based PEF currently available as it encompasses all of the halide ions in the series and provides a complete description of the properties of water for small clusters up to the condensed phase. The MB-nrg PEF is based on the many-body expansion of the interaction energy, and the terms for the 2-body halide ion−water interactions were explicitly fit to high-level CCSD(T) data using permutationally invariant polynomials. The fitting procedure captures the short-range interactions, which are partially due to the quantum effects, such as charge transfer, charge penetration, and Pauli repulsion, with all of the long-range interactions described by the i-TTM model.91 MBnrg uses the MB-pol PEF for the water−water interactions as MB-pol is also ab initiobased and predicated on the many-body expansion with the 2- and 3-body terms fit to CCSD(T) data with permutationally invariant polynomials, and with the remaining N-body terms described by classical induction.92,93 MB-pol has successfully predicted theoretical and experimental data for the structural,94−98 thermodynamic,96,98−102 dynamical,96,99,100,103 and spectroscopic94,98,104−111 properties of small water clusters through to bulk liquid water and ice. Thus, by extension, the MB-nrg PEF is the ideal choice for our studies of small (H2O)NX− systems and their isotopologues. Note that the MB-nrg PEF has recently been used to obtain complete, anharmonic vibrational spectra for the H2OX− and D2OX− dimers (X− = F−, Cl−, Br−, I−) with frequencies and tunnelling splittings that are in good agreement with experimental results.112 We start with employing the diffusion Monte Carlo (DMC) method to probe the nuclear quantum effects in the smallest, (H2O)1−2X−, halide ion−water clusters and their isotopologues by permuting the positions of deuteriums (or hydrogens) to create different isotopomers. Furthermore, the principle aim of this work is to demonstrate the remarkable agreement between DMC and the harmonic approximation (HA) for the ground state solvation energies in (H2O)NX− clusters and their isotopologues. We can afford to run DMC simulations only for relatively small cluster sizes (i.e., N = 1−5), as the numerical effort of DMC increases rapidly with the system size. However, motivated by this agreement, we then extend our study to larger cluster sizes, i.e., N = 6−21, and to the analysis of the finite temperature properties for these clusters using the harmonic superposition method (HSM).113 While the latter approach (which is based essentially on the normal-mode analysis involving multiple energy minima on the potential energy surface (PES)) may appear pedestrian at first glance, the results that we obtain are supported by the established accuracy of HA when compared to the numerically exact DMC results.

NW

Ψ(r; τ ) =

∑ δ[r − rj(τ)] (1)

j=1

(In the DMC version that we employ, each random walker is defined by its position, rj(τ), in the configuration space, and carries the same unit weight.) At each time step Δτ, all random walkers are shifted randomly, rj(τ + Δτ ) = rj(τ ) + zj

(2)

according to a multivariate Gaussian distribution, defined by the value of Δτ and the particle masses. In addition, a branching procedure is implemented such that some walkers are replicated and some are killed. This stabilizes the variance in the instantaneous estimate of the ground state energy Eref = Eref(τ). The latter is updated at each time step to reflect changes in the random walker population, thus, maintaining a stationary distribution of the total random walker population. At sufficiently long projection times, both NW and Eref fluctuate in a stationary manner, and the time average of the latter provides an estimate for the ground state energy: E0 ≈

1 τmax

∫0

τmax

Eref (τ ) dτ

(3)

The ground state wave function is then estimated as a time average over the evolving random walker population, Ψ0(r) ≈

1 τmax

∫0

τmax

Ψ(r; τ ) dτ

(4)

The exact value of E0 and wave function Ψ0(r) are obtained only in the limit of NW → ∞, Δτ → 0, and τmax → ∞. A more detailed description of the DMC algorithm and its implementation that we adapted here is available in ref 116. There is a well-known problem in the DMC approach associated with computation of physical observables that do not commute with the system Hamiltonian. (Note that essentially any interesting physical observable, other than the ground state energy, would face this problem.) For example, consider a coordinate-dependent function O(r), and its expectation value ⟨O⟩ ≔

∫ O(r)|Ψ0(r)|2 dr ∫ Ψ0(r)|2 dr

(5)

Apparently, given the random walker population satisfying eq 4, there is no numerically straightforward procedure to evaluate the above expression, while a related Ψ-averaged quantity (albeit not a true physical observable), B

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A ⟨O⟩1 ≔

∫ O(r)Ψ0(r) dr ∫ Ψ0(r) dr

Z HSM(T ) =

⟨O⟩1 ≈

τmax

∫0

τmax

⎧ ⎪ 1 dτ ⎨ ⎪N ⎩ W

For each energy minimum α, we use the rigid rotor-harmonic oscillator approximation to evaluate the corresponding partition function, i.e., including only vibrational and rotational contributions,

⎫ ⎪ O [ r ( )] τ ∑ j ⎪⎬ ⎭ j=1 NW

Zα(T ) = (7)

1 τmax

⎧ ∑NW O[r (τ )]w (τ ) ⎫ τmax ⎪ j=1 ⎪ j j ⎬ dτ ⎨ NW ⎪ ⎪ 0 ∑ j = 1 wj(τ ) ⎩ ⎭

wα(T ) = Zα(T )/Z HSM(T )

∑ i=1

ℏωiα 2

(12)

Another quantity of interest is the ion solvation free energy ΔFsolvation(T). At zero temperature it is unambiguously defined as the ground state energy difference ΔEsolvation = E[(H 2O)N X−] − E[(H 2O)N ]

(13)

where the energy of an isolated X− is set to zero on the MB-nrg PEF. This expression assumes that the true ground state energies of both (H2O)NX− and (H2O)N can be obtained, which is the case with DMC. On the other hand, when using HA one is required to specify a particular energy minimum for each system. As will be seen, in most cases considered in this study the DMC calculation results in a ground state wave function that can be definitively assigned to a single minimum. However, we also encounter situations (see ref 94) where the ground state wave function is delocalized over more than one minimum with all of the minima having well-defined isomer fractions, or where there exist two or more nearly isoenergetic minima separated by high energy barriers. In the former case we use the minimum corresponding to the highest isomer fraction, while in the latter case of quasi-degenerate states (i.e., for different isotopologues with various permutational isotopomers) the most natural choice corresponds to the minimum with the lowest harmonic energy. Within the HSM framework we can also define the temperature dependent Helmholtz free energy using

(8)

We have found that observables which correspond to slow, intermolecular degrees of freedom (such as certain spatial correlation functions) may converge very slowly with τDW, while others that report on, e.g., fast, intramolecular degrees of freedom converge rapidly. As such, we cannot guarantee that all of the spatial correlation functions and isomer fractions computed using DW in this work are correct from a quantitative standpoint, but their qualitative accuracy is still superior to that of the corresponding Ψ-averaged results. In addition to the DMC calculations for cluster sizes N = 1− 5, the energetic and structural properties of (H2O)NX− clusters with N = 1−21 and their corresponding (D2O)N counterparts were treated in the framework of HA. Namely, for each system in question we performed a long classical Monte Carlo (MC) exploration of the PES at sufficiently high temperature (T ∼ 300 K) quenching the resultant configuration once in every 104 MC steps using the conjugate gradient method. This procedure resulted in a number of minimum energy structures, i.e., isomers, which were then compared (see ref 94 for details), so that a single copy of each isomer would be stored in the isomer library. For each configuration in the isomer library, the 3n − 6 nonzero normal-mode frequencies ωαi were computed, with n defining the number of atoms in the cluster. Consequently, for each isomer α the HA ground state energy was obtained according to 3n − 6

(11)

where (k = 1, 2, 3) are the principal moments of inertia, β = 1/kBT, and σα is the order of the isomer point group, which is equal to unity for nearly all cases encountered. Note that the above expression only includes vibrational and rotational contributions, but completely ignores the effects of pressure and volume. Within HSM the isomer fractions are estimated by



α + E0α = Emin

⎛ I αI αI α ⎞1/2 ⎜ 1 2 3 ⎟ Πi3=n −1 6[1 − exp(−β ℏωiα)] ⎝ σα2β 3ℏ6 ⎠ exp( −βE0α)

Iαk

In principle, expression 5 can be evaluated using the so-called descendant weighting (DW) approach,119,127,128 the primary aim of which is to generate a second copy of the ground state wave function Ψ0(rj(τ)). This is done by counting the number of random walkers that have descended, i.e., progeny from a given walker rj(τ) after waiting a finite projection time τDW into the future. Upon accumulation of the weights wj(τ), the random walker population from DMC provides the first copy of the ground state wave function and the weights provide a second copy, thereby enabling the Ψ2-averaged observable to be computed as follows ⟨O⟩ ≈

(10)

α

(6)

can be readily computed within the standard DMC framework by summing over all the random walkers and averaging over time: 1

∑ Zα(T )

FHSM(T ) = −kBT ln Z HSM(T )

(14)

Thus, from HSM (including only rotational and vibrational contributions), the quantum solvation free energy is computed by ΔFsolvation(T ) = FHSM[(H 2O)N X−] − FHSM[(H 2O)N ]



(15)

RESULTS AND DISCUSSION Ground State Energies and Structures of (H2O)1−2X− and their Deuterated Isotopologues. In this subsection, we explore the properties of the deuterated isotopologues for the smallest, (H2O)1−2X− (X− = F−, Cl−, Br− and I−), clusters by replacing one, two or three hydrogens with deuteriums. Independent DMC simulations were all run starting from different initial conditions by placing the H or D atom in different positions of the MB-nrg PEF global minimum

(9)

In addition, we implemented the harmonic superposition method (HSM)113,129 to estimate the free energies and isomer fractions, i.e., the isomer populations at finite temperatures. Within HSM, the total canonical ensemble partition function is defined as a sum of the isomer partition functions taken over all of the relevant potential energy minima: C

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A configuration to obtain the DMC ground state energies E0, the corresponding binding energies D0 ≔ N × E[(H 2O)] − E[(H 2O)N X−]

Table 2. Favored Isotope Substitution Positions (cf. Figures 1 and 2), DMC Ground State Energies E0, Binding Energies D0, and Isomer Fractions for Different (H2O)2X− Isotopologues (X− = F−, Cl−, Br−, I−)a

(16)

(with analogueus expressions used for the other isotopologues and where E[X−] = 0) and, if possible, the isomer fractions. The chemical formulas of the different isotopologues are all listed in Tables 1 and 2, and the various permutational isotopomers, i.e., the different D (or H) atom positions, are color coded and shown in Figures 1 and 2.

E0

D0

HOD−H2O−F− HOD−H2O−F− HOD−H2O−Cl− HOD−H2O−Br− HOD−H2O−I−

1, 3, 1, 1, 1,

2 4 2, 3 2, 3 2, 3

−20.34 −20.21 −2.03 1.11 4.65

45.07 44.94 26.76 23.62 20.09

0.52, 0.48 0.45, 0.55 − − −

1, 3, 3, 1, 1,

2 4 4 3, 4 3, 4

−24.04 −24.17 −6.01 −2.86 0.69

45.21 45.34 27.18 24.03 20.48

0.47, 0.55, 0.10, 0.01, 0.06,

0.53 0.45 0.90 0.06, 0.93 0.03, 0.91

2 2 2 2

−22.19 −4.01 −0.87 2.65

45.12 26.94 23.80 20.28

0.59, 0.10, 0.07, 0.03,

0.41 0.90 0.93 0.97

System

Table 1. Favored Isotope Substitution Positions (cf. Figures 1 and 2), DMC Ground State Energies E0, Binding Energies D0, and Isomer Fractions for Different H2OX− Isotopologues (X− = F−, Cl−, Br−, I−)a

Structure

Isomer Fraction

System

Structure

E0

D0

Isomer Fraction

HOD−D2O−F− HOD−D2O−F− HOD−D2O−Cl− HOD−D2O−Br− HOD−D2O−I−

HOD−F− HOD−F− HOD−Cl− HOD−Br− HOD−I−

1 2 1, 2 1, 2 1, 2

−15.16 −15.09 −2.37 −0.44 1.63

26.65 26.58 13.86 11.93 9.86

1.00 1.00 0.99, 0.01 0.97, 0.03 0.85, 0.15

H2O−D2O−F− H2O−D2O−Cl− H2O−D2O−Br− H2O−D2O−I−

1, 1, 1, 1,

H2OF− H2OCl− H2OBr− H2OI−

− − − −

−13.27 −0.41 1.50 3.52

26.51 13.65 11.74 9.73

1.00 1.00 1.00 1.00

(H2O)2F− (H2O)2Cl− (H2O)2Br− (H2O)2I−

− − − −

−18.37 −0.07 3.08 6.61

44.86 26.55 23.40 19.88

1.00 1.00 1.00 1.00

D2OF− D2OCl− D2OBr− D2OI−

− − − −

−16.98 −4.19 −2.28 −0.22

26.67 13.88 11.96 9.90

1.00 1.00 1.00 1.00

(D2O)2F− (D2O)2Cl− (D2O)2Br− (D2O)2I−

− − − −

−26.02 −7.90 −4.76 −1.21

45.40 27.27 24.13 20.58

1.00 1.00 1.00 1.00

a

The nomenclature, DMC parameters, and statistical errors are the same as in Table 1 except that NW = 3 × 104. For the case of the HOD−H2O−X− clusters (X = Cl−, Br−, I−), three practically isoenergetic configurations (1, 2, 3) are identified, but the isomer fractions are not converged (see text).

The “−” label for the “Structure” column indicates that a single minimum structure exists for the corresponding system. The ground state energies of the isolated HOD, H2O, and D2O monomers used to compute D0 are 11.49, 13.24, 9.69 kcal/mol, respectively, and all energies are reported in kcal/mol. The DMC simulations were run with NW = 2 × 104 random walkers, a time step of Δτ = 10 au, and a projection time of τmax = 2 × 106 au. The statistical errors for the E0 and D0 values are all on the order of 10−3 kcal/mol or better. HOD−I− is the only case for which the ground state wavefunction is substantially delocalized over the two minima (1 and 2). The isomer fractions were computed using DW with τDW = 500 au. a

The H2OX− dimer and its isotopologues correspond to a simple case with only one classical energy minimum, such that one of the H/D atoms points toward the halide ion (position 1), while the other H/D atom is oriented away from the halide ion in a dangling (free) position (position 2). The DMC ground state energies for the different isotopologues are listed in Table 1 along with their corresponding binding energies D0 and the isomer fractions of the various permutational isotopomers. For HOD−F−, the barrier to the rocking motion for flipping between positions 1 and 2 and vice versa is high such that the H/D atom remains either hydrogen-bonded to the F− ion or in the dangling position. However, the barrier height slowly decreases as the polarizability of the halide ion increases with the isomer fraction of the dangling position growing to the order of ∼0.01 for HOD−Cl− and HOD−Br−. Upon reaching HOD−I−, the ground state wave function is effectively delocalized between positions 1 and 2 having isomer fractions of 0.85 and 0.15, respectively. This delocalization of the ground state wave function for HOD−I− proves that the

Figure 1. Images of the permutational isotopomers of the (H2O)1−2X− (X− = Cl−, Br−, I−) clusters studied in this work. The color coding and numbering scheme corresponds to the numbered positions of the substituted isotope listed in Tables 1 and 2.

two configurations are separated by a relatively low energy barrier, which is consistent with tunnelling splitting results for D

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

flipping a free H atom into a hydrogen-bonded position declines in the HOD−D2O−X− isotopologue as the polarizability of the X− ion increases. The situation is significantly more nontrivial for the HOD− H2O−X− (X− = Cl−, Br−, I−) isotopologue. There are three configurations (1, 2, 3) that are separated by relatively high energy barriers, but appear to be practically isoenergetic (ΔE ∼ 0.01 kcal/mol). Apparently, the DMC ground state energy and the associated binding energy are poor criteria for determining the convergence of the simulation for the ground state wave function as the random walker population is always incorrectly distributed between the quasi-degenerate states pursuant to the choice of the initial condition. Under these circumstances, DMC cannot be trusted to estimate, in particular, the isomer fractions. In such a case a more adequate structural analysis of these systems could be carried out using finite temperature approaches. However, note that the use of, in principle, numerically exact path integral-based methods would likely be problematic as well, especially at sufficiently low temperatures for the same reasons that DMC encounters sampling problems. Apparently, HSM provides a meaningful and computationally affordable (albeit manifestly approximate) means for exploring the interplay between quantum and thermal fluctuations in these small clusters, and assuming that all of the relevant potential energy minima have been found, its implementation has no convergence problems. To this end, Figure 3 shows the

Figure 2. Same as Figure 1, but for the (H2O)1−2F− clusters.

H2OI− and D2OI− obtained in ref 112. Note that the preference of the D atom to participate in hydrogen bonding with the halide ion (albeit at finite temperatures) has also been examined in path-integral molecular dynamics (PIMD) simulations using thermodynamic integration.64,65 On moving to the (H2O)2X− trimer and its isotopologues (see Table 2), the situation is simple for the cases corresponding to a single structure as marked in the “Structure” column by “−”. There is a high barrier to the rocking motion that prevents hydrogen-bonded and free H/D atoms from switching positions in all of the F− isotopologues as the water molecules are strongly attracted to the F− ion. However, a very small barrier exists for a concerted change of the orientations of the two water molecules such that significant delocalization of the ground state wave function is observed for the two isotopomers of H2O−D2O−F−. For this isotopologue, the isomer fractions for the two states are 0.59 and 0.41, respectively, with the D2O molecule preferentially residing at the position with the longer hydrogen bond. The ground state wave function is also delocalized in a similar fashion between the two states corresponding to the hydrogen-bonded (1,2) and free (3,4) H/D positions in the HOD−H2O−F− and HOD−D2O−F− isotopologues (cf. Figure 2). Upon considering the more polarizable halide ions (X− = Cl−, Br−, I−), another hydrogen bond forms between the two water molecules in the classical minimum structure (cf. Figure 1) such that the two water molecules can no longer easily change positions as they can for F−. Instead, the two water molecules remain effectively tethered to each other and to halide ion with the D2O molecule overwhelmingly occupying position 2 in H2O−D2O−X− with isomer fractions ≥0.90. Likewise, the D atoms are invariably present in all of the possible hydrogenbonded positions (1, 2, 3) in the HOD−D2O−X− isotopologues (again with isomer fractions ≥0.90) leaving the lone H atom in position 4, i.e., in the only dangling or free position. The isomer fractions of the first water molecule in H2O−D2O− X− decline, whereas configurations 1 and 3 become more accessible in HOD−D2O−X− upon going from Cl− to I−. This indicates that the X−-water hydrogen bonds become gradually weaker, but at the same time, the water−water hydrogen bond grows stronger as one moves down the halide series. Consequently, the barrier for a concerted change between the two water molecules increases for the H 2 O−D 2 O−X − isotopologue, while the barrier to the rocking motion for

Figure 3. H or D atom fractions computed by HSM (with a rotational correction, cf. eq 11) for some selected isotopologues of the (H2O)1−2I− clusters from Tables 1 and 2. The color coding denoting the positions of the unique H or D atoms for the permutational isotopomers is the same as in Figure 1.

isomer fractions as a function of temperature computed by HSM (with a rotational correction, cf. eq 11) for selected isotopologues of the (H2O)1−2I− clusters. These HSM results are consistent with those in Tables 1 and 2, and thus, they complement each other. For example, the D atom has a preference in position 1 for the HOD−I− cluster at low temperatures with its isomer fraction approaching the asymptotic value of 0.5 at T ∼ 300 K. Furthermore, the D atom in the HOD−H2O−I− cluster starts to be almost equally distributed between positions 1, 2, and 3 (which are listed as quasi-degenerate in Table 2 at relatively low temperatures (T ∼ 50 K), etc. The quantum isotope effect is apparent at low E

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A temperatures, but gets washed out by the thermal fluctuations as the temperature is increased. Similar results for the other small clusters (X− = F−, Cl−, Br−) are reported in the Supporting Information and are in complete agreement with the DMC results. DMC versus HA: Solvation Energies and Structural Properties of (H2O)1−5X− Clusters. Here we consider the (H2O)1−5X− clusters and compare the solvation energies obtained from DMC with those from HA. For a given cluster a typical DMC calculation was initialized by distributing the random walkers equally among all the local minima on the PES that had been found prior to performing the simulation. A single long simulation would then be performed for time τmax, with statistical errors estimated by breaking the simulation into 10 blocks of the same length. The DMC simulation parameters, i.e., the random walker populations NW and the total projection times τmax for each cluster size N are shown in Table 3. A time

can be found in the Supporting Information). Upon comparison of the two types of results, the most striking observation is that HA provides practically exact values (i.e., with accuracy on the order of ΔE ∼ 0.1 kcal/mol or better) for the solvation energies, assuming that the DMC results are numerically exact. This remarkable accuracy of HA is surprising because, for example, the absolute ground state energies alone are associated with significant anharmonic shifts (cf. the third and fourth columns of Table 4). In particular, the contribution to the latter from only one OH stretch mode is about 0.6 kcal/ mol. Clearly, to achieve a total discrepancy on the order of 0.1 kcal/mol, massive error cancellations must take place. At the same time, the mechanism for the error cancellations is not obvious because, for example, the addition of a halide ion does distort the hydrogen bond network in the water cluster substantially. This effect is demonstrated in Figure 4, where the most pronounced distortions occur in the (H2O)5I− cluster. In addition, some minor structural changes, such as flipping of free hydrogens do occur in the smaller (H2O)3−4 and (H2O)3−4I− clusters, but the hydrogen bonding pattern for the water molecules still retains its nearly equilateral triangular or square planar shape, respectively. Despite the apparently trivial disruption of the water−water motif upon addition of the halide ion in the smaller clusters, we do observe a significant increase in the O−O distances even for (H2O)3I− versus (H2O)3. In what follows we report on our results for the structural analysis of (H2O)3X− clusters using various spatial correlation functions. As already mentioned, the true physical observables, eq 5, other than the ground state energies, cannot be computed in a straightforward manner within the basic DMC framework. Instead, either the ψ-averaged quantities, eq 6, are reported, or methods such as DW are employed. Our numerical tests using DW led us to conclude that the ψ2-averaged spatial correlation functions, which only depend on fast (typically intramolecular) degrees of freedom, converge well, while those depending on slow (often intermolecular) degrees of freedom may converge poorly or may not converge at all as it is necessary to use large values for τDW, leading in turn to large statistical errors. This behavior with respect to τDW is expected for the intermolecular degrees of freedom because they evolve on comparatively slower time scales during the DMC simulations, and thus, they require longer correlation times that in turn necessitate larger τDW values. Namely, for the X−−H, X−−O, and O−O distance distribution functions, we were not able to obtain reliable results using DW. Consequently, for these observables we report only the ψ-averaged quantities, while for the O−H distance and H−O−H angle distribution functions we do use the DW method. Note, however, that although the ψ-averaged quantities cannot be compared directly to the true physical observables, they do often provide meaningful physical information. In particular, the ψ-averaged distribution functions are used here to estimate the mean X−−H, X−−O, and O−O interatomic distances. Since for most cases the DMC ground state wave function could be associated with a particular potential energy minimum, and hence with a particular harmonic ground state wave function, the HA estimates (either ψ- or ψ2-averaged) are also provided together with the corresponding DMC results. To this end, Figures 5 and 6 show, respectively, the O−H distance and H−O−H angle distributions for (H2O)3 and (H2O)3X− clusters (X− = F−, Cl−, Br−, I−), with the O and H atoms belonging to the same monomer. Apparently, the O−H

Table 3. Random Walker Populations NW and Projection Times τmax Used in the DMC Simulations for the (H2O)1−5X− Clusters (X− = F−, Cl−, Br−, I−) and the Corresponding (H2O)1−5 Pure Water Clustersa System −

τmax (au)

NW 4

2 × 10 1.96 × 104

2 × 106 2 × 106

(H2O)2X− (H2O)2

3 × 104 2 × 104

2 × 106 2 × 106

(H2O)3X− (H2O)3

7 × 104 6 × 104

6.3 × 105 4 × 106

(H2O)4X− (H2O)4

9 × 104 8 × 104

2.7 × 105 2 × 106

(H2O)5X− (H2O)5

2.1 × 105 2 × 105

6 × 104 2 × 106

H2OX H2O

All of the DMC simulations were performed with a time step of Δτ = 10 au. For the N = 3−5 systems, the random walker population was initialized in the global minimum and among the various local minima on the PES.

a

step of Δτ = 10 au was used for all of the DMC simulations such that the systematic error in Δτ is on the order of 10−2 kcal/mol (cf. the Supporting Information for details). While the τmax values reported in the table are rather short, particularly for the larger (H2O)3−5X− clusters, we have used large random walker populations NW to compensate for the shorter projection times as the statistical error scales according to ∼(NWτmax)−1/2 for a fixed Δτ. This strategy of employing sizable NW values improves the convergence of the DMC simulations by reducing the systematic error (bias) in NW, while simultaneously keeping the magnitude of the statistical error for all of the cluster sizes on the order of 10−2 kcal/mol or better. Note that both the systematic and statistical errors for our DMC calculations are much smaller than the reported accuracy of the MB-nrg PEF.88 The last two columns of Table 4 show the ground state solvation energies computed with DMC and HA (i.e., ΔEDMC and ΔEHA, respectively) for the (H2O)1−5 (X− = F−, Cl−, Br−, I−) systems. (Similar results for the (D2O)1−5X− isotopologues F

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Table 4. Classical Minimum Energies (Emin) and Absolute Ground State Energies from DMC (EDMC), HA (EHA), and Scaled HA (ESHA) along with Solvation Energies from DMC (ΔEDMC) and HA (ΔEHA) for the (H2O)1−5X− Clusters (X− = F−, Cl−, Br−, I−) and the Pure (H2O)1−6 Clustersa EDMC

ΔEHA

ESHA

−27.25 −47.82 −66.28 −82.33 −97.09

−13.27 −18.37 −20.18 −20.12 −18.80

± ± ± ± ±

0.0005 0.002 0.003 0.02 0.03

−13.02 −17.41 −18.74 −18.45 −16.74

−13.37 −18.15 −19.89 −20.00 −18.68

−26.51 −41.71 −49.06 −53.53 −58.86

± ± ± ± ±

0.0007 0.002 0.003 0.02 0.03

−26.50 −41.45 −48.88 −53.28 −58.65

H2OCl− (H2O)2Cl− (H2O)3Cl− (H2O)4Cl− (H2O)5Cl−

−14.74 −29.85 −46.25 −59.90 −73.51

−0.41 −0.07 −0.20 1.84 4.09

± ± ± ± ±

0.0009 0.002 0.005 0.02 0.05

−0.11 0.67 0.96 3.38 5.94

−0.46 −0.07 −0.19 1.85 4.02

−13.65 −23.41 −29.08 −31.57 −35.96

± ± ± ± ±

0.001 0.002 0.005 0.02 0.05

−13.58 −23.37 −29.19 −31.45 −35.97

H2OBr− (H2O)2Br− (H2O)3Br− (H2O)4Br− (H2O)5Br−

−12.76 −26.62 −42.37 −55.75 −68.67

1.50 3.08 3.62 6.00 8.78

± ± ± ± ±

0.0009 0.002 0.004 0.02 0.03

1.80 3.80 4.77 7.47 10.60

1.45 3.06 3.63 5.94 8.67

−11.74 −20.26 −25.26 −27.41 −31.28

± ± ± ± ±

0.001 0.002 0.004 0.02 0.03

−11.67 −20.24 −25.37 −27.36 −31.32

H2OI− (H2O)2I− (H2O)3I− (H2O)4I− (H2O)5I−

−10.60 −22.98 −37.86 −50.87 −62.75

3.52 6.61 8.06 10.82 14.46

± ± ± ± ±

0.0007 0.002 0.004 0.01 0.02

3.85 7.31 9.20 12.27 16.26

3.50 6.58 8.06 10.74 14.35

−9.73 −16.73 −20.82 −22.59 −25.59

± ± ± ± ±

0.0008 0.003 0.004 0.01 0.02

−9.62 −16.73 −20.95 −22.56 −25.65

H2O (H2O)2 (H2O)3 (H2O)4 (H2O)5 (H2O)6

8.00 × 10−5 −4.96 −15.69 −27.18 −35.55 −45.94

13.24 23.34 28.88 33.41 40.06 46.45

± ± ± ± ± ±

0.0005 0.001 0.0009 0.005 0.005 0.004

13.47 24.04 30.14 34.84 41.91 48.73

13.15 23.34 29.03 33.33 40.04 46.44

− − − − − −

System

Emin

ΔEDMC

EHA

H2OF− (H2O)2F− (H2O)3F− (H2O)4F− (H2O)5F−

− − − − − −

a A single scaling factor (λ = 0.976) was applied to the harmonic frequencies of all the systems shown above in order to map the absolute HA energies onto the numerically exact DMC energies (see text). All energies are reported in kcal/mol. The statistical errors for the DMC simulations are on the order of ∼10−2 kcal/mol or better.

DMC versus HA: Absolute Ground State Energies and Scaling Factors. We now shift our discussion to the absolute ground state energies as it is instructive to entertain the possibility of using scaling factors to bring the HA ground state energies into agreement with those from DMC. In the spirit of the common practice exploited, for example, by Temelso et al.,130,131 particularly for water clusters, the scaling procedure here was used to map the HA ground state energies onto the numerically exact DMC energies for the (H2O)1−5X− (X− = F−, Cl−, Br−, I−) and for the pure (H2O)1−6 clusters. In order to define K scaling factors λ1, ..., λK we split the whole frequency range into K intervals by choosing the values h0 = 0 < h1 < ... < hK = ∞. Each [hk−1, hk) interval then corresponds to scaling factor λk with the frequencies in this interval scaled as

distribution function is only slightly perturbed by the addition of a halide ion, whereas the H−O−H angle distribution function shifts to noticeably smaller bond angles upon adding a halide ion. Figures 7 and 8 show, respectively, the X−−H and X−−O distance distributions, and these two distributions depend strongly on the halide ion X−. Note also that, in general, HA accounts only partially for nuclear quantum effects as, in particular, it does not allow for possible changes of internuclear distances caused by the asymmetries in the potential energy well around the minimum. Yet, for the four cases mentioned above the HA and DMC results agree very well with each other, and the anharmonic shifts associated with the centroid positions of the peaks in the distribution functions are very small. However, the behavior of the O−O distance distribution (see Figure 9) is very different. Not only does it depend significantly on the halide ion, but also the anharmonic and nuclear quantum effects manifested in the DMC O−O distribution function cause a substantial shift of the centroid peak positions to longer O−O distances. The effects described above are also demonstrated directly in Figure 10, where the corresponding interatomic distances and angles are plotted as a function of the ion.

ωĩ = λk ωi , hk − 1 ≤ ωi < hk

(17)

The following linear least-squares problem is then solved for the unknowns λk: 3n − 6

∑ α

G

α Emin

+

∑ i=1

ℏωĩ α α − E DMC 2

2

→ min

λ1,..., λK

(18)

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 6. ψ2-averaged H−O−H angle distribution function from HA and DMC for the (H2O)3X− clusters.

Figure 4. Images of the lowest classical minimum structures for pure (H2O)3−5 and (H2O)3−5I− clusters, which all happened to correspond to their ground states. The most dramatic structural change upon addition of the halide ion occurs for N = 5.

Figure 7. ψ-averaged X−-H distance distribution function from HA and DMC for the (H2O)3X− clusters.

Figure 5. ψ2-averaged O−H distance distribution function from HA and DMC (τDW = 500 au) for the (H2O)3X− clusters.

with α running over all or a selected subset of configurations, for which the DMC data exists. We first considered a single scaling factor (i.e., K = 1), resulting in λ = 0.976. The corresponding scaled harmonic energies ESHA and the DMC energies EDMC agreed within ΔE ∼ 0.1 kcal/mol. Moreover, nearly the same value for λ was obtained using only part of the data, which did not include the largest clusters, confirming that the scaling procedure can be used to predict accurate values for the ground state energies only using HA. We then considered two scaling factors, one for the intra- and one for the intermolecular degrees of freedom (h1 = 1000 cm−1), and achieved a similar agreement between the scaled HA and DMC energies, i.e., also within ΔE ∼ 0.1

Figure 8. ψ-averaged X−-O distance distribution function from HA and DMC for the (H2O)3X− clusters.

H

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

obtained using the scaled and unscaled frequencies, we found that they are virtually the same, and hence in our following analysis, we only report the thermodynamic properties computed using the unscaled HA frequencies. Solvation Free Energies and Structural Properties of (H2O)1−21X− Clusters Using the Harmonic Superposition Method. In this subsection, we apply HSM to clusters with sizes up to N = 21. Yet, we do not naively assume that the aforementioned remarkable agreement between the numerically exact DMC results and HA for the small (H2O)1−5X− clusters extends to the full gamut of cluster sizes, as it cannot be guaranteed that deviations would not start to occur for larger N. However, the results for the smallest clusters do indicate a promising trend that could be carried on for the larger cluster sizes studied here. To this end, Figure 11 shows the solvation Table 5. Bulk System Gibbs Solvation Free Energies ΔG°solvation for All of the Halide Ions (X−=F−, Cl−, Br−, I−)a

Figure 9. ψ-averaged O−O distance distribution function from HA and DMC for the (H2O)3X− clusters.

System

ΔG°solvation

ΔFsolvation



−102.49 −72.71 −66.30 −57.36

−79.75 −54.27 −49.12 −41.55

F Cl− Br− I−

The ΔG°solvation values were obtained from ref 135. The ° superscript denotes standard temperature (T = 298.15 K) and pressure (P = 1 atm). Also shown are the present HSM results for the Helmholtz free energy of solvation ΔFsolvation (see text) at T = 300 K for the largest clusters, (H2O)21X−, studied in this work. All energies are reported in kcal/mol, and the uncertainty in the Gibbs free energies is ±1.9 kcal/ mol. a

free energies computed by HSM (cf. eqs 10−15) for N = 1−21 at temperatures of T = 0 and 300 K. In reality, all of the finite cluster sizes considered in this work are extremely small relative to that of the corresponding condensed phase systems, such that we do not contend to have reached the thermodynamic limit. Ignoring the PΔV work term as it is small for ion Figure 10. Mean interatomic distances as obtained from the corresponding distribution functions from HA (filled circles) and DMC (X’s) for the (H2O)3 and (H2O)3X− clusters. The mean H−O− H angles are reported in radians.

kcal/mol. In addition, the two scaling factors turned out to be both very close to the above value of the single scaling factor λ = 0.976. We also tried to use three scaling factors (K = 3), i.e., distinguishing between the intermolecular modes, bends, and stretches (h1 = 1000 cm−1, h2 = 2000 cm−1), but because of the Fermi resonance between the intramolecular bending and stretching modes, the results for λ2 and λ3 proved to be unstable. Consequently, in Table 4 we only report the scaled harmonic energies ESHA calculated with a single scaling factor. (Results for the corresponding D2O isotopologues are available in the Supporting Information.) Although seemingly “ad hoc”, the present scaling procedure, based on mapping the HA energies to exact data, can be viewed as a possible recipe for constructing an optimal harmonic approximation, i.e., in the spirit of such approaches as the selfconsistent phonons132,133 or principal mode analysis.134 Thus, following the general idea one could use the scaled frequencies in order to compute properties, other than the ground state energies, in particular, the thermodynamic properties of all the (H2O)NX− clusters. However, by comparing the HSM results

Figure 11. Solvation free energies ΔFsolvation computed by HA as a function of cluster size N for the quantum (solid lines, filled symbols) and classical (dashed lines, open symbols) (H2O)NX− clusters (X− = F−, Cl−, Br−, I−) at T = 0 and 300 K. Note that the results for T = 300 K correspond to vibrational and rotational contributions to the free energies (cf. eq 11). The dotted horizontal lines correspond to the bulk Gibbs solvation free energies as described in Table 5. I

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

Figure 12. Isomer fractions/populations as a function of temperature from HSM for the (H2O)20X− (X− = F−, Cl−, Br−, I−) systems studied in this work. Images of some highly populated isomers are also displayed.

for the (H2O)NX− and pure (H2O)N systems in the same way, such that they cancel out upon evaluating the free energy difference. This is especially true at high temperatures, and thus, the notable difference in ΔFsolvation for F− stems from the changes in the corresponding potential energy function that directly reflect the enhanced proton affinity (i.e., basic properties) of the F− ion. Finally, in order to gain more insight into the equilibrium properties of halide ions solvated in water, we apply HSM to quantum (H2O)20X− (X− = F−, Cl−, Br−, I−) clusters to compute the corresponding isomer populations as a function of temperature (see Figure 12). Images of the most highly populated isomers are also shown in the figure. At low temperatures, i.e., in the T → 0 limit, only a single low-lying isomer dominates, but several other isomers start to become appreciably populated around T ∼ 100 K, while that of the most populated isomer invariably decreases. As the temperature is increased to T ∼ 300 K, the populations of the low-lying isomers are almost completely depleted and become distributed over a large number of configurations. The relevant plots for the corresponding classical clusters and the quantum (D2O)20X− isotopologue are included in the Supporting Information. As the quantum and isotope effects are small, these plots are qualitatively similar to those shown in Figure 12.

solvation, such that ΔGsolvation ° ≈ ΔFsolvation, the bulk solvation free energies (see Table 5) are all ∼20 kcal/mol lower than our own at T = 300 K upon comparison to the ΔFsolvation value for the largest clusters considered in this work, i.e., (H2O)21X−. The visible bumps at various cluster sizes on the solvation free energy curves in Figure 11 for T = 0 K are due to especially stable isomeric configurations of the pure water clusters used to determine the solvation energies. These features are smoothed out for the corresponding ΔFsolvation results at T = 300 K, as numerous isomers make significant contributions to the partition function at high temperatures. Note that the classical and quantum ΔFsolvation curves exhibit very similar behavior with cluster size at both T = 0 and 300 K regardless of the halide ion type, where the “classical” result is obtained by using eqs 10−15 with ℏ → 0. The minor deviations noticeable at T = 0 K are completely nonexistent at T = 300 K, which indicates that the quantum effects for these systems are small and that the thermal fluctuations at T = 300 K effectively drown out all of the quantum fluctuations. Thus, accurate ΔFsolvation values for these systems at sufficiently high temperatures can in principle be obtained from classical molecular dynamics or Monte Carlo simulations using the same PEFs. The isotope effects turned out to be practically nonexistent for ΔFsolvation, but for completeness we do include the corresponding results for the (D2O)NX− isotopologue in the Supporting Information. We further observe that the magnitude of ΔFsolvation for F− is ∼25 kcal/mol lower than that of Cl− for the largest cluster sizes shown in Figure 11, but the deviation between those for the Cl−, Br− and I− ions never exceed ∼10 kcal/mol. The F− ion causes more appreciable structural disruption and rearrangements of the hydrogen bond network than do the other three ions. There is also a substantially stronger attractive interaction between the F− ion and the water molecules than for the rest of the halide ions, which distorts the pure water structure and makes F− the most thermodynamically favorable to solvate. Furthermore, the nuclear quantum effects are influencing both the intra- and intermolecular motions of the water molecules



CONCLUSION We have explored the nuclear quantum effects for small halide ion−water clusters and their various deuterated isotopologues using the ab initio-based MB-nrg PEF in conjunction with HA and DMC. Upon comparing the HA and DMC ground state solvation energies for small (H2O)1−5X− (X− = F−, Cl−, Br−, I−) systems, we revealed surprisingly good agreement (ΔE ∼ 0.1 kcal/mol) between these quantities from the two methodologies. This massive cancellation of errors was unexpected as the addition of a halide ion leads to signicant structural changes in a pure water cluster. Moreover, the anharmonicities in conjunction with nuclear quantum effects J

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A

simulations of ion-enhanced interfacial chemistry on aqueous NaCl aerosols. Science 2000, 288, 301. (2) Collins, K. D.; Neilson, G. W.; Enderby, J. E. Ions in water: Characterizing the forces that control chemical processes and biological structure. Biophys. Chem. 2007, 128, 95. (3) Marcus, Y. Effect of ions on the structure of water: Structure making and breaking. Chem. Rev. 2009, 109, 1346. (4) Kunz, W. Specific ion effects in colloidal and biological systems. Curr. Opin. Colloid Interface Sci. 2010, 15, 34. (5) Jungwirth, P.; Cremer, P. S. Cremer, Beyond Hofmeister. Nat. Chem. 2014, 6, 261. (6) Cabarcos, O. M.; Weinheimer, C. J.; Lisy, J. M.; Xantheas, S. S. Microscopic hydration of the fluoride anion. J. Chem. Phys. 1999, 110, 5. (7) Corcelli, S.; Kelley, J.; Tully, J.; Johnson, M. Infrared characterization of the icosahedral shell closing in Cl− · H2O · Arn (1≤ n ≤ 13) clusters. J. Phys. Chem. A 2002, 106, 4872. (8) Diken, E. G.; Headrick, J. M.; Roscioli, J. R.; Bopp, J. C.; Johnson, M. A.; McCoy, A. B.; Huang, X.; Carter, S.; Bowman, J. M. Argon predissociation spectroscopy of the OH− · H2O and Cl− · H2O complexes in the 1000−1900 cm−1 region: Intramolecular bending transitions and the search for the shared-proton fundamental in the hydroxide monohydrate. J. Phys. Chem. A 2005, 109, 571. (9) Soper, A. K.; Weckström, K. Ion solvation and water structure in potassium halide aqueous solutions. Biophys. Chem. 2006, 124, 180. (10) Mancinelli, R.; Botti, A.; Bruni, F.; Ricci, M.; Soper, A. Perturbation of water structure due to monovalent ions in solution. Phys. Chem. Chem. Phys. 2007, 9, 2959. (11) Tielrooij, K.; Garcia-Araez, N.; Bonn, M.; Bakker, H. Cooperativity in ion hydration. Science 2010, 328, 1006. (12) Funkner, S.; Niehues, G.; Schmidt, D. A.; Heyden, M.; Schwaab, G.; Callahan, K. M.; Tobias, D. J.; Havenith, M. Watching the lowfrequency motions in aqueous salt solutions: The terahertz bibrational signatures of hydrated ions. J. Am. Chem. Soc. 2012, 134, 1030. (13) Li, R.-Z.; Liu, C.-W.; Gao, Y. Q.; Jiang, H.; Xu, H.-G.; Zheng, W.-J. Microsolvation of LiI and CsI in water: Anion photoelectron spectroscopy and ab initio calculations. J. Am. Chem. Soc. 2013, 135, 5190. (14) Wolke, C. T.; Menges, F. S.; Tötsch, N.; Gorlova, O.; Fournier, J. A.; Weddle, G. H.; Johnson, M. A.; Heine, N.; Esser, T. K.; Knorke, H.; et al. Thermodynamics of water dimer dissociation in the primary hydration shell of the iodide ion with temperature-dependent vibrational predissociation spectroscopy. J. Phys. Chem. A 2015, 119, 1859. (15) Antalek, M.; Pace, E.; Hedman, B.; Hodgson, K. O.; Chillemi, G.; Benfatto, M.; Sarangi, R.; Frank, P. Solvation structure of the halides from x-ray absorption spectroscopy. J. Chem. Phys. 2016, 145, 044318. (16) Wang, S.; Fang, W.; Li, T.; Li, F.; Sun, C.; Li, Z.; Huang, Y.; Men, Z. An insight into liquid water networks through hydrogen bonding halide anion: Stimulated Raman scattering. J. Appl. Phys. 2016, 119, 163104. (17) Chakrabarty, S.; Williams, E. R. The effect of halide and iodate anions on the hydrogen-bonding network of water in aqueous nanodrops. Phys. Chem. Chem. Phys. 2016, 18, 25483. (18) Tobias, D. J.; Jungwirth, P.; Parrinello, M. Surface solvation of halogen anions in water clusters: An ab initio molecular dynamics study of the Cl−(H2O)6 complex. J. Chem. Phys. 2001, 114, 7036. (19) Huang, X.; Habershon, S.; Bowman, J. M. Comparison of quantum, classical, and ring-polymer molecular dynamics infra-red spectra of Cl−(H2O) and H+(H2O)2. Chem. Phys. Lett. 2008, 450, 253. (20) Jahangiri, S.; Dolgonos, G.; Frauenheim, T.; Peslherbe, G. H. Parameterization of halogens for the density-functional tight-binding description of halide hydration. J. Chem. Theory Comput. 2013, 9, 3321. (21) Xu, J.-J.; Yi, H.-B.; Li, H.-J.; Chen, Y. Ionic solvation and association in LiCl aqueous solution: A density functional theory, polarised continuum model and molecular dynamics investigation. Mol. Phys. 2014, 112, 1710.

are shown to substantially affect certain structural properties in these systems. In particular, they weaken the hydrogen bond network between the water molecules leading to a significant increase of the O−O distances, albeit without appreciably affecting the other structural properties. The striking agreement for the smallest N = 1−5 clusters motivated a computationally affordable study of the solvation energies and thermodynamic properties, i.e., isomer fractions of larger (H2O)NX− (N = 6−21) clusters at finite temperatures using HSM. While we acknowledge that the agreement between HA and DMC was demonstrated only for the smallest clusters studied here, it is still plausible that this very promising progression does extend to larger cluster sizes. Of course, we cannot guarantee nor do we intend to imply that the large error cancellations observed in this work will occur for systems besides halide ion−water clusters as the system properties that lead to the aforementioned error cancellations are difficult to predict. Thus, ample caution must be exercised before applying HA to study the ground state and/or thermodynamic properties of other systems as, e.g., has been demonstrated in a recent publication on solvated Li+ cations.136



ASSOCIATED CONTENT

S Supporting Information *

The Supporting Information is available free of charge on the ACS Publications website at DOI: 10.1021/acs.jpca.8b00917. DMC convergence studies for the H2OX− and D2OX− dimers (X− = F−, Cl−, Br−, I−). Isomer fractions from HSM as a function of temperature for the F−, Cl−, and Br− trimer isotopologues. Solvation energies in Table 4 from DMC and HA represented graphically as a function of the cluster size N for the H2O1−5X− and D2O1−5X− systems. Solvation free energies as a function of the cluster size N for the (D2O)1−21X− isotopologues. Isomer fractions/populations as a function of temperature for the (D2O)20X− isotopologue and for the classical (H 2 O) 20 X − system. Tables containing the DMC parameters used for the (D2O)1−5X− systems, the mean interatomic distances and bond angles in Figure 10 and the data shown in Table 4 for the corresponding (D2O)1−5X− systems (PDF).



AUTHOR INFORMATION

Corresponding Author

*[email protected]. ORCID

Joel D. Mallory: 0000-0002-0251-5724 Vladimir A. Mandelshtam: 0000-0001-6395-1397 Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS This work was supported by the National Science Foundation (NSF) Grant No. CHE-1566334. We thank Francesco Paesani and Pushp Bajaj for many useful discussions and for sharing the MB-nrg PEF with us. We also thank Shane Flynn and Bridgett Kohno for helpful suggestions.



REFERENCES

(1) Knipping, E.; Lakin, M.; Foster, K.; Jungwirth, P.; Tobias, D.; Gerber, R.; Dabdub, D.; Finlayson-Pitts, B. Experiments and K

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A (22) Singh, M. B.; Dalvi, V. H.; Gaikar, V. G. Investigations of clustering of ions and diffusivity in concentrated aqueous solutions of lithium chloride by molecular dynamic simulations. RSC Adv. 2015, 5, 15328. (23) Liu, H.; Wang, Y.; Bowman, J. M. Quantum calculations of the IR spectrum of liquid water using ab initio and model potential and dipole moment surfaces and comparison with experiment. J. Chem. Phys. 2015, 142, 194502. (24) Choi, J.-H.; Kuwata, K. T.; Cao, Y.-B.; Okumura, M. Vibrational spectroscopy of the Cl−(H2O)n anionic clusters, n = 1−5. J. Phys. Chem. A 1998, 102, 503. (25) Ayotte, P.; Weddle, G. H.; Kim, J.; Johnson, M. A. Vibrational spectroscopy of the ionic hydrogen bond: Fermi resonances and ion− molecule stretching frequencies in the binary X−· H2O (X = Cl, Br, I) complexes via argon predissociation spectroscopy. J. Am. Chem. Soc. 1998, 120, 12361. (26) Ayotte, P.; Weddle, G. H.; Kim, J.; Johnson, M. A. Mass-selected “matrix isolation” infrared spectroscopy of the I−·(H2O)2 complex: making and breaking the inter-water hydrogen-bond. Chem. Phys. 1998, 239, 485. (27) Ayotte, P.; Bailey, C. G.; Weddle, G. H.; Johnson, M. A. Vibrational spectroscopy of small Br−· (H2O)n and I− ·(H2O)n clusters: Infrared characterization of the ionic hydrogen bond. J. Phys. Chem. A 1998, 102, 3067. (28) Ayotte, P.; Nielsen, S. B.; Weddle, G. H.; Johnson, M. A.; Xantheas, S. S. Spectroscopic observation of ion-induced water dimer dissociation in the X−·(H2O)2 (X = F, Cl, Br, I) clusters. J. Phys. Chem. A 1999, 103, 10665. (29) Ayotte, P.; Weddle, G. H.; Kim, J.; Kelley, J.; Johnson, M. A. A cluster study of anionic hydration: Spectroscopic characterization of the I−·Wn, 1 ≤ n ≤ 3, supramolecular complexes at the primary steps of solvation. J. Phys. Chem. A 1999, 103, 443. (30) Ayotte, P.; Weddle, G. H.; Johnson, M. A. An infrared study of the competition between hydrogen-bond networking and ionic solvation: Halide-dependent distortions of the water trimer in the X−·(H2O)3, (X = Cl, Br, I) systems. J. Chem. Phys. 1999, 110, 7129. (31) Ayotte, P.; Kelley, J. A.; Nielsen, S. B.; Johnson, M. A. Vibrational spectroscopy of the F− ·H2O complex via argon predissociation: photoinduced, intracluster proton transfer? Chem. Phys. Lett. 2000, 316, 455. (32) Kelley, J.; Weber, J.; Lisle, K.; Robertson, W.; Ayotte, P.; Johnson, M. The infrared predissociation spectra of Cl−·H2O·Arn (n = 1−5): experimental determination of the influence of Ar solvent atoms. Chem. Phys. Lett. 2000, 327, 1. (33) Roscioli, J. R.; Diken, E. G.; Johnson, M. A.; Horvath, S.; McCoy, A. B. Prying apart a water molecule with anionic h-bonding: A comparative spectroscopic study of the X−· H2O (X = OH, O, F, Cl, and Br) binary complexes in the 600−3800 cm−1 region. J. Phys. Chem. A 2006, 110, 4943. (34) 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. (35) Chaban, G. M.; Xantheas, S. S.; Gerber, R. B. Anharmonic vibrational spectroscopy of the F− (H2O)n complexes, n= 1,2. J. Phys. Chem. A 2003, 107, 4952. (36) Yoo, S.; Lei, Y.; Zeng, X. C. Effect of polarizability of halide anions on the ionic solvation in water clusters. J. Chem. Phys. 2003, 119, 6083. (37) Trumm, M.; Martínez, Y. O. G.; Réal, F.; Masella, M.; Vallet, V.; Schimmelpfennig, B. Modeling the hydration of mono-atomic anions from the gas phase to the bulk phase: The case of the halide ions F−, Cl−, and Br−. J. Chem. Phys. 2012, 136, 044509. (38) Vlcek, L.; Chialvo, A. A.; Simonson, J. M. Correspondence between cluster-ion and bulk solution thermodynamic properties: On the validity of the cluster-pair-based approximation. J. Phys. Chem. A 2013, 117, 11328.

(39) Neogi, S. G.; Chaudhury, P. Structure and spectroscopic aspects of water-halide ion clusters: A study based on a conjunction of stochastic and quantum chemical methods. J. Comput. Chem. 2013, 34, 471. (40) Migliorati, V.; Sessa, F.; Aquilanti, G.; D’Angelo, P. Unraveling halide hydration: A high dilution approach. J. Chem. Phys. 2014, 141, 044509. (41) Vlcek, L.; Uhlik, F.; Moucka, F.; Nezbeda, I.; Chialvo, A. A. Thermodynamics of small alkali metal halide cluster ions: Comparison of classical molecular simulations with experiment and quantum chemistry. J. Phys. Chem. A 2015, 119, 488. (42) Réal, F.; Severo Pereira Gomes, A.; Guerrero Martínez, Y. O.; Ayed, T.; Galland, N.; Masella, M.; Vallet, V. Structural, dynamical, and transport properties of the hydrated halides: How do At− bulk properties compare with those of the other halides, from F− to I−? J. Chem. Phys. 2016, 144, 124513. (43) Wen, H.; Huang, T.; Liu, Y.-R.; Jiang, S.; Peng, X.-Q.; Miao, S.K.; Wang, C.-Y.; Hong, Y.; Huang, W. Structure, temperature effect and bonding order analysis of hydrated bromide clusters. Chem. Phys. 2016, 479, 129. (44) Chiou, M.-F.; Sheu, W.-S. Charge-transfer-to-solvent absorption spectra of I−(H2O)3−5 at a finite temperature via simulation. Int. J. Quantum Chem. 2017, 117, e25404. (45) Li, J.; Wang, F. Accurate prediction of the hydration free energies of 20 salts through adaptive force matching and the proper comparison with experimental references. J. Phys. Chem. B 2017, 121, 6637. (46) Mucha, M.; Frigato, T.; Levering, L. M.; Allen, H. C.; Tobias, D. J.; Dang, L. X.; Jungwirth, P. Unified molecular picture of the surfaces of aqueous acid, base, and salt solutions. J. Phys. Chem. B 2005, 109, 7617. (47) Jungwirth, P.; Tobias, D. J. Specific ion effects at the air/water interface. Chem. Rev. 2006, 106, 1259. (48) Brown, M. A.; D’Auria, R.; Kuo, I.-F. W.; Krisch, M. J.; Starr, D. E.; Bluhm, H.; Tobias, D. J.; Hemminger, J. C. Ion spatial distributions at the liquid-vapor interface of aqueous potassium fluoride solutions. Phys. Chem. Chem. Phys. 2008, 10, 4778. (49) Tobias, D. J.; Stern, A. C.; Baer, M. D.; Levin, Y.; Mundy, C. J. Simulation and theory of ions at atmospherically relevant aqueous liquid-air interfaces. Annu. Rev. Phys. Chem. 2013, 64, 339. (50) Jungwirth, P.; Tobias, D. J. Molecular structure of salt solutions: A new view of the interface with implications for heterogeneous atmospheric chemistry. J. Phys. Chem. B 2001, 105, 10468. (51) Jungwirth, P.; Tobias, D. J. Ions at the air/water interface. J. Phys. Chem. B 2002, 106, 6361. (52) Jungwirth, P.; Tobias, D. J. Chloride anion on aqueous clusters, at the air-water interface, and in liquid water: Solvent effects on Cl− polarizability. J. Phys. Chem. A 2002, 106, 379. (53) Piatkowski, L.; Zhang, Z.; Backus, E. H.; Bakker, H. J.; Bonn, M. Extreme surface propensity of halide ions in water. Nat. Commun. 2014, 5, 4083. (54) Bowman, J. M.; Xantheas, S. S. ”Morphing” of ab initio-based interaction potentials to spectroscopic accuracy: Application to Cl−(H2O). Pure Appl. Chem. 2004, 76, 29. (55) Rheinecker, J.; Bowman, J. M. The calculated infrared spectrum of Cl−H2O using a new full dimensional ab initio potential surface and dipole moment surface. J. Chem. Phys. 2006, 124, 131102. (56) Rheinecker, J.; Bowman, J. M. The calculated infrared spectrum of Cl−H2O using a new full dimensional ab initio potential surface and dipole moment surface. J. Chem. Phys. 2006, 125, 133206. (57) Kamarchik, E.; Bowman, J. M. Quantum vibrational analysis of hydrated ions using an ab initio potential. J. Phys. Chem. A 2010, 114, 12945. (58) Toffoli, D.; Sparta, M.; Christiansen, O. Vibrational spectroscopy of hydrogen-bonded systems: Six-dimensional simulation of the IR spectrum of F−(H2O) complex. Chem. Phys. Lett. 2011, 510, 36. (59) Kamarchik, E.; Bowman, J. M. Coupling of low-and highfrequency vibrational modes: Broadening in the infrared spectrum of F−(H2O)2. J. Phys. Chem. Lett. 2013, 4, 2964. L

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A (60) Kamarchik, E.; Toffoli, D.; Christiansen, O.; Bowman, J. M. Ab initio potential energy and dipole moment surfaces of the F−(H2O) complex. Spectrochim. Acta, Part A 2014, 119, 59. (61) Diken, E. G.; Shin, J.-W.; Price, E. A.; Johnson, M. A. Isotopic fractionation and zero-point effects in anionic H-bonded complexes: a comparison of the I−· HDO and F−· HDO ion−molecule clusters. Chem. Phys. Lett. 2004, 387, 17. (62) Horvath, S.; McCoy, A. B.; Roscioli, J. R.; Johnson, M. A. Vibrationally induced proton transfer in F−(H2O) and F−(D2O). J. Phys. Chem. A 2008, 112, 12337. (63) Horvath, S.; McCoy, A. B.; Elliott, B. M.; Weddle, G. H.; Roscioli, J. R.; Johnson, M. A. Anharmonicities and isotopic effects in the vibrational spectra of X− · H2O, · HDO, and · D2O [X = Cl, Br, and I] binary complexes. J. Phys. Chem. A 2010, 114, 1556. (64) Videla, P. E.; Rossky, P. J.; Laria, D. Isotopic preferential solvation of I− in low-temperature water nanoclusters. J. Phys. Chem. B 2015, 119, 11783. (65) Videla, P. E.; Rossky, P. J.; Laria, D. Isotope effects in aqueous solvation of simple halides. J. Chem. Phys. 2018, 148, 102306. (66) Litman, Y. E.; Videla, P. E.; Rodriguez, J.; Laria, D. Positional isotope exchange in HX·(H2O)n (X = F, I) clusters at low temperatures. J. Phys. Chem. A 2016, 120, 7213. (67) Yang, N.; Duong, C. H.; Kelleher, P. J.; Johnson, M. A.; McCoy, A. B. Isolation of site-specific anharmonicities of individual water molecules in the I− ·(H2O)2 complex using tag-free, isotopomer selective IR-IR double resonance. Chem. Phys. Lett. 2017, 690, 159. (68) Xantheas, S. S.; Dunning, T. H., Jr Structures and energetics of F−(H2O)n, n = 1−3 Clusters from ab Initio calculations. J. Phys. Chem. 1994, 98, 13489. (69) Xantheas, S. S. Quantitative description of hydrogen bonding in chloride-water clusters. J. Phys. Chem. 1996, 100, 9703. (70) Baik, J.; Kim, J.; Majumdar, D.; Kim, K. S. Structures, energetics, and spectra of fluoride-water clusters F−(H2O)n, n = 1−6: Ab initio study. J. Chem. Phys. 1999, 110, 9116. (71) Kim, J.; Lee, H. M.; Suh, S. B.; Majumdar, D.; Kim, K. S. Comparative ab initio study of the structures, energetics and spectra of X−·(H2O)n=1−4 [X = F, Cl, Br, I] clusters. J. Chem. Phys. 2000, 113, 5259. (72) Irle, S.; Bowman, J. M. Direct ab initio variational calculation of vibrational energies of the H2O···Cl− complex and resolution of experimental differences. J. Chem. Phys. 2000, 113, 8401. (73) Lee, H. M.; Kim, D.; Kim, K. S. Structures, spectra, and electronic properties of halide-water pentamers and hexamers, X−(H2O)5,6 (X = F, Cl, Br, I): Ab initio study. J. Chem. Phys. 2002, 116, 5509. (74) Xantheas, S. S.; Dang, L. X. Critical study of fluoride-water interactions. J. Phys. Chem. 1996, 100, 3989. (75) Dorsett, H. E.; Watts, R. O.; Xantheas, S. S. Probing temperature effects on the hydrogen bonding network of the Cl−(H2O)2 cluster. J. Phys. Chem. A 1999, 103, 3351. (76) Kiss, P. T.; Baranyai, A. A new polarizable force field for alkali and halide ions. J. Chem. Phys. 2014, 141, 114501. (77) Werhahn, J. C.; Akase, D.; Xantheas, S. S. Universal scaling of potential energy functions describing intermolecular interactions. II. The halide-water and alkali metal-water interactions. J. Chem. Phys. 2014, 141, 064118. (78) 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. (79) Kamarchik, E.; Wang, Y.; Bowman, J. M. Quantum vibrational analysis and infrared spectra of microhydrated sodium ions using an ab initio potential. J. Chem. Phys. 2011, 134, 114311. (80) Wang, Y.; Bowman, J. M.; Kamarchik, E. Five ab initio potential energy and dipole moment surfaces for hydrated NaCl and NaF. I. Two-body interactions. J. Chem. Phys. 2016, 144, 114311. (81) 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.

(82) 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. (83) Liu, H.; Wang, Y.; Bowman, J. M. Vibrational analysis of an ice Ih model from 0 to 4000 cm−1 using the ab initio WHBB potential energy surface. J. Phys. Chem. B 2013, 117, 10046. (84) Liu, H.; Wang, Y.; Bowman, J. M. Local-monomer calculations of the intramolecular IR spectra of the cage and prism isomers of HOD(D2O)5 and HOD and D2O ice Ih. J. Phys. Chem. B 2014, 118, 14124−14131. (85) 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. (86) Liu, H.; Wang, Y.; Bowman, J. M. Quantum local monomer IR spectrum of liquid D2O at 300 K from 0 to 4000 cm−1 is in nearquantitative agreement with experiment. J. Phys. Chem. B 2016, 120, 2824. (87) Wang, Y.; Bowman, J. M. Calculations of the IR spectra of bend fundamentals of (H2O)n=3, 4, 5 using the WHBB_2 potential and dipole moment surfaces. Phys. Chem. Chem. Phys. 2016, 18, 24057. (88) Bajaj, P.; Götz, A. W.; Paesani, F. Toward chemical accuracy in the description of ion-water interactions through many-body representations. I. Halide-water dimer potential energy surfaces. J. Chem. Theory Comput. 2016, 12, 2698. (89) Riera, M.; Götz, A. W.; Paesani, F. The i-TTM model for ab initio-based ion-water interaction potentials. II. Alkali metal ion-water potential energy functions. Phys. Chem. Chem. Phys. 2016, 18, 30334. (90) Riera, M.; Mardirossian, N.; Bajaj, P.; Götz, A. W.; Paesani, F. Toward chemical accuracy in the description of ion-water interactions through many-body representations. Alkali-water dimer potential energy surfaces. J. Chem. Phys. 2017, 147, 161715. (91) Arismendi-Arrieta, D. J.; Riera, M.; Bajaj, P.; Prosmiti, R.; Paesani, F. i-TTM model for ab initio-based ion-water interaction potentials. 1. Halide-water potential energy functions. J. Phys. Chem. B 2016, 120, 1822. (92) Babin, V.; Leforestier, C.; Paesani, F. Development of a “first principles” water potential with flexible monomers: Dimer potential energy surface, VRT spectrum, and second virial coefficient. J. Chem. Theory Comput. 2013, 9, 5395. (93) Babin, V.; Medders, G. R.; Paesani, F. Development of a “first principles” water potential with flexible monomers. II: Trimer potential energy surface, third virial coefficient, and small clusters. J. Chem. Theory Comput. 2014, 10, 1599. (94) Mallory, J. D.; Mandelshtam, V. A. Diffusion Monte Carlo studies of MB-pol (H2O)2−6 and (D2O)2−6 clusters: Structures and binding energies. J. Chem. Phys. 2016, 145, 064308. (95) 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. (96) Medders, G. R.; Babin, V.; Paesani, F. Development of a “firstprinciples” water potential with flexible monomers. III. Liquid phase properties. J. Chem. Theory Comput. 2014, 10, 2906. (97) Pham, C. H.; Reddy, S. K.; Chen, K.; Knight, C.; Paesani, F. Many-body interactions in ice. J. Chem. Theory Comput. 2017, 13, 1778. (98) Brown, S. E.; Gö tz, A. W.; Cheng, X.; Steele, R. P.; Mandelshtam, V. A.; Paesani, F. Monitoring water clusters “melt” through vibrational spectroscopy. J. Am. Chem. Soc. 2017, 139, 7082. (99) Reddy, S. K.; Straight, S. C.; Bajaj, P.; Huy Pham, C.; Riera, M.; Moberg, D. R.; Morales, M. A.; Knight, C.; Götz, A. W.; Paesani, F. On the accuracy of the MB-pol many-body potential for water: Interaction energies, vibrational frequencies, and classical thermodynamic and dynamical properties from clusters to liquid water and ice. J. Chem. Phys. 2016, 145, 194504. (100) Videla, P. E.; Rossky, P. J.; Laria, D. Communication: Isotopic effects on tunneling motions in the water trimer. J. Chem. Phys. 2016, 144, 061101. M

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX

Article

The Journal of Physical Chemistry A (101) Bishop, K. P.; Roy, P.-N. Quantum mechanical free energy profiles with post-quantization restraints: Binding free energy of the water dimer over a broad range of temperatures. J. Chem. Phys. 2018, 148, 102303. (102) Gaiduk, A.; Pham, T. A.; Govoni, M.; Paesani, F.; Galli, G. Electron affinity of liquid water. Nat. Commun. 2018, 9, 247. (103) Richardson, J. O.; Pérez, C.; Lobsiger, S.; Reid, A. A.; Temelso, B.; Shields, G. C.; Kisiel, Z.; Wales, D. J.; Pate, B. H.; Althorpe, S. C. Concerted hydrogen-bond breaking by quantum tunneling in the water hexamer prism. Science 2016, 351, 1310. (104) Mallory, J. D.; Mandelshtam, V. A. Binding energies from diffusion Monte Carlo for the MB-pol H2O and D2O dimer: A comparison to experimental values. J. Chem. Phys. 2015, 143, 144303. (105) Medders, G. R.; Paesani, F. Infrared and Raman spectroscopy of liquid water through “first-principles” many-body molecular dynamics. J. Chem. Theory Comput. 2015, 11, 1145. (106) Medders, G. R.; Paesani, F. On the interplay of the potential energy and dipole moment surfaces in controlling the infrared activity of liquid water. J. Chem. Phys. 2015, 142, 212411. (107) Medders, G. R.; Götz, A. W.; Morales, M. A.; Bajaj, P.; Paesani, F. On the representation of many-body interactions in water. J. Chem. Phys. 2015, 143, 104102. (108) Medders, G. R.; Paesani, F. Dissecting the molecular structure of the air/water interface from quantum simulations of the sumfrequency generation spectrum. J. Am. Chem. Soc. 2016, 138, 3912. (109) Straight, S. C.; Paesani, F. Exploring electrostatic effects on the hydrogen bond network of liquid water through many-body molecular dynamics. J. Phys. Chem. B 2016, 120, 8539. (110) Moberg, D. R.; Straight, S. C.; Knight, C.; Paesani, F. Molecular origin of the vibrational structure of ice Ih. J. Phys. Chem. Lett. 2017, 8, 2579. (111) Reddy, S. K.; Moberg, D. R.; Straight, S. C.; Paesani, F. Temperature-dependent vibrational spectra and structure of liquid water from classical and quantum simulations with the MB-pol potential energy function. J. Chem. Phys. 2017, 147, 244504. (112) Bajaj, P.; Wang, X.-G.; Carrington, T., Jr; Paesani, F. Vibrational spectra of halide-water dimers: Insights on ion hydration from full-dimensional quantum calculations on many-body potential energy surfaces. J. Chem. Phys. 2018, 148, 102321. (113) Wales, D. J. Coexistence in small inert gas clusters. Mol. Phys. 1993, 78, 151. (114) Anderson, J. B. A random-walk simulation of the Schrdinger equation: H+3 . J. Chem. Phys. 1975, 63, 1499. (115) Anderson, J. B. Quantum chemistry by random walk. H 2P, H+3 D3h 1A1′ , H2 3Σ+u , H4 1Σ+g , Be 1S. J. Chem. Phys. 1976, 65, 4121. (116) Mallory, J. D.; Brown, S. E.; Mandelshtam, V. A. Assessing the performance of the diffusion Monte Carlo method as applied to the water monomer, dimer, and hexamer. J. Phys. Chem. A 2015, 119, 6504. (117) Mallory, J. D.; Mandelshtam, V. A. Quantum melting and isotope effects from diffusion Monte Carlo studies of p-H2 clusters. J. Phys. Chem. A 2017, 121, 6341. (118) Johnson, L. M.; McCoy, A. B. Evolution of structure in CH+5 and its deuterated analogues. J. Phys. Chem. A 2006, 110, 8213. (119) McCoy, A. B. Diffusion Monte Carlo approaches for investigating the structure and vibrational spectra of fluxional systems. Int. Rev. Phys. Chem. 2006, 25, 77. (120) Acioli, P. H.; Xie, Z.; Braams, B. J.; Bowman, J. M. Vibrational ground state properties of H+5 and its isotopomers from diffusion Monte Carlo calculations. J. Chem. Phys. 2008, 128, 104318. (121) Barragán, P.; Pérez de Tudela, R.; Qu, C.; Prosmiti, R.; Bowman, J. M. Full-dimensional quantum calculations of the dissociation energy, zero-point, and 10 k properties of H+7 /D+7 clusters using an ab initio potential energy surface. J. Chem. Phys. 2013, 139, 024308. (122) Mancini, J. S.; Samanta, A. K.; Bowman, J. M.; Reisler, H. Experiment and theory elucidate the multichannel predissociation dynamics of the HCl trimer: Breaking up is hard to do. J. Phys. Chem. A 2014, 118, 8402.

(123) 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. (124) Czakó, G.; Wang, Y.; Bowman, J. M. Communication: Quasiclassical trajectory calculations of correlated product-state distributions for the dissociation of (H2O)2 and (D2O)2. J. Chem. Phys. 2011, 135, 151102. (125) 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. (126) 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. (127) Suhm, M. A.; Watts, R. O. Quantum Monte Carlo studies of vibrational states in molecules and clusters. Phys. Rep. 1991, 204, 293. (128) Warren, G. L.; Hinde, R. J. Population size bias in descendantweighted diffusion quantum Monte Carlo simulations. Phys. Rev. E 2006, 73, 056706. (129) Doye, J. P.; Calvo, F. Entropic effects on the structure of Lennard-Jones clusters. J. Chem. Phys. 2002, 116, 8307. (130) Temelso, B.; Shields, G. C. The role of anharmonicity in hydrogen-bonded systems: The case of water clusters. J. Chem. Theory Comput. 2011, 7, 2804. (131) Temelso, B.; Archer, K. A.; Shields, G. C. Benchmark structures and binding energies of small water clusters with anharmonicity corrections. J. Phys. Chem. A 2011, 115, 12034. (132) Gillis, N. S.; Werthamer, N. R.; Koehler, T. R. Properties of crystalline argon and neon in the self-consistent phonon approximation. Phys. Rev. 1968, 165, 951. (133) Koehler, T. R. Theory of the self-consistent harmonic approximation with application to solid neon. Phys. Rev. Lett. 1966, 17, 89. (134) Brooks, B. R.; Janezic, D.; Karplus, M. Harmonic analysis of large systems. I. Methodology. J. Comput. Chem. 1995, 16, 1522. (135) Tissandier, M. D.; Cowen, K. A.; Feng, W. Y.; Gundlach, E.; Cohen, M. H.; Earhart, A. D.; Coe, J. V.; Tuttle, T. R. The proton’s absolute aqueous enthalpy and Gibbs free energy of solvation from cluster-ion solvation data. J. Phys. Chem. A 1998, 102, 7787. (136) Dupuis, R.; Benoit, M.; Tuckerman, M. E.; Méheut, M. Importance of a fully anharmonic treatment of equilibrium isotope fractionation properties of dissolved ionic species as evidenced by Li+(aq). Acc. Chem. Res. 2017, 50, 1597.

N

DOI: 10.1021/acs.jpca.8b00917 J. Phys. Chem. A XXXX, XXX, XXX−XXX