Competition between proton transfer and intermolecular Coulombic decay in water

Intermolecular Coulombic decay (ICD) is a ubiquitous relaxation channel of electronically excited states in weakly bound systems, ranging from dimers to liquids. As it is driven by electron correlation, it was assumed that it will dominate over more established energy loss mechanisms, for example fluorescence. Here, we use electron–electron coincidence spectroscopy to determine the efficiency of the ICD process after 2a1 ionization in water clusters. We show that this efficiency is surprisingly low for small water clusters and that it gradually increases to 40–50% for clusters with hundreds of water units. Ab initio molecular dynamics simulations reveal that proton transfer between neighboring water molecules proceeds on the same timescale as ICD and leads to a configuration in which the ICD channel is closed. This conclusion is further supported by experimental results from deuterated water. Combining experiment and theory, we infer an intrinsic ICD lifetime of 12–52 fs for small water clusters.


Supplementary Figures
Electron-electron coincidence data recorded from a jet of N = 76 water clusters at hν = 62 eV. Intensity is shown as events/pixel of 92 × 67 meV 2 width, as a function of kinetic energy E 1 of the faster and kinetic energy E 2 of the slower electron. A linear color scale is used. The E 1 interval marked by two vertical bars is assigned to events initiated by 2a 1 photoionization. The E 2 interval marked by the green box was used to obtain the coincident electron spectra, which were used to determine the ICD efficiency. The total number of coincident events acquired for each cluster size and photon energy was roughly between 4 × 10 5 and 2 × 10 6 , within an acquisition time of 140-900 s. Further discussion of this type of data in general can be found in Supplementary Ref. 1 Supplementary Figure 2: Outer valence spectrum used for estimating the degree of condensation Outer valence photoelectron spectrum of a jet of N = 76 water clusters at hν = 62 eV.
Components of a least squares fit to disentangle cluster and monomer intensity, and the residual of the fit (top trace) are also shown. The components 'g5' and 'g6' are assigned to the cluster HOMO component, components 'g3' and 'g4' are assigned to the monomer 1b 1 , 'g3' tentatively to the v = 2 vibrational component. A number of constraints were applied to produce consistent and well-defined fit results over the whole series of cluster sizes.  Supplementary Table 1). Small round symbols are from spectra recorded before and after an ICD efficiency measurement; in some cases the two data points are overlapping. An arithmetrical average of the two data points at equal N was used as an estimate for c. Error bars show the standard deviation. Additionally, larger round symbols show the result of neglecting the small auxiliary peak 'g3' when calculating c, obviously the effect is small. See Supplementary Methods for details. Supplementary Figure 4: Peak-to-background separation in the experimental data used to determine the efficiency of ICD Coincident and non-coincident electron spectra recorded after photoionization of N = 76 water clusters at hν = 62 eV: Panels (a) and (c) show p(E ph ), the intensity of photoelectrons detected in the region of 2a 1 binding energies; panels (b) and (d) show the integral of P (E ph , E ICD ) over an interval of [0, 5] eV for E ICD (region marked by a green rectangle in Supplementary Fig. 1). In other words, the latter panels show the intensity of photoelectrons detected in coincidence with an ICD electron. Experimental data in the left and right hand side panels are identical, but different choices for peak-to-background separation are indicated.

Supplementary Tables
Supplementary Table 1: Parameters of cluster source and calculated mean cluster size N Expansion conditions for cluster production, and calculated mean cluster size N . All experiments used a conical nozzle with 80 µm smallest diameter and a 15 • half opening angle. The stagnation pressure p was derived from the vapour pressure of water at the reservoir temperature T r . The nozzle temperature T n was set independently from T r . The H 2 O/D 2 O comparison was measured in a separate beamtime. See Supplementary Methods for a detailed description.   Supplementary Note 4. Potential energy curves in polarizable continuum model To estimate the energetics of the ICD process in liquid water, we performed a scan along the proton transfer coordinate using water dimer embedded in a dielectric continuum as a model. The technical details of these calculations are described below. The results are shown in Supplementary Fig. 5.
The curves for the both initial and final ICD states are qualitatively similar to the gas phase, but both curves are shifted to lower energies as the solvent stabilizes the extra positive charge. This stabilization is more pronounced for the doubly ionized final state, resulting in the ICD electrons with higher kinetic energies. In fact, our results suggest that the ICD channel never closes in liquid water. Note that this is a rather crude model and experimental verification is needed for this conjecture. Nevertheless, the widening of the gap between initial and final state is clearly apparent in the measured data as shown in Fig. 6 of the main article.
Let us now comment on the technical aspects of these calculations. Here, we focused only on the initial singly ionized 2a −1 1 state and the lowest triplet doubly ionized state, which can both be reasonably described by single reference methods as shown in Table 2. We used the PBE0 functional with the 6-311++G** basis set, but the gas phase energies were shifted to match the reference MP4 values. The same approach was used for Fig. 3 in the article.
There are two aspects of these calculations that deserve further discussion: 1. the use of single reference DFT method to obtain the energies of highly excited single ionized state and 2. the choice of the polarizable continuum model.
As explained in the Methods section, highly excited states can be obtained with the help of the Maximum Overlap Method (MOM) 3,4 . However, straightforward use of the MOM leads to energy discontinuities along the proton transfer coordinate. We have also tested a recently published Improved MOM (IMOM) method 5 , but with no improvement for this problem. In the end, we were able to make use of the fact that the inner valence orbital of the donor water molecule is energetically well separated from the other orbitals. Therefore, instead of using the standard MOM criterion, we simply populated the desired orbitals based on the energy criterion. We implemented this approach in the development version of the TeraChem package for both MOM and IMOM methods 6,7 .
To model the solvent effects, we used the non-equilibrium formulation of the C-PCM model as implemented in the TeraChem package 8 . Although, the non-equilibrium solvation takes into account that the ionization is a vertical process, the overall ICD process that we are trying to describe here is not vertical and the surrounding water molecules have certain (albeit short) time to relax around the ionized molecule. The dynamical aspect of this process is hard to describe by the PCM methodology, but these effects would not changes the qualitative conclusions.

Supplementary Methods
Cluster production and size estimation: In our work, water clusters were formed via a supersonic expansion of pure water vapour into vacuum (no seeding gas was used). The jet produced such contains clusters with a broad distribution of sizes N , including monomers.
The mean of this distribution N can be estimated from the parameters of the expansion.
For this estimate, we use an empirical relation due to Bobbert et al. 9 Details, and a critical account of this approach to cluster size estimation, have been given by some of the authors Expression for the ICD efficiency: In the following, we give some details on the derivation of the expression that connects the experimental data with the ICD efficiency.
Following that, we analyze the errors assigned to its individual factors.
The process we consider starts by photoionization, creating a photoelectron of kinetic energy E ph . After that, the target is left in an excited state and emits another electron (the ICD electron) with kinetic energy E ICD . By tuning the photon energy we can always achieve E ph = E ICD , i.e. photoelectron and secondary electron can be distinguished experimentally.
In an ideal experiment, the branching ratio α ICD of the autoionization process (α ICD ∈ [0, 1]) can then be derived from experimentally measurable data as: where P is the rate for detection of electron pairs with energies E ph and E ICD , and p the rate for detection of photoelectrons with energy E ph (with or without a subsequent secondary electron). Rigorously, both photoelectrons and autoionization electrons, are emitted within some intervals of kinetic energy. In the following, E ph and E ICD rather designate the central values of the respective intervals E ph , E ICD , and P (E ph , E ICD ) ≡ P (E ph , E ICD ), p(E ph ) ≡ p(E ph ). In our experiments, we always had E ph > E ICD , therefore we designate the photoand ICD electron the 'first' and 'second' electron to arrive at the detector, with kinetic energies E 1 and E 2 .
In our actual experiment, several corrections of the measured data must be taken into account, leading to a more complicated version of supplementary equation (1). These are: 1. the less-than-unity detection efficiency of the spectrometer, 2. the background of uncondensed monomers in the cluster jet, 3. losses by intracluster inelastic electron scattering, 4. a difference in the outer valence photoionization cross section of monomers and clusters.
The first two points lead to a correction factor of cγ(E ICD ) −1 in the expression for α ICD .
Here, c is the degree of condensation (ratio of molecules being a part of a cluster to total number of molecules in the interaction region) and γ(E ICD ) the detection efficiency of the spectrometer for ICD electrons. The detection efficiency for photoelectrons γ(E ph ) cancels from the expression.
Further complications occur due to the overlap of the inner valence monomer and cluster photoelectron lines for molecular clusters. We consider the observed rate of photoelectrons p(E ph ) = γ(E ph )r ph . It contains a rate p(E ph ) >p(E ph ) = γ(E ph )r ph c from clusters, which is not directly observable due to the overlap. We would also like to take into account the effect of inelastic intracluster photoelectron scattering, which diminishesp(E ph ). (Densities in our jet are such that inelastic scattering at other molecules or clusters can be neglected.) We describe the inelastic losses by a factor f := 1−(lost fraction). With that, the observed rate of photoelectrons becomes because the correction by f must only be applied to the cluster part of the intensity. We assume that no correction for intercluster scattering is needed for the ICD electrons, because their energy is low and low cross sections for inelastic electron scattering in amorphous ice were found at these energies 11,12 . For the ratio with the coincident events we now have The degree of condensation c in our experiments is determined by electron spectroscopy as well, and results from a comparison of the features from condensed molecules vs. monomers in a region of the spectrum were they can be distinguished. Practically, the HOMO levels of molecules and clusters were compared (see below). We use: where c and m refer to clusters and monomers, respectively. This implicitly supposes that the measured photoelectron intensities are proportional to the respective numbers of condensed and uncondensed molecules in the interaction region, which is the case if their outer valence photoionization cross sections are equal-an assumption often made in PE studies of clusters, which, however, has rarely been checked rigorously 13 . We would like to allow for a correction factor taking into account possible differences between the water HOMO cross section for the monomer (σ m , the 1b 1 cross section) and for a single water molecule within a cluster, σ c . Defining x := σ m /σ c we can write a cross-section corrected degree of condensation as The expected case is σ c < σ m , because the orbitals are more diffuse in a cluster, and because in larger clusters molecules at the side of the cluster that is exposed to the photon beam may shadow others. In this case x > 1, and the actual degree of condensationc is larger than the calculation of c given in supplementary equation (4). The connection to the measured c is given byc Usingc instead of c in supplementary equation (3) after some algebra yields our final result We would like to note that most likely f is a number slightly smaller than one, while x should be slightly larger than one, as argued above. In this case, the approximate neglect of the last factor in supplementary equation (7) seems justified. This approximation was used in earlier work on rare gas clusters by some of the authors 14,15 and will again be used here.
Additional justification arising from details of our experimental procedures is given below.
In the following, we describe how the individual factors making up supplementary equation (7) were determined.
Spectrometer efficiency: Experiments were carried out with a magnetic bottle time-offlight spectrometer for electrons 16 . The transmission function γ(E) of this instrument is the product of accepted solid angle (as a fraction of 4π) times probability to register a charged particle on the detector, and may weakly depend on kinetic energy E. In supplementary equation (7), we need the value of γ for the kinetic energy interval in which ICD may occur. We determined this property of the spectrometer from measurements of Xe NOO photoelectron-Auger electron coincidences, as described in Supplementary ref. 16. Values for kinetic energies of approx. 0.5 to 5.5 eV were found inbetween 0.5 and 0.6, with no significant energy dependence. In the following, we use their average of γ(E ICD ) = 0.58 (4).
For the data in the comparison of normal to deuterated water, a value of γ(E ICD ) = 0.39 was used. This is however less well bounded than the one for the pure water data. The significance of the normal to deuterated water comparison however is not influenced by that.
Degree of condensation: In supplementary equation (7), a correction of the degree of condensation c in the cluster jet is applied, because the 2a 1 photolines from monomers and clusters cannot be separated. To experimentally determine c, we measured the outer valence photoelectron spectrum of our cluster jet, as it was shown that cluster and monomer signals for the HOMO can be distinguished. 10 We determined c from the area ratio of the two features. In monomers, this line results from ionization of the lone-pair orbital and is dominated by a sharp v = 0 vibrational component. In clusters, it develops a broad peak at a lower binding energy than the monomer. Outer valence spectra were measured at hν = 62 eV before and after the spectra used for ICD efficiency determination. A retardation voltage of V ret = −15 V was used to improve the energy resolution. Nevertheless, peak fitting had to be used to disentangle the cluster from the monomer features ( Supplementary Fig. 2). Results for the degree of condensation c (supplementary equation (4)) are shown in Supplementary   Fig. 3. Here, symbols of different color are shown for an analysis using either only the fit component 'g4', or both 'g4' and the small auxiliary component 'g3' as the monomer area.
Differences between these two approaches and differences between spectra acquired before and after the ICD efficiency measurement amount to a compound estimated error for c of ±0.05 (standard deviation).
The comparison of normal to deuterated water was measured in a separate beamtime.
Here, a similar approach was used to yield c = 0.75 (3)  give an order of magnitude, f ∼ 0.8. Here we would like to take a different approach, however. Recalling that xσ c = σ m , we observe that inelastic losses will also influence the area of cluster outer valence photoelectron lines used in the determination of the degree of condensation. Therefore, if we ignore the slight dependence of inelastic scattering cross section on kinetic energy, 1/f will be one of the components making up x, the other being 'intrinsic' changes in the photoionization cross section due to changes of the orbital shape upon aggregation. The change in orbital shape might be substantial 10,18 , but whether that affects the photoionization cross section is currently unknown. In the absence of further information, we think it is a fair approximation to propose cancellation of the product f x, which will yield unity for the whole round bracketed factor in supplementary equation (7).
Modelling the possible influence of this factor when relaxing the former approximation, we find that it mostly plays a role when the degree of condensation is low. For c = 0.8, a ±20% intrinsic difference in the photoionization cross sections will have less than 5% influence on the result for α ICD .
Extraction of spectral intensities from the experimental spectra: Finally, values of the coincident signal intensity P and non-coincident signal intensity p in supplementary equation (7) were determined from the number of events with electrons registered in some interval of kinetic energies. Typical coincident and non-coincident electron spectra are shown in Supplementary Fig. 4, together with the two background models we have applied.