Simulations of Thin Water Layers on Reduced TiO - American

Apr 3, 2012 - Compared to the defect-free surface, the simulations highlight how the interstitial ... DFT+U method in the other.26 This detailed compa...
0 downloads 0 Views 4MB Size
Article pubs.acs.org/JPCC

DFT-GGA and DFT+U Simulations of Thin Water Layers on Reduced TiO2 Anatase Antonio Tilocca*,† and Annabella Selloni‡ †

Department of Chemistry, University College London, London WC1H 0AJ, U.K. Department of Chemistry, Princeton University, Princeton, New Jersey 08544, United States



ABSTRACT: We investigated the effect of a reducing subsurface Ti interstitial defect on the structure and reactivity of thin water layers adsorbed on the majority anatase(101) surface of titanium oxide using ab initio molecular dynamics simulations. We find that standard DFT-GGA and the DFT+U method predict similar energetics and dissociation barrier for a single water molecule adsorbed on the reduced surface; moreover, the two approaches also lead to very similar structural features and reactivity for an adsorbed water monolayer (ML) on the same surface. This allows us to model 1, 2, and 3 water layers on the reduced surface through Car−Parrinello molecular dynamics simulations up to 20 ps long. Compared to the defect-free surface, the simulations highlight how the interstitial defect alters the stability of surface adsorption sites, substantially enhancing the surface reactivity and leading to a markedly different structure of the first water layers adsorbed on the reduced surface.



INTRODUCTION The interaction of water with titania surfaces continues to attract considerable attention for both its fundamental scientific interest and its important role in technological applications.1−5 After many years of intense experimental and theoretical research, some aspects of this interaction, such as water adsorption in the dilute limit, have been elucidated to a large extent, but other features remain poorly understood. Among the latter, there are the properties of adsorbed thin water layers and how these are affected by the underlying surface structure. Knowledge of these properties is important for a better understanding of the mechanisms of TiO2 photocatalysis and would be very useful in other fields as well. For instance, immobilizing peptides, protein, and DNA strands on titanium oxide surfaces is of key importance for biomolecular recognition to develop biosensors and improve the performances of Ti-based biomaterials.6 The interaction of peptides with TiO2 is mediated by structured water adsorbed at the surface: proteins, collagen, and DNA bases bind to a thin water film in direct contact with the surface, rather than directly to the surface itself.7,8 The properties of the first water layers in close contact with the surface of titanium oxide and other biomaterials can thus control their biological response.9,10 Recently, there have been several computational studies of water on the majority (101) surface of anatase, the TiO2 polymorph most relevant for photocatalytic applications. These studies have shown that the stability of various adsorbates on anatase TiO2 is significantly influenced by coadsorbed water.11−13 Ab-initio simulations have also shown that the main structural features of the first water layers in direct contact with the surface14 are maintained at the interface with bulk liquid water.15 While these results highlight the interest of investigating the properties of thin water layers, they generally © 2012 American Chemical Society

refer to the stoichiometric, defect-free surface, or the reduced surface with a surface oxygen vacancy.17 On anatase (101), however, subsurface defects are usually present.18−20 In this paper we thus use density functional theory (DFT) calculations and first principles molecular dynamics (FPMD) simulations to examine how the structure of up to three monolayers thick water films is affected by a subsurface reducing defect. The structures of a water monolayer (ML), bilayer (BL), and trilayer (TL) on defect-free anatase (101) have been studied by us previously.14 We also investigated the structures of a water ML and a BL on anatase (101) in the presence of a surface oxygen vacancy.17 Although informative, those studies were somewhat limited by the use of a rather thin slab model and the relatively short duration of the simulations. Moreover, the simulations on the reduced surface17 were based on standard DFT in the generalized gradient approximation (GGA), which is now known to provide a rather unsatisfactory description of the electronic structure of intrinsic defects in TiO2 and other metal oxide materials.21−23 In the present study, we overcome these limitations by using thicker slab models and performing significantly longer FPMD simulations. To model the reduced surface, we consider a subsurface Ti interstitial, a frequent intrinsic reducing defect in TiO2.18,20,24,25 We assess the performance of DFT-GGA by comparing the structure and dynamics of an adsorbed ML obtained from two parallel simulations performed using DFT-GGA in one case and the DFT+U method in the other.26 This detailed comparison shows that the ML structure and dynamics from DFT-GGA are actually quite similar to those obtained using the DFT+U Received: February 18, 2012 Revised: March 31, 2012 Published: April 3, 2012 9114

dx.doi.org/10.1021/jp301624v | J. Phys. Chem. C 2012, 116, 9114−9121

The Journal of Physical Chemistry C



ASSESSING THE PERFORMANCE OF DFT-GGA VS DFT+U Single-Molecule Adsorption and Dynamics. We considered two possible stable locations for the subsurface Ti interstitial, previously labeled as T4 and T5.20 As shown in Figure 1, the Ti interstitial is closer to the surface in T4, which

method. On the basis of these results, we then use DFT-GGA for extensive simulations of water ML, BL, and TL on the stoichiometric and reduced surfaces. The latter simulations show that the interstitial defect significantly perturbs the adsorbed water layers, affecting the adsorption strength and the reactivity of the surface sites nearest to it and altering the vertical and in-plane ordering of the water layers observed on the defect-free surface.



Article

COMPUTATIONAL PROCEDURE

All calculations were performed using density functional theory (DFT) within the plane wave-pseudopotential framework as implemented in the Quantum-ESPRESSO (QE) package.27 Both the generalized gradient approximation (GGA) and the DFT+U method26 were employed, and in both cases the PBE28 exchange-correlation functional was adopted. The DFT-GGA calculations were spin-restricted, whereas spin polarization was included in the DFT+U calculations. For the latter, U = 3.5 eV was used, as determined previously.16 Core−valence interactions were represented through ultrasoft pseudopotentials,29 with plane-wave basis set cutoffs of 25 and 200 Ry for the smooth part of the wave functions and the augmented charge, respectively. Because of the relatively large supercell size, kpoint sampling was restricted to the Γ point only. The anatase(101) surface was represented as a periodically repeated slab, using previously calculated bulk lattice parameters20 (a = 3.77 Å and c = 9.54 Å). We used three-layers thick slab, with 1 × 3 (10.26 × 11.31 Å2) surface supercells, corresponding to 36 TiO2 units per simulation cell for the stoichiometric surface; for the reduced surface, an additional Ti atom was included. A vacuum gap of 14.5 Å separated the periodic images of the slab along the c direction.14 In all simulations, the Ti and O atoms in the bottom layer were fixed at their bulk positions. This computational setup has been thoroughly validated in our previous studies of water interacting with TiO2 surfaces.11,14,16 DFT-GGA-based Car−Parrinello molecular dynamics (CPMD)30 simulations were carried out with a fictitious electronic mass of 700 atomic units and a time step of 0.17 fs, using the deuterium mass for hydrogen atoms. The ionic temperature was controlled through a Nosé thermostat in all simulations; for the runs on the reduced surface, an additional Nosé thermostat on the electronic degrees of freedom was needed to avoid drifts in the fictitious kinetic energy and preserve adiabaticity.32 Spin-polarized DFT+U-based simulations were performed using Born−Oppenheimer molecular dynamics (BOMD)33 with a time step of 1 fs and a Berendsen thermostat to control the ionic temperature. Despite the longer time step, the spin-polarized DFT+U BOMD runs are obviously much more time-consuming34 than the spinrestricted DFT-GGA CPMD runs. The total electronic magnetization was unconstrained in the DFT+U BOMD simulations and converged to S = 4.0, which corresponds to a fully ferromagnetic arrangement of the four excess electrons associated with the Ti interstitial. The nudged elastic band (NEB) method35 was used to determine the DFT+U energy barrier for water dissociation on the reduced surface. Nine replicas were employed in the calculation; as in the BOMD runs, the unconstrained total magnetization converged to S = 4.0 for all replicas.

Figure 1. Side (left) and top (right) views of the anatase(101) surface with T4 (top panels) and T5 (bottom panels) titanium interstitials. The Ti5c sites exposed on the two surfaces are numbered as in the text. Ti and O atoms are represented as silver and red spheres, respectively, with larger spheres highlighting the Ti interstitials.

is deemed to result in a stronger effect on the surface reactivity, compared to the deeper and more stable T5.31 After relaxing the bare slabs, we modeled the adsorption of a single water molecule on all the inequivalent 5-coordinated surface Ti atoms (Ti5c) on the T4 and T5 surfaces; these sites are labeled 1, 2(3), 4(5), and 6 in Figure 1 for T4, whereas for T5 the inequivalent sites are those labeled 1(3), 2, 4, and 5(6). All structures were optimized at the DFT+U level and their adsorption energies compared with those previously obtained by standard DFT-GGA16,31 in Table 1. Even though the Table 1. Adsorption Energy of a Single Water Adsorbed at Different Ti5c Sites in the Presence of T4 and T5 Interstitialsa DFT+U surface T4

T5

site 1 2 4 6 1 2 4 5

(3) (5) (3)

(6)

DFT-GGA

Eads

ΔE

Eads16,31

ΔE

0.57 1.39 (1.06) 0.85 1.08 0.73 1.12 (0.92) 0.86 0.89

0.82 0 (0.33) 0.54 0.31 0.39 0 (0.20) 0.26 0.23

0.67 1.10 (0.84) 0.79 0.86 0.71 1.05 0.73 0.74

0.43 0 (0.26) 0.31 0.24 0.34 0 0.32 0.31

a Both absolute energies (Eads, eV) and energies relative to the most stable adsorption mode (ΔE) are reported. Bold (normal) values denote dissociative (molecular) adsorption.

9115

dx.doi.org/10.1021/jp301624v | J. Phys. Chem. C 2012, 116, 9114−9121

The Journal of Physical Chemistry C

Article

absolute values are different, the relative trends of DFT+U and DFT-GGA adsorption energies are exactly the same: the increasing order of stability calculated by the two approaches is 1 < 4 < 6 < 2 both for T4 and T5. Both DFT+U and DFTGGA identify the dissociative state at site 2 (close to the interstitial) as the most stable adsorption mode, whereas adsorption is weaker and molecular on the other Ti5c sites. On T4, dissociation at site 2 is more favorable than molecular adsorption by 0.33 eV (0.26 eV with DFT-GGA): the initial and final states are shown in Figure 2. The main effect of the +U repulsion is to enhance the adsorption energy difference between the different states.

Figure 3. Selected snapshots extracted from the DFT-GGA CPMD simulation of a water monolayer on the T4 surface. The oxygens of the three water molecules on the row atop the interstitial are highlighted in yellow (W1), green (W2), and pink. Figure 2. Molecular (left) and dissociative (right) adsorption states of a single water molecule on site 2 on surface T4. Color codes as in Figure 1, with water oxygen highlighted in yellow and hydrogen white.

Starting from the same initial configuration used for the DFT-GGA CPMD trajectory, a DFT+U BOMD simulation at 160 K showed a dynamical evolution very similar to that observed in the DFT-GGA run: W1 desorbs from its initial site and moves to the same stable configuration where it assists W2 dissociation and coordinates its proton transfer to a BO (see Figure 4). These processes occur however on a slightly shorter time scale compared to the DFT-GGA runs. As in the DFTGGA run, Figure 4 shows that also with DFT+U the W1 molecule leaves the initial site and reaches its stable configuration in ∼1 ps; at this point, however, it immediately gives rise to the proton transfer events leading to W2 dissociation and protonation of a BO, whereas the reaction occurred after an additional ∼7 ps in the DFT-GGA run. The driving force for the observed process likely is the interstitialinduced destabilization of site 1 (see Table 1). This could indeed explain the water desorption from that site and migration to a more favorable configuration, which in turn enables dissociation of an adjacent water molecule. This interpretation is supported by the fact that a CPMD run for surface T5 also led to desorption from site 1, which represents the site with the lowest single-water adsorption energy for the surface with the T5 interstitial as well. We determined the relative stabilities of the intermediate structures observed during the dynamical evolution of the water monolayer by calculating the adsorption energies of the geometry-optimized structures in Figure 5 (see Table 2). The initial flat ML (Figure 5a) corresponds to the most stable structure on the defect-free surface.14 Migration of W1 is driven by the higher stability of structure (b), whose energy is actually very close to that of the corresponding dissociated structure (c). The similar stability of molecular and dissociated structures leads to the frequent proton transfers observed during the dynamics, corresponding to breaking and reassociation of O−H bonds. The trend calculated at the DFT-GGA level is the same, even thoughas in the case of single water adsorptionthe absolute adsorption energies are lower. Vertical Ordering. Figure 6 shows the vertical distance from the surface of each water oxygen during the DFT-GGA and

Previous NEB/DFT-GGA calculations found a low energy barrier (0.24 eV) for water dissociation at site 2 on surface T4.16 With DFT+U, NEB calculations give an even lower barrier, 0.14 eV. The strong tendency of a water molecule to dissociate at site 2 on the T4 surface was confirmed by a 110 K DFT+U BOMD trajectory: the adsorbed molecule dissociated after 0.7 ps and the resulting pair of OH groups (right panel of Figure 2) remained stable for the rest of the trajectory. The dynamic behavior obtained with the DFT+U approach is in close agreement with that observed in DFT-GGA CPMD simulations.16 In summary, we find that DFT-GGA reproduces well both the energetics and the kinetics of single-water adsorption in the presence of Ti interstitials given by the DFT +U calculations. Monolayer. Reactivity. A DFT-GGA CPMD trajectory at 160 K was started after placing a water molecule ∼2.5 Å above each of the six Ti5c exposed on the T4 surface. Figure 3 illustrates the time evolution of the water monolayer. The water molecules in the left row, above the interstitial, show a significantly richer dynamic behavior compared to the molecules in the adjacent row. The water molecule initially adsorbed on site 1 (W1 hereafter) leaves this site and rapidly moves closer to the region above the interstitial, between the water molecule on site 2 (W2 hereafter) and the adjacent row of 2-fold coordinated bridging oxygens (BOs), accepting a hydrogen bond (Hb) from W2 and donating Hbs to two BOs. W1 reaches this stable configuration after 1 ps and maintains it for the rest of the trajectory. In this configuration, the molecule assists the proton transfer from W2 to one of the BOs: this process starts after 8 ps, when W2 dissociates and donates a proton to W1, which rapidly transfers another proton to a BO. Whereas W2 remains dissociated for the rest of the trajectory, frequent proton exchanges occur between W1 and the two BOs, essentially resulting in two different configurations, shown at 12 and 28 ps in Figure 3. The latter configuration is observed more frequently. 9116

dx.doi.org/10.1021/jp301624v | J. Phys. Chem. C 2012, 116, 9114−9121

The Journal of Physical Chemistry C

Article

Figure 6. Distribution (left panels) and trajectory (right) of watersurface vertical distances. Δz = z − z,̅ where z̅ is the average z coordinate of the six Ti5c. Top panels: DFT-GGA; bottom panels: DFT+U.

DFT+U MD trajectories of the water monolayer. In both cases, there is a noticeable difference between the height of the water molecules in the row atop the interstitial and those in the adjacent row, which are closer to the surface. An additional gap is present between the W1 and W2−W3 molecules atop the interstitial: W1, not being coordinated to a Ti5c, is farther away from the surface than W2 and W3. These features are common to both DFT-GGA and DFT+U trajectories, indicating that the perturbation induced by the interstitial on the ML vertical ordering is qualitatively reproduced already at the DFT-GGA level. However, the W4−W6 molecules are slightly (0.15 Å) closer and the W5−W6 slightly (0.1 Å) farther from the surface in the DFT+U simulation compared to the vertical distances in the DFT-GGA CPMD run. This results in a more pronounced vertical split between W4−W6 and W2−W3 molecules, as shown by the vertical distance (Δz) distributions in the left panels of Figure 6. This difference could partially reflect the fact that the proton transfer events take place earlier in the DFT+U runs, whereas, as mentioned above, no dissociation occurs in the DFT-GGA runs within the time frame covered in Figure 6. Water−Water and Water−Surface Interactions. The radial distribution functions for pair interactions involving water molecules in the ML are displayed in Figure 7. We can see that DFT-GGA and DFT+U yield quite similar local structures for the water ML. The characteristic peaks for the water−surface interaction are shown in the bottom panel of the figure. The main Ti−Ow peak (where Ow is a water oxygen) is split in two due the contributions from the Ti−OH bond (around 1.9 Å) formed upon W2 dissociation, much shorter than the Ti−OH2 distance of ∼2.3 Å. The Ow−Os (where Os is an oxygen belonging to the TiO2 slab) main peak at 2.8 Å arises from hydrogen bonds formed between W2 and the two adjacent bridging oxygens; this peak has a short-distance shoulder around 2.5 Å, which reflects a shorter Ti−OsH···Ow hydrogen bond formed upon W1 dissociation. The water−water radial distribution functions are displayed in the top panel of Figure 7. Whereas the Ow−Ow peak at 3.8 Å reflects water molecules not hydrogen-bonded and in registry with the lattice of Ti5c sites, as in the unperturbed ML on the defect-free surface (see ref 14 and also the discussion below), new O−H and O−O peaks at 1.5 and 2.5−2.7 Å reveal

Figure 4. Time evolution of characteristic interatomic distances during the DFT+U and DFT-GGA MD runs of water ML on the T4 surface. The top panels (a) illustrate the migration process of the W1 molecule leaving its initial site and settling into the final configuration; panels (b) highlight the proton transfers to the bridging oxygens; panels (c) describe the dissociation of molecule W2; and panels (d) describe the proton transfers involving W1. The vertical dashed line identifies the start of the dissociation processes in each case.

Figure 5. DFT+U-optimized structures of a water monolayer on the T4 surface. Color codes as in Figure 1, with water oxygen highlighted in yellow and hydrogen white.

Table 2. Adsorption Energy (eV/H2O) of a Water Monolayer on T4 and T5 Surfaces (Structures Are Labeled as in Figure 5) Eads surface

structure

DFT+U

DFT-GGA

T4

(a) (b) (c) (a) (b) (c)

0.80 0.88 0.91 0.78 0.82 0.81

0.77 0.81 0.82 0.73 0.76 0.76

T5

9117

dx.doi.org/10.1021/jp301624v | J. Phys. Chem. C 2012, 116, 9114−9121

The Journal of Physical Chemistry C

Article

Figure 7. Water−water (top) and water−surface (bottom) radial distribution functions for a ML on surface T4, calculated from the DFT-GGA and DFT+U MD simulations. The DFT+U curves are shifted vertically for clarity. Figure 8. Distribution (left panels) and trajectory (right) of water− surface vertical distances for mono-, bi-, and trilayer adsorbed on the defect-free surface. Δz = z − z,̅ where z̅ is the average z coordinate of the six Ti5c.

additional hydrogen bonds within the ML involving almost exclusively W1 and W2. These features characterize both the DFT+U and DFT-GGA description of the ML−surface and intra-ML interactions. The most significant difference is the Ow−Ow H-bond distance, which is 0.2 Å shorter in the DFTGGA model.



INFLUENCE OF THE TI INTERSTITIAL ON THE STRUCTURE OF 1−3 WATER LAYERS From the results in the previous section, it is reasonable to conclude that the DFT+U and DFT-GGA descriptions of structural and dynamical features of a water monolayer on the reduced surface are sufficiently close to justify the use of the standard DFT-GGA approach to investigate systems which would be too cumbersome to simulate with spin-polarized DFT +U-based BOMD. Therefore, we now use DFT-GGA CPMD to study the effects of the interstitial defect on the structure of up to three layers of adsorbed water. For this purpose, we use as a reference the corresponding structures on the defect-free surfaces, obtained with the same computational settings. Thin Water Layers on the Stoichiometric Surface. Vertical Ordering. DFT-GGA CPMD simulations of one, two, and three water layers on the defect-free surface were carried out at 160 K for 18 ps. The vertical distances in Figure 8 show that the ML and BL quickly settle into stable configurations, characterized by well-defined distances from the surface for each molecule and corresponding sharp peaks in the distance distributions at Δz = 2 and 3 Å. The left and middle panel in Figure 9 shows that water molecules in the first layer (Δz = 2 Å) are coordinated to Ti5c, whereas the second layer (Δz = 3 Å) is formed by strongly oriented molecules, all arranged with one hydrogen pointing upward and the other toward the surface and forming a strong H-bond with a BO. The vertical ordering of the first two layers in contact with the surface is maintained also in the TL: the position of the first two peaks in the z distribution is the same already observed for the ML and BL. However the presence of the additional water molecules

Figure 9. Structure of a water mono-, bi-, and trilayer on the defectfree surface. In each case, the last configurations of the corresponding MD runs at 160 K is shown. H2O molecules bonded to Ti5c and Os are highlighted in yellow and blue, respectively, whereas the remaining molecules are green.

introduces some disorder in the first layer, with a molecule becoming slightly higher than the others. Moreover, only three molecules are now in the second layer, H-bonded to BOs, while a third layer is present at Δz = 4.5 Å, formed by molecules Hbonded to the first two water layers (Figure 9, right panel). The corresponding peak in the height distribution has a tail extending up to 6 Å, suggesting a higher mobility and a more disordered arrangement within the top layer. Intra- and Interlayer Structure. H2O−H2O Interactions. The radial distribution functions (rdfs) in Figure 10 allow us to further characterize the structure of the thin layers on the defect-free surface. The top panel shows that Ti5c-coordinated water molecules in the ML do not significantly interact with each other, due to the large Ti5c−Ti5c distance.14 Clear signatures of a H-bonded system (intermolecular Ow−Hw and Ow−Ow peaks at 1.7 and 2.7 Å, respectively) emerge for the BL, and are maintained for the TL, with a generally more disordered arrangement evident from the broader peaks. H2O−Surface Interactions. The persistent Ti−Ow peak around 2.3 Å indicates that water molecules remain strongly 9118

dx.doi.org/10.1021/jp301624v | J. Phys. Chem. C 2012, 116, 9114−9121

The Journal of Physical Chemistry C

Article

to an intermediate position between the Ti5c and BO rows (Figure 11, left panel). The molecules directly above the interstitial tend to stay farther away from the surface, leading to a secondary peak in the z-distribution (Figure 12, bottom panels).

Figure 10. Radial distribution functions for water−water (top) and water−surface (bottom) interactions for mono-, bi-, and trilayers on the defect-free surface.

coordinated to Ti5c at the different coverages, and the strength of this interaction increases at coverages higher than one ML. Water molecules in the ML are anchored to the Ti5c sites and can only form weak Hbs with the surface BO (bottom panel, dashed red and cyan curves). On the other hand, the BL and TL also include molecules that are not coordinated to Ti5c and are thus able to form stronger Hbs with the surface. In particular, water−surface Hb interactions are strongest in the BL due to the contribution from strongly oriented molecules pointing a hydrogen down (Figure 9, middle panel), which leads to sharp Hw−Os and Ow−Os peaks at 1.5 and 2.5 Å. Even though the weight of molecules in this specific configuration is lower in the TL, for which the sharp peaks noted above essentially disappear, water−surface interactions are still significant for the TL. Thin Water Layers on the Reduced Surface. Vertical Ordering. DFT-GGA simulations of 1−3 water layers on the T4 surface were carried out for 20 ps at 160 K. The final configurations are shown in Figure 11. As already pointed out, the presence of the interstitial affects substantially the ML structure compared to the defect-free surface. The ordered arrangement of Ti5c-coordinated molecules observed on the stoichiometric surface changes into a partially dissociated layer with an empty Ti5c site, due to the migration of the W1 water

Figure 12. Distribution (left panels) and trajectory (right) of water− surface vertical distances for mono-, bi-, and trilayer adsorbed on the surface with T4 interstitial. Δz = z − z,̅ where z̅ is the average z coordinate of the six Ti5c.

The Ti interstitial has an even more marked influence on the BL structure, whose vertical arrangement changes from the ordered double layer with sharp peaks seen on the defect-free surface to a multipeak distribution reaching up to 8 Å from the T4 surface. Even though the peaks at 2 and 3 Å characteristic of the BL on the defect-free surface are still present, the structure of these layers is significantly perturbed, with a dissociated molecule in the first layer and only half of the molecules left in the second layer, H-bonded to BOs (Figure 11, middle panel). The remaining molecules produce an additional peak at z = 4.2 Å and further extend up to 8 Å away. For the TL, the structure of the first two layers in contact with the surface (Figure 11, right panel) is similar to that just described for the BL, with a dissociated water in the first layer and a depleted second layer of molecules H-bonded to the surface. Above these two layers in direct contact with the surface sites, the remaining molecules in the BL and TL are significantly spread for about 3 Å perpendicularly to the surface; this upper region is uniformly filled only at TL coverage, whereas low-density gaps are evident for the BL (Figures 11 and 12). Intra- and Interlayer Structure. Water−Water Interactions. The radial distribution functions for water on the T4 surface (Figure 13) further detail the substantial influence of the interstitial defect on the structure of adsorbed water. At variance with the defect-free surface, signatures of hydrogen bonds between water molecules are observed already in the ML. In this case, Hbs always involve OH groups produced by

Figure 11. Structure of a water mono-, bi-, and trilayer on the T4 interstitial surface. In each case, the last configurations of the corresponding CPMD runs at 160 K is shown. Color codes as in Figure 9. 9119

dx.doi.org/10.1021/jp301624v | J. Phys. Chem. C 2012, 116, 9114−9121

The Journal of Physical Chemistry C

Article

BOMD simulations for an adsorbed monolayer on the defected surface, we have shown that DFT-GGA can adequately describe the structure and dynamics of thin water films on reduced TiO2, even though this approach provides an unsatisfactory description of the electronic structure of localized defect states in the material.21−23 In this context, a particularly important finding is that DFT-GGA and DFT+U predict the same order of stability for single-water adsorption on the surface with a subsurface interstitial defect, analogously to what previously observed for the surface with a subsurface O vacancy.16 The improved description of the localized electronic states associated with the interstitial defect, obtained with the +U correction, does not affect significantly the main structural features of the adsorbed water layers. The only small differences concern the time scale of water dissociation and the Ow−Ow distance between two H-bonded water molecules in a ML. The shorter time for the onset of dissociation in the GGA+U approach is consistent with the smaller computed water dissociation barrier. Moreover, taking into account the wellknown GGA underestimation of the O−O distance in liquid water,38 the longer O−O distance between two H-bonded adsorbed water molecules obtained with the +U correction suggests that DFT+U can indirectly improve the description of H-bonds in adsorbed water, similarly to the effect of hybrid functionals for liquid water.39 Experimental techniques such as infrared reflection absorption spectroscopy (IRRAS)40,41 could probe the H-bond structure of the adsorbed layers and complement our findings. On the stoichiometric anatase(101) surface, the main factors which control multilayer water adsorption are the balance between water−surface and water−water interactions (the latter being weaker14) and the surface specific geometry which does not favor interactions between water molecules coordinated to the surface sites. This leads to substantial vertical and in-plane ordering of the first two layers in direct contact with the surface, which is maintained even at higher coverages, although surface-induced ordering is weaker above the second layer. On the reduced (101) surface, on the other hand, the structure of water overlayers is affected also by the defect-induced stress, which destabilizes H2O−Ti5c bonds and breaks the ordered arrangements formed on the stoichiometric surface. We previously observed a similar effect caused by a surface oxygen vacancy,17 which also leads to a partially dissociated and highly disrupted water ML, with H2O−Ti5c transformed into H2O−Os, as well as water−water and water− surface Hbs not present on the clean surface. Intermolecular water interactions at higher coverages have previously been found to favor water dissociation in the presence of both surface17 and subsurface16 oxygen vacancies. This effect is less noticeable for an interstitial due to the higher reactivity of the defect which strongly favors dissociation already at low coverage. The mechanism of O vacancy- and Ti interstitial-induced dissociation may also differ: whereas a surface oxygen vacancy can directly participate to the dissociation mechanism at low and high coverages (the driving force being the filling of the vacancy by an OH group),17,36 a subsurface interstitial affects the process indirectly, through the destabilization and consequent mobilization of water molecules which enables water dissociation by a concerted mechanism. In any case, the present simulations confirm that defects, either surface or more likely subsurface,18 are required to enable water dissociation on anatase(101) at any coverage.

Figure 13. Radial distribution functions for water−water (top) and water−surface (bottom) interactions for mono-, bi-, and trilayers on the T4 surface.

water dissociation and are stronger than Hbs between undissociated water molecules, which become dominant in the bi- and trilayer: Ow−Ow and Ow−Hw peaks lie at 2.5 and 1.45 Å in the ML, compared to values of 2.7 and 1.7 Å in the BL and TL. The water−water hydrogen bond patterns in the BL and TL are very similar, with distances approaching those of a liquid water environment. Water−Surface Interactions. The water−surface rdfs are shown in the bottom panel of Figure 13. The formation of a Ti−OH bond results in the short-distance peak in the Ti−Ow rdf, whose weight is obviously higher in the ML case. The Ti− Ow rdf remains structured for the ML and increasingly less going to BL and TL. At variance with the defect-free surface, hydrogen bonds between water and the reduced surface are already present at ML coverage (compare the Hw−Os rdfs between 1.5 and 2 Å), mainly involving the W1 molecule. The more prominent water−surface Hbs peaks in the BL than in TL case are somewhat a remnant of the highly ordered pattern for the BL on the defect-free surface.



DISCUSSION AND FINAL REMARKS Recent scanning tunneling microscopy (STM) and temperature-programmed desorption (TPD) studies of water on reduced anatase (101) surfaces have provided evidence of a higher water desorption temperature on highly reduced surfaces as well as a clear preference for water adsorption close to subsurface defects.16 These findings were supported by DFT calculations, which showed higher adsorption energies at sites above subsurface oxygen vacancies and Ti interstitials. Lowcoverage water dissociation turned out to be more favorable on Ti interstitials than on O vacancies, indicating a higher reactivity of interstitial defects.16 Motivated by these results, in this paper we have focused on Ti interstitials, and in particular we have studied the influence of a subsurface Ti interstitial on the structure and dynamics of thin water layers on the anatase (101) surface. By performing a detailed comparison of DFT-GGA CPMD and DFT+U 9120

dx.doi.org/10.1021/jp301624v | J. Phys. Chem. C 2012, 116, 9114−9121

The Journal of Physical Chemistry C

Article

(3) Fujishima, A.; Zhang, X. T.; Tryk, D. A. Surf. Sci. Rep. 2008, 63, 515. (4) Diebold, U. Surf. Sci. Rep. 2003, 48, 53. (5) Sun, C. H.; Liu, L.; Selloni, A.; Lu, G.; Smith, S. C. J. Mater. Chem. 2010, 20, 10319. (6) Cohavi, O.; Corni, S.; De Rienzo, F.; Di Felice, R.; Gottschalk, K. E.; Hoefling, M.; Kokh, D.; Molinari, E.; Schreiber, G.; Vaskevich, A.; Wade, R. C. J. Mol. Recognit. 2010, 23, 259. (7) Monti, S. J. Phys. Chem. C 2007, 111, 16962. (8) Monti, S.; Walsh, T. R. J. Phys. Chem. C 2011, 105, 24238. (9) Skelton, A. A.; Liang, T.; Walsh, T. R. ACS Appl. Mater. Interfaces 2009, 1, 1482. (10) Tilocca, A.; Cormack, A. N. ACS Appl. Mater. Interfaces 2009, 1, 1324. (11) He, Y.; Tilocca, A.; Dulub, O.; Selloni, A.; Diebold, U. Nat. Mater. 2009, 8, 585. (12) Miller, K. L.; Musgrave, C. B.; Falconer, J. L.; Medlin, J. W. J. Phys. Chem. C 2011, 115, 2738. (13) Sanchez, V. M.; de la LLave, E.; Scherlis, D. A. Langmuir 2011, 27, 2411. (14) Tilocca, A.; Selloni, A. Langmuir 2004, 20, 8379. (15) Sumita, M.; Hu, C.; Tateyama, Y. J. Phys. Chem. C 2010, 114, 18529. (16) Aschauer, U.; He, Y.; Cheng, H.; Li, S.; Diebold, U.; Selloni, A. J. Phys. Chem. C 2010, 114, 1278. (17) Tilocca, A.; Selloni, A. J. Phys. Chem. B 2004, 8, 4743. (18) He, Y.; Dulub, O.; Cheng, H.; Selloni, A.; Diebold, U. Phys. Rev. Lett. 2009, 102, 106015. (19) Cheng, H.; Selloni, A. Phys. Rev. B 2009, 79, 092101. (20) Cheng, H.; Selloni, A. J. Chem. Phys. 2009, 131, 054703. (21) Di Valentin, C.; Pacchioni, G.; Selloni, A. Phys. Rev. Lett. 2006, 97, 166803. (22) Finazzi, E.; Di Valentin, C.; Pacchioni, G.; Selloni, A. J. Chem. Phys. 2008, 129, 154113. (23) Ganduglia-Pirovano, V. M.; Hofmann, A.; Sauer, J. Surf. Sci. Rep. 2007, 62, 219. (24) Wendt, S.; Sprunger, P. T.; Lira, E.; Madsen, G. K.; Li, Z.; Hansen, J.; Matthiesen, J.; Blekinge-Rasmussen, A.; Laegsgaard, E.; Hammer, B.; Besenbacher, F. Science 2008, 320, 1755. (25) Henderson, M. A. Surf. Sci. 1999, 419, 174. (26) Anisimov, V. I.; Zaanen, J.; Andersen, O. K. Phys. Rev. B 1991, 44, 943. (27) Giannozzi, P.; et al. J. Phys.: Condens. Matter 2009, 21, 395502. (28) Perdew, J. P.; Burke, K.; Ernzerhof, M. Phys. Rev. Lett. 1996, 77, 3865. (29) Vanderbilt, D. Phys. Rev. B 1990, 41, 7892. (30) Car, R.; Parrinello, M. Phys. Rev. Lett. 1985, 55, 2471. (31) Aschauer, U.; Selloni, A. Proc. SPIE 2010, 7758, 77580B. (32) Blöchl, P. E.; Parrinello, M. Phys. Rev. B 1992, 45, 9413. (33) Marx, D.; Hutter, J. In Modern Methods and Algorithms of Quantum Chemistry; Grotendorst, J., Ed.; NIC Series; FZ Julich: Germany, 2000; Vol. 1. (34) Using 32 cores connected by Infiniband on the Della Linux cluster (http://www.princeton.edu/researchcomputing/), a trajectory of 1 ps required 300−400 wall clock hours using spin-polarized DFT +U BOMD and 16 h with spin-restricted DFT-GGA CPMD. (35) Henkelman, G.; Uberuaga, B. P.; Jónsson, H. J. Chem. Phys. 2000, 113, 9901. (36) Tilocca, A.; Selloni, A. J. Chem. Phys. 2003, 119, 7445. (37) Finazzi, E.; Di Valentin, G.; Pacchioni, G. J. Phys. Chem. C 2009, 113, 3382. (38) Sit, P.H.-L.; Marzari, N. J. Chem. Phys. 2005, 122, 204510. (39) Zhang, C.; Donadio, D.; Gygi, F.; Galli, G. J. Chem. Theory Comput. 2011, 7, 1443. (40) Xu, M.; Noei, H.; Buchholz, M.; Muhler, M.; Wöll, C.; Wang, Y. Catal. Today 2012, 182, 12. (41) Kimmel, G. A.; Baer, M. A.; Petrik, N. G.; VandeVondele, J. A.; Rousseau, R.; Mundy, C. J. J. Phys. Chem. Lett. 2012, 3, 778. (42) Cheng, H.; Selloni, A. Langmuir 2010, 26, 11518.

The structure of thin water layers on the defect-free surface that we have found in this study is very close to that previously obtained with less demanding simulations involving thinner slabs and shorter MD runs.14,17 However, among the different BL configurations observed in our previous simulations,14 the present results strongly favor the highly ordered arrangement where each molecule in the second layer is oriented with a hydrogen atom pointing upward. This arrangement denotes unfavorable lateral interactions between molecules H-bonded to a surface BO (Os). Alternative BL configurations can only be attained by removing a molecule from either a Ti5c site or an Os site,14 but we have not observed these alternative and more disordered BL configurations in the present simulations, suggesting that their energetic cost is too high when a more realistic surface model is adopted. Therefore, the most stable structures of a ML and a BL on the defect-free surface contain water molecules in registry with all the Ti5c and Os exposed on the surface; only when a third water layer is adsorbed, a less ordered arrangement appears also in the first two layers, which do not longer match all the surface sites. The Ti interstitial strongly perturbs the ordered adsorption structure present on the stoichiometric surface. Water molecules in the row atop the defect are much more perturbed than those in the adjacent row; moreover, adsorption just above the interstitial is strengthened whereas adsorption at the adjacent sites is weakened, ultimately leading to W1 desorption and W2 dissociation. The disrupted structure of the first two contact layers at BL and TL coverage represents a less favorable landscape for the attachment of additional molecules, which could at least in part explain the significant vertical spread on the reduced compared to the defect-free surface (compare Figures 8 and 12). These effects could also reflect the excess (negative) charge introduced on the surface by the interstitial.37 In particular, the higher spread and also the larger distance from the surface of molecules directly above the interstitial suggest that reducing defects provide some hydrophobic character to the surface. Additional simulations of bulk liquid water on the reduced surface would be needed to assess the long-range consequences of this effect: the results of this work show that this task can be effectively accomplished through standard DFT-GGA CPMD simulations.15,42



AUTHOR INFORMATION

Corresponding Author

*E-mail: [email protected]. Notes

The authors declare no competing financial interest.



ACKNOWLEDGMENTS A.T. acknowledges financial support from the U.K.’s Royal Society (University Research Fellowship). A.S. thanks the support of DoE-BES, Chemical Sciences, Geosciences and Biosciences Division, Contract DE-FG02-12ER16286. We used resources of the TIGRESS high performance computer center at Princeton University, which is jointly supported by the Princeton Institute for Computational Science and Engineering and the Princeton University Office of Information Technology.



REFERENCES

(1) Henderson, M. A. Surf. Sci. Rep. 2002, 46, 1. (2) Henderson, M. A. Surf. Sci. Rep. 2011, 66, 185. 9121

dx.doi.org/10.1021/jp301624v | J. Phys. Chem. C 2012, 116, 9114−9121