Unravelling the influence of quantum proton delocalization on electronic charge transfer through the hydrogen bond

TL;DR

This paper investigates how nuclear quantum effects influence electronic charge transfer in hydrogen bonds, revealing that quantum nuclear motion significantly alters charge transfer magnitude and temperature dependence, crucial for accurate modeling.

Contribution

It demonstrates that nuclear quantum effects are essential for understanding charge transfer in hydrogen bonds and identifies a geometric parameter that governs this transfer.

Findings

Quantum nuclear effects cause large changes in charge transfer magnitude.

Charge transfer exhibits significant temperature dependence due to quantum effects.

A single geometric parameter determines charge transfer in hydrogen bonds.

Abstract

Upon hydrogen bond formation, electronic charge density is transferred between the donor and acceptor, impacting processes ranging from hydration to spectroscopy. Here we use ab initio path integral simulations to elucidate the role of nuclear quantum effects in determining the charge transfer in a range of hydrogen bonded species in the gas and liquid phase. We show that the quantization of the nuclei gives rise to large changes in the magnitude of the charge transfer as well as its temperature dependence. We then explain how a single geometric parameter determines the charge transfer through the hydrogen bond. These results thus demonstrate that nuclear quantum effects are vital for the accurate description of charge transfer and offer a physically transparent way to understand how hydrogen bonding gives rise to it.

Tables

Table 1. Table 1: Parameters of the non-linear least square fit of charge transfer as a function of the proton sharing coordinate τ 𝜏 \tau for NPA, Bader’s QTAIM and the Hirshfeld method assuming the functional form given by Eq.  3 .

	NPA	QTAIM	Hirshfeld
$a$	$0.6254$	$0.2107$	$0.4872$
$b$	$-0.5932$	$-0.6464$	$-0.6526$
$c$	$2.5937$	$1.8677$	$1.4376$

Equations

$$\displaystyle\Delta q(\tau)=\left\{\begin{array}[]{ll}a\cdot\left(\tau-b\right)^{c}&\quad\text{if}\quad\tau>b\\ 0&\quad\text{if}\quad\tau\leq b\end{array}\right.,$$

$$q_i=\sum _{n=1}^{N_{\text{D}}}\Delta q({\tau }_n)-\sum _{m=1}^{N_{\text{A}}}\Delta q({\tau }_m),$$

Full text

charge transfer, nuclear quantum effects, hydrogen bonding

I Introduction

The transfer of electronic charge density through a hydrogen bond manifests in effects ranging from the stabilization and spectroscopy of hydrated ions and aqueous solutions Zhao, Rogers, and Beck (2010); Ben-Amotz (2011); Thompson and Hynes (2000); Torii (2014) to the charging of polar residues in proteins Ufimtsev, Luehr, and Martinez (2011), and has even been implicated in explaining the observed zeta potential at hydrophobic interfaces Takahashi (2005); Vácha et al. (2011a); Ben-Amotz (2011). However, defining the amount of charge transfer (CT) remains a subject of significant debate, since partitioning the electron density and assigning it to particular atoms or molecules is not uniquely defined. This has led to the introduction of a number of approaches to calculate the charge transfer and charge transfer energy, with different approaches exhibiting both different quantitative and qualitative behavior Mulliken (1955); Hirshfeld (1977); Reed, Weinstock, and Weinhold (1985); Stevens and Fink (1987); Breneman and Wiberg (1990); Bader (1991); Stone (1993); Lillestolen and Wheatley (2009); Azar and Head-Gordon (2012). For example, methods based on natural bond orbitals suggest that without CT the 0 K water dimer would adopt a different structure Reed and Weinhold (1983); Glendening (2005). However, studies using other methods suggest a much less substantial role of CT in stabilizing the water dimer hydrogen bond Khaliullin, Bell, and Head-Gordon (2009); Ronca, Belpassi, and Tarantelli (2014); Stone and Misquitta (2009); Stone (2017).

Previous CT studies of hydrogen bonded systems have analyzed minimum energy structures Glendening (2005); Khaliullin, Bell, and Head-Gordon (2009); Ronca, Belpassi, and Tarantelli (2014) or those obtained from simulations that treat the nuclei classically Lee and Rick (2011); Vácha et al. (2011a). However, nuclear quantum effects (NQEs), such as zero-point energy and tunnelling, are known to significantly change the geometry of hydrogen bonded systems Tuckerman et al. (1997); Marx et al. (1999); Morrone and Car (2008); Zeidler et al. (2012); Ceriotti et al. (2013); Wang, Ceriotti, and Markland (2014). The zero-point energy along the hydrogen bond increases the ability of the proton to delocalize between the donor and acceptor, resulting in shorter and stronger hydrogen bonds. On the other hand, the distortion of the hydrogen bond due to quantization of the librational modes perpendicular to the hydrogen bond leads to weaker hydrogen bonding. The relative importance of these two competing quantum effects Habershon, Markland, and Manolopoulos (2009) leads to cases where NQEs strengthen Li, Walker, and Michaelides (2011) or weaken McKenzie et al. (2014) hydrogen bonds, depending on the geometry and chemical environment. This competition has been used to elucidate the seemingly anomalously small NQEs and resulting isotope effects observed in many systems Ceriotti et al. (2016).

Here we investigate the role of NQEs in modulating the CT characteristics of a range of hydrogen bonded systems involving water and ions in both the gas and condensed phase and show how these can be understood based on the changes in the hydrogen bond geometry. We demonstrate that NQEs give rise to large changes in both the amount of CT and its contribution to the donor-acceptor interaction energy, independent of the method used to study charge transfer. In particular, including NQEs increases CT in some cases, while decreasing it in others. It also leads to a much milder temperature dependence of CT in all cases. To uncover the reasons for this, we show that the amount of charge transfer is almost exclusively determined by a scaled proton sharing coordinate. We use this observation to demonstrate that the charge of a water molecule in liquid water can be accurately estimated from the values of this coordinate for each hydrogen bond the molecule forms, thus provides a straightforward approach to modelling CT in solution.

II Simulation Details

To evaluate NQEs on electronic properties in the gas phase, $20\text{\,}\mathrm{ps}$ ab initio path integral molecular dynamics (AI-PIMD) simulations at $25\text{\,}\mathrm{K}$ and $300\text{\,}\mathrm{K}$ of the water dimer as well as the water-fluoride and water-iodide complex were performed. These simulations were performed with a time step of $0.5\text{\,}\mathrm{fs}$ with 256 and 32 replicas, at the two respective temperatures. The B3LYP Becke (1988); Lee, Yang, and Parr (1988) exchange correlation functional and the 6-311++G** Pople basis set Krishnan et al. (1980); Clark et al. (1983) for the water dimer and fluoride-water complex as well as the def2-SVPD basis set Peterson et al. (2003); Rappoport and Furche (2010) for the iodide-water complex were used to calculate electronic interactions on the fly during the simulations. An additional $1\text{\,}\mathrm{ps}$ at the beginning of each simulation was discarded as equilibration after starting from the optimized minimum energy configuration. Langevin thermostats Leimkuhler and Matthews (2013) with optimal coupling to each staging mode Ceriotti et al. (2010) were used to sample the canonical distribution. Two independent $21\text{\,}\mathrm{ps}$ classical ab initio molecular dynamics (AIMD) simulations of each system were carried out using the same electronic structure method, time step and temperature, together with massive Nosé Hoover chains thermostatting Martyna et al. (1996) with a chain length of 6. The simulations were performed with the in-house Python program PyMD Marsalek et al. (2014) that uses IPython parallelization Pérez and Granger (2007) and the Atomic Simulation Environment Bahn and Jacobsen (2002) to link to Gaussian 09 Frisch et al. (2009) for the electronic structure calculation. Our choice of basis set and number of replicas gives good convergence for geometrical properties (see SI). Since it is known that density functional theory (DFT) can cause spurious charge delocalization effects Mori-Sánchez, Cohen, and Yang (2008), we validated our choice of methods by also conducting MP2 simulations of the water dimer, as detailed in the SI. The geometries were further validated by comparing to the very accurate MBPol force field Babin, Leforestier, and Paesani (2013); a detailed analysis can be found in the SI.

Population analyses of the gas-phase configurations were performed with the B3LYP density functional and the augmented cc-pVTZ Dunning basis Kendall, Jr, and Harrison (1992) (water dimer and fluoride-water complex) or the def2-QZVPPD basis Peterson et al. (2003); Rappoport and Furche (2010) (iodide-water complex) using Gaussian 09 Frisch et al. (2009) where the basis sets were chosen to give basis set converged results (see SI). For the AIMD trajectories, configurations every $5\text{\,}\mathrm{fs}$ were analyzed, whereas for the AI-PIMD simulations at $25\text{\,}\mathrm{K}$ configurations every $25\text{\,}\mathrm{fs}$ of every 10${}^{\text{th}}$ bead were used. Configurations of the AI-PIMD simulations at $300\text{\,}\mathrm{K}$ were selected from every 5${}^{\text{th}}$ bead every $5\text{\,}\mathrm{fs}$. Note that we applied a standard hydrogen bond criterion (cutoff of 3.5 Å on the oxygen-oxygen distance and 130 on the O-H-O angle) for the water dimer at $300\text{\,}\mathrm{K}$ to discard non-bonded configurations caused by thermal fluctuations.

Bader’s quantum theory of atoms in molecules (QTAIM) charges Bader (1991) were calculated with the AIMAll program Keith from the wavefunction provided by Gaussian 09, Hirshfeld charges were calculated directly with Gaussian 09 and natural population analysis (NPA) charges with calls to the NBO 6.0 program Glendening, Landis, and Weinhold (2013). SAPT(2+3) interaction energies, including the extension for charge transfer Hohenstein and Sherrill (2010); Stone and Misquitta (2009), were obtained using Psi4 Turney et al. (2012) with the aug-cc-pVTZ basis set and the corresponding auxiliary basis set for density fitting.

To analyze electronic properties in liquid water including an explicit treatment of nuclear quantum effects, AI-PIMD and AIMD trajectories obtained in our previous study Wang, Ceriotti, and Markland (2014) were used. The electronic structure was evaluated with the Quickstep module of the CP2K package VandeVondele et al. (2005) using the BLYP-D3 functional Becke (1988); Lee, Yang, and Parr (1988); Grimme et al. (2010) and a DZVP basis. The i-PI program Ceriotti, More, and Manolopoulos (2014) was used and in the AI-PIMD simulations each atom was represented by 6 beads using the PIGLET algorithm Ceriotti and Manolopoulos (2012).

For the analysis of the fluoride ion solvated in water, AI-PIMD and AIMD simulations of the ion with 64 water molecules in periodic boundary conditions with a box size of 12.414 Å were run for $27\text{\,}\mathrm{ps}$ and $40\text{\,}\mathrm{ps}$, respectively, with a plane wave cutoff of 280 Ry and settings otherwise identical to those for the liquid water simulations.

To obtain the molecular charges in the condensed phase (both neat water and hydrated fluoride), the full electronic structure of the box was calculated for snapshots taken from our trajectories. Configurations of the classical AIMD trajectory and of one bead of the AI-PIMD trajectory were analyzed every $12.5\text{\,}\mathrm{fs}$ using the Quickstep module of the CP2K package VandeVondele et al. (2005). Stricter electronic structure settings were used, based on thorough testing of the convergence of the charges (see SI). Specifically, the BLYP-D3 functional was kept and the SCF convergence was increased to $\epsilon_{\rm{SCF}}$=$1.0\text{\times}{10}^{-10}$, while the cutoff for the plane-wave representation of the charge density was set to 900 Ry. The m-TZV2P basis set optimized for molecular systems (molopt) VandeVondele and Hutter (2007) was employed and the three population analyses also used for the gas-phase systems were performed in the following way. To obtain Bader charges, the electronic density on a grid was analyzed with the Bader Charge Analysis program of the Henkelman group Tang, Sanville, and Henkelman (2009). Hirshfeld charges were determined directly with CP2K. NPA charges are not available directly from CP2K, but were obtained by writing out the atomic orbital overlap and density matrix of every configuration and preparing NBO 6.0 Glendening, Landis, and Weinhold (2013) input.

To enable the decomposition of charge populations on molecules in liquid water, individual charge transfer contributions were calculated for hydrogen bonded dimers extracted from the bulk geometries. For this, hydrogen bonded configurations were detected by applying a cutoff on the oxygen-oxygen distance $r_{\text{OO}}$ of 3.5 Å, which corresponds to the first minimum of the radial distribution function, and a cutoff on the hydrogen bond angle $\measuredangle_{\text{OHO}}$ of 130. Overall, $20\,000$ dimer configurations were randomly extracted from the quantum ensemble, since it features a wider range of structures than the classical ensemble.

Population analyses for extracted dimer configurations were performed using the BLYP-D3 functional with otherwise similar specifications as for the gas-phase dimers. We repeated the population analyses also using MP2 to ensure that the CT is not an artifact of DFT, obtaining results comparable to those from DFT (see Fig. S12 in the SI for details).

III Results

In what follows we investigate CT in several hydrogen-bonded systems over a range of temperatures in both the gas and condensed phase. The water dimer was chosen as it is the prototypical hydrogen bond, the iodide-water complex was chosen as an intermediate strength hydrogen bond that is dominated by the charge of the diffuse and polarizable iodide and the fluoride-water complex was studied as an example with a very strong hydrogen bond. First we consider these systems in the gas phase at the low temperature of 25 K as well as at ambient temperature to isolate the ground state contribution to CT. Second, to understand molecular charge populations in the condensed phase, we investigate CT in liquid water and around the hydrated fluoride anion. Finally, we show that it is possible to construct a simple model of CT based on proton sharing in hydrogen bonds.

For clarity, in the main text we report only NPA charge populations, but results obtained with QTAIM and the Hirshfeld method are reported in the SI. It should be noted that although the absolute amount of charge varies with the different population analyses, the trends and conclusions we demonstrate are universal.

III.1 Gas-phase systems

We start by examining the amount of charge transfer in the three gas phase complexes — water-water, water-iodide, and water-fluoride — at temperatures of $25\text{\,}\mathrm{K}$ and $300\text{\,}\mathrm{K}$. At $25\text{\,}\mathrm{K}$, the properties of all the systems are dominated by the vibrational ground state. Hence, the influence of NQEs when minimal thermal fluctuations are present can be assessed. The results of the NPA CT analysis are shown in the top 6 panels of Fig. 1, which compares the classical and quantum ensembles as well as the minimum energy structure. The CT is largest for the strongly hydrogen-bonded fluoride-water complex, less substantial for the iodide-water complex, which is of intermediate strength, and smallest for the neutral water dimer. At $25\text{\,}\mathrm{K}$, the AIMD simulations, which include only classical thermal fluctuations, show a narrow symmetric distribution of CT around an average which is very close to the value obtained for the minimum energy structure. In contrast, when NQEs are included in the description of the system using AI-PIMD simulations, the distributions become substantially broader and more asymmetric with a long tail of high CT values. This shifts the average CT away from that obtained for the minimum energy structure. However, the direction of this shift differs for the three complexes, with the water dimer average CT decreased by 27 %, while the fluoride-water complex CT is increased by 16 %. These results thus show that, due to the asymmetry and the shift of the average in the quantum simulations, the minimum energy structure does not provide an adequate approximation to the quantum mechanical ground state.

Upon increasing the temperature to 300 K, the difference between the classical and quantum CT distributions markedly decreases. This is largely due to changes in the classical distributions, while the quantum distributions are essentially unchanged from those at 25 K, suggesting that even at 300 K, CT is dominated by the ground state contribution. While for the water dimer the classical and quantum distributions look essentially identical, in the case of the two halide-water complexes, there are small residual differences, with the quantum average matching the minimum structure in the case of iodide and the classical average matching it for fluoride. It is important to note, however, that given the shape and spread of the distributions, these matches should be considered entirely coincidental and should not be relied on in general.

To assess the structural changes which modulate the CT, we consider the elliptical coordinates defined in Fig. 2. These describe the position of the proton in the hydrogen bond in terms of a proton sharing coordinate $\tau$, a dimensionless version of another commonly used proton sharing coordinate Lobaugh and Voth (1996); Marx et al. (1999); Ceriotti et al. (2013), and an orthogonal coordinate $\sigma$ which describes the deviation from linearity. A higher value of $\tau$ corresponds to a hydrogen bond with the proton closer to the acceptor heavy atom, with $\tau$=0 representing a proton being shared equally between the two heavy atoms. As shown in Fig. S6, this coordinate strongly correlates with the amount of charge transfer. This arises since, when the proton is more shared, the electron density shifts more towards the acceptor molecule, resulting in an increase in charge transfer. As depicted in Fig. 1, the changes of the distribution of $\tau$ (bottom 6 panels) are thus mirrored in the behaviour of the CT distribution (top 6 panels). At the lower temperature, the distribution of $\tau$ is spread out considerably upon including NQEs for all systems. In addition, for the water-dimer, NQEs decrease the average value of $\tau$, while for iodide-water it remains unchanged and for fluoride-water the average of $\tau$ is increased. At ambient temperature, the classical and quantum distributions become almost identical for these gas phase systems, again matching the behavior of the magnitude of CT.

Thus far, our analysis was based on the partitioning of the electronic density to determine the amount of electrons partially transferred through the hydrogen bond. To provide a measure of how this CT affects the binding energy of the complex, we performed calculations for the water dimer and fluoride-water complexes using the symmetry adapted perturbation theory (SAPT) method with an extension to partition the induction energy into charge transfer and polarization contributions. Fig. 3 shows that the charge transfer energy is strongly correlated with the amount of CT. Overall, the CT energy is in the range of $\sim$0–4 kJ mol*-1* for the water dimer, which corresponds to 15 % of the total binding energy. As already implied by the much larger amount of CT, the range of CT energies for the fluoride-water complex is over one order of magnitude larger, $\sim$20–100 kJ mol*-1*, and makes up roughly 30 % of the total binding energy. Owing to this correlation, the CT energies are influenced by NQEs in a way analogous to the amount of CT and the proton sharing coordinate.

III.2 Condensed phase

We now consider molecular CT in the condensed phase for ambient temperature liquid water as well as the hydrated fluoride anion. The distributions of molecular NPA charge populations for these two systems with a classical and quantum description of the nuclei are depicted in Fig. 4, while results for the other population analyses are shown in Fig. S9. In the case of the fluoride anion, CT only occurs in one direction — electronic density is transferred from the fluoride anion to the surrounding water molecules through the hydrogen bonds that they donate to the ion. Any deviation of charge population on the fluoride anion from the formal charge of -1 is thus the amount of CT due to all hydrogen bonds it receives in the liquid. The hydrated fluoride anion in the liquid phase has an average CT of 0.19 e in the classical case. This is only $\sim$50 % higher than for the monohydrated anion in the gas phase, despite the increase in the number of hydrogen bonds it receives from 1 to an average of $\sim$4.5–5.0. NQEs increase the average CT by 11 % (slightly less than the 18 % observed in the gas phase) and markedly increase the range of values observed towards higher charge transfer. Water molecules, unlike the fluoride anion, can both accept and donate hydrogen bonds. Due to the symmetry of the hydrogen bonding in the bulk liquid, the distribution of molecular charges is symmetric and has an average of zero and thus only the fluctuations can vary. Our simulations show that NQEs increase these fluctuations around zero charge — the quantum distribution in the top panel of Fig. 4 is considerably wider than the classical distribution.

III.3 Analysis of charge transfer in liquid water

Based on our observations above, we now demonstrate how one can provide a simple but accurate model of CT through the hydrogen bond. Previous investigations of charge transfer between water molecules in the liquid phase have modelled the molecular charge distribution by using a procedure in which formation of a hydrogen bond contributes a constant amount of charge transfer Lee and Rick (2011); Vácha et al. (2011b, a); Rick (2016). In these studies, each molecule was assigned a hydrogen bond donor-acceptor balance $\Delta_{\text{DA}}$, calculated as the difference between the total number of donor $N_{\text{D}}$ and the total number of acceptor hydrogen bonds $N_{\text{A}}$ of a given molecule, $\Delta_{\text{DA}}=N_{\text{D}}-N_{\text{A}}$. This direct charge transfer model then gives a linear relationship between $\Delta_{\text{DA}}$ and the charge population on the molecule Lee and Rick (2011). However, as shown in Ref. 7 and in Fig. S10 for our liquid water trajectories, while the presence of a hydrogen bond shifts the charge population on a water molecule, that alone is not sufficient to explain the CT, as seen from the width and overlap of the distributions in Fig. S10. Indeed, using this criterion there would be essentially no change expected in the CT upon including NQEs because the number of hydrogen bonds is only minimally perturbed. This suggests that the CT is determined by other parameters not captured by the binary criterion of hydrogen bond existence. Motivated by the strong dependence of charge transfer on the proton sharing coordinate $\tau$ in our gas-phase results (Fig. 1 and Fig. S6) we show how incorporating this coordinate in a condensed-phase model of CT yields much improved predictions.

To parametrize such a model for CT we begin by extracting $20\,000$ dimers from our quantum simulations, which were chosen as they sample a wider range of CT values (see top panel of Fig. 4). We then calculated the CT for each of these dimers separately to quantify the CT between that pair of water molecules. The resulting CT values obtained are plotted against $\tau$ in the top panel of Fig. 5. As for the dimer configurations obtained from the gas phase simulations, these show a strong correlation of the CT with $\tau$. On the other hand, the CT shows no dependence on the coordinate $\sigma$ orthogonal to $\tau$, as shown in the bottom panel of Fig. 5, suggesting that including this coordinate in the CT model would offer little benefit. This relationship can be used to determine the expected amount of CT given a value of $\tau$. To fit the dependence we used the following simple functional form that uses only three parameters

[TABLE]

where the values of the parameters $a$, $b$ and $c$ listed in Tab. 1 were determined by a least squares fit to the data. When this function is used to predict the CT of the extracted dimer set, a correlation of the predicted charge with the charge from the population analysis of 96.4 % is observed. This high correlation coefficient shows that the other degrees of freedom account only for under 4 % of the CT and that a simple functional form is thus able to describe the dependence appropriately.

With this model, the total charge on each molecule in the bulk liquid can be predicted by identifying the hydrogen bond donors and acceptors of each molecule and their $\tau$ values. One can then calculate the charge $q_{i}$ of water molecule $i$ by summing the incoming and outgoing CT contributions,

[TABLE]

where $\Delta q(\tau)$ is obtained from Eq. 3. Fig. 6 shows the predictions from this simple model compared to the molecular charges obtained from a full electronic structure calculation for the quantum liquid water ensemble. Excellent agreement is obtained for the NPA data shown, with the model giving slightly smaller molecular charges than the full electronic structure calculations. For NPA, the prediction works the best with a correlation coefficient of 95.2 %, while QTAIM and the Hirshfeld method yield correlation coefficients of 86.2 % and 83.6 %. The prediction for the classical ensemble shows comparably high correlation coefficients of 93 %, 80 % and 81 %, respectively. This shows that the prediction of molecular charges based solely on the proton sharing coordinate $\tau$ in hydrogen bonds is a good approximation for the total charge population in the liquid phase.

Since we have demonstrated that CT and the resulting molecular charge in liquid water are almost entirely dictated by a single geometric parameter of the hydrogen bond, we can now elucidate the difference between the quantum and classical ensembles. Fig. 7 shows the distributions of this geometric parameter, the proton sharing coordinate $\tau$ for liquid water. Although the classical and quantum ensembles match closely at more negative values of the coordinate, as one approaches zero, corresponding to a highly shared proton, the quantum distribution exhibits significantly higher probability density. This increase in proton sharing, observed in a number of previous studies  Ceriotti et al. (2013); Wang, Ceriotti, and Markland (2014); Ceriotti et al. (2016); Marsalek and Markland (2016, 2017), is vital in this case, as it solely determines the amount of CT. This is in stark contrast to many other properties which exhibit a cancellation between quantum effects along $\tau$ and the orthogonal coordinate $\sigma$ Habershon, Markland, and Manolopoulos (2009); McKenzie et al. (2014); Li, Walker, and Michaelides (2011).

IV Conclusions

In conclusion, by comparing quantum and classical simulations over a wide range of hydrogen bonded structures in both the gas and condensed phases we have shown that NQEs can make large qualitative and quantitative changes to the magnitude of CT. We have demonstrated that these results are robust with respect to the choice of electronic structure method and population analysis approach. In particular, the insensitivity of CT to the coordinates orthogonal to the hydrogen bond direction means that the higher CT due to increased proton sharing is not compensated by decreases due to other degrees of freedom. Due to the large amount of zero point energy in the $\tau$ coordinate in the quantum case, the additional thermal energy at the higher temperature only changes the position distribution by a small amount, in contrast to the classical case. This leads to a very small temperature dependence of the quantum CT distribution compared to that observed classically. As such, one can easily be misled by considering classical configurations or minimum energy structures alone when assessing the amount of CT.

Further, we were able to show that it is possible to accurately model the charge population on a water molecule in the condensed phase by considering only the sharing of protons with its immediate hydrogen bonding neighbors. This approach is able to explain 80–95 % of the charge population on a given molecule even though it lacks higher body effects, showing that the wider environment plays only a minor role. Indeed, recent work has shown that 2-body effects are sufficient to model the dipole moment surface of liquid water Liu, Wang, and Bowman (2015), further suggesting that the electrostatics of water are well captured at this level. Our analysis thus provides physical insight into the origin of charge transfer in hydrogen bonded systems as well as a simple basis for its accurate modelling.

Acknowledgements.

This material is based upon work supported by the National Science Foundation under Grant No. CHE-1652960. T.E.M also acknowledges support from a Cottrell Scholarship from the Research Corporation for Science Advancement. We would like to thank Stanford University and the Stanford Research Computing Center for providing computational resources and support that have contributed to these research results. C.S. gratefully acknowledges financial support from a scholarship from the German Academic Scholarship Foundation (Studienstiftung des Deutschen Volkes).