The water/ceria(111) interface: Computational overview and new structures

Thin film structures of water on the CeO 2 (111) surface for coverages between 0.5 and 2.0 water monolayers have been optimized and analyzed using density functional theory (optPBE-vdW functional). We present a new 1.0 ML structure that is both the lowest in energy published and features a hydrogen-bond network extending the surface in one-dimension, contrary to what has been found in the literature, and contrary to what has been expected due to the large bulk ceria cell dimension. The adsorption energies for the monolayer and multilayered water structures agree well with experimental temperature programmed desorption results from the literature, and we discuss the stability window of CeO 2 (111) surfaces covered with 0.5–2.0 ML of water.


I. INTRODUCTION
The interaction of water with solid oxide surfaces is of critical importance in many environmentally and technologically relevant processes, including electrochemistry, corrosion, catalysis, atmospheric chemistry, and photocatalytic water splitting (see, for example, Refs. 1 and 2). One such oxide is CeO 2 (ceria) relevant in all the aforementioned applications as well as in a range of emerging ones such as solid oxide electrolysis cells, 3 hydrophobic coatings, 4 and bio-medicine. 5 Scientifically, the mechanisms underpinning these applications are all very complex and particularly so by virtue of the structure and properties of the oxide material itself, its surface terminations, its acid-base properties, and the possible presence of defects, all of which influence adsorption of water. Molecular-level insight concerning where and how the water molecules bind to the ceria surface is a prerequisite to building mechanistic models that can handle such scenarios. In this paper, we aim to review previous findings and elaborate on the interface structures at low water coverages on ceria(111), in search of extended hydrogen-bond patterns using density functional theory (DFT).
In the last two decades, considerable efforts have been made to explore thin water films (up to a few monolayers) on ceria using a range of experimental methods such as X-ray photoelectron spectroscopy (XPS), [6][7][8][9] temperature programmed desorption (TPD), [6][7][8][9][10] and scanning force microscopy (SFM). 11,12 The experimental evidence suggests that water molecules in thin water films on ceria(111) partially dissociate, i.e., there is coexistence of H 2 O and OH. [7][8][9]11 Moreover, it is quite well established that water molecules (intact or dissociated) adsorb above exposed cerium cations as shown in Ref. 11, but the details of the interface structure have not been possible to decipher by experimental methods. Here, computational studies can provide valuable information.
Indeed, a significant number of computational studies have been published for water on ceria(111) with the aim to elucidate the interface structure and quantify the water-water and watersurface interactions. 4,[13][14][15][16][17][18][19][20][21][22][23][24][25] All periodic quantum-mechanical (QM) studies that we are aware of for water/ceria(111) are listed in Table I, except for Ref. 23, which is an ab initio MD simulation of a water film on stoichiometric ceria(111) at room temperature. It should be noted that Refs. 19 and 21 contained a large number of monomer calculations; all their smaller cell results were not included in our table. In summary, the rightmost column in Table I shows that, to date, all stable structures found by total energy minimization at coverages up to and including 1.0 ML coverage (monolayer is defined in Sec. II) consist of individual water molecules or finite clusters of water molecules, adsorbed on the surface. For higher water TABLE I. Theoretical work on adsorption of water on the ceria (CeO 2 ) (111) surface, ordered by surface coverage. First column: Reference. Second column: Water coverage (see definition in text). Third column: Q is the fraction of dissociated water molecules out of the total number of adsorbed molecules. Fourth column: Surface supercell size used. Fifth column: Reported adsorption energies. For some entries, two energy values are given; the first refers to the most stable structure obtained and the second to the second most stable structure (when reported with the same coverage, cell, Q, and method). For some other entries, additional values in brackets, [⋅], are given, namely, in the cases where we have started from the reported structure and reoptimized using optPBE-vdW with our computational protocol; in all those cases, we reproduced the published structure as far as we could judge. Sixth column: Program and method. Seventh column: A short description of our interpretation of the hydrogen-bond motif described in each reference. General: The last three entries are three new structures found in this work.

References
Coverage  coverages, starting at 1.75 ML, infinite hydrogen-bonded networks have been found, spanning the two surface dimensions, but for lower coverages, 1.0 ML or less, this is not the case.
In this paper, we continue the exploration of thin water films on ceria with special focus on the 1.0 ML case for the fully oxidized ceria(111) surface with the aim to find if a single monolayer can form a continuous hydrogen bond network extending the surface. The answer is yes and we will present a 1.0 ML structure, which not only is (slightly) more stable than all other monolayer water/ceria(111) structures presented in the literature so far but also forms a periodically repeated H-bond network in one dimension.
Another focus of this paper is the consequences-or not-of water dissociation. Considering the extremes for monolayer coverage, i.e., zero and full dissociation, the energy difference is about 0.4 eV/H 2 O, 18 in contrast to the tiny energy difference observed between an intact isolated water molecule on the ceria(111) surface and its dissociated counterpart (see Table I). Here, we present monolayer structures that span the full dissociation range: 0%, 25%, 50%, 75%, and 100%. Using DFT, we have optimized 13 unique monolayer structures, including the three previously reported in the literature. 16,18,21 Although these 13 represent just a subset of all the existing energy minima on the 1.0 ML potential energy surface, they allow us to discuss stability vs dissociation degree and structural motif since these structures display a range of H-bond patterns. This is the second aim of this paper. For the 1.0 ML coverage, we find that partial dissociation is, in fact, more stable than the fully intact layer. Furthermore, we attempt to validate our results by simulating the TPD desorption temperatures. This paper is organized as follows: Section II describes the electronic structure calculations, geometric structure optimization, and water adsorption starting structures. Section III A addresses our first aim and provides a condensed review of the existing literature for different coverages and a presentation of the most stable structures that we have found. Here, we also present predicted scanning tunneling microscopy (STM) images primarily as a good way to highlight differences between structural patterns. Section III B addresses the second aim, water dissociation. In Sec. III C, we relate the adsorption energies to TPD experiments from the literature, and calculations of the stability windows for the different coverages are presented in Sec. III D.  subsequently, they were all reoptimized by us using the optPBE-vdW method and consistent settings (as described above and in Sec. II B) concerning the number of atomic slab layers, the vacuum gap size, the cell parameters, and our computational protocol. The other 21 structures were constructed by placing water at plausible orientations in mixes of intact and dissociated water molecules according to 0%, 25%, 50%, 75%, and 100% dissociation. The degree of dissociation is henceforth abbreviated as Q-N, where N is the fraction of dissociated water molecules, i.e., Q-0.0, Q-0.25, Q-0.5, Q-0.75, and Q-1.0.
After full geometry optimization, only 17 unique structures remained. They are listed in Table II and span coverages from 0.5 to 2.0 ML. As mentioned, 13 of these pertain to the 1.0 ML coverage. In the five cases where we started with structures from the literature, our optimized structures essentially gave the same structure as the published ones. The adsorption energies for these cases are reported in brackets in the fifth column in Table I, and all our adsorption energies are listed in Table II. A reference calculation for the gas-phase water molecule [H 2 O(g)] was simulated in a 13 Å × 13 Å × 13 Å box.

B. Electronic structure calculations
All calculations were performed using the Vienna ab initio Simulation Package (VASP version 5.4.4) 27-30 using the optPBE-vdW functional, 31,32 found to quite accurately describe adsorption of water on oxides 33 and in condensed systems containing H 2 O and ions. 34 The valence states, Ce(4f, 5s, 5p, 5d, 6s), O(2s, 2p), and H(1s), were treated explicitly, and valence-core interactions were treated with the projector-augmented wave (PAW) method. 35,36 The calculations were not spin-polarized or dipole corrected. The plane-wave basis-set was truncated at 1000 eV, and the electronic structure convergence criterion in the self-consistent-field (SCF) cycle was set to 10 −7 eV.
The slab system's Brillouin zones were sampled using a Γ-centered 5 × 5 × 1 Monkhorst-Pack grid and Blöchl tetrahedron TABLE II. Results from our own optPBE-vdW structure optimizations. The upper part lists all 13 1.0 ML structures in order of decreasing stability. The structure label in the first column refers to the number displayed in Fig. 3. The lower part lists the stable structures that we examined for each coverage. A, 1, and B were constructed and optimized by us (and are also listed on the last three rows in Table I), while C and D were started from published structures (as is clarified in the last column).  Atomic positions were optimized in the fixed supercells using the conjugate-gradient algorithm converged to forces less than 0.01 eV/Å on any and all atoms (individually and summed together). The adsorption energy was converged (to 0.005 eV) with respect to energy cutoff, k-points, and number of layers.

Structure label
A Hubbard U-type correction, 37 which is often used in computational studies of ceria systems, was purposely not selected for this work as reduced cerium ions or oxygen vacancies were neither expected nor observed.

C. Properties
The water adsorption energy was calculated from the total energy of the optimized slab with water molecules adsorbed E[H 2 O/CeO 2 ], the optimized clean slab E[CeO 2 ], and the optimized gas-phase molecule E[H 2 O(g)] according to the following formula: where n is the number of adsorbed water molecules per supercell (intact or not), meaning that the adsorption energy is calculated per water molecule (intact or not).
In addition, we have simulated scanning tunneling microscopy (STM) images using the Tersoff-Hamann 38 approximation at a height of 2.0 Å above the highest located oxygen with a bias voltage of −3 eV corresponding to filled-state imaging for visualization purposes.

III. RESULTS AND DISCUSSION
A. Finite and extended water structures on ceria(111) One of the aims of this paper, as mentioned above, is to investigate whether stable extended hydrogen-bond networks can be found for water adsorbed on the ceria(111) surface. But first, let us place this question in the context of what has been published so far.

Single water molecules
Adsorption of a single water molecule has been studied extensively by computational methods in Refs. 14,15,17,18,20,21,24, and 25 nicely summarized and somewhat definitively answered in Ref. 19. Two types of surface water species have been found, intact and dissociated water molecules [cf. Fig. 1(b)]. The dissociated water molecule is formed by transferring a H (from H 2 O) to a neighboring surface oxygen, resulting in a "free" hydroxide OH f , and a surface hydroxide OsH (note that the negative charge is not written out in this notation). All these studies are included in Table I, labeled "monomer" and "dissociated monomer," respectively, in the rightmost column. From these works, it has become clear that the adsorption energy difference between the intact and the dissociated water monomer is small and that adsorbed water can hydrogen bond to surface oxygen.

0.5-1.0 monolayer
With the formation of a thin water film, the adsorbed species (intact or dissociated water) can not only interact with the surface but also with other adsorbed water molecules. The 0.5 and 1.0 ML coverages were first explored by Kumar and Schelling, 16 showing a marginal gain in adsorption energy per water molecule with increasing water coverage and small interactions between water molecules, but not observing any hydrogen bonds between them. Again we refer to Table I for quantitative information. A later study, by Fronzi et al. 17 of coverages up to 0.5 ML, also observed that adsorption energies were not coverage dependent up to 0.5 ML. Molinari et al. 18 studied 0.17-1.0 ML coverages, finding no hydrogen bonds between adsorbed water in their 1.0 ML structure, and attributed it to the largely spaced sublattice of the ceria(111) surface, an effect that had earlier been hypothesized by Henderson et al. 6 In Ref. 18, each water molecule was estimated to occupy only 10 Å 2 of the available 13 Å 2 surface area, as calculated from the 3.88 Å distance between water adsorption sites, supporting Henderson's hypothesis. The fully dissociated layer has also been investigated but was found to be far less stable than the fully intact layer. Kropp et al. 21 proposed a new structure for the monolayer where the adsorbed water molecules form hydrogen bonds to each other in a water tetramer, with partial dissociation. Despite the increase in the total number of H-bonds in the system, the system of Kropp et al. was about as stable as those previously investigated. This structure included something not observed previously, namely, hydrogen bonding between adsorbed water, but not extending over to the periodically repeated image.
In this study, we found a water monolayer structure that is somewhat more stable than other published monolayer structures on ceria(111) [ Fig. 2(a), E ads = 0.773 eV/H 2 O]. It consists of an H-bond network of alternating intact and dissociated water molecules. Each OH f ion accepts two hydrogen bonds from two water molecules (R(H⋯O) = 1.79 Å) but does not  Table II). (b) 1.0 ML coverage, Q-0.33 (Structure A in Table II). For atom colors, see Fig. 1. , similar to what has been observed in calculations for the tetragonal-ZrO 2 (101) surface 39 and in a combined experimental-theoretical study for the RuO 2 (110) surface. 40 To ensure that this, our most stable, 1 ML structure was not constrained by the use of the (1 × 2) supercell, we also investigated a cell that was twice as large along the chain direction. The resulting structure was the same with respect to motifs and bond lengths.

The Journal of Chemical Physics
Furthermore, we performed a computer experiment to examine the interchain interaction in this structure. The shortest oxygenoxygen interchain distance is 5.0 Å between OH f and H 2 O in adjacent chains. Every second chain in Fig. 2(a) was removed, resulting in an oxygen-oxygen interchain distance of 11.5 Å after optimization and an E ads value of 0.740 eV/H 2 O. This is 0.03 eV/H 2 O less stable than the original 1.0 ML chain structure, indicating only a weak interaction between 1D-chains at the 1.0 ML coverage. The 0.5 ML structure so obtained is listed in its own right at the bottom of Table I. For coverages above 1.0 ML, we only studied a few structures and here, we also present a new stable structure for 1.5 ML coverage (C in II); the structure consists of a bottom layer, identical to the stable 1.0 ML structure, with an additional water molecule which bridges the chains (E ads = 0.729 eV/H 2 O). This new water molecule accepts a very long OH f ⋯OH 2 hydrogen-bond with R(H⋯ O) = 2.13 Å, and it donates two long H-bonds to the firstlayer water molecules with R(H⋯O) = 2.09 Å, as visualized in Fig. 2(b), thereby bridging the separate chains discussed for the monolayer. The resulting structure has a 2D-network of hydrogen bonds extending the surface.
The two stable 1.0 and 1.5 ML structures in Fig. 2 are each visualized together with a simulated STM image. These simulations serve two purposes. First, they provide an excellent way of illustrating patterns formed by the adsorbed water molecules. Second, they could provide a reference for future STM studies, although it may prove difficult to perform such experiments not least due to the low conductivity of CeO 2 . With the STM settings chosen (a negative bias of −3 V), it is the highest occupied electronic states of the adsorbed water molecules and OH f hydroxide groups that contribute to the contrast. For the structures presented in the current investigation, the brightness in the image is closely related to the topography of the structure.
In summary, we have presented two new extended H-bonded water structures on ceria(111): an infinite 1D chain structure for 1.0 ML coverage and an infinite network extending over two dimensions for 1.5 ML. A third new structure, our 0.5 ML structure, is also a stable, low-energy structure. It was obtained by a simple (at least computationally) modification of the monolayer. All three are listed at the end of Table I and constitute the A, 1 and B cases in Table II.  Figure 3 displays the adsorption energy per water molecule with respect to the degree of dissociation, Q, for our various variants of the monolayer coverage. The lowest-energy structures are interpolated with a line to highlight the convex hull. The inset shows the adsorption energy with respect to coverage. In it, the 1.0 ML is more stable than the structures found at both lower and higher coverages. All the data in the figure are also given in Table II. At 1.0 ML coverage, in the system of fully intact water molecules (Q-0), two stable structures were found (labeled 4 and 5 in Fig. 3). Although the adsorption energy difference is negligible (0.003 eV, i.e., the same in the convergence criterion) they show different bonding motifs. In structure 5 all water molecules bind to the surface but not to each other. Structure 4 is stabilized by forming H-bonded dimers but this will simultaneously destabilize the water-surface bonds. The net result is that 4 and 5 have virtually the same adsorption energy. Our attempts to extend the H-bond network were unsuccessful due to the large distances between adsorption sites. Upon geometry optimization of tested water chains, they reverted back to isolated water molecules or dimers. This limitation is overcome by partial dissociation of the surface water species. The OH f ions allow for an increased flexibility of H-bond geometries. A mix of intact and dissociated water molecules can then form H-bonded networks involving more than two water  Table II.

ARTICLE
scitation.org/journal/jcp molecules, and we find H-bonded tetramers at Q-0.25 and the 1Dchain described above at Q-0.5. Dissociation at Q-0.5 yields the most stable structure. From previous studies, we believe that hydrogen bonds strengths are generally ordered as HOH⋯OH f > HOH⋯OH 2 > OH f ⋯OH 2 (see Fig. 5 in Ref. 41). At Q-0.5, the number of strong hydrogen bonds can be maximized, as in the case of our 1.0 ML chain structure.
With increasing dissociation, destabilization from increasing hydroxide interactions occur but is initially stabilized by favorable hydrogen bonds in the adsorbed layer. As the dissociation surpasses Q-0.5, the possibilities for hydrogen bonding diminish and the destabilization will dominate. Note, however, that the energy differences between many structures at 1.0 ML coverage are small, as exemplified in Fig. 3, which is attributed to different hydrogen-bond arrangements. In reality, many of these structures may come into play simultaneously, even at low temperatures.

C. Comparison with experimental TPD
At least four experimental TPD studies of water desorption from the stoichiometric CeO 2 (111) surface have been published in the literature, as mentioned in the Introduction. Supported by XPS investigations, several spectral features have been identified and discussed, namely, a broad water peak, roughly centered around 200 K, and shifting toward lower temperatures for higher initial water coverage, 6-9 and a smaller peak with a maximum in the range of 270-290 K. [6][7][8] The latter peak has been interpreted as water molecules resulting from a re-combination of OH f and OsH hydroxide ions in an activated process, making it less connected to the adsorption energy and it will not be further discussed here. The broad large peak centered around 200 K has been interpreted as consisting of several contributions, namely, (i) a large peak below 200 K (maximum in the range 150-200 K in different studies [6][7][8][9], originating from multilayer water molecules not in contact with the oxide surface, and (ii) a shoulder believed to originate from the first-layer water molecules and peaking at about 200 K 6,7,9 and for very small coverages even at higher temperatures. 6 We have estimated the maximum position of the temperature programmed desorption peaks (simply denoted desorption temperatures here) based on our calculated structures and adsorption energies. A thermal correction was applied by adding the zero point energy contributions to the adsorption energy from the gas-phase water molecule and surface-bound water (decreasing the magnitude of the adsorption energy). This contribution was estimated from harmonic vibrational normal mode calculations that we performed using the phonopy package. 42 We obtained values 0.093 and 0.089 eV/H 2 O for structures 1 and 4, which are 1 ML structures with Q-0.5 and Q-0.0, respectively (cf. Table II; the values in the table are the uncorrected values). Rather than explicitly calculating the contribution for each monolayer structure, a representative value of 0.090 eV/H 2 O was applied to all those adsorption energies in the TPD analysis. The value 0.090 eV/H 2 O is consistent with Hofmann and Sauer's findings in Ref. 39 and, as said, decreases the magnitude of the adsorption energy.
The Redhead equation 43 was then used to relate the adsorption energy with the desorption temperature using an experimentally determined prefactor, ν, of 10 14.6 Hz, as given by Ref. 10, and a heating rate, β, of 2 K/s. First, we estimated a desorption temperature from our calculated E ads values in Table II for structures at coverages 1.0 ML (structure 1) and 2.0 ML (structure E) with the added thermal correction, subtracting one from the other according to 2 ⋅ E ads,2.0ML − E ads,1.0ML . This approach, one could say, yields a crude estimate of the "multilayer" feature (i) in the experimental TPD spectrum. Our desorption temperature obtained this way is 177 K, i.e., in the range of experimental observations (150-200 K).
Second, the desorption temperature for the water molecules in monolayer structure 1 was addressed. Our E ads value here does not refer only to desorption of intact water molecules as it is an average over the 50% intact water and 50% hydroxide species. Accepting this conceptual mismatch (whose effect may be modest), our 1.0 ML E ads value with the added thermal correction gives a desorption temperature of 227 K, slightly higher than the peak value of the comparable broad experimental peak (200 K).
In summary, our estimated desorption temperature values based on optPBE-vdW adsorption energies agree quite well with the experimental findings. Nevertheless, it should be noted that it is not straightforward to relate experimental desorption temperatures to adsorption energies (see, for example, the discussion by Campbell et. al. 10 ). Furthermore, the Redhead equation depends on the assumed parameter settings and will yield a 15 K shift in desorption temperature for an order of magnitude difference in the assumed parameters ν or β.

D. Surface energy
In order to compare the stabilities of surfaces with different water coverages, we have constructed a surface phase diagram (Fig. 4). The key quantity in the construction of such a diagram is the surface free energy, γ, as a function of the chemical potential of where n is the number of water molecules in the system, A is the area, and m is the number of CeO 2 formula units in the slab supercell. We also assume that volume changes in the slab can be neglected. The chemical potential of gas-phase water is referenced to the electronic (DFT) energy and includes the rotational, translational, and vibrational contributions for an ideal gas and is a function of the temperature and pressure of the gas. In this work, we have calculated this dependence from thermodynamic tables 44 in accordance with Refs. 39 and 45. As water is adsorbed on one of the two exposed surfaces of the CeO 2 slab, the γopposite term needs to be included; it is the surface free energy contribution of the clean surface opposite to the surface of adsorbed water. Figure 4 shows the resulting surface energy γ as a function of μ H 2 O . It is seen that, at low chemical potentials of water, up to −0.68 eV, the clean surface is more stable than the water covered surfaces. At chemical potentials in the window from −0.68 to −0.55 eV, the 1.0 ML coverage dominates. For a higher chemical potential, the 1.5 ML coverage is stabilized in a narrow window around −0.55 eV after which 1.75 ML dominates. Above −0.46 eV, the 2.0 ML coverage is stabilized.
In conclusion, our calculations suggest the monolayer could dominate under ambient conditions of 300 K and 10 3 Pa, where the latter value was selected to represent a typical water pressure at ambient conditions. At room temperature and ultrahigh vacuum, a clean surface would be expected. It is further noted that lowering the temperature or increasing the pressure results in a gradual buildup of thicker layers.

IV. CONCLUSION
In previous quantum-mechanical studies of thin film water adsorption on ceria(111), up to and including 1.0 ML coverage, much emphasis has been placed on the largely spaced ceria sublattice that tends to hinder efficient hydrogen bonding between the water molecules, making them essentially isolated, regardless of whether they are dissociated or not. One exception is the finite surface clusters observed in the QM calculations of Ref. 21. In the present study, we have investigated thirteen 1.0 monolayer structures with the aim to search for extended hydrogen-bonded networks and explore the influence of water dissociation on the interface stability. We report a new structure for the 1.0 ML coverage, consisting of a onedimensional ⋯OH f ⋯HOH⋯OH f ⋯HOH⋯ chain with an infinite hydrogen-bonded network in one dimension. This is the most stable structure reported for 1.0 ML on ceria(111). At 1.5 ML, we found another hydrogen-bond network, this time spanning the two surface dimensions.
We also related the adsorption energies for the first layer and for multilayer adsorption to experimental TPD measurements, finding a close agreement between our zero-point energy corrected E ads values calculated at the optPBE-vdW level, and experimental results. Finally, we calculated the stability windows for our 1.0 and 1.5 ML structures plus for two highercoverage structures at 1.75 and 2.0 ML, taken from the literature. With our computational method, the structure at 1.0 ML is the stable phase under ambient conditions (T = 300 K, p H 2 O = 10 3 Pa).