Yiwei Zhang

Guillaume Stirnemann

James T Hynes

Damien Laage

Water dynamics at electrified graphene interfaces: a jump model perspective

The reorientation dynamics of water at electrified graphene interfaces was recently shown [

Introduction

Interfaces between an aqueous electrolyte and an electrified material play a key role in fundamental electrochemical processes involved in, e.g., energy storage in batteries 1,2 and photocatalytic energy conversion 3,4 . The great importance of these interfaces has recently prompted a considerable effort, both experimental and theoretical (see, e.g. Refs 5-13 and references therein), to gain a better understanding of their unique properties.

While most studies so far have focused on structural aspects, the much less frequently examined dynamics of interfacial water molecules is essential for the transport of redox species to and away from these interfaces, and for the kinetics of charge transfer reactions taking place next to the electrode 14 . Recent simulations by some of us 15 showed that the reorientational and translational dynamics of water molecules in contact with an electrified graphene interface are affected in dramatically different ways by electrode potentials with opposite signs: while water reorientation dynamics is increasingly slowed down as the magnitude of an applied positive electrode potential increases, it changes in a non-monotonic fashion with the potential at negatively charged electrodes. In the bulk, water reorientation essentially proceeds via exchanges of hydrogen-bond (H-bond) partners, which occur through large amplitude jumps 16,17 . Some of us have previously shown that in the vicinity of solutes or interfaces, jumps generally remain the main reorientation pathway, and the impact of solutes and interfaces on water reorientation dynamics mostly arises from their effect on the jump dynamics 17 . In an interfacial context, the Extended Jump Model was shown to provide a molecular description of water dynamics at a broad range of interfaces, ranging from hydrophobic 18 and hydrophilic 19 mineral interfaces, to biomolecular [20][21][22] hydration layers.

Using the extended jump framework, we elucidate here the molecular origin of the radically different water dynamical responses at the positive and negative electrodes. We show that the concepts that we have previously developed in other contexts still apply at charged interfaces, and are able to quantitatively predict the surprising non-monotonous behavior at negatively charged interfaces, as well as the effect of positively charged electrodes.

In the remainder of this contribution, we first describe our simulation methodology and characterize the geometric arrangement of interfacial water molecules. We then provide an analysis of the water H-bond jump dynamics at negative and positive electrode interfaces, and show that the extended jump model can successfully describe the electrode potential's effect on water reorientation dynamics. Next we characterize the correlations between the charges on the graphene electrode atoms, and show that interfacial water dynamics does not exhibit significant heterogeneities. Finally, we offer some concluding remarks.

Methodology

Our simulation methodology follows the procedure described in our recent work 15 . We use classical molecular dynamics to simulate a slab of liquid water between two parallel single sheet graphene plates. The two rigid graphene sheets are constructed from an ideal hexagonal lattice with a 1.42 Å distance between neighboring carbon atoms; each plate contains 240 atoms and the two sheets are 6 nm apart along the z axis. The box contains 1158 water molecules; this number was determined iteratively via an insertion/deletion procedure such that the density in the middle of the box equals the experimental liquid water density. Water is modelled by the widely employed SPC/E potential, which has been shown 23 to provide an excellent description of water dynamics at ambient temperature. Graphene carbon Lennard-Jones parameters are taken from ref. 24 where they were optimized to reproduce DFT-based molecular dynamics simulation results. Long-range electrostatic interactions are described via a two-dimensional Ewald summation method 25 with 10 4 force accuracy.

The graphene carbon atomic charges are not fixed; they are allowed to fluctuate in order to maintain the electrostatic potential at each electrode atom fixed. This approach is based on pioneering work by Siepmann and Sprik 26 , who introduced an extended Lagrangian method in which the electrode atomic charges are treated as dynamical variables and subjected to the constraint that the average electrode potential should be constant. Our simulations use the constant potential method subsequently developed by Madden and coworkers 27 , where the constraint of constant electrode potential is enforced on each electrode atom at every simulation step. Graphene carbon charges are described by Gaussian distributions whose amplitudes are determined by solving a set of linear equations at every step, using the constant potential implementation by Laird and coworkers 28 . The Gaussian inverse width is fixed to be 1.979 Å 1 , following ref. 27 (narrower distributions were suggested 29 for carbon nanotubes but we verified that this did not affect our water dynamics results). The constant potential method and related techniques have been successfully used for a broad range of electrode-electrolyte interfaces (see, e.g., refs. 27,28,30-34). (Other constant electrode potential simulation techniques including, e.g., constant Fermi-level simulations 35,36 , are being developed but require much more computationally demanding explicit electronic structure calculations; recently developed 37 classical polarizable potentials for metals based on the embedded atom method can also be combined with the Siepmann and Sprik extended Lagrangian and would be a promising route to describe constant potential electrodes.) These classical simulations were shown 15 to provide descriptions of interfacial water structure 24 and water sum-frequency generation spectra 38 in very good agreement with DFT-based simulations, while allowing the propagation of long trajectories required for our dynamical analysis. Polarizability is most important for graphene 39 , and is described here by the constant potential approach. Concerning water, the SPC/E model is a fixed-charge non-polarizable model, but we will show in Section 5.4 that, for the relatively weak electrode fields and small electrode charges studied here, this model's approximations are valid. There is a further point noteworthy in connection with our treatment's inclusion of the electrode polarizability. For water dynamics in aqueous salt solutions, where the ions provide an additional electric field source, an effective description of electronic polarizability via scaling of the solute charge without modifying the SPC/E water model was shown to provide an excellent description 40 . In the present case of charged electrodes, such an approach would simply shift the potential-(and thus charge-) dependence of water dynamics, but -importantly -would not change the overall behavior.

Molecular dynamics simulations are run with the LAMMPS program 41 , at a series of electrode potentials ranging from 0 to ±2.5 V. This range of potentials goes beyond the stability domain of water, but recent experiments 42 on water-graphene interfaces have shown that ±1.5 V potentials are accessible, owing to the important overpotential 43 , and the largest ±2.5 V potentials in our simulations are only used here to highlight the trends. To provide a connection between the fixed electrode potential used in our simulations and experimentally accessible potential measurements, we determined the potential of zero electrode charge (pzc). As shown in the ESI †, the pzc value from our simulations is '-0.02 V. In the following, the reported potential values can thus be approximately considered as the potential differences with respect to the pzc.

Each system is first equilibrated for 1 ns with a thermostat at 298 K. Trajectories are then propagated with a 1 fs time step in the microcanonical ensemble for 4 ns.

Interfacial water reorientation dynamics

To characterize the effect of graphene electrodes on water dynamics, we focus on the reorientation time of water OH groups at increasing distances from the electrode. We calculate the reorientation time correlation function

C OH 2 (t) = hP 2 [u OH (t) • u OH (0)]i (1)
where u OH (t) is the water OH group direction, P 2 is the second order Legendre polynomial and the h...i ensemble average is performed over all water OH groups in a given region at t = 0. The t OH 2 reorientation time can be measured in the bulk via a range of experimental techniques including NMR and ultrafast vibrational spectroscopies 17,44 ; measurements of t OH 2 at interfaces are more challenging but recent developments in surface-specific ultrafast vibrational spectroscopy showed that it is accessible at aqueous interfaces 45 . Here t OH 2 is determined as the long-time decay of this reorientation time correlation function, determined via an exponential fit on the 2-10 ps interval.

Our results in Figure 1a show that, as recently reported in ref. 15, the water reorientation times in the first layer in direct contact with the graphene sheets dramatically change with the electrode potential. Positive electrode potentials cause a marked slowdown in water dynamics, while increasingly negative potentials first induce a surprising acceleration of water dynamics, followed by a retardation for very negative potentials.

But this strong effect is short-ranged, since the reorientation times of water molecules initially in the second layer are very close to the values determined in the middle of the box (Figure 1b): only a small residual retardation is found in the second layer next to positively charged electrodes, and the dynamics is already essentially bulk-like within the negative electrodes' second hydration layer. In contrast, much longer-ranged effects would be anticipated from traditional continuum electrostatic models of the water dipole -electrode-induced field interaction. (We note that this difficulty is apparently circumvented -in heavily screened ionic liquids -by new mean-field theories 46,47 accounting for the water dipoles' rearrangement). However, a common -and critical -limitation of these theories is the point dipole treatment of water molecules. This description is not correct at interfaces where the water-electrode distance is of the same order as the water molecular dimensions; a key direct consequence relevant here is that these simple models would not account for the asymmetry between the behaviors at the positive and negative electrodes. A molecular picture explaining both the magnitude and the range of the electrodes' impact on water dynamics is thus needed.

Arrangement of water molecules at electrified graphene interface

In order to determine the molecular origin of the just-described dramatic effects caused by electrified interfaces on neighboring water dynamics, we first characterize the arrangement of water molecules at these interfaces.

Figure 2 addresses the water molecules' orientations by showing the probability distributions of the angle q between each water OH group and the interface normal vector for increasing distance Dz to the electrode. Next to the 0 V graphene interface, three main orientations are identified within the first hydration layer (Figure 2a). The first of these corresponds to a small number of OH groups ('5% fraction) pointing towards the electrode (q .40 ). Since the graphene interface is approximately apolar, these OH groups are not engaged in any strong interaction; indeed they are analogous to the dangling OH groups found at hydrophobic interfaces, including, e.g., the air-water interface 45,[49][50][51] . This population of OH groups pointing towards the electrode surface is noted S in the following. The second main orientation arises from OH groups tangent to the interface, i.e. q '90 ; these OH groups are typically donating a H-bond to other interfacial water molecules, and is noted W T . Finally, a third OH group orientation is found to be due to first-layer OH groups pointing towards the solution, and forming a H-bond with second-layer (almost bulk-like) water molecules; it is noted W B in the following. Figure 2 shows that in the second hydration layer, preferred OH orientations are much less pronounced than in the first layer. The water structure and these main H-bond orientations are schematically represented in Figure 3a.

We pause to note that the W B population implies the presence of H-bonds between the first and second hydration layers next to graphene; this therefore markedly contrasts with the situation found at Pt 100 and Pt 111 interfaces with water 31 , where simulations suggested that the adlayer interacts strongly and exclusively with the metal, but not with the second layer, thus creating a hydrophobic interface for the subsequent water layers.

The application of a negative potential on the graphene electrode induces a pronounced restructuring of the electrode first hydration layer. First, molecules in the first hydration layer move closer to the interface 15 . But most importantly, Figure 2b shows that at -1 V, the S population is strongly enhanced, while instead the W T and W B populations exhibit a small decrease. The attractive electrostatic interaction between the negative graphene charges and the positive charge on water hydrogen atoms favors the reorientation of interfacial water molecules which point one of their OH groups toward the electrode, while same molecule's other OH group is approximately tangent to the interface (Figure 3b). Such an arrangement has been found in a series of prior water-electrified interface simulations, both with classical forcefield and DFT-based molecular dynamics 30,52 ; it is also consistent with the configuration of water molecules next to anions in aqueous solutions obtained from simulations 53 and neutron-scattering experiments 54 . This configuration contrasts with a suggested 55 different arrangement of water molecules at negatively charged interfaces, where water dipoles would align with the surfaceinduced electric field, and both OH groups would point toward the interface. In our simulations, such configurations were not found (no peaks are found around q ' 50 in figure 2), even for unrealistically large negative electrode potentials (< 5 V).

In contrast to the significant rearrangement induced by negative potentials, our simulations show that positively charged electrodes less dramatically affect the structure of the interfacial layer 1-11 | 3 (Figure 3c). Water oxygen atoms are attracted by the electrode positive charges, resulting in shorter distances to the surface. Figure 2c shows that at +1 V, the W T and W B populations are slightly enhanced with respect to the 0 V case, and that the small S population has vanished, due to the unfavorable electrostatic water OH group -electrode positive charges interaction.

Water jump dynamics

Some of us have shown 16,17 that water reorientation essentially proceeds through sudden large-amplitude jumps when a water OH group switches H-bond acceptors. In a series of prior works, [17][18][19]21 we demonstrated that the difference of the bulk water reorientation dynamics and those at an interface (with an extended surface or with a small solute) is mostly due to the interface's impact on the rate constant of these H-bond jumps. Accordingly, we will employ this jump perspective to elucidate how electrified interfaces affect water dynamics. One key novel fea- ture here -and a main difference -compared to previously investigated systems lies in the fluctuating character of the electrode charges; prior systems included solutes and surfaces with fixed charges.

V = 0 V S W B W T V = -1 V S W B W T V = +1 V W B W T

Extended jump model

We first quickly summarize the key features of the extended jump model 16,17 . The main reorientation pathway for water molecules in the liquid is via large-amplitude angular jumps when OH groups exchange H-bond acceptors 16,17 (Figure 4). These hydrogen-bond acceptor exchanges of average amplitude Dq can be fruitfully viewed as a chemical reaction whose inverse rate constant is the jump time t jump . The water reorientation time associated with the reorientation time correlation function eqn 1 was shown 16,17 to be well described by an analytic extended jump model, which combines the jump reorientation with the slower diffusive tumbling time t frame of the coordinate frame for a water OH engaged in an intact H-bond between jump events,

1 t reor = 1 t jump  1 1 5 sin (5Dq /2) sin (Dq /2) + 1 t frame (2)
Studies of water reorientation next to a wide range of solutes and interfaces have shown 16,17 that the changes in water reorientation dynamics mainly result from a change in t jump . Jump times are determined 16,17 from the cross time correlation function between the initial (I) and final (F) states

1 hp I (0)p F (t)i = exp t/t jump (3)
where p I,F (t) is 1 if the system is in state I (resp. F) and 0 otherwise. States I and F are defined within the Stable States Picture 56 with strict geometric criteria to remove contributions from fast, unproductive, barrier crossing, and absorbing boundary conditions in the product state ensure that the forward rate constant is calculated (In the following, S states are defined by Dz O <4.5 Å and q < 30 , W T,B by Dz O <4.5 Å and the presence of a water oxygen H-bond acceptor O a next to the considered OH group such that OO a <3.1 Å , HO a <2.1 Å and \ HOO a < 20 ). The mechanisms of two main types of jumps (W T,B ! W T,B and S ! W T are schematically shown in Figure 4, highlighting the important molecular motions that will be extensively discussed in the following.

W W S W b a

Fig. 4 Schematic illustrations of the H-bond jump exchange mechanism for (a) W T,B ! W T,B and (b) S ! W T jump events. Importantly, this process can be fruitfully viewed as a chemical reaction that proceeds from an initial reactant state (left) to a final product state (right) through a transition state (middle). Key molecular motions are indicated by black arrows. 16,17 The slow motion frame reorientation time t frame is determined by calculating the reorientation time correlation function eqn 1 during time intervals when the H-bond is intact.

Negative graphene interface

The striking non-monotonic change in water reorientation time with increasingly negative electrode potential reported in Figure 1 is reminiscent of the behavior observed 19,57 for the dynamics of water next to a silica surface whose hydrophilicity was changed by progressively switching the surface's atomic charges from the apolar limit of zero up to their full values. For that system, the non-monotonic character of the reorientation time was shown 19 to arise from a change in the OH group populations with different H-bond acceptors combined with the surface charge effect on the jump times of these different OH group populations. We will follow an approach similar to that of ref 19 here.

We first examine how the increasingly negative electrode po-tential affects the fractions of interfacial groups with two different H-bond types: first the OH groups pointing towards the electrode (S state, see Figure 2) and second the interfacial OH groups donating H-bonds to other water molecules; the latter includes the W T and W B populations which in a first approximation are now described as a single entity W. Determining the two state populations from geometric probability distributions is often ambiguous because of the necessarily arbitrary choice of state boundaries; we therefore prefer to determine the state populations via their equilibrium constant obtained from the less ambiguous jump times eqn 3 between the two states. Figure 5a shows that increasingly negative potentials lead to an increase in the S fraction, which is stabilized by its constituent water molecules' favorable interaction with the growing negative charges on graphene. We now consider the jump dynamics of the above two OH group populations at the graphene interface. Due to their different orientations, these two OH group subsets are affected in very different ways by the negative graphene charges. The W population is either tangent to the interface or pointing away from it, and it is engaged in a H-bond with a water molecule. Figure 5b shows that the H-bond jump time of OH groups initially in the W state is very little sensitive to the electrode potential: the jump time is slightly retarded with respect to the 3.2 ps value determined in the middle of the cell. This slowdown was shown in prior work 19,58 to be caused by an entropic excluded volume effect induced by the interface which hinders the approach of a new H-bond acceptor and thus slows the H-bond jumps (as shown in ESI †, while the W jump time is almost unaffected by the electrode potential, the slower frame reorientation does display a small change).

In contrast, the jump time of OH groups initially in the S state dramatically increases with increasingly negative electrode potential. This is explained by the transition-state H-bond strength effect 19,59 which affects H-bond jumps starting from H-bond acceptors of varying strength. In particular, the jump mechanism involves the concerted breaking and making of H-bonds (black arrows in Figure 4), and reaching the jump transition state thus requires weakening of the initial H-bond. In the present case, the increasingly negative electrode potentials lead to larger negative charges on the graphene carbon atoms, and thus stronger interactions with water OH groups in the S state. This thereby increases the free energy cost to stretch the H-bond with initial acceptor and reach the jump transition state, thus slowing down the jumps.

Both effects can be combined to reconstruct the overall extended jump model reorientation time eqn 2 (see details in ESI †). Figure 5c shows that the model provides an excellent description of the reorientation time non-monotonic behavior with the electrode potential, and thus offers a molecular interpretation of this behavior, as follows. We begin with the start from the apolar electrode at 0 V, with the increasingly negative potential. This increase first leads to an increase in the S fraction, this fraction's whose reorientation is initially much faster than in the W state, because in the latter the initial H-bond with water has to be stretched while in the former there is no initial H-bond. The net effect is that the average interfacial water reorientation speeds

1-11 | 5 -2.5 -2 -1.5 -1 -0.5 0
Electrode potential (V) up (decreasing time) when the negative electrode potential increases in magnitude. But this does not continue: when the negative potential magnitude increases further, the electrode's interaction with S state OH groups is strengthened, and the jump time increases and reorientation from the S state slows down. This ultimately leads to a slowdown in the average interfacial reorientation.

The combined increase in the S population and the progressive slowdown of S-state reorientation dynamics cause an extremum in the average interfacial reorientation time. Our simulations suggest that this extremum is located in the vicinity of -1.5 V. While the water reorientation time at this point is almost equal to its bulk value, the dynamics of interfacial water is nonetheless not bulk-like: it is much more anisotropic than in the bulk, with the average value concealing the presence of two subpopulations with very different dynamics. The extremum is determined by the competition between the population transfer between the W and S states, and the pronounced change in the S-state dynamics (together with a limited change in the W-state dynamics). The former's reorientation time is essentially determined by an entropic, temperature-independent, excluded-volume factor induced by the presence of the interface hindering the approach of new acceptors. The latter's reorientation time changes with the strength of the mostly enthalpic interaction with the electrode charges. This competition suggests that the location and amplitude of this extremum should be sensitive to temperature, and to the presence of ions. In addition, the different values of the surface charge densities at the extrema in the electrified graphene and scaled-charge silica cases (approximately -0.3 e/nm 2 for graphene and -0.9 e/nm 2 for the silica silanols 19 ) suggest that the surface sites configuration also plays a role.

Positive graphene interface

We now turn to the change in interfacial water reorientation dynamics at the positive graphene electrode. The simulation results in Figure 1a show a pronounced slowdown when the electrode potential increases. In contrast with the situation described at the negative interface, the slowdown at the positive electrode cannot be explained by a transfer between H-bond populations; Figure 2c shows that -in contrast with the large change in populations of OH groups with distinct orientations observed at the negative interface -the OH group populations do not significantly change in the positive electrode's neighborhood, except for the disappearance of the S state which in any case represented only a very small fraction at 0 V. As we will now see, the explanation of the water reorientational behavior with increasingly positive electrode potential is to be found in a marked slowdown in the jump dynamics, caused by interfacial water structural changes. To begin, we see that Figure 6a shows that the radial distribution function of water oxygen atoms around central water oxygen atoms lying within the first hydration layer exhibits increasingly pronounced peaks and minima for increasing positive potentials. While the locations of these extrema are not affected, their amplitudes are, as can also be seen in Figure 2b. Remarkably, despite its apparently quite modest character, this enhanced structuring has a strong impact on the kinetics of the water structural rearrangements, explained via the jump model. The next step in the argument is to recall that the free energy barrier of water hydrogen-bond jump exchanges results from both elongation of the initial H-bond and approach of the future Hbond partner 48,60 . And some of us recently showed -in connection with a quantum nuclear treatment of the jumps 61 -that the jump model can quantitatively connect changes in the oxygenoxygen radial distribution function -and thus in the free energy costs to stretch and elongate O-O distances -to changes in the jump free energy barrier. In order to calculate the impact of Fig-

0 0.5 1 1.5 2 2.5
Potential (V) Radial distribution function ure 2b's structural enhancements on the interfacial water's reorientation for the positive electrode, we will need to first summarize the key model features. The jump time is the inverse of the H-bond jump exchange rate constant, determined via Transition State Theory by the jump free energy barrier; the latter can be well approximated as the sum of two independent contributions, from the initial H-bond elongation and from the compression of the distance to the future H-bond partner

0 V 0.25 V 0.5 V 1 V 1.5 V 2.
DG ‡ ' DG ‡ elong + DG ‡ compr (4)
whose physical aspects are the following. Elongation takes the initial acceptor from its stable H-bonded distance within the first shell to the transition-state geometry where it lies half-way between first and second shells (cf Figure 6a); in the compression, the future partner has to go at the same time from the second shell to that same transition state distance. This identification allows the connection to the radial distribution in Figure 6a via eqn 5: each free energy barrier component can be approximated 61 using the potential of mean force W (r) corresponding to the O-O radial distribution function

W (r) = k B T ln [g OO (r)], DG ‡ ' [W (r min ) W (r max1 )] + [W (r min ) W (r max2 )] (5)
where r min,max1,max2 are respectively the OO distances at the radial distribution function's first minimum, first maximum, and second maximum. Finally, the ratio between the water jump time next to an electrode at potential F and its value at 0 V is thus given by t jump (F) t jump (0) = g (r max1 ; F) g (r max2 ; F) g (r min ;0) 2 g (r max1 ;0) g (r max2 ;0) g (r min ; F) 2

Figure 6b shows that this model eqn 6 provides an excellent description of the changes in water jump times induced by increasingly positive electrode potentials. The jump model again offers a molecular interpretation of the electrodes' main impact on water reorientation dynamics.

It is noteworthy that this enhanced structuring effect was not observed at the negatively charged electrode. Our calculated O-O radial distribution functions at this interface actually show a slight de-structuring with increasingly negative potentials (see ESI †), consistent with the change in density profiles reported in ref. 15. This probably arises from the reduction in water-water Hbonds due to the growing interaction with the negatively charged electrode.

Spatial range of electrode effects

We stress that the molecular insight provided by the jump model also explains why positively charged electrodes have a longer ranged effect on water reorientation times than negatively charged electrodes (Figure 1b). The enhanced structuring caused by the positive potential can propagate further to the second shell, while the H-bond population change induced by the negative potential involves the specific S population which is only present in the first layer, and this effect does not significantly affect the subsequent layers.

Here, we examine why the dynamics of water molecules beyond the two interfacial layers are not dramatically affected by the electrode potential (Figure 1b). We have calculated the mean electrical potential profile across the cell, by solving the Poisson equation with a finite-difference approach and boundary conditions fixed by the constant potential at the electrodes 62 . The potential profiles shown in Figure 7a display pronounced oscillations at the interfaces, due to the successive layers of water molecules, in agreement with the molecular density profiles 15 . The results reveal that the first two water layers at each interface screen the electrode potential very efficiently. This screening is remarkable at the positive electrode, where all water molecules have their negatively charged oxygen site facing the interface; consequently, the residual potential after the double layer is practically independent of the electrode potential (Figure 7a). At the negative electrode, the potential drop arises from the fraction of water OH groups pointing towards the electrode (S state). Because of the competition between interactions with the negative electrode and with neighboring water H-bond acceptor sites, only a fraction of interfacial groups are in the S state. The resulting screening is thus partial, as is manifest below -0.5 V in Figure 7a.

We now turn to the behavior of the electric field along the z axis at increasing distances from the electrodes, and we start by considering some useful magnitudes for typical fields in bulk water. In bulk liquid water in absence of any external field, typi-cal local fields induced by neighboring water molecules are on the order of 2 V/Å i.e. 20,000 MV/m 45 . When static external fields are applied, field strengths above 100 MV/m are typically considered to lead to dielectric saturation in liquid water 63 , and recent simulations 64,65 found that electronic polarization effects become important above 500 MV/m field strength. Our simulations show that in the vicinity of graphene electrodes, the fast potential spatial fluctuations (Figure 7a) lead to electric fields whose magnitude never exceeds 2 V/Å , thus remaining smaller than or comparable to typical bulk local fields. In the middle of the cell, the residual potential beyond the interfacial layer leads to a modest but constant electric field; Figure 7b shows that this static field's magnitude increases with the electrode potential but remains always much below typical local fields and -except for the largest ±2.5 V electrode potentials -below the dielectric saturation threshold. This residual static field therefore has a very moderate impact on the H-bond dynamics of water molecules far from the electrodes. However, despite this dynamical insensitivity and despite the comparatively small amplitude, this static field can be important in some contexts; for example such weak fields were found 15,66 to play a major role in sum-frequency generation spectra by breaking the local symmetry.

Correlations and heterogeneities

As pointed out in Section 5, an important difference between the electrified interfaces addressed here and prior studies of water jump dynamics next to a range of solutes and interfaces lies in the fluctuating character of the electrode charges. Ideally polarizable metallic (or semi-metallic for graphene) electrodes could support collective polarization effects, leading to correlated motions of interfacial water molecules and dynamical heterogeneities. Such correlated motions with dynamical heterogeneities relaxing on very slow ns timescales were recently found 31,67 in simulations of water-Pt interfaces. We therefore now examine the possibility of such collective effects at the graphene-water interface.

Electrode charges

We first study the distribution of electrode charges. As detailed below, our analysis reveals that neighboring electrode charges are correlated, and that the fluctuating charge with constant potential method used here is essential for the proper description of interfacial water dynamics at electrodes.

First, figure 8a shows that at ±1 V the total charge on each electrode essentially follows a Gaussian distribution, with the notable difference that probabilities to find very charged electrodes are slightly enhanced with respect to the normal distribution. The latter shows that atomic charges on the graphene carbon sites are not strictly independent, and suggests that the water arrangements which polarize the electrode induce correlations between neighboring sites. Second, we now examine the local charges on graphene carbon atoms. Within the constant potential method 27 used in our simulations, these charges are not fixed and are determined at every step such that the potential is constant on each carbon atom. Figure 8b shows that the local charge distributions are typically bimodal, with a smaller band at more negative charges, whose amplitude grows with increasingly negative potentials. These two bands arise from the two main orientations of interfacial water OH groups, either pointing towards another water molecule or pointing towards the electrode. The latter orientation corresponds to the S state in Figures 2 and5: the proximity of the positively charged water hydrogen atom polarizes the electrode negatively. The significant width of the local charge distributions and their bimodal character stress the importance of the constant potential description with respect to a fixed charge approach. The latter would not allow the proper stabilization of the S population, which would have critical consequences for interfacial water dynamics given the important role of this geometry described in Section 5.2. An approximate description of an electrode held at constant potential with a set of fixed charges would thus miss some key features; in addition to the limitations already identified in prior work 28,30,31,68 , our results show that this approximation would fail to provide a correct description of interfacial water dynamics.

The spatial correlations between graphene carbon charges in Figure 9a shows that neighboring charges are correlated over a distance of '3 Å. This is approximately the size of a water molecule, which -together with the feature that the correlation lengthscale does not exhibit any strong change with the electrode potential -suggests that correlations found in the total charge distribution arise from the simultaneous polarization of several carbon sites by the same water molecule. Fig. 8 Semilog plots of charge distributions for a) the total electrode charge at -1 V (green dashes) and +1 V (solid green), together with a Gaussian fit (blue), along the reduced d Q/s Q coordinate, where d Q = Q hQi is the instantaneous deviation from the average total charge and s Q = p hd Q 2 i is the distribution's second moment, and b) the atomic carbon charges at the positive and negative graphene interfaces for growing potentials, shown both as an effective local surface charge density (lower axis) and as an atomic charge on the carbon atoms (upper axis).

Figure 9b shows that the graphene carbon charges fluctuate over a fast sub-ps timescale, responding to the rapid rearrangements of the hydration layer structure. This relaxation time is not strongly affected by the electrode potential, and it is similar at the positive and negative electrodes, except for the sub-100 fs initial decay shown in Figure 9b inset. This initial decay arises from water librational motions, and the disparity in ultrafast decay timescales is caused by the different types of H-bonded populations: at the negative electrode, the graphene charges are determined by the S state librations, which take place on a comparably slower timescale because of the weak water-electrode interaction; in contrast, at the positive electrode, the dominant population is the W T configuration whose librations are faster due to the stronger H-bonds donated to water acceptors.

Absence of dynamical heterogeneity within interfacial water layer

Finally, we study the impact of these local charge correlations on the dynamics of the interfacial water molecules. Here it is useful to consider the H-bond jump time correlation functions eqn 3. In the case of homogeneous dynamics with a single jump rate constant, a single exponential decay is expected for the correlation function. For water molecules initially within the first layer and donating a H-bond to another water, our results in Figure 10 show that the jump time correlation functions eqn 3 exhibit single exponential ps decays at every electrode potential, with no visible sign of slower ns components. This absence of dynamical heterogeneity therefore markedly contrasts with the situation previously reported 32,67 at the Pt-water interface. This suggests that compared to the strongly interacting Pt adlayer, the weaker interaction of water with graphene facilitates the water structural rearrangements, which accordingly require no very cooperative motions.

Concluding remarks

The present study has shown that the jump picture for water reorientation both qualitatively and quantitatively characterizes in detail how interfacial water dynamics are affected by electrified graphene surfaces. Among these characterizations, the jump model reveals the molecular origin of the asymmetry between the impacts of positive and negative electrodes on water reorientation dynamics. This study has also shown that the graphene interfacial water layers do not exhibit significant dynamical heterogeneity. Further efforts will address the effects of electrolytes on interfacial water dynamics at charged electrodes, and how ions affect the water jump dynamics 53 at these interfaces. We consider that

a

  PASTEUR, Department of Chemistry, École normale supérieure, PSL University, Sorbonne Université, CNRS, 75005 Paris, France; E-mail: damien.laage@ens.fr

Fig. 1

 1 Fig. 1 Water reorientation times t OH 2 a) for water OH groups initially in the first layer (|z O z C | <4.5 Å), and b) for water OH groups initially in the first layer (red), in the second layer (green, 4.5< |z O z C | <8.5 Å) and beyond (black), together with the bulk liquid water value 48 (dots).

Fig. 2

 2 Fig.2Normalized probability distributions for the angle q between a water OH and the surface normal vector along the distance Dz to the surface plane at a) 0 V, b) -1 V, c) + 1 V. See inset in panel a) for the q and Dz definitions. (The approximate population fractions in the S/W T /W B states are respectively 0.05/0.80/0.15 at 0 V, 0.15/0.75/0.10 at -1 V and 0.00/0.82/0.18 at +1 V.)

Fig. 3

 3 Fig. 3 Schematic illustrations of water structure and H-bond populations at (a) neutral, (b) negatively charged and (c) positively charged electrodes. Second hydration layer water molecules are shaded, while water molecules in further layers are not shown.

Fig. 5

 5 Fig. 5 Negative electrode interface. a) HB populations in the S (orange) and W (blue) states within the first hydration layer; b) HB jump times for an OH group initially in the S (orange) and W (blue) states, and in the middle of the cell (black circles) 48 , towards any other state; c) water reorientation times from simulations (black) and from the extended jump model (red).

Fig. 6

 6 Fig. 6 Positive electrode interface. a) Oxygen radial distribution function around water oxygen atoms within the positive interface first hydration layer; b) H-bond jump times for an OH group initially at the positive interface and switching between W or bulk stable H-bond acceptors, obtained from simulations (black) and from our model eqn 6 (red) for increasing electrode potentials.

Fig. 7 a

 7 Fig. 7 a) Mean electrical potential profile across the cell; b) average electric field along z in the middle of the cell (-10 Å< z < +10 Å).

Fig. 9 a

 9 Fig. 9 a) Normalized atomic charge spatial correlation functions at the positive (crosses) and negative (circles) interfaces; b) normalized atomic charge time-correlation functions at the positive (solid lines) and negative (dashes) interfaces.

[1][2][3][4][5][6][7][8][9][10][11]

Jump time correlation function the molecular picture that has been established here and in ref. 15 will be essential for future descriptions of chemical reaction dynamics and ion transport at electrode interfaces.

Conflicts of interest

There are no conflicts to declare.

Notes and references