J. Phys. Chem. B 1998, 102, 2017-2026
2017
Adsorption of Methanol on TiO2(110): A First-Principles Investigation S. P. Bates*,† and M. J. Gillan Physics Department, Keele UniVersity, Staffordshire ST5 5BG, U.K.
G. Kresse Institut fu¨ r Theoretische Physik, Technische UniVersita¨ t Wien, A-1040 Wien, Austria ReceiVed: December 11, 1997
We have performed first-principles static and dynamic calculations based on density functional theory and the pseudopotential method to investigate the adsorption and deprotonation of methanol on the stoichiometric (110) surface of TiO2. Static calculations, employing full relaxation of adsorbate and substrate atom positions, are performed. In the high-coverage limit (θ ) 1), we find that there are several structures of approximately equal stability. In two of these, the methanol molecule is dissociated, resulting from scission of the O-H or C-O bonds. In the third, methanol is molecularly adsorbed. Other structures of approximately equivalent energy contain 1:1 mixtures of these conformations. At lower coverage (θ ) 1/2), we find that the two dissociative modes of adsorption found at θ ) 1 are favored over molecular adsorption by 19 kJ/mol (O-H scission) and 7 kJ/mol (C-O scission). The adsorption energy of the most stable θ ) 1/2 conformation changes by approximately (5% as the coverage is reduced to θ ) 1/3 and θ ) 1/4. Intermolecular attractions and repulsions are found to play a crucial role in determining the stability of different conformations at different coverages. Conversion of the metastable θ ) 1/2 molecularly adsorbed complex via O-H scission to a dissociated complex is predicted to be barrierless. First-principles molecular dynamics calculations on this system in which the methanol molecule approaches the surface predict spontaneous dissociation by rupture of the O-H bond and also that C-O bond breaking is likely to be an activated process. Further dynamical simulations indicate that the probability of finding conformations other than that obtained after O-H bond rupture is small.
1. Introduction In the past few years, first-principles calculations based on density functional theory (DFT) and the pseudopotential method have made increasingly significant contributions to our understanding of the nature of oxide surfaces1-9 and the interaction of these surfaces with small adsorbates.10-14 These methods are able to predict both thermodynamic and kinetic aspects of simple surface reactions, and the quality of the predictions that are obtained is now becoming comparable to that of experimental results, without being subject to the same real-world limitations. In this paper, we use these calculation methods to investigate the nature of the adsorption of methanol on the stoichiometric (110) surface of rutile (TiO2). TiO2 has received a considerable amount of attention in recent years, both experimental and theoretical. This is driven by a fundamental interest in this material as a model transition metal oxide and by the practical uses of TiO2, e.g., as a pigment in the paint industry or as a catalyst support. The structures of the low-index surfaces, of which (110) is the most stable, are now well understood.2,3,15 The (110) surface contains two different types of Ti (five- and six-coordinate) and O atoms (in-plane and bridging), as illustrated in Figure 1. (The labels for these atoms in Figure 1 are used throughout the following text.) These atoms are displaced from their bulk positions; inplane oxygen and six-coordinate titanium atoms are displaced † Present address: Department of Chemistry, West Mains Road, Edinburgh EH9 3JJ, U.K.
Figure 1. The TiO2(110) stoichiometric surface. Smaller light gray circles are titanium, larger dark gray circles are oxygens. Optimized bond lengths for the clean surface are given in angstroms.
outward, whereas bridging oxygens and five-coordinate titaniums are displaced inward. There is now semiquantitative agreement between the surface displacements predicted by theoretical methods and those determined experimentally.2,3,15 A more challenging problem is to understand the nature of the interaction of adsorbates with these surfaces. First-principles calculations that employ periodic boundary conditions are an effective way to model surfaces and the
S1089-5647(98)00499-4 CCC: $15.00 © 1998 American Chemical Society Published on Web 02/18/1998
2018 J. Phys. Chem. B, Vol. 102, No. 11, 1998 processes that occur on them. They have recently been applied to study the nature of the interaction of water with TiO2(110).1,10 Dissociative adsorption was found to be favored at low coverage (θ ) 1/2), with a structure containing equal numbers of molecules in dissociated and molecular conformations most stable at θ ) 1. The nature of the adsorbed water species was predicted to be a delicate balance of several factors, including coverage and temperature. The OH functional group of methanol may be expected to show modes of adsorption similar to that for water, but the methanol molecule can also adsorb in other ways, as we show later. One of the overall aims of our calculations on TiO2-adsorbate systems is to investigate the nature of adsorption as a function of acidity of the adsorbate molecule. The OH bond of methanol is substantially weaker than that of water,16 making it an ideal candidate for investigation. To our knowledge, no theoretical study of methanol on TiO2 has been published. There have, however, been a number of experimental investigations of this system. Onishi et al.17 used ultraviolet photoelectron spectroscopy (UPS) to ascertain the nature of methanol on the stoichiometric (110) and (441) surfaces of a single crystal at approximately θ )0.5. They concluded that methanol binds in a molecular state on both surfaces. In contrast, infrared (IR) measurements of Ramis et al.18 conclude the adsorption is dissociative, as indicated by an absence of hydroxyl stretching (ν) and torsional (δ) modes in their spectrum, obtained from a self-supporting disk of powdered rutile. A similar result is obtained by Suda et al.19 A temperature-programmed desorption (TPD) and UPS study by Kim and Barteau20 on (001) single-crystal surfaces suggests that a dissociative mode of adsorption predominates on the stoichiometric surface, whereas both dissociated and molecular adsorbates are present on the reduced surface. The coexistence of two such forms has also been suggested on the basis of studies on polycrystalline powder samples.21 More recent studies, such as the X-ray photoelectron spectroscopy (XPS) and TPD study of Aas et al.22 and the IR measurements of Hussein et al.,23 have considered the subsequent reaction of methanol on TiO2 and have identified methane to be the dominant product. In both cases, the methane originates from dissociatively adsorbed methanol in the form of a methoxide surface species. The lack of consistency in these findings clearly demonstrates the need for a thorough investigation. In particular, two key questions are raised by the experimental information available for this system. The first concerns the nature of the adsorbed methanol, prior to any reaction, and how this is influenced by factors such as coverage etc. The second question is the nature of the dissociatively adsorbed species. In most experimental studies, the O-H bond of methanol is found to break, generating a methoxy anion and a proton. However, scission of the C-O bond of methanol is also possible, yielding a methyl cation and hydroxyl anion as products. In the gas phase, this process is far easier than breaking the O-H bond.16 With the exception of Carrizosa et al.,24 who have suggested that surface CH3 groups are responsible for methane production by subsequent reaction, few authors have considered this mode of adsorption. One of the main aims of this paper is to provide answers to these questions. In a number of these experimental studies, the (110) surface is found to contain oxygen vacancies. Such a reduced surface will bind oxygen-containing adsorbates more strongly than the stoichiometric one. Our calculations are only concerned with the stoichiometric surface and can therefore be considered to be the lower limit of adsorption energies of methanol on “real” TiO2 surfaces. A similar argument follows for surface defects
Bates et al. such as steps or kinks; the perfect surface is likely to bind adsorbates more weakly than structural defect sites. An investigation into the behavior of the stoichiometric, defectfree surface forms the starting point for calculations that include vacancies or defects. (We mention in passing that firstprinciples calculations have recently been applied to study both electronic and structural aspects of defect sites on oxide surfaces1,11,25.) The remainder of the paper is organized as follows: the next section contains a thorough description of our calculation methodology before we present our results, discussion, and conclusions. 2. Calculation Methodology 2.1. General Details. We have performed density functional theory calculations that use a plane-wave basis set. The use of such a basis set necessitates the representation of tightly bound core electrons by a pseudopotential, while the valence electrons are treated explicitly. The general principles of the method are well-documented elsewhere26-31 and the approach has been successfully applied to a wide variety of oxide surfaces1-9 and molecular adsorption processes.10-14 We have used the VASP code32-34 to perform both static calculations (with full relaxation of all atoms) and ab initio molecular dynamics calculations. Ultrasoft pseudopotentials35 are used, allowing us to use a considerably smaller basis set than would be needed with standard norm-conserving pseudopotentials. The size of the basis set, i.e., the number of plane waves used, is limited by inclusion of only those plane waves whose energy is less than a certain cutoff value (referred to as the plane-wave cutoff energy in the rest of the paper). The Ti and O pseudopotentials used in this work have been shown to predict zone-center vibrational modes for bulk TiO236 in very good agreement with experiment.37 These pseudopotentials are identical to those used in our previous study of the surface structure of TiO2(110).2 Bulk TiO2 lattice parameters (a and c) are also well reproduced using these pseudopotentials.2 In all our calculations, we have used the generalized gradient approximation (GGA) of Perdew and Wang,38 known as GGAII, as this improves the accuracy of the energetics of molecular dissociation when compared with the local density approximation (LDA). The oxide surface is represented as a slab of finite thickness to which full periodic boundary conditions are applied, thus generating an infinite stack of quasi-two-dimensional slabs, each separated from its neighbors by a certain vacuum width. We have previously shown2 that the separation of neighboring slabs must be at least 5 Å and the thickness ideally five or more TiO2 layers (where a “layer” is defined as a (110) plane that contains both Ti and O atoms). All our calculations employ vacuum widths of at least 5 Å or larger, and our static calculations use slabs of five layers. As a tradeoff against computational cost, our molecular dynamics simulations employ three-layer slabs. Our experience with calculations on the bare (110) surface2 and water adsorbed on the (110) surface39 suggests that although use of a three-layer slab may introduce significant quantitative uncertainties, qualitative features are preserved. This is confirmed by molecular dynamics tests described later in this section. Bulk TiO2 is tetragonal and is described by two lattice parameters, a and c; thus, the bulk unit cell has dimensions a × a × c. The (110) surface unit cell lies in the plane formed by the (1h10) and (001) directions, as illustrated in Figure 1. The dimensions of this surface unit cell are x2a × c (along
Adsorption of Methanol on TiO2(110)
J. Phys. Chem. B, Vol. 102, No. 11, 1998 2019
TABLE 1: Comparison of Calculated and Experimental Bond Lengths, r, Bond Angle, ∠COH, Dipole Moment, G, Homolytic O-H Bond Dissociation Energy, E(hom), and Heterolytic O-H Bond Dissociation Energy, E(het), for Methanola r (OH) (Å) r (CO) (Å) r (CH) (Å) ∠COH (deg) F (D) E(hom) (kJ/mol) E(het) (kJ/mol)
calc
expt16,46,47
0.9738 (+1.1) 1.4320 (+0.8) 1.0998 (+0.5) 107.8 (-0.2) 1.756 (+3.2) 424 (-2.9) 1582 (-0.9)
0.963 1.421 1.101 108.0 1.70 437 1597
a Figures in parentheses are the percentage deviation from the experimental value.
(1h10) and (001), respectively), where a and c are the same bulk lattice vectors mentioned above. Most of our calculations use surface unit cells of x2a × nc, which we refer to as multiples of the smallest surface unit cell, i.e., 1 × n. This corresponds to a surface coverage of methanol of θ ) 1/n. In addition, we have used a surface unit cell of 2x2a × 2c ( i.e., 2 × 2) with a surface coverage of methanol of θ ) 1/4. 2.2. Confidence Tests. We have previously reported extensive tests of the DFT-pseudopotential method applied to the surface energy and structure of the stoichiometric (110) surface.2 We have also investigated the structural and energetic properties of the other component of the system, i.e., the isolated methanol molecule. To ensure comparability, we have used exactly the same DFT-pseudopotential method. The calculation details were identical to those used for all the calculations of the bare surface, and each methanol molecule was placed at the origin of a cubic cell of side 7 Å. (This procedure has recently been used to evaluate the performance of the DFTpseudopotential method against an atom-centered Gaussian basis set for a number of small molecules.40) Further technical details of our calculations will be described elsewhere,39 but Table 1 shows the calculated structural and energetic properties of the gas-phase methanol molecule. We find that the experimental values are well reproduced by the calculations, in particular the heterolytic splitting of the O-H bond. 2.3. Static Calculations. The static relaxation calculations used a plane-wave cutoff energy of 400 eV and the full-coverage calculations used the lowest-order Monkhorst-Pack41 set of two k-points, with the component parallel to the surface normal set to zero, i.e., fractional coordinates of ((1/4, 1/4, 0) in the coordinate system of the surface unit cell. Both cutoff energy and k-points were the same as we have previously found to accurately describe the surface energy and structure of stoichiometric (110).2 At the lower coverage of θ ) 1/2, doubling the unit cell size along the (001) direction allowed us to halve the number of special k-points in this direction. Thus, we used one k-point with reciprocal coordinates (1/4, 1/2, 0) for coverages of θ ) 1/2 or less. Iterative relaxation of atomic positions was stopped when the change in total energy between successive steps was less than 0.001 eV. With this criterion, forces on the atoms were generally less than 0.1 eV/Å. No symmetry was used in any of the calculations. Adsorption energies were calculated as the difference between the surface-methanol system energy and the sum of the bare surface and gas-phase methanol molecule. A positiVe adsorption energy is favorable and indicates stabilization. 2.4. Barrier to Proton Transfer. We have calculated the energy barrier for transfer of the methanol proton to/from the
(110) surface. To do this we have used the nudged elastic band method (NEB) of Jo´nsson and Mills.42,43 The NEB method represents a general way of locating transition-state structures between two minima. A series of structures is generated between initial and final states on either side of the barrier, and these are then connected by a series of “elastic bands”. The ionic positions of all the images are iteratively relaxed simultaneously until the forces on these ions (which include an elastic band force resulting from the position of ions in adjacent images) are sufficiently small to satisfy the optimization criteria. A transformation of the force acting on each image in such a way as to zero out the perpendicular component of the force between images ensures that the minimum energy path between initial and final states is found. We have used a system of four images, created by linear interpolation between the coordinates of initial and final states. Minimization criteria and calculation parameters (plane-wave cutoff energy etc.) were the same as for the respective static calculations from which we obtained the initial and final conformations. 2.5. Dynamic Simulations. We have also used the VASP package to perform ab initio molecular dynamics (MD) calculations. We have performed two different types of calculations in which the molecule is either sent in from the gas phase or adsorbed on the surface in thermal equilibrium. These two types of calculation provide information on the dynamics of the adsorption process and the probability of finding different conformations of the molecule on the surface. In an MD calculation, the forces on the ions are determined and used in the integration of Newton’s equations of motion. In our ab initio MD calculations, the forces are found from the evaluation of the exact electronic ground state at each time-step and we use a Verlet44 algorithm to move the ions. More details can be found elsewhere.32-34 The efficiency of this method is comparable to a conventional Car-Parrinello-based MD simulation,45 in which electrons and ions are updated at the same time. Our dynamics calculations are intended to provide semiquantitative insight into the dynamic behavior of the methanolTiO2 system. We have reduced the expense of these simulations by lowering the plane-wave energy cutoff from 400 to 300 eV and using a slab thickness of three rather than five layers. We have performed extensive tests to ascertain the effects of these changes. We find that a reduction in the plane-wave energy cutoff has a very small effect on the results of the dynamics calculations. A thinner oxide slab has a more pronounced effectsthe distortions and displacements of the TiO2 atoms upon adsorption of methanol tend to be more pronounced on a threelayer slab than on one with five layers. This is in accordance with what we have previously found for the displacements of surface atoms on the bare surface as a function of slab thickness.2 However, we find that the same dynamical phenomena are shown by both slab thicknesses. We have used a time-step of 1 fs and increased the mass of the hydrogen atoms to 3 amu. This increase of the hydrogen mass has also been used in similar MD calculations of water on TiO2.10,1 Time-dependent properties such as vibrational frequencies will of course be changed, as will the attempt rate of dynamical processes involving the hydrogen atom (e.g., proton transfer). However, time-independent properties will remain unchanged. With this combination of parameters, we are able to perform simulations of the order of a few ps for coverages of θ ) 1 and 1/2. Typically, 2 ps is long enough for the system to explore the majority of the phase space. For dynamical simulations at a coverage of θ ) 1, the same set of two k-points was used as in the static calculations. Using
2020 J. Phys. Chem. B, Vol. 102, No. 11, 1998
Bates et al.
Figure 2. The different bonding modes of methanol on the (110) surface. The series of diagrams on the left are side-on views of the surface, those on the right are top views. The labels M1, M2, etc. are defined in the text. The shading of Ti and O atoms is identical to Figure 1, while the methanol carbon and hydrogen atoms are represented by mid-gray and white circles, respectively.
only one k-point produced far larger forces on the ions in this small unit cell. For all lower coverages, the Brillouin zone was sampled only at the Γ-point. For all calculations, unless stated otherwise, we used a temperature of 600 K. 3. Results 3.1. Static Calculations. 3.1.1. Full-Coverage Adsorption (θ ) 1). In this section, we present results for calculations performed with one methanol molecule per surface unit cell of dimensions x2a × c. The relative sizes of the surface unit cell and the methanol molecule mean that this corresponds to full surface coverage. This correspondence between the size of the methanol molecule and the size of the smallest surface unit cell reduces the number of possible conformations of the methanol on the surface. We identified six candidate structures on the basis of favorable electrostatic interactions, i.e., matching up electrostatically negative parts of the adsorbate molecule with electrostatically positive regions of the surface and vice versa. Two of these candidate structures optimized to the same conformation, hence only five different conformations were obtained.
These conformations are illustrated in Figure 2. Two of these involve methanol binding molecularly to the surface in different orientations, which we label M1 and M2. The key difference between these two is that in configuration M1 the hydrogen of the methanol hydroxyl points toward the bridging oxygen, while in M2 it points toward the in-plane oxygen. The dissociated conformations are labeled D1, D2, and D3. Conformation D1 involves scission of the methanol O-H bond. The resulting fragments, H+ and CH3O-, bind to the bridging oxygen and 5-fold titanium atoms, respectively. Conformation D2 involves C-O bond breaking, with the OH- and CH3+ fragments bound to the 5-fold titanium and bridging oxygen, respectively. The final conformation, D3, is a variant of D1 in which the CH3O and H fragments lean toward each other, facilitating formation of a hydrogen bond. The optimized bond lengths and adsorption energies are given in Table 2. The molecular conformations M1 and M2 are similar in geometry and adsorption energy. In conformation M1, the methanol acidic proton forms a hydrogen bond with the surface bridging oxygen, while in M2 this hydrogen is closest to the in-plane oxygen of TiO2. This twist or rotation
Adsorption of Methanol on TiO2(110)
J. Phys. Chem. B, Vol. 102, No. 11, 1998 2021
TABLE 2: Optimized Bond Lengths (in Å) and Adsorption Energies (Eads, in kJ/mol) for the Different Conformations of Methanol on TiO2(110)a θ)1
θ ) 1/ 2
Cm-Om Om-Hm Hm-Ob Om-Ti5 HMe-Ob
1.434 0.997 1.905 2.189 2.714
1.438 1.003 1.772 2.132 3.072
Cm-Om Om-Hm Hm-Ob Hm-Oip Om-Ti5
1.454 0.987 2.759 2.349 2.342
Cm-Om Om-Ti5 Hm-Ob Ob-Ti6
1.377 1.821 0.969 2.042
Cm-Ob Ob-Ti6 Om-Hm Om-Ti5
1.428 2.022 0.991 1.905
1.423 2.020 0.957 1.812
Cm-Om Om-Hm Hm-Ob Om-Ti5 HMe-Ob
1.408 1.778 1.006 1.884 2.973
1.413 1.837 0.998 1.850 2.986
θ)1
θ ) 1/2
M1 Ob-Ti6 Ti6-Oip Oip-Ti5 Eads
1.874 2.013 1.933 70.4
1.881 2.036 1.933 100.8
M2 HMe-Ob Ob-Ti6 Ob-Ti6 Oip-Ti5 Eads
2.397 1.864 2.022 1.945 64.6
D1 Ti6-Oip Oip-Ti5 Eads
1.962 1.981 44.4
D2 Ti6-Oip Oip-Ti5 Eads
1.984 1.951 78.2
2.015 1.977 107.8
D3 Ob-Ti6 Ti6-Oip Oip-Ti5 Eads
2.032 1.991 1.960 74.3
2.021 2.022 1.969 119.1
a The conformations M1, M2, etc. are described in the text and shown in Figure 2. The labels for the surface TiO2 atoms are defined in Figure 1. Labels Hm, HMe refer to the methanol acidic proton and the nearest proton of the CH3 group to the surface, respectively.
that forms M2 from M1 weakens the interaction of the methanol oxygen with the surface (as shown by the bond lengths in Table 2). This is the most likely reason for the slightly reduced stability of M2 compared to M1. In both conformations, the molecule is physisorbed on the surface and this is reflected in the relatively small perturbations of Ti-O bond lengths between the surface atoms. The high coverage results in neighboring methanol molecules that are close to one another, resulting in repulsive interactions between methyl hydrogens of adjacent molecules, which can approach as close as 1.8 Å. The range of adsorption energies found for the dissociated conformations is somewhat larger. Conformation D1, which corresponds to scission of the O-H bond, yields the smallest adsorption energy. The separation of the fragments is too large to facilitate any intermolecular bonding, and there are repulsions between neighboring methyl groups. The perturbation of the surface is greater than for the molecular conformations; the bond between the 6-fold titanium and the bridging oxygen lengthens by approximately 10% on adsorption, with neighboring bonds between the 6-fold titanium and the in-plane oxygen contracting. The 5-fold titanium is strongly displaced out of the surface in all three dissociated conformations. The tilted counterpart to D1, D3, is more strongly bound by approximately 30 kJ/mol. This stabilization by tilting symmetrically arranged dissociation products has been predicted previously from similar static calculations of water on TiO2 by Goniakowski and Gillan13 and the MD calculations of Lindan et al.10 The reason for the stabilization is the formation of a strong hydrogen bond (1.78 Å) between the methanol hydrogen fragment and the CH3O fragment. The O-H and C-O bonds lengthen slightly, and Figure 2 shows significant buckling of
the surface. The uppermost surface atoms move in such a way as to facilitate the close approach of the two fragments. The similarity between conformations D3 and M1 is striking, and we return to this later when we consider the process of proton transfer. The most stable conformation of all five studied is D2, and there are two reasons for this. The first is that the C-O bond is weaker than the O-H bond (heterolytic dissociation energies are 1160 and 1596 kJ/mol in the gas phase, respectively16). A second reason for the larger stabilization is the intermolecular bonding between neighboring O-H fragments. The hydrogen of one O-H fragment is separated from the oxygen of its neighbor by only 2 Å, resulting in a chain of O-H fragments parallel to the (001) direction. We have also performed calculations on a unit cell doubled in length in the (001) direction (x2a × 2c) with two molecules of methanol in different conformations to each other. Thus, a coverage of θ ) 1.0 is maintained, but we investigate the stability of “mixed” structures, which similar calculations have predicted to be most stable for water on TiO2(110).1 We have restricted our considerations to combinations of the three most stable conformations, i.e., M1, D2, and D3. We find that all three combinations (M1/D3, M1/D2, and D2/D3) are of approximately equal stability (with adsorption energies between 74 and 77 kJ/mol per methanol molecule). In the calculations involving conformation D2, the favorable intermolecular hydrogen-bonding chain along the (001) direction is destroyed by every other molecule being in a different conformation and the O-H bond vector of the D2 fragment is now normal to the surface. This destabilization by removal of hydrogen bonding must be balanced by a less strained or perturbed surface. These calculations on mixtures further confirm that there is little difference energetically between the three structures M1, D2, and D3 at full-coverage. 3.1.2. Half-Coverage Adsorption (θ ) 1/2). The conformations described in the previous section are all subject to a significant amount of intermolecular repulsion due to close proximity of neighboring methanol molecules. We have performed a series of calculations with half-coverage (cell dimensions x2a × 2c and 1 methanol molecule per unit cell), restricted to the three most favorable conformations described above, i.e., M1, D2, and D3. The optimized bond lengths for these half-coverage calculations are given in Table 2. In general, these are very similar to what was found for the full-coverage case. The half-coverage adsorption energies are larger than those at full-coverage; the reasons for this are the removal of the unfavorable intermolecular repulsion present at full-coverage and the less perturbed surface. We now find that the D3 conformation is the most stable and the difference in energy between this and the most stable molecular conformation increases to nearly 20 kJ/mol. The D2 half-coverage calculation yields two interesting features worthy of discussion. The first is that the favorable intermolecular interactions are lost when half the methanol molecules are removed. The O-H bond vector is found to be perpendicular to the surface, as in the mixed-full-coverage calculations described above. This loss of favorable interaction destabilizes the complex somewhat, and accounts for the fact that the D2 conformation is found to be less stable than D3 at half-coverage. The second feature concerns the location of the dissociated fragments of the D2 conformation. For a dissociated conformation at a coverage of θ ) 1/2, the dissociated fragments may be separated by either a/x2 (in the same 1 × 1 unit cell) or
xa2/2+c2 (in adjacent 1 × 1 unit cells).
We find that the
2022 J. Phys. Chem. B, Vol. 102, No. 11, 1998
Figure 3. A schematic view of the pattern of surface atomic displacement caused by half-coverage adsorption of methanol in the D2 conformation. Dashed lines represent the surface unit cell, while dotted lines show bonds between atoms. The atom labels (5, 6, b, and ip) are defined in Figure 1, and Ti and O atoms with bold edges indicate the locations of the OH and CH3 fragments, respectively (these are not shown for clarity). Arrows indicate the direction of displacement of the surface atoms, and the values on the arrowheads are the displacement in the (110) direction (in angstroms), relative to the calculated positions of the atoms for the clean, relaxed surface. Part a represents the case where the two dissociated fragments of methanol are placed in the same 1 × 1 surface unit cell, while in part b they are in neighboring cells.
larger separation gives the more favorable adsorption energy by almost 10 kJ/mol, but not because of repulsions between fragments. The reason for the enhanced stability at the greater separation is due to the TiO2 surface being less strained. At a fragment separation of a/x2, rows of Ti5-Ob-Ti5 atoms in the (1h10) direction are alternately highly strained (due to methanol fragments displacing these surface atoms upward) and relaxed (no adsorbate). At the larger fragment separation, each row in the (1h10) direction contains both strained and relaxed atoms. An alternate way of expressing this is that the two fragment separations generate a surface atom displacement pattern that is either striped (smaller separation) or checked (larger separation). The latter is more stable, and it is this value that is included in Table 2. The pattern of surface atomic displacements is shown schematically for both cases in Figure 3. This choice of separation of the fragments is not relevant for the other dissociated conformation we have studied at a coverage of θ ) 1/2 (D3), since separating the fragments removes the hydrogen bond between them, and the structure optimizes to a separated version of D1. Half-coverage is also achieved by adsorbing one molecule of methanol in a unit cell of dimensions 2x2a × c. This yields alternately filled and empty rows of adsorbates along (001). It is unlikely that this arrangement will produce structures that are more stable than those we have described above; the
Bates et al. main source of intermolecular repulsions due to close approach of neighboring methyl groups is still present. We note in passing that coverages lower than θ ) 1 allow the possibility of a number of other conformations of the methanol molecule in which the C-O bond vector is aligned approximately parallel to the (001) direction (rather than along (1h10) or perpendicular to the surface, as in the five conformations in Figure 2). We have tried to locate other adsorption minima, but have found that the structures are either significantly less stable than D3, D2, and M1 or collapse into one of these structures during optimization. 3.1.3. Lower Coverage Adsorption (θ ) 1/3,1/4). We have performed lower-coverage calculations, to investigate the convergence of the adsorption energy to the low-coverage limit. We have restricted ourselves to the most stable conformation at θ ) 1/2, D3. Using 1 × 3 and 1 × 4 unit cells, we have obtained adsorption energies of 115 and 125 kJ/mol, respectively. The adsorption energy of the D3 conformation as a function of surface unit cell size (and hence surface coverage) is shown graphically in Figure 4. The adsorption energy is converged to within (5% at a coverage of θ ) 1/2. We have also investigated θ ) 1/4 coverage using a 2 × 2 surface unit cell, occupied by a single methanol molecule. We found an adsorption energy for the D3 conformation of 105 kJ/mol. The reason for the smaller adsorption energy for θ ) 1/ using the 2 × 2 unit cell, compared to the 1 × 4 unit cell, 4 is that the hydrogen bond between dissociated fragments is longer and therefore weaker. The 2 × 2 surface is far less puckered than the 1 × 4 surface as a result of this longer hydrogen bond. Thus, the systems appears to show a thermodynamic preference toward the formation of adsorbate rows along (1h10). We could reasonably assume that the same trends in the adsorption energy as a function of coverage would be found for conformation M1. Conformation D2, however, might be somewhat different because of the choice of position of dissociated fragments. However, it is unlikely even if this is the case that the overall energetic picture of low-coverage adsorption will change significantly. 3.1.4. Barriers to Proton Transfer. In addition to the thermodynamic information we have presented in the previous sections, we have also investigated the barriers to proton transfer to/from the methanol molecule to gain insight into the kinetics of the methanol-TiO2 system. As stated earlier, we have used the NEB method, which requires initial and final conformations on either side of the barrier to be known, and the M1 and D3 structures, which are shown in Figure 2 to be very similar with the exception of the location of the acidic methanol hydrogen. We have investigated the barrier to proton transfer at two coverages, θ ) 1 and θ ) 1/2. For the full-coverage case, we used a slab comprising three layers as a preliminary test. The half-coverage barrier calculation used the larger five-layer slab. In both cases, we found that the barrier to proton transfer was negligible (less than 5 kJ/mol). 3.2. Dynamic Simulations. In addition to the static relaxation of ionic positions, we have also performed ab initio MD simulations for methanol coverages of θ ) 1 and θ ) 1/2. As stated earlier, we have used MD simulations to semiquantitatively investigate both the dynamics of the adsorption process (which is done by sending the molecule in from the gas phase above different regions of the surface) and the probability of finding different conformations (which is done by a simulation of the methanol molecule in thermal equilibrium with the surface).
Adsorption of Methanol on TiO2(110) 3.2.1. Full-Coverage MD Simulations (θ ) 1.0). To make a preliminary investigation into probabilities of different conformations, we have performed full-coverage MD simulations. The full-coverage calculations were intended as an initial test before considering lower coverages but there are features of the MD trajectories from these introductory calculations that merit discussion. A simulation at 600 K of duration 2 ps was performed, with the initial conformation of the methanol molecule chosen to be D3. During the course of the simulation, we found that the surface exhibits a considerable amount of distortion, however this is not unexpected given the simulations were performed on a three-layer slab. We did not see the reformation of the methanol molecule via donation of the proton from the surface to the methoxy fragment. The Ob-Hm bond was lengthened significantly at various times, but was never broken. In these three-layer calculations, the energy difference between M1 and D3 conformations was 15 kJ/mol, substantially larger than the 4 kJ/mol found for the five-layer calculations (Table 2). Therefore, at a temperature of 600 K, the M1 state is far less likely to be populated compared to the D3 state, and it is not surprising that we see no proton transfer back to the methanol molecule during our simulations. The fact that the energy gap between these two conformations is far larger for a three-layer slab illustrates the quantitative inaccuracies introduced with thin slabs. We have observed similar behavior for water adsorption on TiO2.39 Examination of the simulations shows that the potential energy undergoes large fluctuations. The structures associated with configurations of highest potential energy (i.e., least stable) are characterized by unfavorable interactions between methyl hydrogens on neighboring methanol molecules. In certain case, hydrogens atoms approached as close as 1.8 Å, resulting in a structure that was approximately 50 kJ/mol less stable than the average energy of the system in a D3-like conformation. A subsequent simulation at 800 K exhibited broadly the same features as that at 600 K. 3.2.2. Half-Coverage MD Simulations (θ ) 1/2). More comprehensive simulations were carried out at a methanol coverage of θ ) 1/2. Two different sets of simulations were considered: the methanol molecule approaching different regions of the surface from the gas phase, and a simulation of the molecule dissociated on the surface in the D3 configuration (predicted to be the most stable from our half-coverage static calculations described earlier). In simulating the approach of the methanol molecule to the surface, we performed four separate calculations. These were: the oxygen of the methanol molecule directly over the 5-fold titanium of the surface (to bias the conformation to that of M1 or D3); the oxygen of the methanol molecule directly over the in-plane oxygen; the methanol oxygen directly over the bridging oxygen and finally the C-O bond directly over the bridging oxygen (to bias the conformation to that of D2). The first three of these differ only in the starting position of the methanol molecule along the (1h10) direction, and they are represented by the labels I, II, and III in Figure 5a. In these three simulations, the methanol molecule had the same orientation (C-O bond vector pointing toward surface) but was positioned above different surface atoms. In the fourth simulation, the orientation was altered so that the C-O bond vector was parallel to the surface and the midpoint was directly above the bridging oxygen atom. In all cases, the methanol molecule started at a distance of approximately 3 Å above the surface, and the methanol atoms were not given initial velocities to direct the molecule toward the surface but were allowed to “find” the
J. Phys. Chem. B, Vol. 102, No. 11, 1998 2023
Figure 4. The adsorption energy of the D3 conformation as a function of surface unit cell size. One molecule of methanol is placed in each unit cell, thus the abscissa represents a surface coverage of 1/n.
Figure 5. Snapshots of the half-coverage MD simulation in which the methanol molecule is sent in from the gas phase with its oxygen atom above the 5-fold titanium atom of the surface. In part a, a hydrogen bond forms between acidic methanol hydrogen and surface bridging oxygen. In part b, the methanol oxygen interacts with the 5-fold titanium. In part c, the methanol proton is transferred from the methanol molecule to the bridging oxygen (part d). Atom colors are as defined in Figures 1 and 2. The labels I, II, and III indicate the positions of the methanol molecule at the start of simulations.
surface during the course of the simulations. The reason for this was our aim to avoid introduction of a bias toward a particular conformation by choice of initial velocities. For the simulation in which the methanol was positioned directly above Ti5, the molecule moves toward the surface and a hydrogen bond develops between methanol hydrogen and bridging oxygen within 0.2 ps. The methanol oxygen interacts with the Ti5 atom, and scission of the O-H bond follows. Snapshots of this process are shown in Figure 5. When the methanol molecule is started above the in-plane oxygen (a displacement of a/2x2 along (1h10) compared to the previous simulation) the methanol molecule takes longer to find the surface. It first drifts back toward the Ti5 position halfway between neighboring bridging oxygen atoms and rotates around the C-O bond. A similar trajectory to that described previously then evolves, with the acidic proton first hydrogen bonding to a bridging oxygen, followed rapidly by dissociation of the O-H bond. The third simulation of this series started with the methanol oxygen directly above the bridging oxygen of the surface. During the course of the simulation, the molecule is repelled by the surface and “bounces” away into the vacuum gap. (In actual fact, it crosses the vacuum gap between slabs
2024 J. Phys. Chem. B, Vol. 102, No. 11, 1998
Bates et al. in both conformations. As in the case of the full-coverage simulations, we do not observe proton transfer back to the methoxy fragment. For the majority of the simulation time, the methanol molecule is in a conformation which is D3-like (Figure 6d).
Figure 6. Snapshots of the half-coverage MD simulation in which the methanol molecule is dissociated in the D3 conformation. Figure 6a illustrates the conformational freedom of the OCH3 groupsthe C-O bond vector lies completely along (001). Parts b and c of Figure 6 illustrate the “flipping” of the methanol H+ fragment. Figure 6d shows a typical D3-like conformation, in which the system spends the majority of the time. Atom colors are as defined in Figures 1 and 2.
and finds the lower surface of the next slab and dissociates on this surface by O-H bond scission.) The final simulation of this series started from a conformation of the methanol biased toward the formation of D2. We found that the methanol molecule did not dissociate via C-O bond breaking. Instead, the molecule once again migrated to a position approximately midway between neighboring rows of bridging oxygens and attached to the surface by interaction of the methanol oxygen with the 5-fold titanium atom. Dissociation to the D3 conformation was then observed. This series of MD simulations starting from different initial positions illustrates that the most probable structure of the methanol molecule on the surface at half-coverage was D3 (as we also found from the static calculations). To show that the D3 dissociated conformation really is the most probable and does not transform into another conformation, we have performed a 2 ps calculation with the methanol molecule starting in conformation D3. We do find that the system remains in this conformation for the duration of the simulation. The most notable features of the simulation are shown in a series of snapshots in Figure 6. One of the most striking features of the animated trajectory of the atoms is the freedom of movement enjoyed by the OCH3 group. As can be seen from Figure 6, the C-O bond vector has a large freedom of movement. In the full-coverage calculations the O atom positions were confined to an arc of a circle in the (1h10)-(110) plane, whereas at half-coverage this atom is able to move in the (001) direction as well, and the locus of positions falls on a disk. Figure 6a shows a conformation in which the C-O bond is leaning strongly toward the (001) direction, while Figure 6d shows it leaning strongly in the (1h10) direction. In addition to this motion, the CH3 group rotates very freely about the C-O bond, as can be seen from Figure 6. We have previously stated that in the full-coverage MD simulations, repulsive interactions between hydrogen atoms on neighboring molecules cause a large steric hindrance. With half the methanol molecules removed, this hindrance disappears completely and this gives the larger freedom of movement that we observe. Another striking feature of the snapshots is the mobility of the H+ fragment. Several times during the course of the simulations the methanol H+ fragment “flips”, from one conformation that is D3-like to another D3-like conformation in which the O-H bond vector along (1h10) points in the opposite direction (Figure 6b,c). The hydrogen bond between the H+ and OCH3- fragments is present
4. Discussion Our static calculations predict that the most stable conformation of methanol on TiO2 is sensitive to the surface coverage and also that no single conformation is unambiguously preferred over others. At full monolayer coverage, there are four different single-conformation structures within an energy window of 13 kJ/mol. The two most stable conformations are separated by only 4 kJ/mol, and in both cases the methanol molecule is predicted to be in a dissociated state. In our confidence tests on the dissociation energies of the gas-phase methanol molecule, we were able to reproduce experimental values within 15 kJ/ mol. If we assume approximately the same magnitude of uncertainty in our calculations of methanol adsorbed on the (110) surface, we find that we cannot conclusively state which conformation will prevail. Further confirmation of how complicated the picture is at full-coverage is provided by the mixed conformation calculations, which comprise 1:1 mixtures of the three most stable single conformations. All three of these possible conformations are well within the original energy window for the four single-component calculations. At a coverage of θ ) 1/2, the most stable conformation corresponds to dissociation of the O-H bond, with the two dissociated fragments leaning toward each other to generate a hydrogen bond. This is predicted to be 12 kJ/mol more stable than another dissociated conformation, in which the C-O bond is broken, which is in turn 7 kJ/mol more stable than the most stable molecular conformation. Once again, the energy differences are rather small and comparable to the uncertainty of our calculations. However, our findings are in agreement with the majority of experimental data that report the methanol molecule to be dissociated on the surface. Experimental evidence that both molecularly and dissociatively adsorbed complexes can coexist on the surface is supported by the small energy differences we calculate and also the effectively barrierless transfer of the methanol proton to/from the surface. There are several factors that govern the relative stability of different conformations at different coverages. An obvious factor is intermolecular interaction, both attractive and repulsive. At full-coverage, all conformations are sterically hindered by the close approach of neighboring methyl hydrogens. This is demonstrated by the large increase in adsorption energy for all conformations when the coverage is reduced from θ ) 1 to θ ) 1/2. Full-coverage intermolecular interactions are not all bad, however, as the D2 and D3 conformations show. Hydrogen bonding between dissociated fragments in the same molecule (D3) or in neighboring molecules (D2) contributes significantly to the overall stabilization of these conformations. A second factor that controls overall stability is the displacement of surface atoms upon adsorption of methanol. At a coverage of θ ) 1/2, we found that different arrangements of the dissociated fragments of conformation D2 produced adsorption energies that differed by 8 kJ/mol, and the origin of this difference was the displacement of the surface atoms. A final factor governing stability is the strength of bond in the methanol molecule that is being broken upon dissociation. We have previously stated that there is a large difference in the gas-phase heterolytic dissociation energy of the C-O and O-H bonds. Different combinations and magnitudes of these factors are present for different conformations.
Adsorption of Methanol on TiO2(110) The overall stability of any one conformation depends on entropic terms as well as enthalpic. At higher temperatures, the entropy contributions become more important. When comparing the M1 and D3 conformations, the atoms of the methanol molecule in the M1 structure are far more confined than those in the D3 structure. The freedom of movement of the D3 fragments, especially at half-coverage, suggest this conformation is likely to be favored on entropic as well as enthalpic grounds. The results from our molecular dynamics calculations are entirely consistent with those from the static relaxations. Stable structures predicted from our static calculations are observed in the snapshots of the MD trajectory. In addition, finite temperature effects are observed in the dynamics simulations that are not accessible to our static calculations. Examples of these effects include the fluxional behavior of the methanol proton shown in Figure 6. We observed the methanol molecule dissociating spontaneously on the surface via O-H bond breaking. This corroborates the results of our proton-transfer calculations which showed this process to be nonactivated. Our static calculations suggest that there is little difference energetically between C-O and O-H bond scission, yet we do not observe the C-O bond breaking in the dynamics calculations, despite biasing a simulation toward this conformation. Thus, it would appear that despite being energetically favorable, the scission of the C-O bond is activated on the (110) surface. This is a possible explanation as to why this mode of adsorption is not observed experimentally. There are strong similarities between our calculation results for the adsorption of methanol on TiO2 and previous calculations of water adsorption.10,1 In both cases, a balance between molecular and dissociated modes of adsorption exists. We are currently studying the adsorption of other ROH molecules on TiO2, including hydrogen peroxide (HOOH) and formic acid (HCOOH). The results of these calculations will be reported in a future publication.39 These calculations represent only the initial stages of methanol sorption on TiO2(110); conversion of methanol to formate and then methane is observed experimentally. The calculations presented here simulate surface coverages down to θ ) 1/4, which is well within the range of surface coverages usually investigated experimentally. The present pace of development of computer hardware and software will make more complex surface reactions at lower coverages far more accessible within a few years. 5. Conclusions We have performed first-principles static and dynamic calculations to investigate the adsorption of methanol on the stoichiometric surface of TiO2(110). The conclusions we draw from these calculations are summarized below. There are several different adsorption conformationssboth molecular and dissociativesthat are within a narrow energy window for surface coverages of both θ ) 1 and θ ) 1/2. In both cases, we find that dissociative adsorption is slightly favored over molecular. We find that the transfer of the methanol proton to the (110) surface is barrierless at both fulland half-coverage. These findings may explain the conflicting assignments of the nature of adsorbed methanol from experimental studies. We find that dissociation via C-O bond scission is equally favorable on energetic grounds as O-H bond scission, despite large differences in the gas-phase heterolytic bond dissociation energies. Nevertheless, we do not observe C-O bond scission
J. Phys. Chem. B, Vol. 102, No. 11, 1998 2025 in our MD simulations, and this suggests that this is likely to be an activated process on the (110) surface. Our full-coverage MD simulations show that the methanol molecules are severely sterically hindered and highly unfavorable conformations are formed transiently by methyl hydrogens in neighboring methanol molecules approaching within 2 Å or less. Half-coverage MD simulations show that the methanol molecule will dissociate spontaneously on the (110) surface for a variety of different starting conformations. Only when the methanol molecule approaches the surface from directly above the bridging oxygen atom is it repelled away from the surface. As has been seen previously for water adsorption on the (110) surface,10 intermolecular hydrogen bonding between sorbate molecules or fragments is critical in determining the energetically most favorable structure. Acknowledgment. The work of S.P.B. was supported by EPSRC grant GR/J34842. The project was performed within the framework of the UK Car-Parrinello Consortium. We are grateful to the High Performance Computing Initiative (HPCI) for an allocation of time on the Cray T3D at EPCC. The authors acknowledge the exploratory calculations on this system performed by Dr. S. Pugh at Keele in 1995-6. References and Notes (1) Lindan, P.; Harrison, N.; Gillan, M.; White, J. A. Phys. ReV. B. 1997, 55, 15919-15927. (2) Bates, S.; Kresse, G.; Gillan, M. Surf. Sci. 1997, 385, 386-394. (3) Ramamoorthy, M.; King-Smith, R.; Vanderbilt, D. Phys. ReV. B 1994, 49, 16721. (4) Manassidis, I.; Vita, A. D.; Gillan, M. Surf. Sci. Lett. 1993, 285, L517. (5) Manassidis, I.; Gillan, M. J. Am. Ceram. Soc. 1994, 77, 335. (6) Pugh, S.; Gillan, M. Surf. Sci. 1994, 320, 331. (7) Goniakowski, J.; Holender, J.; Kantorovich, L.; Gillan, M.; White, J. Phys. ReV. B 1996, 53, 957. (8) Manassidis, I.; Gillan, M. Surf. ReV. Lett. 1994, 1, 491. (9) Manassidis, I.; Goniakowski, J.; Kantorovich, L.; Gillan, M. Surf. Sci. 1995, 339, 258. (10) Lindan, P.; Harrison, N.; Holender, J.; Gillan, M. Chem. Phys. Lett. 1996, 261, 246-252. (11) Lindan, P.; Muscat, J.; Bates, S. P.; Harrison, N.; Gillan, M. Faraday Discuss. 1997, 106, 135-154. (12) Kantorovich, L.; Gillan, M. Surf. Sci. 1997, 374, 373-386. (13) Goniakowski, J.; Gillan, M. Surf. Sci. 1996, 350, 145. (14) Szyman´ski, M.; Gillan, M. Surf. Sci. 1996, 367, 135-148. (15) Charlton, G.; Howes, P.; Nicklin, C.; Steadman, P.; Taylor, J.; Muryn, C.; Harte, S.; Mercer, J.; McGrath, R.; Norman, D.; Turner, T.; Thornton, G. Phys. ReV. Lett. 1997, 78, 495-498. (16) Lide, D. Handbook of Chemistry and Physics, 74th ed.; Boca Raton, FL, CRC Press: 1994. (17) Onishi, H.; Aruga, T.; Egawa, C.; Iwasawa, Y. Surf. Sci. 1988, 193, 33. (18) Ramis, G.; Busca, G.; Lorenzelli, V. J. Chem. Soc., Faraday Trans. 1 1987, 83, 1591-1599. (19) Suda, Y.; Morimoto, T.; Nagao, M. Langmuir 1987, 3, 99. (20) Kim, K. S.; Barteau, M. A. Surf. Sci. 1989, 223, 13-32. (21) Kim, K. S.; Barteau, M. A.; Farneth, W. Langmuir 1989, 4, 533. (22) Aas, N.; Pringle, T. J.; Bowker, M. J. Chem. Soc., Faraday Trans. 1994, 90, 1015-1022. (23) Hussein, G. A. M.; Sheppard, N.; Zaki, M. I.; Fahim, R. B. J. Chem. Soc., Faraday Trans. 1991, 87, 2655-2659. (24) Carrizosa, I.; Castaner, S.; Munuera, G. J. Catal. 1977, 49, 265. (25) Kantorovich, L.; Holender, J.; Gillan, M. Surf. Sci. 1995, 343, 221239. (26) Hohenberg, P.; Kohn, W. Phys. ReV. 1964, 136, 864. (27) Kohn, W.; Sham, L. Phys. ReV. 1965, 140, 1133. (28) Jones, R.; Gunnarsson, O. ReV. Mod. Phys. 1989, 61, 689. (29) Gillan, M. In Proc. NATO ASI on Computer Simulation in Materials Science, Aussois; Mayer, M., Pontikis, V., Eds.; Kluwer: Dordrecht, 1991; p 257. (30) Payne, M.; Teter, M.; Allan, D.; Arias, T.; Joannopoulos, J. ReV. Mod. Phys. 1992, 64, 1045.
2026 J. Phys. Chem. B, Vol. 102, No. 11, 1998 (31) Gillan, M. Contemporary Phys. 1997, 38, 115. (32) Kresse, G.; Hafner, J. Phys. ReV. B 1993, 47, RC558. (33) Kresse, G.; Furthmu¨ller, J. Comput. Mater. Sci. 1996, 6, 15-50. (34) Kresse, G.; Furthmu¨ller, J. Phys. ReV. B 1996, 54, 11169-11182. (35) Vanderbilt, D. Phys. ReV. B 1990, 41, 7892. (36) Kresse, G. Unpublished results, 1996. (37) Traylor, J.; Smith, H.; Nicklow, R.; Wilkinson, M. Phys. ReV. B 1971, 3, 3457. (38) Perdew, J.; Chevary, J.; Vosko, S.; Jackson, K.; Pederson, M.; Singh, D.; Fiolhais, C. Phys. ReV. B 1992, 46, 6671. (39) Bates, S.; Kresse, G.; Gillan, M. Surf. Sci., submitted.
Bates et al. (40) Andrews, S.; Burton, N.; Hillier, I.; Holender, J.; Gillan, M. Chem. Phys. Lett. 1996, 261, 521-526. (41) Monkhorst, H.; Pack, J. Phys. ReV. B 1976, 13, 5188. (42) Mills, G.; Jo´nsson, H.; Schenter, G. K. Surf. Sci. 1995, 324, 305. (43) Jo´nsson, H.; Mills, G. Chem. Phys., to be published. (44) Verlet, L. Phys. ReV. 1967, 159, 989. (45) Car, R.; Parrinello, M. Phys. ReV. Lett. 1985, 55, 2471-2474. (46) Pople, J. A.; Hehre, W. J.; Radom, L.; Schleyer, P. v. R. Ab initio Molecular Orbital Theory, 5th ed.; Wiley: New York, 1986. (47) Berkowitz, J.; Ellison, G. B.; D. G. J. Phys. Chem. 1994, 98, 27442765.