Deep vacancy induced low-density fluxional interfacial water

Interfacial water on transition metal oxides such as TiO$_2$ has been widely studied because of their structural complexity and scientific relevance in e.g. photocatalysis and ice growth. Using ab initio molecular dynamics, we find that interfacial water on anatase (101) surface features an unconventional fluxional structure with reduced contact layer density. The density reduction and flexibility of interfacial water are induced by oxygen vacancy defects deeply located in the subsurface. Our study proposes a fresh perspective of the anatase/water interface, raising the importance of non-trivial long-range effects caused by deep defects. These often-neglected effects highlight the necessity and challenges of the state-of-the-art simulation and experimental probing of solid liquid interfaces.

tional theory (DFT) with the Vienna ab initio simulation package(VASP) [32]. The Van der Waals inclusive optB86-vdW exchange correlation functional was employed [33,34]. Additional simulations with SCAN-rvv10 exchange correlation functional [35,36] and with Hubbard corrections (U=3.9 eV) [37] were performed to test the validity of conclusions. Electron-ion interactions were described using the projector augmented wave method, with Ti(3d, 4s), O(2s, 2p), and H(1s) electrons treated explicitly as valence electrons [38]. Wave functions were expanded in plane-waves up to a kinetic energy cutoff of 500 eV. K-point sampling was restricted to the Γ point only. The deuterium mass (2 amu) was used for hydrogen to allow for a larger time step of 1 fs. Temperature was controlled using Nosé-Hoover thermostat at 330 K within NVT ensemble. Ten different configurations were considered for anatase(101), including a stoichiometric substrate and nine substrates each containing one oxygen vacancy (V O ) located at different depths. The oxygen vacancy concentration in our supercell is 1 120 , which is within the typical range in reduced TiO 2 . Taking the stoichiometric configuration as an example, a periodically repeated 1 × 3 supercell of size 10.21 × 11.328 × 31.58 was used in our simulation, corresponding to a composition of ten Ti-layers and twenty O-layers (Ti 60 O 120 ) and a (101) surface area of 10.21 × 11.328 Å 2 . Water slabs including 52 water molecules were constructed with the bulk water density (1.02 g/cm 3 ).
After reaching equilibrium in e.g. 10 ps, trajectories were collected for up to 120 ps. The pressure is within 50 MPa throughout the simulations. We carry out additional simulations with 15 Å vacuum above 76 water molecules and a substrate with 36 Ti atoms and 72/71 oxygen atoms senting two ordered layers of water. This structure of water over-layer is also reproduced in many previous simulations [39][40][41][42]. For liquid water-anatase interface at room temperature, however, the structure is not clear [43]. [101] [101] [101] [010] [101] [010] [101]   Prior to a systematic discussion of the impact of defects on the interfacial water structure, we first discuss another typical interfacial water structure identified in our simulations. In this structure ( Fig.1c-d), water molecules on the O 2c row mostly form no more than one hydrogen bond with O 2c , in contrast to two hydrogen bonds in the high-density model. In addition, water molecules on Ti 5c may directly form hydrogen bond with the O 2c rather than with adsorbed water molecules (orange molecules in Fig.1a,c). In this unconventional structure, fewer water molecules feature in the contact layer, thus we refer it as the low-density interfacial water [45] Fig. 1e further shows that the low-density interfacial water has fewer hydrogen bonds within the contact layer, and fewer hydrogen bonds between water and the surface. Such differences are also shown in the radial distribution function of water hydrogen (H W ) around O 2c (Fig.1f). In the high-density structure, more H atoms prefer to stay bonded to O 2c s, while the density of H atoms is obviously lower in the low-density model.
Another measure to distinguish the low-density and the high-density interfacial water is the g(z) profile. Fig.2b,c shows all g(z) profiles of ten different simulations, including a pristine substrate and nine substrates with oxygen vacancy, which are labeled from V O 1 to V O 9 according to the depth of the oxygen vacancy. The relative stability of these 9 vacancy sites is yet to be established, therefore as a first step it is important to understand how each vacancy site impacts the interfacial water structure. In reality, vacancies may migrate across the substrates and have a dynamical influence on the interfacial water structure. Fig.2a shows a sideview of the pristine model, and the colored circles highlight the positions of oxygen vacancy in other nine substrates. Overall from g(z) profiles, interfacial water in these ten simulations can be categorized into 4 high-density ( Fig.2b) and 6 low-density (Fig.2c) structures. For both the high-density and the low-density interfacial water, a sharp peak at about 2 Å corresponds to water adsorbed on Ti 5c . In the high-density model, water molecules bonded to O 2c s give rise to a second peak at about 2.6 Å, while the peak is missing for the low-density model (grey shaded area in Fig.2b,c). An immediate consequence of such structural changes is the density reduction of interfacial water, and a quantitative analysis will be discussed later.
Besides the density reduction, the low-density interfacial water is accompanied by a fluxional feature, which are reflected in Fig.3a-b by projected trajectories of O W in x-y plane. There are six black trajectories indicating six water molecules on Ti 5c . Because water bonding energy to Ti 5c is significantly higher than hydrogen bonding energy to O 2c [46,47], these molecules remain adsorbed in all simulations. Difference appears for water molecules on O 2c s (viridis lines), whose trajectories are concentrated in the high-density model, while they are extended in the low-density model.
Another feature indicated by the trajectory projection is that the high-density interfacial water layer is ordered, whereas the fluxional interfacial water layer of the low-density model is disordered. To quantify the ordering, we define an order parameter σ based on the lateral angular distribution of water. The order parameter is calculated by σ = 1 N ∑ N i cos(3θ ) in polar coordinate system, where N is the total number of water molecules within a radius of 3 Å from O 2c and θ is the polar angle of the molecule (on Ti 5c )-molecule (on O 2c ) line with respect to the x-axis (Fig.   S4a). Orderly adsorbed water molecules on Ti 5c and O 2c sites feature a distribution with θ = 0, θ = 2 3 π and θ = 4 3 π, where σ gets a maximum value of 1. In contrast, σ has a minimum value of -1 in directions around θ = 1 3 π, θ = π, and θ = 5 3 π. As a result, for relatively ordered distribution of water in the high-density model, σ approaches 1, whereas the fluxional distribution in the low-density model results in a smaller, possibly negative, value of σ . Furthermore, Fig.3c indicates there is a relation between the density of interfacial water and the order parameter σ . The number density of water molecules in Fig.3c is defined by the number of water molecules per unit cell within the region of 3 Å above the Ti 5c surface. Overall, the higher the density is the larger the order parameter σ is, meaning the structure is more ordered. However, when the density is lower the structure is less ordered and more fluxional. The disordering of the contact layer is also reflected on the orientation of water molecules. Fig.3d-e plots the distribution of molecular orientation φ of water above O 2c rows, which are hydrogen bonded with O 2c s directly or with the water adsorbed on Ti 5c . For the ordered high-density model, there is a clear preference of water molecules aligning in a double-leg orientation, which peaks at about 160 degree with respect to the surface normal (Fig.3d). In contrast, a typical fluxional interfacial model would break such a preference (Fig.3e). This difference in orientation indicates the fluxional molecules alternately form donating hydrogen bonds with O 2c and the bulk liquid. The distribution profiles of all ten systems are shown in the SI (Fig. S5, Fig. S6), supporting the analyses discussed.
In addition, the fluxional behaviors are accompanied by dynamic exchanges between the contact layer and the bulk liquid (Fig. S8).
From above results, a particularly interesting finding is that the density reduction and the fluxional interfacial water layer is observed on substrates with deep vacancies. Even for the V O 7 system, showing a similar second peak on g(z) to the conventional high-density interfacial water model, other analyses of order parameter (Fig.3c), water orientation (Fig. S5) and hydrogen bonded water number (Fig. S6) clarify its fluxional nature. V O 9, located ∼ 9 Å under the interfacial water layer, is the deepest vacancy considered in our simulations, and our results suggest that deeper vacancies would have similar effects. This means that only the pristine substrate and the substrates with surface vacancies (V O 1 and V O 2) maintain the high-density interfacial water.
Given the fact that in TiO 2 oxygen vacancies are ubiquitous at different depths across the substrate, low-density fluxional interfacial water may appear generally. The density reduction appeared in our simulation is in line with the experimental results of I.M.Nadeem et al., in which they found that a 50% coverage of water on O 2c fits best with measurements [43].
It is also worth noting that our results do not contradict with previous studies where water dissociation was discussed. Recently, Calegari Andrade et al. investigated the spontaneous dissociation of water on TiO 2 surface using machine learning force fields and suggested that spontaneous water dissociation requires nanoseconds for equilibrium [22]. Here we focus on the interfacial structure of water, characterized by the distribution of oxygen positions, whose equilibrium is often reached within a few tens of picoseconds. In fact, water dissociation induced by surface V O is also observed in our simulations with V O 1 and V O 2. Dissociation, however, does not lead to the density reduction of interfacial water. The fact that surface vacancies (V O 1 and V O 2) have minor effects on the interfacial water structure is rather counter-intuitive because surface oxygen vacancies are reactive to water, and lead to the formation of hydroxide groups on the surface.
Our counter-intuitive results suggest there is very delicate balance of competing interactions at the interface associated with the details of substrate. To provide deeper physical insights, we further compared the position and charge distribution of the V O 8 system and the pristine substrate. Fig.4a-b plots the Bader charge of oxygen and titanium atoms on each layer across the pristine and V O 8 substrates. In comparison to the pristine system, oxygen vacancies result in excess electrons in the substrate. We note that although the difference is a few percent, this magnitude of change is large enough to induce an electrostatic energy difference of tens of meV simply based on electrostatics, which is comparable to the 330K temperature (28 meV). Such a difference is enough to tip the balance of the competition between hydrogen bonding network and adsorption on bridging oxygen, thus leading to a new interfacial water structure on TiO 2 . For example, previous studies have shown a competition between two water molecules for hydrogen bonding to the bridging oxygen that raises up the adsorption energy by about 33 meV [48]. Slight electron deficit on bridging oxygen (Fig.4a) in defected system leads to weaker hydrogen bonding, hence pushing one of the competing water molecules away from the bridging oxygen.
Previous studies have shown that excess electrons can migrate from subsurface to the surface, affecting hydrogen bonding, adsorption and dissociation of water [13,29,30]. We note that the behavior of excess electrons is quite sensitive to DFT methods used. For example, excess electrons may localize and show polaronic behavior when DFT+U methods are used. We have also carried out simulations with SCAN-rvv10 and PBE+U, which show consistent results that density reduction can be induced by deep oxygen vacancies (Fig. S3, Fig. S7-8, Fig. S10). These additional tests suggest the low-density fluxional interfacial water is not an artifact and not sensitive to the description of localized electrons and the extent of water dissociation.  (Fig. S12). Overall, here we propose that the long-range effects to the interfacial water structure may entangle with minuscule atomic position and charge distribution charges of the substrate, which should be paid attention to and further investigated in the future.
Before concluding, we note that the delicate nature indicates that experimental identification of the interfacial water structure on anatase(101)-water interface is of great challenge, requiring not only advanced characterization techniques but also clean and well controlled substrate samples.
The interfacial water density change reflects to some extent on the hydrogen bonding change at the interface, which can be detected experimentally using vibrational spectroscopies. and D1 between pristine and defected systems leads to the difference in the interfacial structure.
Specifically, oxygen vacancy shifts the preference of local hydrogen bonding structure from D1 to C1, which is also highlighted in the dipole orientation profiles (Fig. 3d-e).
Finally, to examine whether our findings can extend to other systems, we have considered a few other substrates, namely rutile(110), anatase(001), and silica(0001) (Fig. S13, Fig. S14). These substrates cover three different kinds: (i) on rutile(110) water remains intact; (ii) on anatase(001) water dissociates; and on silica(0001) the surface is readily hydroxylated. In the end, on these substrates tested we do not see strong influences of deep oxygen vacancies on the interfacial water structure. Therefore, whether there are other substrates that possess similar interesting behaviors identified on anatase(101) remains an open question.
In summary, we have reported ab initio simulations of water-anatase(101) interface, where we have identified the low-density fluxional interfacial water due to oxygen vacancies in the substrate.
We show that the influence of oxygen vacancies can be rather long-range and indirect, which lead to the finding that even the deepest vacancy considered in our simulations can strongly impact the interfacial water structure. We also reveal the density reduction and induced flexibility is closely related, and the behaviors are likely to be physically amplified by minor changes in the substrate due to the presence of oxygen vacancy. The non-trivial long range effect identified in our study is another highlight of the complexity of interfacial water structure, and it is also potentially an opportunity for materials engineering to manipulate interfacial water properties by controlling the inner part of materials.