Understanding water behaviour on 2D material interfaces through single-molecule motion on h-BN and graphene

Domain: physics
Published: 2025-11-25
Editorial outcome: minor_revision

Abstract
Understanding water behaviour on 2D materials is crucial for applications in sensing, microfluidics, and tribology. While graphene-water interactions are well studied, water on hexagonal boron nitride (h-BN) remains largely unexplored. Despite its structural similarity to graphene, h-BN possesses polar B-N bonds that give rise to distinct electronic and chemical properties. Most previous studies have also focused on multilayer water, leaving single-molecule dynamics poorly understood. Here we show how individual water molecules diffuse on h-BN compared to graphene using helium spin-echo spectroscopy and ab initio calculations. On h-BN/Ni, water exhibits coupled rotational-translational motion, in contrast to the discrete hopping observed on graphene. Water molecules rotate freely around their centre of mass, and although binding energies are similar on both materials, the activation energy for water dynamics on h-BN is 2.5 times lower than on graphene. These dynamics, which classical models fail to capture, highlight the fundamentally different nature of water transport on polar 2D surfaces. We further demonstrate that the supporting substrate strongly influences water friction, with h-BN/Ni showing markedly lower friction than graphene/Ni, opposite to the behaviour of free-standing layers. These findings challenge assumptions and offer insights for designing microfluidic devices requiring precise control of water mobility.

Introduction
Graphene and hexagonal boron nitride (h-BN) are among the most important and well-investigated 2D materials and find applications in fields as varied as coatings, constructions, nanotechnology, sensors, biomedics and microfluidics1–3. Graphene and h-BN both have a honeycomb structure with very similar lattice constants but completely different electronic properties. Isolated and pristine graphene is a zero-bandgap material, while h-BN is an electric insulator with a band gap close to 6eV4. Both materials have extraordinary mechanical resistance and thermal conductivity, which make them desirable for chemical applications such as catalysis and filtration1,5,6. Water-surface interactions are ubiquitous in these processes, determining structural, dynamic, and chemical properties in applications, such as water purification2,7,8, drug delivery in aqueous media9–11, and hydrogen storage12,13.

Understanding the molecular basis of water-surface interactions is thus crucial for both technological advancements and fundamental physics14,15, providing valuable insight into the nature of water/2D material interfaces, and opportunities for the design of advanced materials by tuning their nanoscopic properties16–19. However, while properties such as wettability and ice nucleation at the macroscopic level are well understood20, measurements at the atomic or molecular level are scarce and focused on metal substrates21–24 or metal-supported thin ionic films such as NaCl25. In addition, previous ab initio computational studies on the interaction of water with 2D material surfaces have primarily examined the adsorption of water, with dynamic behaviour often modelled using classical molecular dynamics simulations26. On the one hand, 2D materials can be used as protective coatings of metal and semiconducting surfaces, to reduce corrosion and inhibit interdiffusion, wetting and icing8,27–29. On the other hand, defect-free 2D material synthesis is best achieved via low-pressure chemical vapour deposition (CVD), and here we have specifically chosen Ni(111) as a substrate in the CVD growth because of the perfect lattice match with both graphene and h-BN.

In the context of water interactions with 2D materials, friction plays a crucial role and understanding the relationship between the energy landscape and surface vibrations (i.e., phononic or “mechanical” friction) remains a central challenge. Recent theory suggests a quantum contribution to friction at the solid/liquid interface, explaining the anomalously high friction of water on graphite compared to graphene30, a phenomenon also linked to electron cooling at the water/graphene interface31. Since h-BN is an insulator, unlike graphene, water dynamics on h-BN are especially interesting due to potentially different electronic couplings. Tocci et al. predicted higher macroscopic friction on h-BN via classical and ab initio simulations32,33, and recent work shows water’s vibrational modes influence interfacial friction with graphene34. However, many such studies rely on classical forcefields35. Here, we instead use an ab initio description of a single H2O molecule, together with experimental data, to go beyond classical models.

Experimental investigations of water/2D interfaces are limited by the challenges of imaging picosecond and sub-angstrom dynamics36. The high proton mobility and delocalisation in H2O, even at low temperatures, hinder atomic-level characterisation37. Furthermore, water is easily damaged by electrons and high-energy photons22,38,39. Helium spin-echo (HeSE), with its picosecond resolution and use of non-destructive, low-energy He atoms, circumvents these issues40–43. A previously reported HeSE study of water motion on graphene44 illustrated the stark contrast between water diffusion on weakly interacting 2D materials such as graphene26,45 and metal substrates, where water is bound much stronger and further influenced by hydrogen bonding between the molecules23,39,46. To establish whether this is a general trend on 2D materials, further studies extending to, e.g., h-BN, are required. Moreover, the mentioned graphene study did neither consider the role of additional molecular degrees of freedom, nor the influence of single-molecular friction on water mobility for a supported 2D material26.

In the present study, we demonstrate that for H2O diffusion on epitaxial graphene and h-BN, the single-molecule perspective is highly complex and incompatible with both a single-point particle diffusion model and with liquid water behaviour. Using the HeSE technique, illustrated in Fig. 1a, we measure surface correlations in the water motion. In contrast to water motion over graphene44, single molecules can easily reorient on h-BN, and their motion can no longer be treated as a series of jump-like motions by point-like particles occupying equivalent adsorption sites. By analysing the dephasing rates in correlation measurements, we find that despite a water adsorption energy similar to that on graphene (within 16%), the activation energy on h-BN is much smaller than on graphene. The dynamics involve a fast rearrangement of the water molecule orientation during diffusion, in which the O-H bonds precess around the water centre of mass during translational motion. We rationalise these findings in the context of the multidimensional potential energy surface (PES) during diffusion and molecular friction on h-BN/Ni(111). With ab initio calculations, including a detailed quantum chemical study using density functional theory (DFT) and ab initio molecular dynamics (AIMD) simulations, we show how water behaves differently on h-BN than on graphene and establish the influence of the supporting metal substrate on the friction coefficient. The concept of friction as a single parameter/average over a microscopic ensemble is typically used at the macroscopic level, e.g., for liquid/solid interfaces, as mentioned above32,33. However, we illustrate that atomistic details such as the potential energy corrugation and the changes in the electronic structure induced by the supporting substrate are important, e.g., the behaviour of water on h-BN and graphene is reversed if the support is explicitly included in the theoretical model.

Results
Multidimensional energy landscape of H2O on h-BN

We first address the deficiency of single-molecule studies for water on h-BN and other supported 2D materials by establishing the energy landscape for water binding and adsorption on h-BN, including the complete system with the corresponding metal substrate. With the exception of ice structure studies on metal-supported h-BN forming a Moiré superlattice47, where similarities to water on a graphene Moiré superlattice48,49 and on functionalised graphite50 were found, no single-molecule studies on h-BN/Ni(111) are available. Thus, to get a full picture of the adsorption behaviour of water on h-BN/Ni(111), we investigated the system both experimentally by conducting extensive adsorption and desorption measurements and theoretically using detailed vdW-corrected DFT calculations.

In the experiments, the Ni(111) substrate was first prepared under ultrahigh vacuum (UHV) conditions47, and h-BN was grown following a chemical vapour deposition (CVD) process according to Auwärter et al.51 (see Methods - Experiment and sample preparation and Supplementary Note 1.1). To verify the surface quality of the h-BN overlayer on top, we performed diffraction scans in the high-symmetry ΓM orientation before and immediately after the CVD growth (see Fig. 1). The position of the first order diffraction peak at ΔK ≈ 2.9 Å−1 clearly shows that the h-BN layer exhibit the same periodicity as the Ni(111) surface and that a clean, ordered overlayer of h-BN has been formed. We then started to determine the conditions under which individual water molecules stick on the surface and remain mobile, by conducting extensive adsorption and desorption experiments on the h-BN/Ni surface.

Continuous dosing with water below 120K (Supplementary Note 1.3 and Supplementary Fig. 2) results in the formation of an amorphous layer of solid water covering the entire surface44, while increasing the temperature to 135K prevents full surface coverage, indicating hydrophobic behaviour and island formation on h-BN/Ni, similar to water adsorption on graphene/Ni44. As such it does not exclude the formation of clusters or islands; however, once small aggregates such as dimers have formed, they are unlikely to dissociate, and, as discussed below, exhibits significantly reduced diffusivity. The elastic component, from the dynamics measurements shown in the bottom panel of Fig. 1b, confirms that no ordered H2O structure is present for the experimental conditions of the dynamics measurements52. By performing thermal desorption measurements (Supplementary Fig. 4), in which the He reflectivity is measured during successive surface heating, we estimate the desorption energy to ≈(0.53±0.06)eV, which is very similar to the experimentally determined 0.52 eV for water on graphene/Ni(111)44.

To obtain a better understanding of the water-surface interaction of a single molecule, we performed a detailed quantum chemical study of the geometries and adsorption energies (Eads) of H2O on h-BN/Ni(111) (see Methods section Theoretical methods). Initially, seven high-symmetry sites were chosen (Supplementary Fig. 6) to construct a PES as a function of the oxygen position. The oxygen coordinates were fixed at a given position, and the hydrogen atoms were allowed to relax. The molecular orientation was optimised to find the most favourable conformation for each site (the more negative Eads, the stronger the adsorption). Further optimisations in which the oxygen atom was allowed to fully relax were included to find the global minimum (the most stable configuration of the adsorbate-surface system). The energy differences between different adsorption sites indicate a weakly corrugated PES (Fig. 2b), where minor changes in the adsorption height of the water molecule result in minimal energy variation (Fig. 2c).

As shown in Fig. 2a, H2O physisorbs near a boron atom and exhibits an Eads of −0.25eV (compared to an Eads of −0.21eV for graphene/Ni, as shown in Supplementary Fig. 6). This physisorption site is favoured due to the interaction between the oxygen lone pair and the partially positively charged surface boron atom, forming a weak intermolecular bond. Additionally, one of the hydrogen atoms in the water molecule and the nearest surface nitrogen atom engage in hydrogen bonding, which is characterised by Hirshfeld charges of +0.17e and −0.16e on the boron and nitrogen atoms, respectively, further enhancing the stability of this configuration. The high electron density and localised bond population of each B-N bond render the “bridge" sites of h-BN unfavourable for H2O adsorption. Consequently, the oxygen atom tends to relax slightly away from the site directly above the bond.

To determine the effect of the supporting metal substrate, the adsorption energetics and geometries of H2O on freestanding h-BN were also calculated and analysed. Similar to that on the h-BN/Ni surface, our DFT calculations of water adsorption on a free-standing h-BN layer reveal the formation of a hydrogen bond between the water molecule and a surface nitrogen atom (see Supplementary Fig. 6). The hydrogen bond stabilises this configuration, resulting in an Eads of −0.18eV, which indicates slightly weaker adsorption compared to h-BN/Ni. In the case of free standing h-BN, the distal hydrogen atom is oriented towards another surface nitrogen atom, with a H-N distance of 3.5 Å: This arrangement maximises the long-range hydrogen bonding. The differences between adsorption sites, orientations and bonding interactions between the freestanding h-BN and h-BN/Ni can be attributed to differences in charge localisation and polarity. In the case of freestanding h-BN, the surface nitrogen exhibits a higher Hirshfeld charge of +0.21e and -0.21e on the boron and nitrogen atoms, respectively, indicating an increased hydrogen bonding strength. Thus, the adsorption site that maximises hydrogen bonding is considered to be the most favourable in the context of free-standing h-BN.

A low barrier for water dynamics on h-BN

In the next step, we establish the energies which are necessary for the onset of single-molecular H2O dynamics on h-BN/Ni and compare those to the graphene/Ni system. To study the molecular motion of water on the h-BN/Ni(111) surface with a temporal sensitivity on the picosecond timescale, we performed HeSE experiments in the temperature range of 120–135K at a H2O coverage between 0.12 and 0.20 monolayer (ML) (see Supplementary Note 1.2). At this coverage well below one monolayer, the effects of correlated motion are negligible, allowing us to study the single-molecule dynamics. During a single spin-echo measurement, the decay of the intermediate scattering function (ISF) versus the spin-echo time tSE provides a measure of the surface correlation and thus the timescale of motion at the length scale corresponding to the momentum transfer ΔK=∣ΔK∣. By measuring the ISF over a range of temperatures and ΔK values, we obtain details about the diffusive motion and its length scales in real space. As illustrated in the inset of Fig. 1a, each measurement can be described by a single-exponential function:

I(ΔK,tSE)=Aexp(−α(ΔK)tSE)+C(ΔK). (1)

In Eq. (1), A is the amplitude at t=0, the constant C(ΔK) represents the static level of the surface, i.e., the fraction of the signal that arises from immobile parts on the surface. The prefactor in the exponent α(ΔK) is the so-called dephasing rate. The dephasing rate, in units of s−1, measures how fast diffusion on the surface occurs. Its functional dependence on the momentum transfer ΔK parallel to the surface contains a variety of information and provides a signature of both the rate and mechanism of molecular motion.

We measured the temperature-dependence of α at a fixed coverage and a fixed momentum transfer ΔK=0.6 Å−1 for both the high-symmetry ΓM and ΓK directions over a temperature range of 120–135K (see the Arrhenius plot Fig. 3a). For thermally activated motion, the rate is given by Arrhenius’ law,

α=α0exp(−Ea/(kBT)), (2)

where kB is the Boltzmann constant, T is the surface temperature of the sample and Ea is the activation energy for diffusion. The uncertainties in the data points are the corresponding confidence bounds (1σ) of the exponential fits. From the slope of the linear fit in Fig. 3a, we obtain an activation energy of Ea=(24±8)meV. We note that due to the proximity of the diffusion and desorption regimes, the temperature range is constrained to 120–135K, within which coverage remains stable. To reduce the confidence interval of the linear fit in Fig. 3a, we include data points for both the ΓM (blue) and ΓK (green) azimuths, assuming that diffusion occurs at a similar rate along both directions, as also supported by AIMD calculations.

The experimental value can then be compared to theoretical models of the reaction pathway from DFT calculations, where the minimum energy pathway for water migration connects the global minimum energy site (the N atom) to the nearest equivalent site on h-BN/Ni. Following an evaluation of the transition state pathways, the energetically most favourable pathways are shown in Fig. 3b, c. Fig. 3b illustrates a sequence where the water molecule undergoes a dynamic process involving translation and rotations. We note that such multi-step diffusion pathways can, in principle, alter the Arrhenius prefactor through entropic contributions to the transition state, although such effects are not resolved in the experimental dephasing observable α(ΔK), and in general both prefactor and activation energy may vary with temperature, as shown, e.g., in ref. 53 The molecule first crosses a boron site, then a hollow site, before ultimately rotating into an equivalent position adjacent to another nitrogen atom. The calculated activation energy barriers for traversing the boron and hollow sites are 24meV and 40meV, respectively. Following this motion, a barrierless rotation around the nitrogen atom occurs, returning the water molecule to its original orientation. Furthermore, alternative paths involving translation without molecular rotation were explored, yielding different energy barriers depending on the traversed surface sites Fig. 3c. The hollow site exhibit the lowest activation energy barrier at 31meV, followed by the boron site at 32meV, with the nitrogen site presenting the highest energy barrier at 34meV. In each case, the transition state structure closely aligns with the midpoint between each energy minimum.

The calculated activation energy barriers to rotation across a boron site (24meV) and translation without rotation (32–34meV) are consistent with the experimental measurements with an Ea of 24±8meV. As the energy barriers for these pathways are considerably low, each pathway is expected to be readily accessible and can occur with minimal energy input at low temperatures. Repeated traversals along these pathways likely lead to the isolated water molecule undergoing a sequence of jumping and spinning motions across the surface. Apparently, such behaviour of a single water molecule has not been reported previously, while the terminology for the molecular degrees of freedom follows from thiophene diffusion54.

For completeness, we also examined the possibility of dimer formation, dissociation, and diffusion, as detailed in Section 2.2 of the Supplementary Information and illustrated in Supplementary Figs. 7 and 8. While dimer formation is thermodynamically allowed under low-coverage conditions, the associated diffusion barriers are substantially higher than those for single-molecule motion. Specifically, calculated barriers for dimer translation and rotation range from 79 to 220meV, in contrast to the 24meV barrier observed experimentally for monomer diffusion. Moreover, dimer dissociation involves an energy cost of ≈220meV, indicating that dimers, once formed, are kinetically stable but effectively immobile under the experimental conditions. These findings support the conclusion that the observed dynamics predominantly arise from isolated water molecules and that single-molecule simulations provide the most accurate description of the experimental system.

To sum up, based on both experimental results and first-principles calculations, we conclude that the differences in Eads between various adsorption sites is very small, and the corrugation of the potential energy surface for water motion over h-BN/Ni (measured as the activation energy variation on the surface plane, Fig. 2) is particularly weak. Such a small barrier for the onset of motion could be easily overcome by tip interactions in STM measurements at cryogenic temperatures55, whereas the low-energy He atoms do not interfere with the motion of H2O. Most importantly, compared to single-molecule water diffusion on graphene/Ni, with an activation energy of 60±4 meV44, both the experimental and theoretical values for h-BN/Ni are considerably smaller by a factor of 2.5. These differences in activation energy suggest a significant distinctness of the molecule-surface interaction, including the underlying energy landscape and surface dynamical properties, which are explored in the following sections, including also the Molecular friction of water on graphene and h-BN.

Walking motion of water on h-BN

Despite the similarity of the substrate structure and H2O adsorption energy on h-BN/Ni to those on graphene/Ni, we observe a pronounced difference in the multidimensional energy landscape and a low activation energy for molecular dynamics on h-BN/Ni(111). In the following, we address the details of H2O motion on h-BN/Ni, specifically considering the additional molecular degrees of freedom and the stark contrast in water motion on h-BN/Ni compared to graphene/Ni. The blue points in Fig. 4a show the variation in α(ΔK) for water molecules at T=120K along the high-symmetry ΓM direction in reciprocal space.

For each data point, a single ISF was recorded and fitted to Eq. (1). We used an iterative routine to optimise the range for inclusion, in which the data show the expected deviation from an exponential for short times44. The error bars correspond to the confidence bounds (1σ) of the exponential fit. Our data exhibit several signatures suggesting that a complex interplay among several processes characterises the observed motion. A dip in the dephasing rate appears at ΔK≈2.8 Å−1, which is very close to the position in reciprocal space of the first-order diffraction peak (Fig. 1b).

By fitting the Chudley-Elliot (CE) model, which is normally used to describe hopping of a point-like particle between adsorption sites on a Bravais lattice,44,54,56–60 to the data points with ΔK>1 Å−1, as illustrated by the red dashed lines in Fig. 4, we can extract a set of “effective” mass transport coefficients. With this model, we obtain a residence time τ=290ps, an average jump length 〈l〉=5.5 Å and a diffusion constant D=(0.26±0.5)nm2ns−1 at 120K. However, the simple CE model for jump diffusion cannot reproduce our data in the entire ΔK range, especially for small values of ΔK<1 Å−1. In the region of very low ΔK values, the contribution to the dephasing rate from motion parallel to the surface can be neglected, and we can isolate the effects of perpendicular motion, which makes this region particularly interesting for our analysis, whereas at higher ΔK values, parallel motion becomes evident, and both processes contribute to the signal59,61. We observe that the dephasing rate remains almost constant for ΔK<1.5 Å−1 and, more importantly, remains at a nonzero constant value as ΔK approaches 0. Such a behaviour for ΔK→0 is characteristic of confined motion perpendicular or parallel to the surface26,54,59,61. Confined motion can include several processes, such as rotation, flipping or spinning of the molecule with respect to its centre of mass54. The fact that our ISFs always decay to a finite value (see the static component in Fig. 1b) is also a clear indication of confined motion59,61.

To investigate the dynamics in more detail and confirm the nature of the confined motion, we performed AIMD simulations at 150K, which include all degrees of freedom of the adsorbed water molecule. Starting from the global minimum determined by DFT, as shown in Fig. 2, the trajectories and starting velocities were randomised in canonical simulations (NVT), see Methods - Theoretical methods, with an example trajectory plotted in Fig. 5a, while the real-space motion is illustrated by Supplementary Movies obtained from the AIMD simulations. Analysis of the oxygen atom path reveals a preference for trajectories that remain close to the nitrogen atoms, consistent with the PES and transition state energies discussed in the previous section. Furthermore, the z coordinate of the oxygen atom (corresponding essentially to the height of the water molecule) oscillates as the water moves between the nearest neighbour minima (see the inset Fig. 5b).

The observation of such a quasi-continuous, “walking" motion of water on h-BN/Ni(111), characterised by a combination of translational and rotational dynamics, shares conceptual similarities with prior studies of single-molecule water diffusion on NaCl(001) in refs. 25, 62: Water on a NaCl bilayer on Ag(111) undergoes a rocking motion, where translations are coupled with rotation, driven by electrostatic interactions with Na+ ions and hydrogen bonding with Cl− ions25.

The regular oscillations indicate periodic surface features with consistent water-surface interactions and resulting forces. The movement of the hydrogen atoms, which precess around the oxygen trajectory, indicates a fast spinning motion of the molecule perpendicular to its direction of travel, meaning that the position of the O-H bonds with respect to the oxygen atom can easily flip (averaging 1.64 rotations per picosecond). Supplementary Movie 1 illustrates the mentioned trajectory, while Supplementary Movies 2 and 3 illustrate other trajectories, e.g., in the case of Movie 3 with a more confined motion. For the other movies, please refer to Data availability. These trajectories can also be used to analytically calculate the ISFs (see Supplementary Note 2.5). For vanishing ΔK, the analytical ISFs exhibit a decay analogous to our experimental ISFs, the results confirm that the motion along the z−component contributes to the constant offset of α(ΔK). For comparison, Fig. 4b displays α(ΔK) for graphene/Ni, which exhibits a negligible contribution from rotational motion relative to h-BN, and is instead dominated by repulsive intermolecular interactions. For a detailed analysis, see refs. 44, 63 and the supporting information of ref. 44.

To better understand the nature of the observed dynamic motion, the spatial probability distribution of water was calculated on h-BN/Ni and graphene/Ni from AIMD trajectories. A comparison of the probability distributions in Fig. 5c indicates that the water motion in the two systems substantially differs. For the h-BN/Ni system, the probability distribution demonstrates notable uniformity among the sites, including the boron atom and hole sites and their interconnections, with the only exception being the nitrogen atom site. In contrast, the site probability for graphene/Ni as shown in the right panel of Fig. 5c is highly localised, leading to the discrete jumps observed in experiments (see Supplementary Note 2.3). The probabilities for h-BN/Ni coincide with the PES determined through assessment of the water Eads, with higher probabilities at the preferred positions (Fig. 2), which favour close proximity to the boron atom and avoidance of the nitrogen atom.

The more uniform distribution of water on h-BN/Ni is also evident from the present experimental data, while for graphene/Ni, H2O most likely sits on one specific site, which requires the molecule to traverse a larger distance across the surface, resulting in point-like jumps between these sites according to the CE model, as has been previously reported44. A more extensive comparison which includes also the free-standing h-BN and graphene systems can be found in Supplementary Fig. 9. Notably, at a temperature of 150K for h-BN/Ni, the probability of finding the water molecule in the conformation with the minimum energy is 5.8%, while there is a 90% probability of finding the molecule within one of the 48 lowest energy sites (out of 240). These states are distributed across the surface and follow a similar pattern to the AIMD probability distribution, as illustrated in Supplementary Fig. 10. At low temperatures, a variety of water conformations coexist across sites; therefore, water diffusion on h-BN/Ni is more complex than that described by the jump-diffusion model of water on single-crystal metal surfaces.

In conclusion, despite the similarity of the substrate structure and H2O adsorption energy on h-BN/Ni to those on graphene/Ni, both experiments and theory show completely different motions of single water molecules on the two surfaces. The diffusion of water on h-BN is slower than that on graphene (D=0.26nm2ns−1 cf. D=0.32nm2ns−1 extrapolated to 120K), but the activation energy for molecular motion is negligible compared to that on graphene because of the simple rearrangement, i.e., rotational and flipping motion, of the molecule on the substrate. These findings are in line with the underlying PES, which is much more uniform for water on h-BN/Ni than for water on graphene/Ni. In the following, we further rationalise these findings by illustrating the strongly divergent behaviour of the two systems, considering the atomic-scale friction.

Molecular friction of water on graphene and h-BN

In the previous sections, we demonstrated how the nature of the 2D substrate and the additional molecular degrees of freedom significantly impact the mobility of water on 2D materials. We rationalise these findings in term of molecular friction and illustrate the correlation between friction and both vibrational coupling between surface and adsorbate and the corrugation of the PES. Specifically, we show how both changes of the potential energy surface corrugation and the vibrational coupling, upon inclusion of the supporting metal substrate underneath the 2D materials, give rise to a completely different molecular friction when compared to the free-standing 2D materials.

Therefore, we performed a series of AIMD simulations on both h-BN and graphene, with and without the presence of a supporting metal substrate.

We employed the microcanonical ensemble, considering a supercell of the corresponding system (see Methods - Theoretical methods), to determine the friction of water. The friction coefficient λ for each MD simulation was calculated using the Green-Kubo (GK) relationship, which is defined as:32,64

λGK = lim t→∞ λGK(t), (3)

with,

λGK(t) = 1/(2AkBT) ∫0t ⟨F(0)F(t0)⟩ dt0, (4)

where F(t0) is the lateral force acting on the sheet at a time t0 and A is the interfacial lateral area. The λGK(t) up to 2ps were calculated by averaging results from AIMD simulations, each involving a single water molecule on the examined surfaces. Figure 6a shows the λGK(t) for a single water molecule on free-standing graphene and h-BN, as well as a comparison with that for water on metal-supported h-BN and graphene. A single error bar is shown for the graphene/Ni system, which illustrates the statistical uncertainty and the degree of convergence when the simulation is extended to 3ps; further details are provided in Supplementary Note 2.6.

Overall, the λGK values observed for freestanding h-BN (19.9⋅106Nsm−3) are considerably greater than those for freestanding graphene (6.2⋅106Nsm−3), by a factor of ≈3.2. These results are consistent with the MD results reported by Tocci et al.32, where classical force fields were assumed for MD simulations of systems containing ≈400 water molecules and revealed that h-BN has a λ that is a factor of 3.1 greater than that of graphene (with values of 30⋅104Nsm−3 and 9.6⋅104Nsm−3 for h-BN and graphene, respectively). Following the same methodology for a water monomer on the Ni-supported 2D materials, we determined the molecular friction of the systems, with λGK for graphene/Ni (9.5⋅107Nsm−3) being 7.9 times greater than λGK for h-BN/Ni (1.2⋅107Nsm−3) throughout the simulation time of 2ps. Thus, as illustrated in Fig. 6a, by including the substrate (Ni), we observe a major change in friction, and the relative relationship between the λGK of h-BN and graphene is reversed, i.e., the friction of water on supported graphene is much greater than on supported h-BN.

Interestingly, in a previous study of water on NaCl bilayers, the role of the underlying substrate was also highlighted-specifically, the Ag(111) substrate was found to significantly affect water adsorption on NaCl/Ag(111) compared to bulk NaCl25. Although the study did not address friction, it clearly showed that the supporting metal substrate (Ag(111)) modifies the electronic structure and impacts water diffusivity. On the other hand, the diffusion barrier was substantially higher, around 150meV, illustrating that the energetics are clearly different compared to the supported 2D-material systems considered here.

We note that our computational methodology differs from the work of Tocci et al. in the fact that we analysed the motion and friction of single water monomers instead of the average motion of a confined water sheet. Thus, higher λGK values are obtained in our study, which can be attributed to the interaction area between the liquid and the graphene/h-BN sheet, where we define A as the area of a water molecule. By assuming an A value corresponding to the full coverage of water molecules on the freestanding material, the λGK values are of the same order of magnitude as in previous works32,33 (10.5⋅104Nsm−3 and 3.2⋅104Nsm−3 for h-BN and graphene, respectively). Therefore, the frictional changes experienced upon inclusion of the Ni-support are indeed caused by changes of the surface properties and not any coverage-induced effects.

As shown in Fig. 5c, the different probability distributions (and therefore different corrugation of the PES) for h-BN/Ni and graphene/Ni give rise to completely different motions, which in turn are correlated with the different friction coefficients on the studied 2D systems (Fig. 6d). Phononic effects also contribute to friction, and we thus considered the vibrational density of states (DOS) of the 2D material systems. We calculated the difference between the phonon DOS of the freestanding and Ni-supported systems for both graphene and h-BN, as shown in Fig. 6b, c, respectively. The vibrational power spectrum was calculated from the Fourier transform of the velocity autocorrelation function of the water molecule, with further details described in Supplementary Note 2.7.

The vibrational spectra in Fig. 6b-c show that the coupling between the 2D overlayers (graphene and h-BN) and water changes when the 2D material interacts with the Ni substrate, leading to differences in the distribution of vibrational states between the stretching modes (at ≈3800cm−1) and bending modes (at ≈1900cm−1) - see Supplementary Note 2.7 for the assignment and accuracy of the vibrational states65,66. For h-BN, the frequencies of the vibrational bands do not change considerably, while for graphene, the frequency of the 3000cm−1 band decreases upon inclusion of the Ni substrate. Furthermore, for h-BN, as shown in Fig. 6b, only a small dampening of the bending mode at ≈1400cm−1 is observed, while for graphene in Fig. 6c, the bending mode is significantly damped, with a strong redistribution to the stretch mode at 3900cm−1, upon inclusion of the Ni substrate. Comparison of the corresponding peak areas shows that the stretch peak on graphene/Ni is 15 times greater compared to free-standing graphene. Hence, in the case of water on graphene/Ni(111), the substantial increase in the population of the H-O stretching band indicates a stronger molecule-surface coupling with a broader range of available states, with a greater λGK(t) experienced by the molecule. Further analysis of the vibrational power spectrum can be found in Supplementary Note 2.7. As shown in Supplementary Fig. 16, upon inclusion of the Ni(111) substrate, both h-BN and graphene exhibit similar electronic DOS at and around the Fermi level to those without the substrate; thus, we anticipate that the change in friction is not caused by changes in the electronic structure of the material. However, there is a strong correlation between the corrugation of the PES and the corresponding friction, as described above.

Finally, to establish the relationship between the properties of the different 2D systems and λGK, we conducted a regression analysis based on properties derived from the time-evolved spatial coordinates of the AIMD trajectories, with statistical significance of the differences being assessed using one-way ANOVA67, where p<0.01 is considered significant. The strong correlation between the corrugation of the PES and the corresponding friction is demonstrated in Fig. 6d: λGK follows a linear increase as a function of the mean free energy variation ΔG, which describes the overall energy landscape experienced by the water molecule on each of the surfaces. I.e., water molecules situated above the graphene/Ni system encounter a much larger corrugation of the underlying PES according to the mean ΔG, resulting in much larger diffusion process and is thus related to the anomalously high friction on graphene/Ni. For all systems, the curves represent an average of the friction coefficients and phonon DOS estimated from AIMD simulations up to 2ps. Panel d illustrates that the increasing friction of the systems is correlated with both the surface topography (the corrugation of the PES) and substrate-induced vibrational changes experienced by the molecules. The linear regression plots, with error bars representing the margin of error (corresponding to the 98% confidence interval to the mean) from 10 independent NVE simulations per system (see Theoretical methods), depict the mean free energy (ΔG), average H-H distance, and average H-O distance for freestanding h-BN and graphene, as well as for h-BN/Ni and graphene/Ni.

Moreover, the distance between the two hydrogen atoms in H2O is correlated with λGK, as shown in Fig. 6d (R2=0.98), with an increase in the mean H-H distance resulting in an increase in λGK. This correlation is linked to the vibrational bending modes of the water molecule. Furthermore, λGK is correlated with the O-H bond length (R2=0.99), as shown in Fig. 6d, i.e., an increase in the O-H bond length corresponds to an increase in λGK (see also Supplementary Fig. 17). The mean O-H bond length of graphene/Ni is significantly greater than that of freestanding graphene, whereas the mean O-H bond lengths for freestanding h-BN and h-BN/Ni are relatively similar. The observed increase in the bond length is associated with increased energy vibrations due to the O-H bond, accompanied by more extensive stretching vibrations. As both the stretching and bending motions of the water molecule increase, the translational motion is consequently slowed down, resulting in a higher λGK. Essentially, the atomic-scale friction λGK is thus clearly dependent on both the energy landscape of the surface and the vibrations of the water molecule. Most importantly, considering a dynamic 2D material together with its supporting metal substrate is essential, and our results demonstrate the crucial, yet largely unexplored, role of the supporting substrate in studies of 2D material systems.

Methods
Experiment and sample preparation

We performed all of our experiments on the Cambridge helium-3 spin echo (HeSE) apparatus, which generates a nearly monochromatic polarised 3He beam with an incident energy of 8meV, that is scattered off the sample in a fixed 44.4∘ source-sample-detector geometry. Essentially, the HeSE method is based on the manipulation of the nuclear spin of 3He-atoms in a magnetic field. A detailed description of the setup can be found in refs. 61, 81 The schematic principle of HeSE is shown in Fig. 1a. After passing through a magnetic field, the incident helium beam, which is polarised in the x-direction, is split into two wave packets of opposite nuclear spins, which are temporally separated by the spin-echo time tSE. The scattered wave packets are recombined in a second magnetic field. As a result of the surface motion that occurs during tSE, the two spin components will differ, resulting in a loss of polarisation of the detected beam. An individual HeSE measurement provides the intermediate scattering function I(ΔK,t), with ΔK the momentum transfer vector. In practice, the analysis is performed for the projection of ΔK along a chosen azimuthal direction, such that the ISF is expressed as a function of the scalar ∣ΔK∣. In the case of molecular surface diffusion, as in our study, the ISF follows an exponential decay (see refs. 26, 61, 82 for more information).

The Ni single crystal has been cleaned by multiple cycles of Ar+ sputtering and annealing to 1050K. A single layer of h-BN was grown by chemical vapour deposition (CVD) according to the procedure given by Auwärter et al. 51, during which the Ni(111) surface was maintained at 1050 K and exposed to the gas-phase precursor borazine (B3H6N3) for a few hours. The hot Ni surface acts as a catalyst, initiating the chemical reactions such as the breaking of the borazine rings, dehydrogenation of borazine, and resulting in the subsequent formation of a complete, non-rotated epitaxial overlayer perfectly matching the Ni(111) lattice due to the small lattice mismatch47,83. Additional details can be found in the Sample preparation section in Supplementary Note 1.1.

Water was deposited onto h-BN with a microcapillary array beam doser, which was brought near the surface. To maintain identical experimental conditions throughout each individual measurement, the partial pressure of the water was kept constant by using an automatic leak valve. All dynamic measurements were performed at approximately the same coverage (between 0.12 ML and 0.20 ML, corresponding to an attenuation of the reflectivity by a factor of 4 (see Supplementary Note 1.2). This factor was regularly checked to ensure reproducibility.

Theoretical methods

Spin-polarised electronic structure calculations were carried out using CASTEP84, and the Perdew-Burke-Ernzerhof (PBE) functional85 was used to parameterise the exchange-correlation potential in combination with the Tkatchenko and Scheffler dispersion correction method86. All the calculations were performed using Vanderbilt ultrasoft pseudopotentials87, with a cutoff energy of 400eV for the plane wave basis set. The Ni(111) substrate was modelled using a 5-layer slab with the lowest two layers fixed to represent both the bulk and surface structures. Calculations were performed on a (7×7)R19.1∘ unit cell for the h-BN/Ni(111) and graphene/Ni(111) surfaces and a (3×3) unit cell for freestanding h-BN and graphene monolayers to eliminate spurious interactions with periodic images. A vacuum region of 15 Å was introduced to separate the periodically repeated images and avoid spurious interactions. For all calculations, a Monkhorst-Pack88 (4×4×1) grid was used for k-point sampling. The self-consistent field energy tolerance was set to 1⋅10−7eV. In geometry optimisation using the BFGS minimiser89, structures were relaxed until the maximum force on each atom was less than 0.025eV Å−1. Adsorption energies (Eads) were calculated using the standard formula:

Eads = Exy − Ex − Ey, (5)

where Exy is the energy of the adsorbed species, Ex and Ey are the energy of the dissociated species. Transition states were identified using the linear-quadratic-synchronous transit (LST/QST) algorithm90.

The energy barriers reported in Fig. 3, differing by a few meV from the experimental results, lie at the threshold of precision typically attainable with density functional theory (DFT) using generalised gradient approximation (GGA) functionals, such as the Perdew-Burke-Ernzerhof (PBE) functional employed herein. Such minimal energy differences are susceptible to variations arising from numerical convergence criteria, basis set completeness, and the treatment of dispersion interactions, potentially introducing uncertainties on the order of 10−20meV. Nevertheless, relative energy differences rather than absolute energies, as those between distinct adsorption configurations or diffusion pathways-can often be computed with meV-level accuracy within a consistent DFT framework, provided that stringent convergence parameters are enforced. The close agreement between our calculated energy barriers and the experimental diffusion dynamics observed via HeSE spectroscopy indicates that our DFT approach, augmented with Tkatchenko-Scheffler dispersion corrections, reliably captures the essential π-H bonding and van der Waals interactions governing water physisorption on h-BN/Ni(111) and graphene/Ni(111). While such quantitative agreement is not routinely expected from GGA-based DFT calculations, analogous accuracies have been reported for related systems dominated by weak intermolecular interactions44,54,58,91–95. For the series of AIMD simulations, used to investigate the dynamic behaviour of water molecules at the surface, the velocity Verlet algorithm96, as implemented by CASTEP84 was employed. The system was initialised with the coordinates set to the global minimum position, and initial velocities were randomised using single iteration AIMD simulations within the canonical ensemble (NVT). The AIMD trajectories were extended for a minimum duration of 2ps (to ensure λGK(t) had converged) with a timestep of 1fs, and the temperature converged at ≈150K within the microcanonical ensemble (NVE). To ensure comprehensive analysis, AIMD simulations that extended beyond 2ps were segmented into 2ps intervals. A minimum of 10 AIMD trajectories were generated for each system. (7×7)R19.1° unit cells of h-BN/Ni(111) and graphene/Ni(111) were considered, as well as the (3×3) expanded, free-standing layers of graphene and h-BN. In each timestep of the AIMD trajectories, the centre of mass of the water molecule was determined, and the closest surface site was calculated, accounting for the dynamic movement of the underlying surface atoms. By aggregating data from all AIMD simulations conducted for a specific system, we are able to derive the probability of the water being located at distinct surface sites. For comparison, to evaluate the probability distribution produced using all relaxed geometries of h-BN/Ni(111) during the construction of the PES (Fig. 2), the representing configuration space was sampled at different temperatures according to Boltzmann population analysis (Supplementary Fig. 10). The Supplementary Movies show the full length of the AIMD simulation, with every frame representing a 1fs timestep, played at 120 frames per second.

Peer review rounds
[{"round": 1, "reviewer_comments": {"reviewer_1": "In „Single-Molecule Water Motion on h-BN and Graphene: A paradigm shift in Understanding the Behaviour of Water on 2D Material Interfaces“, P. Seiler et al. use He-spin echo measurements and ab initio calculations to compare the diffusivity of monomeric water on four related surfaces, Ni supported and non-supported graphene and boron nitride and its influence on friction. The article is well-written, and the data is sound and thoroughly analyzed. The advantage of He-spin echo versus another common method of investigating diffusivity, scanning tunnelling microscopy, is the much higher experimental temperature and the statistical data obtained in this involved and challenging technique. The topic is timely, and the results deserve publication in Nature Communications. However, the claim that \"such a behaviour of a single water molecule had not been reported previously\" is incorrect. A conceptionally similar adsorption, including hydrogen bonding in addition to surface interaction with the oxygen atom of water, and a motion, including a combination of rotation and translation, has been reported for water adsorption and diffusion on NaCl on Ag(111) (P. Cabrera-Sanfelix and co-authors, J. Phys. Chem. B 2006, 110, 24559, ACS Nano 9, 3572–3578 (2015)). Also, the influence of the metal support on the diffusivity of a bilayer of NaCl was discussed and compared to bulk NaCl. BN and graphene exist free-standing, in contrast to NaCl, and the previous studies did not relate the single molecule diffusivity to friction. Nonetheless, the previous studies need to be cited, discussed, and related to the present study. In this respect, the title claiming a ‘paradigm shift’ seems exaggerated. Even without this claim, the study remains relevant because He-Spin echo probes a largely different temperature than STM and because of the relevance to ‘true’ 2D material and relation to friction.\nFurther major points to address:\nThe desorption energy was determined based on the Redhead formula to be approximately (0.53 ± 0.02) eV. First, the Redhead formula is so approximate that an error bar of 20 meV seems inadequate. More importantly, the formula holds for 1st order desorption, for which the desorption peak maximum is independent of initial coverage. Water desorbs from a large number of weakly interacting surfaces with zero-order desorption. The authors present their desorption only for one coverage. Please include a desorption experiment for another coverage in the supplementary material to confirm 1st order desorption and, thus, that the method is applicable.\nIsland formation does not exist up to 125 K. This is surprising given the very low diffusion barrier of water monomers on BN, suggesting a classical nucleation-growth scenario according to Venables theory. Is there an additional barrier for water dimer formation or attachment of water to existing clusters? Please comment on this apparent discrepancy.\nP.M. is named in the Author contributions, who conducted the experiments. P.M. is not on the author list. P.M. should be a co-author unless it is a typo.\nMinor points to consider:\nA multi-step motion on BN should alter the Arrhenius prefactor from the single-step motion on graphene. Please comment on it if the data exists. If not, this might be beyond the scope of the present article but is interesting from a fundamental point of view of diffusivity.\nThe AIMD simulation presented in Fig. 5 suggests a directed motion, for which the molecule does not fully equilibrate in its new energetic minimum between single diffusion steps and not a diffusive, i.e., random, motion. It might be an artefact of the used simulation, but should be commented on if all trajectories are directed as the one shown. Or, diffusive trajectories could be added to the supplementary information. If these 'long’ jumps are part of the motion, how is this accounted for in the He-spin-echo analysis? Please clarify this point for the less specialized reader.\nInconsistencies in the Journal naming in the reference list will be corrected for the main text by the Editors. The authors should unify journal abbreviations and capitalization in the supplementary material. Also, ensure a space between the number and the unit in the text of the supplementary material.", "reviewer_2": "The manuscript \"Single-Molecule Water Motion on h-BN and Graphene: A Paradigm Shift in Understanding the Behaviour of Water on 2D Material Interfaces\" by P. Seiler et al. reports about the diffusion of water molecules at sub-monolayer regime on h-BN and graphene. The experimental study is carried out by using the helium-3 spin echo technique and the results show a difference in the behavior of water on the two substrates. For the interpretation of data, the authors performed DFT calculations and ab-initio molecular dynamics simulations for estimating adsorption and activation energies. Their conclusion is: the model in which the hopping of point-like particles between two different lattice sites can reproduce the data on graphene but it is not able to reproduce the data on h-BN for the entire range of momentum exchange probed in the experiment. In particular the authors observe that the trend of experimental points at small momentum exchange is related to a dynamic process involving molecular translation and rotation.\nIn my opinion, the paper is technically sound, the manuscript clearly written with proper details, there are new results, but I have some points that the authors should address.\n1) The authors consider a single sigma to estimate the error bars throughout the manuscript and I think this value is optimistic. Moreover, in figure 3a, the activation energy is estimated in a limited range of temperature (120-135 K) for 3 points along ΓM and two points along ΓK. In SI Figure 4, authors show thermal desorption measurements, water film desorbs around 170 K and seems stable up to 150 K, therefore an extension up to 150 K of the Arrhenius plot looks possible and more points in a extended temperature range would increase the reliability of the estimation even with larger error bars.\n2) In the same figure 3, there is a comparison among calculated energy barriers which differ of a few meV, a value too small to be calculated with the necessary precision (and accuracy) by DFT.\n3) In figure 4, the dephasing rate alpha is plotted for h-BN and graphene. The authors comment about the different trend at small momentum exchange for the two cases (DK < 1 Ang^-1). A careful inspection of the two plots reveals that the h-BN case presents almost a constant value with points between 3 and 5 ns^-1, on the other hand, even in the graphene case the points are not tending to zero for DK going to zero but to a non-vanishing constant value. At this point, based on the experiment, even for graphene a small contribution of the molecular rotation cannot be excluded.\n4) In SI Figure 1 the specular reflection intensity shows a decay which is monotonic and exponential for 120 K, but for 130 K there is a clear initial (apparent) linear decay, followed by an exponential decay. The initial part is related to the reduction of reflectivity since clean surface presents higher reflectivity with respect to regions covered by the adsorbate (producing diffuse scattering) and, during the deposition, the percentage of clean surface is reduced causing a reduction of specular intensity. The adsorbate islands grow until they start to touch and merge and the intensity trend could change. In the present case, the \"linear\" decay extends up to 0.03 ML, a value which is very small compared to other cases present in literature (up to 20% of a ML), moreover, a merging of islands at this low 0.03 ML coverage is not justifiable, even taking into account the effective cross section seen by helium. Therefore some doubts arise about the calibration procedure to estimate the coverage. Maybe the measurement of the pressure in a specific position of the chamber is not real pressure at the sample. The authors should try to estimate a reasonable error bar (or confidence interval) on coverage and, maybe, increase the 0.16 ML value.\n5) In figure 6a, friction coefficient λ(t) are reported, my concern is about the graphene/Ni(111) case that look like its value has not yet reached its asymptotic value at 2 ps.\nMinor correction:\nMain text:\n1) in \"Walking motion of Water on h-BN\": \" in canonical simulations ((NVT)\" remove double parenthesis\n2) In \"Theoretical methods\" the unit cell notation \"(√7 × √7)R19:1°\" has to be corrected to \"(√7×√7)R19.1\"\n3) the meaning of \" periodic images were eliminated. \" is not very clear.\nIn SI\n1) 1.2 Coverage calibration, line 39 \"In Supplementary Supplementary Figure 1\", remove repetition\n2) page 2, line 76: \"ν is the preexponential or frequency factor,\" ν is missing in the equation\n3) page 4, line 165 & line 178 \"hBN\"-> \"h-BN\"\nIn conclusion, the manuscript require changes before to be accepted for publication.", "reviewer_3": "The authors combine helium spin‑echo spectroscopy and ab initio calculations to show that single water molecules on h‑BN/Ni(111) diffuse via a coupled rotational‑translational “walking” mechanism—with an activation energy roughly 2.5× lower than the ∼60 meV barrier on graphene/Ni—despite exhibiting similar binding energies on both substrates. They further demonstrate that metal support dramatically alters friction, with h‑BN/Ni displaying markedly lower water friction than graphene/Ni, challenging classical jump‑between‑site models and underscoring the critical role of substrate polarity in designing 2D‑material‑based microfluidic systems.\nThe work is very interesting, it addresses the lack of data for the H2O on h-BN/Ni system, and it shows marked differences between h-BN and its graphene counterpart case. It is shown that despite the similarity of the substrate and H2O adsorption energy on h-BN/Ni and graphene/Ni, both experiments and theory show completely different motions of single water molecules on the two surfaces. However, there are a couple comments that need to be addressed before publication.\n1. What would be the role of entropy in the energetics described in this work? Is entropy accounted for?\n2. DFT methods used here are known to overestimate adsorption energies, especially compared to coupled‑cluster/QMC methods. How would this overestimation affect calculated diffusion barriers?\n3. The calculated vibrational spectra for the OH stretching is in the order of 3900 cm-1. This is much higher than normally found for experimental OH stretches (typically below 3750 cm-1). Why is this calculated frequency much higher than experimentally found?\n4. In the following sentence in page 2, “A previously reported HeSE study of water motion on graphene illustrated the stark contrast between water diffusion on 2D materials such as graphene and metal substrates, where water is bound much stronger and further influenced by hydrogen bonding between the molecules.” It is not clear in the sentence to which case it refers that water is bound much stronger. Does this refer to the latter (metal surface)?", "reviewer_4": "This is a most interesting paper. It combines state-of-the-art experimentation and simulation to probe the nature of the surface water interactions on h-BN and graphene. The present referee is not an expert in the experimental techniques but they appear to have carefully undertaken. The modelling work is of high quality using state of the art techniques and is clearly reported and reproducible. The results are particularly interesting, especially in the contract between the diffusion mechanisms on the two systems. The proposed mechanism on h-BN is of special interest although I think it could be more clearly explained. The paper is well written and will be of general interest. Publication in Nat Comms is recommended but I suggest two modifications:\n1. Can we have a clearer account of the rotational- translational dynamic motion of water, with perhaps a good figure.\n2. Can we remove \"paradigm shift\" from the title. This is very interesting and novel work but it is not a paradigm shift in understanding water at 2D interfaces"}, "author_response": "We thank the referee for their constructive suggestions and for recognising the timeliness and significance of our work. We are pleased that the results were found to be suitable for publication in Nature Communications, and we have addressed all comments in detail below.\nMajor Points\n1) “However, the claim that \"such a behaviour of a single water molecule had not been reported previously\" is incorrect. A conceptionally similar adsorption, including hydrogen bonding in addition to surface interaction with the oxygen atom of water, and a motion, including a combination of rotation and translation, has been reported for water adsorption and diffusion on NaCl on Ag(111) (P. Cabrera-Sanfelix and co-authors, J. Phys. Chem. B 2006, 110, 24559, ACS Nano 9, 3572–3578 (2015)). Also, the influence of the metal support on the diffusivity of a bilayer of NaCl was discussed and compared to bulk NaCl. BN and graphene exist free-standing, in contrast to NaCl, and the previous studies did not relate the single molecule diffusivity to friction. Nonetheless, the previous studies need to be cited, discussed, and related to the present study. In this respect, the title claiming a ‘paradigm shift’ seems exaggerated. Even without this claim, the study remains relevant because He-Spin echo probes a largely different temperature than STM and because of the relevance to ‘true’ 2D material and relation to friction.”\nAuthor reply:\nWe thank the reviewer for the constructive comment and for drawing attention to the relevant prior work on single water molecule dynamics on NaCl surfaces. In response, we have cited these studies as Ref. [25] and [61] in the revised manuscript and have incorporated a discussion of their relevance, as detailed below.\nCabrera-Sanfelix et al. (J. Phys. Chem. B 2006, 110, 24559) examined the adsorption and diffusion of a single water molecule on NaCl(001) using density functional theory (DFT). Their results revealed a concerted translation–rotation diffusion mechanism, facilitated by electrostatic interactions with Na⁺ and hydrogen bonding with Cl⁻. While the quasi-continuous motion they describe is conceptually analogous to our findings, the underlying substrate is ionic rather than covalent or weakly polar, as in the h-BN and graphene surfaces studied here. Moreover, this study does not address the role of friction or metal support effects, which are central to our investigation.\nIn ACS Nano 2015, 9, 3572–3578, Cabrera-Sanfelix and co-authors employed scanning tunnelling microscopy (STM) and DFT to study water dynamics on a NaCl bilayer supported on Ag(111). This work demonstrates that the metal substrate modifies adsorption energies and diffusion pathways, an effect that aligns qualitatively with our finding of Ni(111)-induced modulation of potential energy surface corrugation. However, the focus remains on STM-accessible dynamics at lower temperatures (42–52 K), whereas our work probes single-molecule dynamics at elevated temperatures (120–140 K) via helium spin-echo (HeSE) spectroscopy.\nOur study extends beyond these prior contributions in several critical aspects:\n1. We employ HeSE to investigate single-molecule water dynamics on h-BN and graphene above 100 K, thus accessing a higher temperature regime and a different time scale (picoseconds) than STM.\n2. We quantify atomic-scale friction via the Green-Kubo formalism, a method not utilised in the NaCl studies.\n3. We demonstrate a substrate-induced inversion of frictional response between h-BN and graphene when supported on Ni(111), revealing a previously unreported metal–2D material coupling effect.\n4. Our system is based on covalently bonded 2D materials rather than ionic substrates, introducing fundamentally different interfacial interactions and electronic environments.\nTo reflect these distinctions, we have revised the title and updated the introduction and discussion sections to acknowledge and appropriately contextualise the prior NaCl studies.\nTitle Revision\nOriginal Title: “Single-Molecule Water Motion on h-BN and Graphene: A Paradigm Shift in Understanding the Behaviour of Water on 2D Material Interfaces”\nRevised Title: “Single-Molecule Water Motion on h-BN and Graphene: Understanding the Behaviour of Water on 2D Material Interfaces”\nRevised text\nPage 1 (Introduction):\nPage 5 (Motion):\nPage 7 (Friction and role of the substrate):\nFurther major points to address:\nThe desorption energy was determined based on the Redhead formula to be approximately (0.53 ± 0.02) eV. First, the Redhead formula is so approximate that an error bar of 20 meV seems inadequate. More importantly, the formula holds for 1st order desorption, for which the desorption peak maximum is independent of initial coverage. Water desorbs from a large number of weakly interacting surfaces with zero-order desorption. The authors present their desorption only for one coverage. Please include a desorption experiment for another coverage in the supplementary material to confirm 1st order desorption and, thus, that the method is applicable.\nAuthor reply:\nWe appreciate the referee’s insightful comment concerning the desorption kinetics and the limitations of the Redhead analysis. We agree that the Redhead method assumes first-order desorption and that its precision is inherently limited, particularly in the context of water desorption on weakly interacting surfaces such as hexagonal boron nitride (h-BN). The main goal of including this analysis was to provide a reference point for the desorption temperature regime, as this represents, to our knowledge, the first such experimental data reported for water on h-BN.\nAlthough the desorption kinetics are not the central focus of our study, and our main conclusions regarding diffusion are unaffected by this discussion, we acknowledge the importance of clarifying the desorption order. Unfortunately, we do not have additional thermal desorption spectra for different coverages. Nevertheless, we have experimentally confirmed the desorption temperature regime through the recovery of the helium reflectivity signal.\nTo address the referee’s concern, we have substantially revised the supplementary section. Specifically, we now discuss the potential for zero-order desorption kinetics, which are more commonly observed for water on weakly interacting surfaces such as graphene, due to the formation of hydrogen-bonded water clusters – and calculate also the energy based on that kinetics. We also increase the uncertainty in the derived desorption energy to reflect the approximations inherent in both the first- and zero-order analyses.\nThe reported value in the main text has been corrected; for the added paragraphs see the revised version of the supplementary information – on page 2.\nIsland formation does not exist up to 125 K. This is surprising given the very low diffusion barrier of water monomers on BN, suggesting a classical nucleation-growth scenario according to Venables theory. Is there an additional barrier for water dimer formation or attachment of water to existing clusters? Please comment on this apparent discrepancy.\nAuthor reply:\nWe clarify that our manuscript does not claim an absence of island formation below 125 K and we have also done additional calculations for dimer formation and diffusion, see below.\nAdded on page 2:\nIn response, we have now expanded the discussion of dimer diffusion and dissociation in a newly added Section 2.2 of the Supplementary Information. As shown in Supplementary Figures 7 and 8, water dimers on both freestanding and Ni-supported h-BN exhibit significantly higher energy barriers for diffusion (∼80–100 meV) and dissociation (220 meV) compared to isolated water molecules. These elevated barriers indicate that dimers are effectively immobile under the experimental conditions and thus do not contribute to the observed dynamics. The comparison confirms that the experimentally measured activation energies are best represented by single-molecule behaviour.\nWe have further added the following text on page 4 of the main manuscript to point the reader to the SI:\nP.M. is named in the Author contributions, who conducted the experiments. P.M. is not on the author list. P.M. should be a co-author unless it is a typo.\nAuthor reply:\nWe would like to apologise for this typo. P. M. is Philipp Maier which at the time of writing was his birth surname - after getting married this has now changed to Philipp Seiler. We have corrected this accordingly in the author's contribution.\nMinor points:\n1) “A multi-step motion on BN should alter the Arrhenius prefactor from the single-step motion on graphene. Please comment on it if the data exists. If not, this might be beyond the scope of the present article but is interesting from a fundamental point of view of diffusivity.\nAuthor reply:\nWe appreciate the reviewer’s observation concerning the possible impact of multi-step motion on the Arrhenius prefactor. The experimentally determined dephasing rate α(ΔK) in helium spin-echo measurements is a reciprocal-space observable that integrates over all dynamical processes in real space. Consequently, it does not allow for a direct decomposition into translational, rotational, or multi-step components. In this framework, the extracted Arrhenius parameters represent effective quantities that encapsulate the ensemble-averaged dynamics projected onto a single timescale.\nWe are, of course, mindful that within the framework of transition state theory, a multi-step diffusion mechanism can lead to an enhanced entropy of the transition state, thereby increasing the Arrhenius prefactor. This arises due to the greater number of accessible microstates and vibrational modes contributing to the rate-limiting barrier crossing. Conversely, enhanced frictional coupling to the substrate or internal degrees of freedom may act to suppress the effective attempt frequency. Furthermore, for strongly corrugated potential energy surfaces (PES), the classical picture of single-barrier crossing becomes less applicable, and the prefactor becomes highly sensitive to the detailed nature of the PES and the coupling to dissipative environments.\nGiven these considerations, we note that a more detailed analysis of prefactor variations arising from multi-step dynamics lies beyond the scope of the present study, whose primary focus is on quantifying the effective activation energy and the overall diffusive behaviour – and we have added the following sentence to the manuscript on page 4:\n2) “The AIMD simulation presented in Fig. 5 suggests a directed motion, for which the molecule does not fully equilibrate in its new energetic minimum between single diffusion steps and not a diffusive, i.e., random, motion. It might be an artefact of the used simulation, but should be commented on if all trajectories are directed as the one shown. Or, diffusive trajectories could be added to the supplementary information. If these 'long’ jumps are part of the motion, how is this accounted for in the He-spin-echo analysis? Please clarify this point for the less specialized reader.”\nAuthor reply:\nWe thank the reviewer for the careful observation. The trajectory presented in Figure 5 is intended as a representative illustration of water mobility on h-BN/Ni(111), rather than as an exemplar of universal behaviour across all simulations. While the specific trajectory shown may appear directed, examining the motion in detail shows it is random, consistent with thermally activated diffusion processes.\nTo facilitate a clearer understanding of the dynamic behaviour, we have provided visualisations of several ab initio molecular dynamics (AIMD) trajectories in the form of supplementary videos. These include three additional trajectories that accompany the revised manuscript and are referenced explicitly in the updated text. The trajectory shown in Figure 5 is designated as Video 1.\nAll videos (10 in total) are also permanently archived in the data repository cited under the \"Data Availability\" section, enabling full access to the underlying motion data.\nWe further note that the mentioned phenomenon is consistent with finite-temperature ab initio molecular dynamics at 150 K: Thermal energy permits the molecule to explore a range of configurations around local minima without necessarily achieving full equilibration within each potential well prior to a subsequent transition. This behaviour stands in contrast to the situation at approximately 0 K, where atomic positions relax strictly to the nearest local minimum in a deterministic manner, governed solely by the potential energy landscape and not by thermal fluctuations.\nThis form of motion is expected for thermally activated dynamics where the thermal energy of a system at constant T is comparable to the energy barrier and has been well documented in both theoretical and experimental studies1,2.\n3) “Inconsistencies in the Journal naming in the reference list will be corrected for the main text by the Editors. The authors should unify journal abbreviations and capitalization in the supplementary material. Also, ensure a space between the number and the unit in the text of the supplementary material.”\nAuthor reply:\nWe have checked and updated the supplementary references, and we fixed the missing space between number and text.\n1 https://doi.org/10.1103/PhysRevLett.96.026101\n2 https://pubs.acs.org/doi/full/10.1021/jacs.5b08001\nReviewer #2 (Remarks to the Author):\nThe manuscript \"Single-Molecule Water Motion on h-BN and Graphene: A Paradigm Shift in Understanding the Behaviour of Water on 2D Material Interfaces\" by P. Seiler et al. reports about the diffusion of water molecules at sub-monolayer regime on h-BN and graphene. The experimental study is carried out by using the helium-3 spin echo technique and the results show a difference in the behavior of water on the two substrates. For the interpretation of data, the authors performed DFT calculations and ab-initio molecular dynamics simulations for estimating adsorption and activation energies. Their conclusion is: the model in which the hopping of point-like particles between two different lattice sites can reproduce the data on graphene but it is not able to reproduce the data on h-BN for the entire range of momentum exchange probed in the experiment. In particular the authors observe that the trend of experimental points at small momentum exchange is related to a dynamic process involving molecular translation and rotation.\nIn my opinion, the paper is technically sound, the manuscript clearly written with proper details, there are new results, but I have some points that the authors should address.\nAuthor reply:\nWe thank the referee for the positive assessment of our work and for acknowledging the technical soundness, clarity, and novelty of our results. We appreciate the thoughtful comments and suggestions, which have helped us to improve the manuscript. Below, we address each point raised in detail.\n1) The authors consider a single sigma to estimate the error bars throughout the manuscript and I think this value is optimistic. Moreover, in figure 3a, the activation energy is estimated in a limited range of temperature (120-135 K) for 3 points along ΓM and two points along ΓK. In SI Figure 4, authors show thermal desorption measurements, water film desorbs around 170 K and seems stable up to 150 K, therefore an extension up to 150 K of the Arrhenius plot looks possible and more points in a extended temperature range would increase the reliability of the estimation even with larger error bars.\nAuthor reply:\nWe thank the referee for raising this important point regarding the uncertainty estimates and the limited temperature range used for the Arrhenius analysis. We agree that the determination of the activation energy from a narrow temperature interval can introduce sensitivity to statistical and systematic errors.\nWhile we did conduct an additional measurement at 140 K, we note that this temperature lies near the onset of water desorption, as confirmed by the thermal desorption profile presented in SI Figure 4. Under these conditions, maintaining a constant surface coverage becomes increasingly difficult without applying a continuous water overpressure, which is experimentally challenging and introduces further uncertainties.\nThe temperature range employed in our analysis (120–135 K) was selected to ensure that diffusion dominates over desorption, thereby maintaining well-defined coverage conditions. It is worth noting that similar or even narrower temperature windows are frequently used in scanning tunnelling microscopy (STM) studies of water diffusion, due to desorption constraints, but also to achieve the necessary temporal resolution.\nWe have already combined two independent datasets to maximise the temperature coverage and statistical reliability of the extracted activation energy. Given the limitations described above, we believe this range represents the most robust basis for Arrhenius analysis under the current experimental conditions. Nevertheless, to reflect the acknowledged uncertainties more conservatively, we have increased the reported error bars on the activation energy estimate in the revised manuscript, and we have also added 2 sentences on page 4:\n2) In the same figure 3, there is a comparison among calculated energy barriers which differ of a few meV, a value too small to be calculated with the necessary precision (and accuracy) by DFT.\nAuthor reply:\nWe appreciate the reviewer’s observation regarding the limitations in absolute energy accuracy inherent to density functional theory (DFT), particularly when using generalised gradient approximation (GGA) functionals such as the Perdew-Burke-Ernzerhof (PBE) functional employed in our study. While absolute quantities - such as adsorption energies - are often subject to systematic errors and uncertainties arising from numerical convergence, basis set limitations, and dispersion treatment, relative energy differences between local minima and transition states within the same DFT framework can be determined with a precision of a few meV, provided that computational parameters are carefully converged. This distinction is particularly relevant when comparing diffusion barriers or competing adsorption configurations within a unified computational approach and the close agreement between our calculated energy barriers and the experimental values obtained from helium spin-echo (HeSE) measurements suggests that our DFT calculations effectively capture the key physical interactions governing water molecule dynamics on h-BN/Ni(111) and graphene/Ni(111).\nIn this context, the close agreement between our calculated energy barriers and the experimental values obtained via helium spin-echo (HeSE) spectroscopy suggests that our DFT calculations accurately capture the essential physics underlying water molecule dynamics on h-BN/Ni(111) and graphene/Ni(111). The governing interactions-namely π–H bonding and weak van der Waals forces-are effectively described by the GGA functional in conjunction with the Tkatchenko-Scheffler dispersion correction. Furthermore, to substantiate the reliability of our methodology, we have cited relevant benchmark studies that demonstrate similar levels of accuracy for weakly interacting surface systems [87,90]. Accordingly, while absolute agreement within a few meV is not generally expected, the consistency observed in our study reflects the robustness of our approach in capturing the relative energetics of water physisorption and diffusion on these two-dimensional substrates as also shown in Refs [43,53,56,88,89,91].\nWe have also added the following paragraph in the manuscript:\n3) In figure 4, the dephasing rate alpha is plotted for h-BN and graphene. The authors comment about the different trend at small momentum exchange for the two cases (DK < 1 Ang^-1). A careful inspection of the two plots reveals that the h-BN case presents almost a constant value with points between 3 and 5 ns^-1, on the other hand, even in the graphene case the points are not tending to zero for DK going to zero but to a non-vanishing constant value. At this point, based on the experiment, even for graphene a small contribution of the molecular rotation cannot be excluded.\nAuthor reply:\nWe thank the reviewer for this careful inspection and insightful observation. Indeed, as the reviewer points out, the dephasing rate α for graphene does not vanish entirely in the low-ΔK limit. While we cannot categorically exclude a minor contribution from rotational degrees of freedom for H O adsorbed on graphene, as discussed in Ref [43] and Supplementary material thereof, the constant offset is mostly an artifact of the correction upon extracting the self-diffusion rather than the collective diffusion on graphene where the latter has a strong signature of repulsive interactions. We have also added the following sentence in the manuscript on page 5:\n4) In SI Figure 1 the specular reflection intensity shows a decay which is monotonic and exponential for 120 K, but for 130 K there is a clear initial (apparent) linear decay, followed by an exponential decay. The initial part is related to the reduction of reflectivity since clean surface presents higher reflectivity with respect to regions covered by the adsorbate (producing diffuse scattering) and, during the deposition, the percentage of clean surface is reduced causing a reduction of specular intensity. The adsorbate islands grow until they start to touch and merge and the intensity trend could change. In the present case, the \"linear\" decay extends up to 0.03 ML, a value which is very small compared to other cases present in literature (up to 20% of a ML), moreover, a merging of islands at this low 0.03 ML coverage is not justifiable, even taking into account the effective cross section seen by helium. Therefore some doubts arise about the calibration procedure to estimate the coverage. Maybe the measurement of the pressure in a specific position of the chamber is not real pressure at the sample. The authors should try to estimate a reasonable error bar (or confidence interval) on coverage and, maybe, increase the 0.16 ML value.\nAuthor reply:\nWe appreciate the reviewer’s detailed analysis and agree that pressure- and coverage-calibration procedures are intrinsically challenging. We follow the methodology outlined in Ref [43] and its SI and have added additional clarifying statements in the current SI. We further note that – as described in the SI of Ref [43] for graphene/Ni we had two independent ways of calibrating the coverage as also added in the revised version of the current SI. We have also added an extensive explanation of the uptake curve behaviour in the SI (starting from the last paragraph on page 1) as indicated by the referee.\nIn addition to these challenges, to ensure reproducibility and comparable experimental conditions across individual ISF measurements, we performed the dynamical measurements at approximately equivalent levels of specular intensity attenuation I/I and despite applying consistent procedures, some variation in this attenuation is unavoidable; for the 120 K dataset, for example, values range between 0.27 and 0.35.\nAccordingly, we have expanded the boundaries of our estimated coverage, which now correspond to a range of approximately 0.12 to 0.20 monolayers.\n5) In figure 6a, friction coefficient λ(t) are reported, my concern is about the graphene/Ni(111) case that look like its value has not yet reached its asymptotic value at 2 ps.\nAuthor reply:\nWe thank the reviewer for this comment. We acknowledge that the friction coefficient λ(t) for graphene/Ni(111) exhibits a slight increase around 2ps. Our choice of 2 ps as a cutoff was guided by the work of Tocci et. al. (reference 29), which demonstrated that this duration achieves sufficiently converged friction coefficients for reliable estimation. To ensure convergence, we extended our analysis by running additional AIMD simulations with identical parameters and a minimum duration of 3 ps. Here, we observed the friction coefficient rises from ~9.5×107 Nsm−3 at 2 ps to ~9.7×107 Nsm−3 at 3 ps (an increase of 2 %). The gradient of λ(t) remains marginally positive but exhibits significant noise. We do not expect the minimal changes in friction coefficient λ(t) to alter our conclusions. Please see section 2.6 in the SI together with Supplementary Figure 13 for the corresponding plot.\nWe have further added a single error bar to Figure 6(a) and the corresponding reference to the supplement in the main text on page 6:\nMinor correction:\nMain text:\n1) in \"Walking motion of Water on h-BN\": \" in canonical simulations ((NVT)\" remove double parenthesis\n2) In \"Theoretical methods\" the unit cell notation \"(√7 × √7)R19:1°\" has to be corrected to \"(√7×√7)R19.1\"\n3) the meaning of \" periodic images were eliminated. \" is not very clear.\nIn SI\n1) 1.2 Coverage calibration, line 39 \"In Supplementary Supplementary Figure 1\", remove repetition\n2) page 2, line 76: \"ν is the preexponential or frequency factor,\" ν is missing in the equation\n3) page 4, line 165 & line 178 \"hBN\"-> \"h-BN\"\nIn conclusion, the manuscript require changes before to be accepted for publication.\nAuthor reply:\nHave all been corrected accordingly\nReviewer #3 (Remarks to the Author):\nThe authors combine helium spin-echo spectroscopy and ab initio calculations to show that single water molecules on h-BN/Ni(111) diffuse via a coupled rotational-translational “walking” mechanism—with an activation energy roughly 2.5× lower than the ∼60 meV barrier on graphene/Ni—despite exhibiting similar binding energies on both substrates. They further demonstrate that metal support dramatically alters friction, with h-BN/Ni displaying markedly lower water friction than graphene/Ni, challenging classical jump-between-site models and underscoring the critical role of substrate polarity in designing 2D-material-based microfluidic systems.\nThe work is very interesting, it addresses the lack of data for the H2O on h-BN/Ni system, and it shows marked differences between h-BN and its graphene counterpart case. It is shown that despite the similarity of the substrate and H2O adsorption energy on h-BN/Ni and graphene/Ni, both experiments and theory show completely different motions of single water molecules on the two surfaces. However, there are a couple comments that need to be addressed before publication.\nAuthor reply:\nWe thank the referee for their positive evaluation and for recognising the significance of our findings regarding the contrasting behaviour of water on h-BN/Ni and graphene/Ni. We appreciate the thoughtful comments and have carefully addressed the points raised below.\n1. What would be the role of entropy in the energetics described in this work? Is entropy accounted for?\nAuthor reply:\nWe thank the reviewer for raising this important point. Our study concentrates on the dynamics of isolated water molecules adsorbed on h-BN and graphene/Ni(111) surfaces, in the absence of chemical reactions or cooperative processes such as self-assembly. This stands in contrast to systems where configurational entropy or collective interactions play a dominant role in modulating energetics.3,4 For single-molecule motion in such weakly bound systems, the relative contributions from vibrational, translational, and rotational entropy are expected to vary minimally across different adsorption sites and during lateral diffusion. Consequently, the entropic term TΔS is anticipated to be small relative to the enthalpic component ΔH and unlikely to meaningfully affect the observed energetic hierarchy [see Ref 55].\nWe also acknowledge the computational constraints inherent in accurately estimating ΔS along the entire simulation trajectory. Methods such as thermodynamic integration, vibrational mode analysis, or harmonic approximations - while valuable - are computationally expensive, especially when applied to systems with high configurational dimensionality or when aiming for spatially resolved entropic mapping. For these reasons, and in alignment with common practice in AIMD studies of similar scope, we have not included explicit entropic corrections at each simulation point.\nNonetheless, we emphasize that our ab initio molecular dynamics (AIMD) simulations naturally incorporate thermal fluctuations and permit access to an effective free energy surface under experimental temperature conditions. This approach captures entropic effects implicitly to the extent that they modulate the system’s thermal behaviour over accessible time and length scales.\nIn future work, it would be interesting to explore entropic contributions more explicitly through thermodynamic integration and harmonic approximation for selected configurations to assess their quantitative influence on site preference and diffusion barriers.[55]\nWe have therefore added the following sentence in the conclusion/outlook section on page 8:\n2. DFT methods used here are known to overestimate adsorption energies, especially compared to coupled-cluster/QMC methods. How would this overestimation affect calculated diffusion barriers?\nAuthor reply:\nWe appreciate the reviewer’s insightful comment. It is indeed well known that standard DFT functionals, including those used in this study, can overestimate absolute adsorption energies relative to high-level wavefunction-based methods such as coupled-cluster or quantum Monte Carlo (QMC). However, these methods remain largely impractical for extended periodic systems involving surfaces due to their high computational cost and scaling with system size, especially when modelling surface-adsorbate interactions with full geometry relaxation and performing extensive AIMS simulations.\nIn our case, the diffusion barriers of interest are derived from relative energy differences along the potential energy surface (PES), rather than absolute adsorption energies. These relative barriers are typically less sensitive to systematic errors in the underlying exchange-correlation functional, especially when the nature of bonding and local electronic structure does not vary drastically along the diffusion path. This is supported by the excellent agreement we observe between our calculated diffusion barriers and the experimental values obtained via HeSE, as further demonstrated by previous work on weakly adsorbed systems (Refs [43,53,56,88,89,91] in the revised manuscript, plus we have also added some references to benchmarking / accuracy for weakly adsorbed systems with Refs [87,90] in the revised manuscript).\nThus, while we acknowledge that DFT-based methods may not yield quantitatively exact adsorption energies, we believe the relative energetics governing diffusion are sufficiently accurate for the present study.\n3. The calculated vibrational spectra for the OH stretching is in the order of 3900 cm-1. This is much higher than normally found for experimental OH stretches (typically below 3750 cm-1). Why is this calculated frequency much higher than experimentally found?\nAuthor reply:\nWe thank the reviewer for this pertinent observation. We have now clarified this point explicitly in the revised manuscript and the Supplementary Information (Section 2.6, see Refs. [64,65]). As shown in Fig. 6(b–c) of the main text and detailed in the Supplementary Figure S6, the OH stretching modes appear near 3900 cm-1, consistent across both interfacial and isolated water. To verify the origin of this elevated frequency, we computed the vibrational power spectrum of an isolated water molecule in vacuum. The stretching frequency remains centred near 3900 cm-1, indicating that the observed redshift is not primarily attributable to molecule–surface interactions but rather reflects known limitations of DFT-based harmonic vibrational analyses. Specifically, the overestimation of OH stretching frequencies by 8–12% is a recognised issue arising from the neglect of anharmonic effects and the limitations of commonly used exchange–correlation functionals. Our findings align with benchmark studies – mentioned as Ref 64 and 65 in the revised manuscript – i.e. Medel and Suhm (Phys. Chem. Chem. Phys., 2021, 23, 5629–5643) and Varandas (J. Chem. Phys., 2022, 157, 174110), which report comparable overestimations of OH stretching frequencies and advocate for the application of empirical scaling factors to improve agreement with experimental data.\nImportantly, our analysis centres on relative frequency shifts - such as the comparison between water on h-BN/Ni(111) and on free-standing h-BN - where systematic errors are expected to largely cancel. This renders the interpretation of substrate-induced effects more reliable, despite the inherent limitations in absolute frequency prediction.\nIn the main text we refer now to the SI via the following text on page 7:\nWhile in the SI in 2.7 we have added a new paragraph which explains the mentioned redshift and the recognised issue from DFT in calculating absolute frequencies.\n4. In the following sentence in page 2, “A previously reported HeSE study of water motion on graphene illustrated the stark contrast between water diffusion on 2D materials such as graphene and metal substrates, where water is bound much stronger and further influenced by hydrogen bonding between the molecules.” It is not clear in the sentence to which case it refers that water is bound much stronger. Does this refer to the latter (metal surface)?\nAuthor reply:\nWe have changed the text to: … weakly interacting 2D materials such as graphene and metal substrates, where water is bound much stronger and further influenced by hydrogen bonding between the molecules. 23,36,43\nReviewer #4 (Remarks to the Author):\nThis is a most interesting paper. It combines state-of-the-art experimentation and simulation to probe the nature of the surface water interactions on h-BN and graphene. The present referee is not an expert in the experimental techniques but they appear to have carefully undertaken. The modelling work is of high quality using state of the art techniques and is clearly reported and reproducible. The results are particularly interesting, especially in the contract between the diffusion mechanisms on the two systems. The proposed mechanism on h-BN is of special interest although I think it could be more clearly explained. The paper is well written and will be of general interest. Publication in Nat Comms is recommended but I suggest two modifications:\nAuthor reply:\nWe thank the referee for highlighting the broad relevance of our work and the interest in the proposed mechanism on h-BN. We acknowledge the positive assessment of the quality of our modelling and the recommendation for publication. Based on the constructive suggestions we have revised the manuscript accordingly to improve clarity and address the specific points raised.\n1. Can we have a clearer account of the rotational- translational dynamic motion of water, with perhaps a good figure.\nAuthor reply:\nDue to the number of trajectories and the complexity of the motion we have produced videos of each of the H O/h-BN/Ni simulations. See also our reply to Ref #1, Comment 2:\n2. Can we remove \"paradigm shift\" from the title. This is very interesting and novel work but it is not a paradigm shift in understanding water at 2D interfaces\nAuthor reply:\nWe agree and since this point was also raised by referee 1 we have now removed “paradigm shift” from the title and instead propose a revised title.\nRevised Title: “Single-Molecule Water Motion on h-BN and Graphene: Understanding the Behaviour of Water on 2D Material Interfaces”"}, {"round": 2, "reviewer_comments": {"reviewer_1": "The author has thoughtfully addressed all my comments, responding to them in a satisfactory manner and revising the manuscript appropriately. I fully endorse the acceptance of this article for publication in Nature Communications in its current form.sincerely considered all my remarks. They have answered them to my full satisfaction and altered the manuscript accordingly.\nI fully support the acceptance of this article for publication in Nature Communications in the present form.", "reviewer_2": "The revised version of the article with changed title \"Single-Molecule Water Motion on h-BN and Graphene: Understanding the Behaviour of Water on 2D Material Interfaces\" now addresses all the points I raised in my previous review and I am satisfied with the authors' changes. In my opinion, this revised version can be accepted for publication.", "reviewer_4": "The authors have responded positively and well to my comments and those of the other referees.The additional information , especially on mechanistic aspects is helpful and the new title is more appropriate.\nI recommend that the paper be accepted for publication"}, "author_response": ""}]

Decision letter