Enhanced conductivity of water at the electrified air–water interface: a DFT-MD characterization

Fabrizio Creazzo*a, Simone Pezzottiab, Sana Bouguerouaa, Alessandra Servac, Jiri Sponerd, Franz Saijae, Giuseppe Cassone*de and Marie-Pierre Gaigeot*a
aLAMBE UMR8587, Univ Evry, Université Paris-Saclay, CNRS, 91025 Evry, France. E-mail: fabrizio.creazzo@univ-evry.fr; mgaigeot@univ-evry.fr
bLehrstuhl für Physikalische Chemie II, Ruhr-Universität Bochum, 44780 Bochum, Germany
cSorbonne Université, CNRS, Physico-chimie des électrolytes et nano-systèmes interfaciaux, PHENIX, 75005, Paris, France
dInstitute of Biophysics of the Czech Academy of Sciences, Královopolská 135, 61265 Brno, Czech Republic. E-mail: giuseppe.cassone@ipcf.cnr.it
eCNR-IPCF, Viale Ferdinando Stagno d'Alcontres 37, 98158 Messina, Italy

Received 26th December 2019 , Accepted 12th February 2020

First published on 12th February 2020

DFT-based molecular dynamics simulations of the electrified air–liquid water interface are presented, where a homogeneous field is applied parallel to the surface plane. We unveil the field intensity for the onset of proton transfer and molecular dissociation; the protonic current/proton conductivity is measured as a function of the field intensity/voltage. The air–water interface is shown to exhibit a proton conductivity twice the one in the liquid water for field intensities below 0.40 V Å−1. We show that this difference arises from the very specific organization of water in the binding interfacial layer (BIL, i.e. the air–water interface region) into a 2D-HBond-network that is maintained and enforced at the electrified interface. Beyond fields of 0.40 V Å−1, water in the BIL and in the bulk liquid are aligned in the same way by the rather intense fields, hence leading to the same proton conductivity in both BIL and bulk water.

1 Introduction

The structure of liquid water at the interface with the air is an essential key to rationalize and characterize chemical and physical phenomena observed at such an interface, among which are proton trapping and hopping along “water wires”,1 charge separation/recombination processes,2,3 changes in acidity/basicity with respect to bulk water,4,5 the atypical Pockels effect,6 and surface tension.7

Hassanali et al.1 reported the high affinity of protons for the interface especially in terms of specific proton hopping pathways at the air–water (AW) interface, with protons exchanged between water molecules belonging to the first interfacial layer, via water wires running parallel to the surface. This result strongly suggests that a certain ordering of the water molecules within the surface plane is present at the AW interface. In our recent paper8 – where we have combined density functional theory-based molecular dynamics simulations (DFT-MD) and non-linear vibrational sum frequency generation (vSFG) spectroscopy – we have shown that such an order consists of a two-dimensional (2D) H-bonded network (denoted hereafter as “2DN”), connecting the vast majority of the interfacial water molecules (on average more than 90%) through water–water H-bonds/wires oriented parallel to the instantaneous water surface.9,10 Furthermore, due to the additional constraint imposed by the preferential H-bond orientation, water molecules in the 2DN have fewer degrees of freedom for rotation and libration, which was shown to result in a slower orientational dynamics of the interfacial water molecules and, at the same time, to more dynamical H-bond breaking/reforming processes than in bulk liquid water.9 The structure and dynamics of the 2DN thus provide a framework for the preferential direction of the above-mentioned proton hopping reported in ref. 1 and 11. Interestingly, a recent MD simulation of the AW interface has shown that the application of an electric field perpendicular to the interface induces a less efficient reorientation of water molecules than a field applied parallel to the surface.12 However, the way in which the local structure of interfacial water changes in response to an external static electric field, and how this can affect proton hopping remain poorly understood both at the molecular and macroscopic levels.

We hence report here for the first time, to the best of our knowledge, an ab initio MD study of the microscopic effects produced by an external static and homogeneous electric field applied at the AW interface and oriented parallel to the water surface (i.e. along the −x direction in the simulation box). We reveal the possible perturbations in the 2DN at the AW interface under the influence of an external electric field and the consequence on proton hopping at the electrified interface. Beyond proton hopping, we also characterize the electric conditions for the protolysis reaction 2H2O ⇌ OH + H3O+ to occur, where formally a proton transfer between two water molecules gives rise to the formation of hydroxide (OH) and hydronium (H3O+) ions. The protolysis reaction is however a rare event, and interestingly, Saitta et al.13,14 have shown that it is possible to stimulate the proton transfer process in bulk (liquid and ice) water – and hence, to investigate protolysis in a more systematic fashion – by applying a static external electric field. Based on ref. 13 for liquid water, the electrostatic coupling of interfacial water with an external electric field is expected to perturb the interfacial H-bonded network, hence possibly affecting proton transfer, water dissociation and protolysis at the AW interface.

Some information that can be readily gained upon applying sufficiently strong electric fields is the effective thresholds associated, respectively, with the onset of proton transfer and with the onset of molecular dissociation. In liquid water, fields of ∼0.25 V Å−1 are needed to induce proton transfers and molecular dissociations of water along the 3D H-bonded network,13,15–17 whereas a field intensity of at least 0.35 V Å−1 has to be applied in order to establish a measurable protonic current.13 A further and correlated consequence of the application of static electric fields to liquid water is the gradual alignment of an increasing fraction of molecular dipole moments along the field direction.18 Moreover, as very recently demonstrated by monitoring the IR and Raman spectra of electrified liquid water via ab initio MD,19 static electric fields of intensities beneath the molecular dissociation threshold induce structural changes in the H-bonded network and in the water tetrahedrality, in that the water structure becomes more ice-like.

Here, we show how the proton conductivity is enhanced by the presence of the specific 2DN at the air–water (AW) interface under external fields. This paper is organized as follows. In Section 2, we present the methodology, and Section 3 reports results on the protonic current density–voltage diagram and the structural effects introduced by the external field applied parallel to the air–water 2DN interfacial network. We provide a detailed analysis of the H-bond network, water dissociation and proton conduction properties under increasing field strengths. Concluding remarks are in Section 4.

2 Computational methods

Density functional theory (DFT)-based molecular dynamics (MD) simulations (DFT-MD) have been carried out with the CP2K package,20,21 consisting of Born–Oppenheimer MD by means of the DFT-BLYP22,23 exchange–correlation functional including the Grimme D2 correction for dispersion interactions,24,25 GTH pseudopotentials26 for the oxygen and hydrogen atoms, and a combined Plane-Wave (400 Ry density cutoff) and TZV2P basis set. The simulation box of 19.734 × 19.734 × 35 Å3 is composed of a liquid phase made of 256 water molecules, periodically repeated in the x and y directions and separated by a vacuum layer of 15 Å from the replica in the vertical z direction. See Fig. 1a for a snapshot.
image file: c9cp06970d-f1.tif
Fig. 1 (a) Snapshot extracted from DFT-MD simulations showing the simulation box composed of 256 water molecules. The Willard and Chandler instantaneous surface27 is shown in grey sheets, the identified water layers (BIL and bulk water) are color-coded and discussed in the text. A large slab of 15.0 Å of vacuum is used in order to separate the simulation box from its replicas along the vertical z-direction. (b) Electrified air–water interface: time averaged water density profiles normalized with respect to bulk liquid water obtained for applied electric fields of 0.25 V Å−1 (5 V potential) in black line and 0.40 V Å−1 (8 V potential) in blue line. The density is plotted as a function of the distance from the instantaneous Willard and Chandler surface.27

The 256 neutral air–water (AW) trajectory is the one presented in ref. 8, while the other trajectories in the presence of an external electric field applied parallel to the −x axis have been generated for the present investigation. The non-zero-field regime was explored in the range [0.05; 0.70] V Å−1, the electric field being gradually increased with a step-increment of about 0.05 V Å−1. The implementation of an external electric field in numerical codes based on DFT can be achieved by exploiting the modern theory of polarization and the Berry phase28 (see e.g. ref. 29 for the technical implementation of a static and homogeneous electric field in ab initio codes and ref. 30 for a review of several methods that allow for the application of external fields in various simulation frameworks). In a nutshell, the difficulty in treating finite electric fields in first principles periodic systems is the non-periodic nature of the position operator. Within the modern theory of polarization31,32 and of the Berry phase,28 one can introduce a variational energy functional29

image file: c9cp06970d-t1.tif(1)
where E0[{ψi}] is the energy functional of the system in the zero-field approach and P[{ψi}] is the polarization along the field image file: c9cp06970d-t2.tif direction, as defined by Resta:31
image file: c9cp06970d-t3.tif(2)
where L is the periodicity of the cell and S[{ψi}] is a matrix composed of the following elements
Si,j = 〈ψi|eix/L|ψi (3)
for doubly occupied wavefunctions ψi. Umari and Pasquarello29 demonstrated that this variational approach is valid for treating finite homogeneous electric fields in first-principles calculations and that it can be extended to provide atomic forces in first-principles MD simulations.

We performed simulations at the nominal temperature of 300 K, kept fixed through the coupling of the system with a Nosé–Hoover thermostat. The molecular systems were kept in an isothermal–isochoric (NVT) ensemble and the classical Newton's equations of motion for the nuclei are integrated through the velocity Verlet algorithm with a time-step of 0.4 fs. For each electric field strength, the dynamics was followed for time lengths up to about 30 ps, extending to about 100 ps in the absence of the field. Hence, we globally cumulated a total simulation time of approximately 400 ps.

Analyses of the DFT-MD trajectories into instantaneous surface and water layers (Fig. 1) follow the derivation, respectively, by Willard and Chandler27 and Pezzotti et al.33 Water–water H-bonds are defined through Galli and coworkers’ criteria:34 O(–H)⋯O ≤ 3.2 Å and O–H⋯O angle in the range [140–220]°.

The identification of the water interfacial layers at the neutral AW interface, namely the BIL (binding interfacial layer) and bulk liquid water, has been done following our methodology described in ref. 33 on the basis of water structural descriptors only. As already validated in our previous works9,11,33,35,36 and confirmed by the present results at the electrified AW interface, the BIL is systematically composed of the topmost water molecules located within 3.5 Å from the instantaneous water surface,27 forming less water–water H-bonds (2.9 H-bonds per mol) and being 1.4 times denser than water in the bulk. These water molecules form H-bonds preferentially oriented parallel to the surface plane, resulting in the formation of a collective and extended 2D-Hbond-Network (2DN for short notation) in the BIL.8 This leads to the breaking of centrosymmetry and consequent SFG activity of the BIL.8,11 Further away than 3.5 Å from the surface, centrosymmetric bulk water is recovered (with hence no SFG activity).

One of the three descriptors used to define the BIL,33 namely the water density profile as a function of the vertical distance from the instantaneous water surface, is shown in Fig. 1b for two electrified air–water interfaces (respectively, homogeneous static electric field intensities of 0.25 V Å−1 (5 V potential) and 0.40 V Å−1 (8 V potential) applied along the −x-axis/surface plane). As will become clear in the discussion in the following Section 3, they correspond to crucial field values for the water conductivity in the BIL and bulk regions. The plots reveal that the water density profile is only slightly affected by the increase in the field, and in particular that the first peak (position and intensity, as well as following minimum position) is maintained in the 0.25–0.40 V Å−1 field-regime. The thickness of the BIL is thus not changed in this regime. One can see the onset of small modifications to the second peak and following bulk region in the density profiles with the increase in the field intensity, showing that the field-induced rearrangement kicks-in in the 3D-HB-network of bulk water before it affects the 2D-HB-network of the BIL interface. This will be discussed in more detail later in this manuscript.

According to Ohm's law, the current density is related to the number of charge carriers ΔN flowing through a section area a2 orthogonal to the direction of the electric field within a time interval Δt. With a being the side of the simulation box orthogonal to the field direction and q being the elementary charge (1.6 × 10−19 C), the current density is:

image file: c9cp06970d-t4.tif(4)
expressed in μA nm−2. The protonic conductivity is then calculated as
image file: c9cp06970d-t5.tif(5)
expressed in S cm−1.

3 Results

Although it is now established that applying static electric fields of the order of 0.30 V Å−1 to liquid water favours molecular dissociations,13,16 theoretical modeling of the microscopic behaviour at the air/liquid water (AW) electrified interface provides fundamental insights on the conductivity properties of the interfacial H-bonded network that, to the best of our knowledge, has not been explored so far. From liquid water modelling, the application of an external static field is known to alter the H-bond environment in the bulk, triggering the cleavage of some oxygen–hydrogen (O–H) covalent bonds and thus promoting the hopping of protons along intermolecular H-bonds, resulting in the formation of new ionic complexes such as hydronium (H3O+) and hydroxide (OH) in liquid water. This is due to at least two cooperative roles played by the field, which (i) aligns the water molecule dipole moment vectors along the field direction18 and (ii) elongates/weakens their covalent O–H bonds.9 In neat bulk water, the lowest field intensity able to induce a measurable net proton flow/current has been quantified theoretically to a value of 0.35 V Å−1[thin space (1/6-em)]13,16 (obtained from DFT-MD using the PBE exchange–correlation functional; note that a change in the functional might induce a slight change in the absolute value of the field threshold), while a lower field strength of 0.25 V Å−1 triggers a series of ordered and correlated proton jumps via the Grotthuss mechanism in electrolytic aqueous solutions.16,18,37

In the case of the AW electrified interface here investigated, the first significant molecular dissociation events have been recorded for field strengths equal to 0.30 V Å−1 applied parallel to the air–water surface plane. Moreover, as shown in the protonic current density–voltage diagram plotted in Fig. 2, such a field intensity, which corresponds to the application of a voltage of 6 V at the edges of the employed simulation box, is not only able to trigger water dissociations but also to give rise to a net proton flow both at the AW interface (i.e., in the BIL) and in the bulk liquid.

image file: c9cp06970d-f2.tif
Fig. 2 Left: Protonic current density–voltage diagram calculated in the BIL (green squares) and in bulk water (blue circles). The corresponding electric field strength is given with the top axis. The dotted red line highlights the conductivity threshold discussed in the text. σBIL and σbulk are the conductivity values calculated in the BIL and in bulk water, respectively. Table: for each electric field strength applied (and the related voltage for a cell side of 19.734 Å), a list of protonic current density values calculated in the BIL and bulk water. Data highlighted in red represent the conductivity (σ) threshold discussed in the text.

Molecular dissociation processes (BIL and bulk alike) already start at 0.25 V Å−1 (corresponding to a voltage of 5 V). However, similar to bulk liquid water,13 these events are rare enough, and the created hydronium and hydroxide ions are short-lived (i.e., their lifetime is ≤20–30 fs), which is not enough to give rise to a measurable protonic current. Once a field intensity of 0.30 V Å−1 is applied, the BIL-AW slab shows Ohmic behaviour, as already observed in ref. 13, 16, 18 and 37 for bulk water and electrolytic aqueous solutions. In order to extract the current density contributions arising separately from the BIL and from the bulk liquid, respectively, these two regions have been systematically identified in the simulations based on the procedure presented in ref. 10, 11 and 33. As discussed in the methods section, the BIL includes all water molecules within a slab having a thickness equal to 3.5 Å from the instantaneous water surface, while all remaining water molecules are assigned to the bulk region, as depicted in Fig. 1, independent of the field strength. Importantly, as will be demonstrated later in the text, the water–water BIL-2DN specific 2-dimensional H-bond network is maintained at the electrified AW interfaces, which is of high relevance for the rationalization of our findings for the protonic current densities presented and discussed below.

As shown in Fig. 2, two conductivity regimes can be identified, one for the BIL and one for the bulk liquid. In particular, as displayed in Fig. 2 and in the table included in this figure, when an electric field strength equal to 0.30 V Å−1 (corresponding to 6 V potential) is applied parallel to the water surface, protons start to flow along the field direction, with a higher current density along the water–water 2D-Hbond-network than that in the bulk. While there is also a protonic flow in the liquid, the protonic current density measured in the BIL (1.02 μA nm−2) is roughly twice larger than that of the bulk (0.54 μA nm−2), up to 0.40 V Å−1 fields. Correspondingly, the protonic conductivity in the BIL (σBIL = 3.67 S cm−1) is twice the one of the bulk (σbulk = 1.76 S cm−1). Thus, in the [0.30–0.40] V Å−1 field intensity range (corresponding to 6–8 V potentials), the electrified AW interface (i.e., the BIL) is a much better protonic conductor than the electrified bulk water.

On the other hand, beyond an electric field strength of 0.40 V Å−1 (corresponding to about 8 V potential), the protonic current densities in the BIL and in the bulk liquid become roughly identical. Under such a high-voltage regime (i.e., ≥8 V), the BIL and the bulk protonic conductivities are equal to an average value of ∼4.8 S cm−1 (Fig. 2, right). The lower absolute bulk protonic conductivity found here in comparison to that of the pioneering work of Saitta et al.13 (i.e., 7.8 S cm−1) is presumably due to a combination of differences in the adopted theoretical frameworks between our works (i.e., Born–Oppenheimer vs. Car–Parrinello MD, dispersion-corrected BLYP XC functional vs. PBE, etc.) and to different statistics (i.e., box sizes and simulation timescales).

The rationale behind the significant difference in the conduction properties of the BIL and of the bulk for fields below 0.40 V Å−1 can be ascribed to the specific organization of the interfacial water molecules in the BIL, creating the already mentioned 2DN that connects more than 90% of the water molecules belonging to the BIL within a unique extended and collective network via H-bonds all oriented parallel to the surface plane8 and surviving the application of a static electric field, as demonstrated now.

With the aim of providing a statistical and quantitative analysis of the 2DN in the BIL, Fig. 3 shows the probability distribution Pn(%) of the number of BIL–water molecules (n, x-axis) inter-connected by H-bonds through a non-interrupted 2-dimensional interfacial network. The probability distribution Pn(%) is presented for the zero-field case in Fig. 3a; it is the reference for the two other probability distributions presented here for electric fields of 0.25 V Å−1 and 0.30 V Å−1 (Fig. 3(b) and (c)), the latter being the electric field threshold able to dissociate water molecules and to establish a protonic current.

image file: c9cp06970d-f3.tif
Fig. 3 Probability distribution Pn(%) of the possible structure of water molecules located in the interfacial layer (BIL), 3.5 Å thickness. n (x-axis) is the number of BIL water molecules connected by non-interrupted H-bonds. (a) Zero-field case. (b) Electric field strength of 0.25 V Å−1 (5 V potential). (c) Electric field strength of 0.30 V Å−1 (6 V potential).

As depicted in Fig. 3(a) (and already discussed in ref. 9 and 10), the vast majority of the water molecules (i.e. more than 90%, as obtained by integration of Pn(%) for values n ≥ 38) located in the BIL (binding interfacial layer) form one single collective and extended H-bond structure – i.e., the 2DN – as described in our previous works.9 Less than 5% of interfacial water molecules are found either isolated (n = 1), or involved in dimers (n = 2) or in other small H-bonded structures (n ≤ 5), on average.

Similar considerations hold at 0.25 V Å−1 and 0.30 V Å−1 (Fig. 3(b) and (c)), where the 2DN is not only still present, but is even enforced by the electric field applied parallel to the surface. One can indeed see that the main peak in the Pn(%) distribution is shifted towards a higher central value of water molecules (n) forming the extended and collective 2DN, while the peak distribution is also less broad than in the zero-field case. At both fields shown here, the minimum size of the continuous 2DN motif is obtained for n ∼ 42–45. Not surprisingly, this reveals that the 2DN collectivity benefits from the alignment of the water dipoles under the influence of the external electric field applied along the direction parallel to the water surface (i.e. parallel to the 2DN H-bond direction). Let us also stress here that the current-density in the BIL (Fig. 2) is entirely due to the 90% water molecules that build the special 2DN network at the interface.

Besides, the 2DN is composed of H-bonded water rings, as already emphasized in ref. 1 and 9. These rings are quantified here, following the same method as in ref. 9 for the non-electrified air–water interface. Fig. 4 hence reports the probability distribution Pn(%) of finding ring structures of given sizes in the interfacial BIL-2DN, in the absence of the electric field (Fig. 4(a)), and in the presence of the 0.25 V Å−1 (Fig. 4(b)) and 0.30 V Å−1 (Fig. 4(c)) fields. As far as the zero-field case is concerned, rings composed of four, five, and six H-bonded water molecules are the most likely structural motifs that build the collective 2DN. The most likely ring sizes are 4, 5 and 6, with decreasing order of probability. There are also probabilities to observe rings composed of up to nine water molecules.

image file: c9cp06970d-f4.tif
Fig. 4 Probability distribution Pn(%) of the size of the ring structures formed by the interfacial water molecules within the 2DN. (a) Zero-field. (b) Electric field strength of 0.25 V Å−1 (5 V potential). (c) Electric field strength of 0.30 V Å−1 (6 V potential).

For the two external fields reported in Fig. 4(b) and (c), one can see that the distribution of ring sizes between 4 and 6 is still the most probable, however with a global distribution that now clearly shifts towards the ring size of 5 as the most probable/favored, especially for the 0.30 V Å−1 field applied. The formation of H-bonded rings in the BIL – water with the H-bonds oriented parallel to the water surface plane is the fingerprint of the 2DN at the air–water interface, maintained and even strengthened once the interface was electrified, as shown here.

For proton transfers to occur, protons have to move from one water molecule to the neighbouring one along H-bonded chains of molecules known as “water wires”.1,38 Because of the reduced number of available spatial configurations in the collective BIL-2DN, water molecules within the 2DN have lower degrees of freedom for rotation and libration, which leads to a slower timescale for the orientational dynamics of the interfacial water molecules.10 Somehow counter-intuitively, however, these interfacial water molecules exhibit a H-bond breaking/reforming dynamics that is faster than for the water molecules in the bulk liquid.10 The BIL-2DN and its rings of H-bonded molecules connected to each other through this network of H-bonds formed parallel to the surface (at both zero-field and at the electrified interfaces) indeed create preferential water wires that can favor proton hopping along these wires. The BIL-2DN furthermore leads to an increase of the residence time of protons at the interface, as already reported in ref. 1 and 11. Moreover, the preferential H-bond orientation naturally present in the 2DN along with the fast H-bond dynamics within the surface plane, as reported in ref. 9, enables easier alignment of the water molecules in the BIL in response to an electric field applied along a direction parallel to the surface plane. All these properties highly favor proton hopping along the BIL-2-dimensional-Hbond-wires, more efficiently than along the 3D H-bonded network of the liquid bulk, and also favor more efficient water dissociation and hence higher protonic flows within the BIL, which is indeed what is observed in this work.

A good illustration of these points can be found in Fig. 5, where 3D-plots report the probability of the combined O–O H-bonded distances and O–O orientation of the water–water H-bonds with respect to the applied [E with combining right harpoon above (vector)] field vector (θ in the inset scheme), comparing the results for the water molecules in the BIL (left) and for water in the liquid bulk (right), for the electric field strengths of 0.25 V Å−1 (Fig. 5a, electric field condition before the onset of detectable water dissociations and protonic currents) and 0.40 V Å−1 (Fig. 5b). The probability coding is given by the scale from blue (lower probability) to red (higher probability). Very interestingly and in line with our discussion above, one can see immediately that the 0.25 V Å−1 field-induced reorientation of the H-bonded water molecules measured through θ is more efficient in the BIL-2DN (see Fig. 5a), where the maximum probability (red spots) is observed for values of cos θ between 0.6 and 1.0, compared to those in the liquid bulk where the red spots are found between 0.4 and 0.9. For a field intensity of 0.40 V Å−1 (8 V potential), both BIL-2D and bulk-3D H-bonded networks become equally oriented by the electrostatic driving force. One can indeed see that the 3D-plots presented in Fig. 5b for the BIL and bulk regions are very similar when such a higher field is applied, with the same final net HB-orientation of the water in the two regions. The only appreciable difference is found in the length of the HBs forming the 2D-HB-network in the BIL, which are slightly longer than those of the HBs formed between the bulk water molecules. This was already found at the non-electrified air–water interface9 or at the lower 0.25 V Å−1 field in Fig. 5a.

image file: c9cp06970d-f5.tif
Fig. 5 3D plots of the H-bond patterns between the water molecules in the interfacial layer-BIL (on the left) and between the water molecules in the bulk water (on the right) when an electric field of 0.25 V Å−1 (5 V potential) – (a) at the top – and 0.40 V Å−1 (8 V potential) – (b) at the bottom – is applied. The x-axis reports the O–O distance (Å) between 2 H-bonded water molecules and the y-axis provides the angle (cosine value) between the O–O vector (from donor to acceptor) and the direction (−x axis) of the applied electric field (see inset scheme). The colors represent the probability (P) to find one O–H group forming one H-bond with a given distance and orientation. The probability increases from blue to red, see the scale.

The water wires in the BIL are consequently found to be more oriented along the field direction than the water wires in the bulk, at least at the 0.25 V Å−1 low-field strength. This can also be seen by the eye in Fig. 6(b), and in Fig. 6(c) at the slightly higher 0.30 V Å−1 field strength. As furthermore highlighted in Fig. 6(c), the water wires in the bulk retain their 3D-structure, resulting in proton motions that explore a larger 3D portion of space in the reoriented bulk than in the reoriented 2DN, as illustrated by the two wires in Fig. 6(c). It follows that in order to move any proton from a position A to a position B under an external field applied parallel to the AW surface, a lower number of proton jumps are required along the more aligned water wires in the BIL than along the more spatially spread water wires in the bulk. This leads to the higher conductivity of the BIL in the low-to-moderate field regime, as reported in Fig. 2. Moreover, as shown in Fig. 5a for the 0.25 V Å−1 field strength, the reorientation of the interfacial water molecules in the BIL along the field direction (−x axis) leads to longer and hence weaker H-bonds than in the liquid, which also favors and enhances the proton conductivity.

image file: c9cp06970d-f6.tif
Fig. 6 Snapshots extracted from the DFT-MD simulations representing the instantaneous surface in grey sheets and the two water layers (BIL and bulk water) identified and discussed in the text. Oxygen atoms are colored in red and hydrogen atoms are in white (hydrogen in yellow only in panel (b)). (a) Zero-field. (b) Oriented water molecules along the field direction (−x axis) by an electric field strength equal to 0.25 V Å−1 (corresponding to 5 V potential). (c) Illustration of proton hopping water wires in the BIL and in bulk water under the action of a field intensity of 0.30 V Å−1 (6 V) along the −x axis.

It is important to note that a few H-bonds present in the BIL are naturally weaker (and thus more dynamical) also under the zero-field condition, as a way to satisfy the finite temperature geometrical constraints on the water–water H-bonds and on the rings that thus maintain the extended 2DN structure. The further increase of the number of such weaker H-bonds with increasing field strength is a direct consequence of the additional 1D constraint imposed by the application of the field along one direction only. These weaker H-bonds have an influence on the lifetime of the water wires formed at the interface, which are hence expected to be shorter-lived than the water wires of the bulk due to the increased H-bond dynamics within the 2DN.10 It is well known that autoionization in water is generated by fluctuations of the water dipole moments and is hence connected to librations and to more dynamical water wires that ultimately favour water dissociation.39 The efficient separation of hydronium and hydroxide ions is also due to short-lived water wires, which in turn also reduce the probability of ionic recombination. All these effects play a role for field strengths slightly higher than the water dissociation threshold (0.25 V Å−1). At larger intensities (≥0.40 V Å−1), however, the limited size of the BIL likely leads to saturation of the 2DN conductivity, which cannot be further enhanced by the action of the field. In other words, any differences in structures that exist between the BIL-2DN and the 3D H-bonded network in the bulk are washed out at higher fields, simply because both networks are then equally and completely oriented by the electrostatic driving force.

4 Conclusions

Based on state-of-the-art ab initio molecular dynamics simulations, we have characterized proton transfer and water dissociation at the air–water interface, triggered by intense static and homogeneous electric fields applied parallel to the air–water surface plane. These results have been directly compared with those measured in the bulk portion of the liquid.

We have found that the onset of water dissociation (i.e., the minimum field intensity capable of ionizing water) is not affected by the specific 2-dimensional Hbond network formed by water at the air–water interface. The first formation of hydronium (H3O+) and hydroxide (OH) ions has been recorded at the binding interfacial layer (BIL) and in the bulk at the same field strength (i.e., 0.30 V Å−1). However, the proton transfer activity at low-to-moderate field regimes (≤0.40 V Å−1) is differently influenced in the two regions of the liquid. The response of the current density–voltage diagrams is Ohmic in both cases (provided that a conduction regime has been achieved), the protonic conductivity of the BIL (σBIL = 3.67 S cm−1) is twice the one recorded in the bulk (σbulk = 1.76 S cm−1). By monitoring the behaviour of the H-bond networks in the BIL and in the bulk liquid, respectively, we showed this difference in conductivity to be due to the specifically organised 2-dimensional Hbond network (2DN) shaping the water at the air–water interface, which was shown to enhance the proton transfer events under low-to-moderate (0.30–0.40 V Å−1) electric field strengths applied along the interface plane (i.e. along the 2DN). The reduced dimensionality of the intermolecular network has a clear influence on the behaviour of the water wires responsible for the proton conduction. The better aligned and shorter-lived water wires, as existing in the BIL, lead to more efficient spatially (and temporaly) correlated proton hoppings than those in the 3D liquid bulk. On the other hand, for more intense fields (≥0.40 V Å−1), both BIL and bulk protonic conductivities converge to the same value (∼4.8 S cm−1), because the 1D direction constraint imposed by the stronger electrostatic field now aligns both BIL and bulk water in a similar way and hence reduces the structural differences between the BIL and the bulk H-bonded networks. The insights gained from this investigation certainly could have more practical implications, typically in relation with water splitting in confined electrified/electrocatalytic solid/water environments. According to the present work, any confined environment exhibiting the 2DN structural arrangement of water at the interface would indeed be favorable for water dissociation/splitting, especially under electrified conditions applied parallel to the BIL-2DN surface.

Conflicts of interest

There are no conflicts of interest to declare.


This work was performed under Grants LABEX CHARMA3T 11-LABEX-0039/ANR-11-IDEX-0003-02 ‘Excellence Laboratory’ program of the University Paris-Saclay and ANR DYNAWIN ANR-14-CE35-0011-01 (ANR Agence Nationale de la Recherche), and using HPC resources from GENCI-France Grant 072484 (CINES/IDRIS/TGCC). Discussions with Dr D. R. Galimberti are greatfully acknowledged.


  1. A. Hassanali, F. Giberti, J. Cuny, T. D. Kühne and M. Parrinello, Proc. Natl. Acad. Sci. U. S. A., 2013, 110, 13723–13728 CrossRef CAS PubMed.
  2. V. Venkateshwaran, S. Vembanur and S. Garde, Proc. Natl. Acad. Sci. U. S. A., 2014, 111, 8729–8734 CrossRef CAS PubMed.
  3. J. A. Kattirtzi, D. T. Limmer and A. P. Willard, Proc. Natl. Acad. Sci. U. S. A., 2017, 114, 13374–13379 CrossRef CAS PubMed.
  4. X. Zhao, S. Subrahmanyan and K. B. Eisenthal, Chem. Phys. Lett., 1990, 171, 558–562 CrossRef CAS.
  5. R. Bianco, S. Wang and J. T. Hynes, J. Phys. Chem. A, 2008, 112, 9467–9476 CrossRef CAS PubMed.
  6. Y. Suzuki, K. Osawa, S. Yukita, T. Kobayashi and E. Tokunaga, Appl. Phys. Lett., 2016, 108, 191103 CrossRef.
  7. R. C. Tolman, J. Chem. Phys., 1949, 17, 333–337 CrossRef CAS.
  8. S. Pezzotti, D. R. Galimberti and M.-P. Gaigeot, J. Phys. Chem. Lett., 2017, 8, 3133–3141 CrossRef CAS PubMed.
  9. A. Serva, S. Pezzotti, S. Bougueroua, D. R. Galimberti and M.-P. Gaigeot, J. Mol. Struct., 2018, 1165, 71–78 CrossRef CAS.
  10. S. Pezzotti, A. Serva and M.-P. Gaigeot, J. Chem. Phys., 2018, 148, 174701 CrossRef PubMed.
  11. S. Pezzotti and M.-P. Gaigeot, Atmosphere, 2018, 9, 396 CrossRef CAS.
  12. M. Nikzad, A. R. Azimian, M. Rezaei and S. Nikzad, J. Chem. Phys., 2017, 147, 204701 CrossRef PubMed.
  13. A. M. Saitta, F. Saija and P. V. Giaquinta, Phys. Rev. Lett., 2012, 108, 207801 CrossRef PubMed.
  14. G. Cassone, P. V. Giaquinta, F. Saija and A. M. Saitta, J. Phys. Chem. B, 2014, 118, 4419–4424 CrossRef CAS PubMed.
  15. E. M. Stuve, Chem. Phys. Lett., 2012, 519, 1–17 CrossRef.
  16. G. Cassone, F. Creazzo, P. V. Giaquinta, J. Sponer and F. Saija, Phys. Chem. Chem. Phys., 2017, 19, 20420–20429 RSC.
  17. Z. Hammadi, M. Descoins, E. Salançon and R. Morin, Appl. Phys. Lett., 2012, 101, 243110 CrossRef.
  18. G. Cassone, F. Creazzo, P. V. Giaquinta, F. Saija and A. M. Saitta, Phys. Chem. Chem. Phys., 2016, 18, 23164–23173 RSC.
  19. G. Cassone, J. Sponer, S. Trusso and F. Saija, Phys. Chem. Chem. Phys., 2019, 21, 21205–21212 RSC.
  20. J. Hutter, M. Iannuzzi, F. Schiffmann and J. VandeVondele, WIREs Comput. Mol. Sci., 2014, 4, 15–25 CrossRef CAS.
  21. J. VandeVondele, M. Krack, F. Mohamed, M. Parrinello, T. Chassaing and J. Hutter, Comput. Phys. Commun., 2005, 167, 103–128 CrossRef CAS.
  22. A. D. Becke, Phys. Rev. A: At., Mol., Opt. Phys., 1988, 38, 3098–3100 CrossRef CAS PubMed.
  23. C. Lee, W. Yang and R. G. Parr, Phys. Rev. B: Condens. Matter Mater. Phys., 1988, 37, 785–789 CrossRef CAS PubMed.
  24. S. Grimme, J. Comput. Chem., 2004, 25, 1463–1473 CrossRef CAS.
  25. S. Grimme, J. Comput. Chem., 2006, 27, 1787–1799 CrossRef CAS.
  26. S. Goedecker, M. Teter and J. Hutter, Phys. Rev. B: Condens. Matter Mater. Phys., 1996, 54, 1703 CrossRef CAS.
  27. A. Willard and D. Chandler, J. Phys. Chem. B, 2010, 114, 1954–1958 CrossRef CAS PubMed.
  28. M. V. Berry, Proc. R. Soc. London, 1984, 392, 45–57 Search PubMed.
  29. P. Umari and A. Pasquarello, Phys. Rev. Lett., 2002, 89, 157602 CrossRef CAS PubMed.
  30. N. J. English and C. J. Waldron, Phys. Chem. Chem. Phys., 2015, 17, 12407–12440 RSC.
  31. R. Resta, Rev. Mod. Phys., 1994, 66, 899–915 CrossRef CAS.
  32. R. D. King-Smith and D. Vanderbilt, Phys. Rev. B: Condens. Matter Mater. Phys., 1993, 47, 1651–1654 CrossRef CAS PubMed.
  33. S. Pezzotti, D. R. Galimberti, Y. R. Shen and M.-P. Gaigeot, Phys. Chem. Chem. Phys., 2018, 20, 5190–5199 RSC.
  34. J. A. White, E. Schwegler, G. Galli and F. Gygi, J. Chem. Phys., 2000, 113, 4668–4673 CrossRef CAS.
  35. S. Pezzotti, D. Galimberti, Y. Shen and M.-P. Gaigeot, Minerals, 2018, 8, 305 CrossRef.
  36. F. Creazzo, D. R. Galimberti, S. Pezzotti and M.-P. Gaigeot, J. Chem. Phys., 2019, 150, 041721 CrossRef PubMed.
  37. G. Cassone, F. Creazzo and F. Saija, Mol. Simul., 2019, 45, 373–380 CrossRef CAS.
  38. T. von Grotthuß, Mémoire sur la décomposition de l'eau et des corps qu'elle tient en dissolution à l'aide de l'électricité galvanique, 1805 Search PubMed.
  39. P. L. Geissler, C. Dellago, D. Chandler, J. Hutter and M. Parrinello, Science, 2001, 291, 2121–2124 CrossRef CAS PubMed.

This journal is © the Owner Societies 2020