Catalysis Today 105 (2005) 44–65 www.elsevier.com/locate/cattod
First principles calculations of the adsorption and diffusion of hydrogen on Fe(1 0 0) surface and in the bulk
Dan C. Sorescu
U.S. Department of Energy, National Energy Technology Laboratory, Pittsburgh, PA 15236, USA Available online 31 May 2005
Abstract First-principles plane wave calculations based on spin-polarized density-functional theory (DFT) and the generalized gradient approximation (GGA) have been used to study the adsorption of hydrogen on Fe(1 0 0) surface and in the bulk. It was found that H2 adsorption takes place dissociatively with a classical activation energy of about 3.5 kcal/mol. In the low coverage regime at u = 0.25, H atom adsorbs at both two-folded and four-folded sites with a slight preference for the four-folded site. In the full coverage regime, there is a clear distinction between two-folded and four-folded adsorption sites with a net preference for adsorption at four-folded site. The dependence of H binding energy on coverage in the range 1.0 u 3.0 was also determined and the corresponding sequence of sites ?lling has been analyzed. After ?lling all four-folded sites, it was found that occupation of two-folded followed by one-folded sites is possible while adsorption at nearby mixed two-folded and one-folded sites leads to H–H recombination. The minimum energy pathways for surface diffusion of atomic H between selected pairs of local minima indicate the existence of small classical barriers with values of about 1.9 kcal/mol. These barriers increase slightly with the increase of coverage. When H diffuses from surface to subsurface sites, the corresponding barriers are larger than on the surface with values in the range 7.5–9.5 kcal/mol. At these subsurface sites, the absorption energy is still exothermic relative to gas phase H2 and increases with coverage. Once H penetrates the ?rst two surface layers, the corresponding diffusion barriers decrease to values close to those obtained in bulk Fe. Absorption of H in bulk bcc Fe is endothermic relative to isolated gas phase H2 and takes place at tetrahedral sites. The most favorable diffusion pathway among tetrahedral sites was found to pass through a trigonal site and has a low barrier of about 1.1 kcal/mol. Published by Elsevier B.V.
Keywords: Density-functional theory periodic calculations; Chemisorption; Diffusion; Iron; Hydrogen
1. Introduction Continuous review over the past decades of the chemisorption properties of hydrogen on iron surfaces with techniques of increasing sophistication has been motivated by the important applications of this system, particularly in the ?eld of heterogeneous catalysis. Examples of such applications are the ammonia synthesis or Fischer–Tropsch synthesis of aliphatic hydrocarbons [1–3]. Additionally, hydrogen is known to play an important role in metallurgical properties of steel as hydrogen embrittlement can lead to an enhanced decohesion of the steel with degradation of its mechanical performances and internal cracking . From a
E-mail address: Dan.Sorescu@netl.doe.gov. 0920-5861/$ – see front matter. Published by Elsevier B.V. doi:10.1016/j.cattod.2005.04.010
fundamental point of view, a central issue in these technological problems is to understand various elementary processes, such as chemisorption and diffusion of hydrogen both on Fe surfaces and in the bulk material. The interaction of hydrogen with iron single crystals has been the subject of several previous experimental studies [5–10]. Here, we will review some representative results related to adsorption on Fe(1 0 0) surface. In an initial investigation, Boszo et al.  have analyzed based on lowenergy diffraction (LEED), thermal desorption spectroscopy (TDS), work function measurements and ultraviolet photoelectron spectroscopy (UPS) the chemisorption of H on low indices single crystal Fe surfaces. It was found that on (1 1 0), (1 0 0) and (1 1 1) surfaces hydrogen adsorbs dissociatively with adsorption energies of 26, 24 and 21 kcal/mol,
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
respectively. On Fe(1 0 0) surface, two desorption states denoted b1 and b2, have been identi?ed as having energies of about 18 and 24 kcal/mol. By increasing the coverage, a continuous shift of the b2 state towards lower temperatures was observed indicating the weakening of the Fe–H bonds. The problem of H adsorption of Fe(1 0 0) surface has also been analyzed in several papers emerging from Madix and co-workers [6–8]. They have con?rmed based on X-ray photoelectron spectroscopy (XPS) and temperature programmed desorption (TPD) that for exposures of less than 10 L (1 L = 10?6 Torr s) adsorption of H2 takes place dissociatively and that desorption spectrum presents two states, in agreement with previous observations reported by Boszo et al. . At low H coverages of less than 3 ? 1014 atoms/cm2, the b2 state has been characterized as having a desorption peak near 400 K and a corresponding activation energy of desorption of 20.6 kcal/mol. This value has been revised to 23 kcal/mol in a subsequent paper . The second desorption state b1 is predominant at higher coverages and has a desorption temperature lower than the one of b2 state. For this b2 state, an activation energy of desorption of 17 kcal/mol has been determined , value which was also revised to 14.1 kcal/mol in subsequent work . Based on the experimental activation energy of desorption of b1-H2 peak, a lower limit to the Fe–H bond energy of 59 kcal/mol has been estimated . Further insight into the nature of the b1 and b2 states has been provided by Merrill and Madix based on highresolution electron energy loss spectroscopy (HREELS) measurements . At low coverages, the high temperature b2 state exhibits a low vibrational frequency loss at 700 cm?1. This state converts at high coverages to the b1 state which presents a frequency loss at 1000 cm?1. Based on the combined TPD and HREELS data, Merrill and Madix have assigned the loss feature at 700 cm?1 to a high ligancy state, namely the four-fold (4F) hollow site. This state converts with the increase in coverage to a state of smaller ligancy, tentatively assigned as a pseudo three-fold site, asymmetric within the 4F hollow site . Despite the relatively large number of experimental techniques used to investigate various aspects of the adsorption of H2 on iron surfaces, there are but a few theoretical studies focused on this subject [11–14]. Among these, the earlier theoretical work [11,12] was based on the use of semiempirical methods. Using a modi?ed MINDO/ SR semiempirical SCF method in conjunction to a 12 atom Fe cluster, Blyholder et al.  have determined the adsorption of H at various symmetric sites on Fe(1 0 0) surface. A stability order bridge > four-fold > on-top has been determined with the most stable con?guration corresponding to a displaced bridge-center site. This displaced bridge-center site has been used by Merrill and Madix  to argue in favor of their proposed pseudo threefold site. However, Blyholder et al.  have also shown that vibrational Fe–H frequencies at bridge center and at 4F site are 1343 and 1151 cm?1, respectively. Both these values are
in sharp contrast to those observed by Merrill and Madix who assigned the 1000 cm?1 loss to the three-fold site and the 700 cm?1 loss to the 4F hollow site. Raeker and DePristo  have also analyzed the chemisorption of H on Fe(1 0 0) surface using semi-empirical corrected effective medium theory and a slab with 2D periodicity. They found that the most stable adsorption site is the 4F hollow site. The bridge site was found to be about 3.7 kcal/mol less favorable than the 4F hollow site, result which contradicts the stability order provided by Blyholder et al. . From the coverage dependence of the binding energy at 4F site a small increase from 62.9 kcal/mol at u = 1/ 4 to 63.4 kcal/mol at u = 1 has been determined . A ?rst attempt to describe the interaction of H atoms with bcc iron using ab initio methods was done by Walch . Using an ab initio effective core potential method and a cluster model, he studied both the adsorption of H on (1 0 0) surface as well as in bulk bcc iron. His results support the fact that adsorption at 4F site is more stable than at the 2F site, in disagreement with results reported by Blyholder et al. . Walch has also analyzed H penetration into subsurface layers and has shown that adsorption at subsurface sites is less favorable than at the surface sites. In particular, it was found that penetration of the surface at 4F site toward a second layer encounters a large barrier while surface penetration at 2F site is energetically more favorable and can lead to absorption at a tetrahedral interior site. Additionally, for the case of H adsorption in bulk Fe, Walch determined that the most stable con?guration corresponds to tetrahedral site while octahedral site was found to correspond to a maximum on the potential energy surface. More recently, some initial results of H adsorption on Fe(1 0 0) surface based on density-functional theory (DFT) and 3D slab models have been reported by Eder et al.  in conjunction with oxidation of (1 0 0) and (1 1 0) surface by water molecules. In the regime of low coverages addressed in this study it was determined that adsorption can take place at both the bridge and hollow sites with almost identical adsorption energies. A more comprehensive study related to hydrogen adsorption, dissolution and diffusion in ferromagnetic bcc iron has been reported very recently by Jiang and Carter . Based on DFT–GGA (PBE) calculations, it was found that there is a preference for H to stay on the Fe surface instead of subsurfaces or in bulk. On the Fe(1 0 0) surface adsorption at four-fold hollow site was found to be energetically the most stable. Additionally, H diffusion to subsurface layers and dissolution into the bulk were determined to be endothermic. In the bulk, the minimum energy pathway for H diffusion was found to correspond to hops among the neighbor tetrahedral sites. In an attempt to provide further insight into the adsorption and diffusion properties of H on Fe(1 0 0) surface and in the bulk, we report here the results of ?rst principles densityfunctional theory calculations within the pseudopotential approximation. The periodic nature of the surface, which is essentially neglected in simple cluster models, has been
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
considered in the present case by the aid of a slab model (supercell model) repeated periodically in all three dimensions. We consider ?rst the process of dissociative chemisorption of H2 on Fe(1 0 0) surface. The chemisorption properties of atomic H on Fe(1 0 0) surface together with the role played by lateral H? ? ?H interactions on the adsorption process have then been considered. Further, we analyze the minimum energy pathway for diffusion of H between different surface sites and for diffusion from surface to subsurface sites. Finally, we consider the problem of H adsorption and diffusion in bulk bcc Fe. Throughout the paper the terms ‘‘one-fold’’ and ‘‘ontop’’ con?guration will be used interchangeably, as will be ‘‘two-fold’’ and ‘‘bridged’’ con?guration. The organization of the paper is as follows. In Section 2, we describe the computational methods. The results of total energy calculations for dissociative chemisorption of H2, adsorption and diffusion of H on surface and in the bulk are presented in Section 3. Finally, we summarize the main conclusions in Section 4.
2. Computational method The calculations performed in this study were done using the Vienna ab initio simulation package (VASP) [16–18]. This program evaluates the total energy of periodically repeating geometries based on density-functional theory and the pseudopotential approximation. In this case, the electron-ion interaction can be described by fully non-local optimized ultrasoft pseudopotentials (USPPs) similar to those introduced by Vanderbilt [19,20]. The USPPs include non-linear core corrections to account for possible valencecore charge interactions due to the 3d valence states overlapping with the 3p semi-core states in iron . For comparison, in majority of cases investigated we checked the accuracy of the USPPs using the projector augmented ¨ wave (PAW) method of Blochl  in the implementation of Kresse and Joubert . This method which is an allelectron DFT technique within the frozen core approximation improves the description of magnetic transition metals . The surface and the adsorbate–surface systems have been simulated using a slab model. Periodic boundary conditions are used, with the one electron pseudo-orbitals expanded over a plane wave basis set. The expansion includes all plane waves whose kinetic energy is less than a predetermined cutoff energy Ecut, i.e., h?2k2/2m < Ecut, where k is the wave vector and m the electronic mass. The k points are obtained from the Monkhorst–Pack scheme . The corresponding kinetic energy cutoffs were 495 eV in the case of USPPs and 400 eV for PAW calculations, respectively. Electron smearing is employed via the Methfessel– Paxton technique , with a smearing width consistent with s = 0.1 eV, in order to minimize the errors in the Hellmann–Feynman forces due to the entropic contribution
to the electronic free energy [17,18]. All energies are extrapolated to T = 0 K. The minimization of the electronic free energy is carried out using an ef?cient iterative matrixdiagonalization routine based on a sequential band-by-band residuum minimization method (RMM) [21,22], or alternatively, based on preconditioned band-by-band conjugategradient (CG) minimization . The optimization of different atomic con?gurations is based upon a conjugategradient minimization of the total energy. The value of Ecut and the k-point grid are chosen to ensure the convergence of energies and structures. The minimum energy paths between different minima were optimized by using the nudged elastic band (NEB) ? method of Jonsson and co-workers . In this approach, the reaction path is ‘‘discretized’’, with the discrete con?gurations, or images, between minima being connected by elastic springs to prevent the images from sliding to the minima in the optimization. The energies of the intermediate states along the reaction path are simultaneously minimized but the atomic motion is restricted to a hyperplane perpendicular to the reaction path. Unless otherwise stated, the calculations have been done using the spin polarized generalized gradient approximation (GGA) of Perdew et al. . In the case of PAW calculations, we have also included the spin interpolation for the correlation energy proposed by Vosko–Wilk–Nusair (VWN) . The effect of using a different exchangecorrelation functional on the adsorption energies has also been investigated by using the revised form of the Perdew, Burke and Ernzerhof (PBE) functional as introduced by Hammer et al.  (RPBE). These tests have been done only for the case of PAW-GGA potentials. In a number of instances, we have evaluated the vibrational frequencies of adsorbed molecular or atomic H species. This has been done under the frozen phonon mode approximation in which the metal atoms are ?xed at the relaxed geometries. The Hessian matrix has been determined based on a ?nite difference approach with a step ? size of 0.02 A for the displacements of H atoms along each Cartesian coordinate. By diagonalization of the massweighted Hessian matrix, the corresponding frequencies and normal modes have been determined.
3. Results and discussion 3.1. Calculations for bulk Fe The accuracy of the computational method used in this study has been tested initially to describe the properties of bulk iron, the bare surface and isolated H2 molecule. In the case of bulk properties of bcc Fe the results predicted using GGA-USPP have already been reported by us in our previous study related to adsorption of CO molecule on Fe(1 0 0) surface . Here, we summarize in Table 1, the corresponding USPP results from Ref.  for the
D.C. Sorescu / Catalysis Today 105 (2005) 44–65 Table 1 Calculated equilibrium lattice constant (a0), bulk modulus (B0), magnetic moment (M0) and saturation magnetization (Msatur) for ferromagnetic bcc Fe using PW91 functional and USPP and PAW potentials, respectively ? B0 (Mbar) M0 (mB) Msatur (G) Method a0 (A) USPP  USPP  PAW (present calc.) PAW  PAW  Exp.  2.8653 2.86 2.8341 2.83 2.834 2.87 1.597 1.55 1.697 1.74 1.74 1.68 2.33 2.32 2.20 2.20 2.21 2.22 1844 – 1792 – – 1750
3.3. Test calculations for the isolated H2 molecule A ?nal test we have performed was related to predictions of equilibrium properties of H2 molecule in gas phase. In this case, an isolated H2 molecule has been optimized in a cubic ? box of 10 A side. Using USPP (PAW) calculations, we ? obtained an equilibrium bond distance re = 0.749 A ? ), a bond dissociation energy De = 105.16 kcal/ (0.749 A mol (104.92 kcal/mol) and a vibrational frequency n = 4395 cm?1 (4315 cm?1). These values agree reasonable ? with experimental data (re = 0.741 A, De = 109.54 kcal/mol, ?1 n = 4401 cm )  as well as to other similar GGA results . The data obtained in the above-presented tests indicate that there is a good agreement of our calculated results for bulk Fe, the bare surface, and for the isolated H2 molecule with data obtained in other theoretical or experimental studies. These characteristics make us con?dent to pursue the next step of our investigations, namely interaction of H atoms and molecules with Fe(1 0 0) surface. 3.4. H2 dissociative adsorption on Fe(1 0 0) surface A ?rst problem we analyzed in our investigations was adsorption of H2 molecule on Fe(1 0 0) surface. As indicated by previous experimental studies [5,6], dissociative chemisorption of H2 can readily take place on Fe surface however there is not a clear evidence about the corresponding pathway followed by H2 during this process. Moreover, it presents interests to analyze if our computational methods are able to predict H2 dissociative chemisorption on Fe surface and evaluate the corresponding activation energy for this process. In this set of calculations, we have considered the adsorption at high symmetry sites 1F, 2F and 4F starting from a con?guration of H2 molecule parallel to the surface. These calculations have been done with a two-step approach. In this ?rst step, the H2 molecule was moved sequentially towards the surface and the atomic positions were allowed to relax in a plane parallel to the surface. It was found that indeed, for close proximity of the molecule to the surface dissociation takes place. Once the ?nal adsorbed con?gurations have been determined we have re?ned our previous calculations using the nudged elastic band method. Representative results related to adsorption at 1F and 2F sites are indicated in Fig. 1. In panels (a) and (c) of this ?gure, we provide the calculated minimum energy path as function of the relative height of the H2 molecule. Panels (b) and (d) indicate the corresponding variation of the H? ? ?H bond during the adsorption processes. When H2 molecule is approaching the surface above an Fe atom (see Fig. 1a) there is a small local minimum corresponding to a weakly bound molecular state at a height ? of about 1.752 A above the surface. This con?guration is represented in the second panel from the right of the top set of molecular structures in Fig. 1. In this con?guration, the
equilibrium lattice constant, bulk modulus, magnetic moment and saturation magnetization of bcc ferromagnetic Fe together with corresponding experimental data and other previous calculations. In Table 1, we also provide the results obtained in this work using GGA-PAW method and a k-point mesh of 15 ? 15 ? 15. By comparing the USPP and PAW results it can be seen that the agreement of geometrical lattice parameters to experimental data is very good with slightly better results for USPP values than PAW results. The corresponding differences of lattice parameters from experimental values are ?0.16% and ?1.25% for USPP and PAW, respectively. However, as indicated in Table 1, PAW potentials provide slightly more accurate data relative to experiment for bulk modulus, the local magnetic moment and the saturation magnetization than USPP does, in agreement with previous ?ndings . 3.2. Slab calculations for Fe(1 0 0) surface In our previous study , we have provided a detailed analysis of the type of relaxations that take place in a Fe(1 0 0) slab as function of the number of layers of the slab. It has been shown that a slab model with at least six layers should be adequate for adsorption studies. Particularly, we have shown based on USPP calculations that interlayer separation between top two layers Dd1–2 is compressed while the separations between the next two layers Dd2–3 expands relative to bulk values. The corresponding relaxations were found to be ?2.6% and 3.1% for a slab with six layers and ?3.8% and 2.7% for a slab with seven layers. These results are similar to those obtained by Eder et al.  of ?3.5% and 2.3% using a computational method similar to that presented here. These calculated values are also in agreement with the low-energy electron diffraction results of Wang et al.  who found values of Dd1–2/ d = ?5.0 ? 2% and Dd2–3/d = 5.0 ? 2%. The USPP results reported in this study were done using 2 ? 2 slab models containing either six or seven layers. As it will be shown in the next sections, there are small differences among these two models. The GGA-PAW calculations reported in this study were done using a 2 ? 2 slab model with seven layers. In this case, the corresponding relaxations of the top two layers were found to be Dd1–2/d = ?3.0% and Dd2–3/d = 2.71%, values consistent to those obtained in the case of USPP calculations.
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Fig. 1. Potential energy surface for dissociative adsorption of H2 molecule at 1F (panel (a)) and 2F (panel (c)) sites as function of the relative height of H2 molecule above the surface. In both cases the molecule is oriented initially parallel to the surface. Zero of energy was taken when H2 molecule is positioned at large separations from the surface. The corresponding variations of the H–H bond distances are provided in panels (b) and (d). For the top set of con?gurations, the second panel from the right indicates the physisorption state of H2 above a Fe surface atom.
? Fe–H distance is about 1.80 A and the H–H bond is stretched ? ? to 0.824 A relative to the gas phase value of 0.749 A. From this local minimum, the H2 molecule can overcome a small barrier of about 3.75 kcal/mol followed by dissociative chemisorption. If zero-point energies corrections are considered this barrier is found to be equal to 2.01 kcal/ mol. Along this process, the H? ? ?H distance further ? increases from 0.825 to 1.158 A at the top of the barrier
? and further to 2.864 A in the ?nal state. The ?nal state of this process corresponds to two adsorbed H atoms at 2F sites. We need to emphasize however that the weakly bounded molecular state described above should be considered carefully as current implementations of DFT do not have the right physics to account for van der Waals interactions and consequently might not provide an accurate description of such weakly bound states. Nevertheless, we have further
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
characterized this local minimum on the potential energy surface based on frequency calculations in the frozen phonon mode approximation. The results obtained indicate that all calculated vibrations are positive de?ned. Additionally, we found that the H–H stretching frequency decreases for this state from the gas phase value of 4395– 3125 cm?1. This variation is consistent to the observed elongation and weakening of the H–H bond. A second case we investigated corresponds to adsorption of H2 molecule above a 2F site, with H–H bond perpendicular to the Fe–Fe bond direction. This process is represented in the bottom set of panels in Fig. 1. As in previous case, upon approaching the surface the H2 molecule dissociates leading to two adsorbed H atoms at 4F hollow sites. From Fig. 1c, it can be seen that upon approaching the surface the molecule overcomes initially a slowly increasing barrier of about 3.48 kcal/mol. In this process, the corresponding H–H bond lengthens continuously. Once the distance between the molecule and surface ? drops below 1.2 A the molecule dissociates leading to formation of two chemisorbed atomic species at 4F sites. A similar dissociation barrier of about 3.78 kcal/mol has been determined in the case when H2 molecule is positioned initially above the middle of the 4F–2F distance, parallel to the surface. In this case, the dissociation pathway is no longer symmetric and the H2 molecule will tilt initially towards the 4F site followed by dissociation and adsorption at two neighbor 4F sites. We have also analyzed the case when H2 adsorption takes place above a 4F hollow site with a concerted dissociation leading to two H atoms adsorbed at 2F sites, symmetric distributed around the 4F site. The corresponding barrier found was signi?cantly larger with a value of 14.1 kcal/mol, indicating that such a symmetric pathway is less ef?cient. Based on the results presented in this section, it can be concluded that H2 dissociation on Fe(1 0 0) surface can take place with very small activation energies, particularly when the impinging molecule is parallel to the surface and positioned above 1F or 2F sites. In these cases, the calculated barriers have values in the range 3.5–3.8 kcal/mol. For the case of dissociation pathway that involves the weakly adsorbed molecular state above the on-top site, the zeropoint-energy corrections lead to a barrier of about 2.0 kcal/ mol. These results agree with previous experimental ?ndings [5,7], which indicate that dissociative chemisorption of H2 can readily take place. In particular, Burke and Madix have estimated an upper limit to the activation barrier for adsorption of about 2 kcal/mol . 3.5. Hydrogen adsorption on Fe(1 0 0) surface Characterization of the chemisorption properties of hydrogen on Fe(1 0 0) surface has been performed using four sets of calculations. In the ?rst set calculations have been done using a 2 ? 2 slab with six layers containing a total of 24 Fe atoms. In this case, adsorption on a single side
of the slab has been considered. This model (denoted herewith as M1) is consistent with the one we have used previously to investigate adsorption of CO on Fe(1 0 0) surface . In order to check if results were not affected by signi?cant dipole–dipole interactions, we have also analyzed the case when adsorption takes place on both sides of a slab model with seven layers containing a total of 28 Fe atoms. This model will be denoted as M2. For both M1 and M2 models, calculations have been done using GGA-USPP method and PW91 functional. The effect of using GGAPAW potentials has been tested in models M3 and M4. In both these last two cases we have considered 2 ? 2 slabs with seven layers in which adsorption takes place on both sides of the slab. In the case of M3 and M4 models, calculations have been done using PW91 and RPBE exchange-correlation functionals, respectively. In all models mentioned above, we have relaxed the coordinates of the H atoms and of the Fe atoms in the ?rst two layers of the slab. For the case of double side adsorption relaxations of Fe atoms in top two layers of each side of the slab have been considered. The remaining atoms of the slab have been frozen at the bulk-optimized positions. The dependence of the adsorption energies on surface coverage has been analyzed by considering several H atoms in the supercell. As each simulation model contains four surface sites of type 1F, 2F and 4F, we will abbreviate adsorption at speci?c sites as nF(m) where the index n denotes the type of sites (n = 1, 2 or 4 for one-fold, two-fold and four-fold) while the index m denotes the number of atoms in the supercell adsorbed at particular types of sites (1 m 4). For the optimized adsorption con?gurations, we report here the corresponding adsorption energies calculated with respect to isolated H2 molecule (Eads1) or to isolated H atom (Eads2). Speci?cally, Eads1 ? and Eads2 ? NEH ? Eslab ? E?H?slab? N (2) 1=2NEH2 ? Eslab ? E?H?slab? N (1)
where EH2 and EH are the energies of the isolated H2 molecule or H atom, N the total number of H atoms in the simulation box, Eslab the total energy of the relaxed slab and E(H + slab) is the total energy of the optimized adsorbate– slab system. A positive Eads1 (Eads2) corresponds to a stable adsorbate–slab system relative to the corresponding reference state. The same Brillouin-zone sampling has been used to calculate the properties of the bare slab and of the adsorbate–slab systems. In the case of H2 molecular adsorption state, the adsorption energy has been determined as: Eads ? NEH2 ? Eslab ? E?H2 ?slab? N (3)
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Table 2 ? Calculated equilibrium distances (A) and adsorption energies (kcal/mol) using PW91-USPP for molecular and atomic H adsorbed on the Fe(1 0 0) surface at different coverages Con?guration u r(FeI–H)a r(FeII–H) h(H-surf)b – 1.752 1.685 1.789 – – – – – – 1.744 1.756 1.774 1.784 1.795 1.799 1.793 1.680 1.691 1.587 1.576 1.571 1.063 1.026 0.952 0.382 0.373 0.357 0.337 ?0.321 ?0.285 ?0.314 ?1.484 ?2.186 – 1.585 1.586 1.561 1.062 1.019 0.959 0.385 0.376 0.357 0.346 D12 (%)c ?3.58 1.76 ?1.33 0.14 ?1.10 ?0.94 0.43 1.98 0.65 2.40 ?4.94 ?3.74 ?1.22 0.94 ?1.80 ?4.10 ?6.40 6.18 ?4.48 ?4.73 ?1.23 ?3.36 ?2.53 1.31 0.48 ?0.54 ?6.12 ?4.72 ?1.29 ?0.38 D23 (%) 2.09 1.82 0.83 0.45 1.74 1.64 1.60 0.74 ?1.01 ?2.89 4.28 3.71 3.05 2.30 4.49 6.96 13.00 6.71 8.19 2.66 2.45 1.98 2.24 1.33 ?0.25 0.88 4.81 4.10 2.06 1.86 Eads1d – 2.44 1.88 0.80 ?4.59 ?7.71 ?10.14 7.39 6.49 5.17 7.22 7.86 8.41 8.74 0.06 0.49 1.24 ?7.57 ?5.37 – ?4.94 ?6.05 ?10.50 7.53 6.13 5.01 8.06 8.19 8.83 8.91 Eads2 – – – – 48.04 44.93 42.49 60.03 59.13 57.82 59.87 60.50 61.05 61.38 52.70 53.13 53.88 45.06 47.26 – 47.69 46.58 42.13 60.17 58.78 57.66 60.71 60.84 61.47 61.56
Model M1, six layers slab model, single side adsorption Slab 0 – – Molecular hydrogen H2-1F(1) H2-1F(2) H2-1F(4) con?gurations 0.25 1.800 0.50 1.732 1.00 1.835
Surface atomic hydrogen con?gurations 0.25 1.587 1F(1)e 1F(2) 0.50 1.576 1F(4) 1.00 1.571 2F(1) 0.25 1.720 2F(2) 0.50 1.712 2F(4) 1.00 1.687 4F(1) 0.25 2.088 4F(2) 0.50 1.944 4F(3) 0.75 1.928 4F(4) 1.00 1.874 Subsurface atomic hydrogen con?gurations S1(1) 0.25 1.689 S1(2) 0.50 1.693 S1(4) 1.00 1.708 S2 0.25 1.723 S3 0.25 –
Model M2, seven layers slab model, double side adsorption Slab 0.0 – – Surface atomic hydrogen con?gurations 0.25 1.585 1F(1)e 1F(2) 0.50 1.586 1F(4) 1.00 1.561 2F(1) 0.50 1.720 2F(2) 0.50 1.710 2F(4) 1.00 1.689 4F(1) 1.00 2.089 4F(2) 1.00 1.955 4F(3) 1.00 1.928 4F(4) 1.00 1.876
– – – – – – 1.731 1.744 1.774 1.775
FeI and FeII denote iron atoms in the ?rst and second layer, respectively. This distance is measured perpendicular to the surface. c ? The relative distances are taken with respect to the bulk-optimized interlayer distance of 1.4326 A. d The adsorption energies Eads1, Eads2 are calculated with respect to energies of isolated H2 and H species. e Notation mF(n) with m = 1, 2, 4 and n = 1, 2, 4 refers to adsorption at sites of type ‘‘m’’ (m = one-fold, two-fold and four-fold, respectively) while the index ‘‘n’’ denotes how many such sites among a maximum of 4 in the supercell are occupied.
where E?H2 ?slab? is the total energy of the H2–slab system. The calculated adsorption energies reported in this study do not include zero point energy corrections. 3.5.1. Molecular adsorption states of H2 on Fe(1 0 0) surface We have indicated in Section 3.4 that a weak adsorption state of undissociated H2 can be found at an on-top position. In this con?guration represented in the second top panel from the right of Fig. 1, the molecule is positioned parallel to the surface, symmetrically bonded to a Fe atom underneath. For this adsorption con?guration, the corresponding binding energy calculated based on Eq. (3) (see entry H2-1F(1) in
Table 2) is small with a value of 2.44 kcal/mol indicating a weakly physisorbed state. By increasing the H2 surface coverage such that all four on-top sites in the simulation box are occupied, we ?nd that H2 adsorption energy decreases to about 0.88 kcal/mol. This weakening of Fe–H2 interaction is also re?ected by the vertical separation between the molecule and the surface, h(H-surf), which increases from ? ? about 1.752 A at u = 0.25 to about 1.789 A at u = 1.0. This weakening was also probed based on vibrational frequencies calculations in the frozen phonon approximation. Indeed, our results indicate a decrease of the H–H vibrational frequency from 3125 cm?1 at u = 0.25 to 3089 cm?1 at u = 1.0.
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Table 3 ? Calculated equilibrium distances (A) and adsorption energies (kcal/mol) using PAW potentials for molecular and atomic H adsorbed on the Fe(1 0 0) surface at different coverages u Con?guration u r(FeI–H)a r(FeII–H) h(H-surf)b D12 (%)c D23 (%) 2.73 2.38 1.61 0.01 2.74 2.84 1.50 ?0.28 ?1.50 4.68 4.72 3.65 2.80 5.00 7.33 12.39 8.09 9.91 Eads1d – 2.88 2.06 1.69 ?4.57 ?10.07 8.37 6.66 5.81 9.27 9.29 9.34 9.54 1.28 1.59 1.86 ?6.49 ?5.06 (?7.45) (?13.04) (4.88) (3.03) (1.70) (5.20) (5.24) (5.36) (5.36) Eads2 – – – – 47.90 42.40 60.84 59.13 58.28 61.74 61.81 61.89 62.01 53.75 54.00 54.34 45.97 47.41 (45.20) (39.60) (57.53) (55.68) (54.35) (57.85) (57.90) (58.01) (58.01)
Models M3 (PW91) and M4 (RPBE): 7 layers slab model, double side adsorption Slab 0 – – – ?3.18 Molecular hydrogen con?gurations 0.25 1.738 H2-1F(1) 0.50 1.711 H2-1F(2) H2-1F(4) 1.00 1.802 Surface atomic hydrogen con?gurations 1F(1) 0.25 1.566 1F(4) 1.00 1.554 2F(1) 0.25 1.706 2F(2) 0.50 1.697 2F(4) 1.00 1.675 4F(1) 0.25 2.073 4F(2) 0.50 2.039 4F(3) 0.75 2.037 4F(4) 1.00 2.034 Subsurface atomic hydroggen S1(1) 0.25 S1(2) 0.50 S1(4) 1.00 S2 0.25 S3 0.25 con?gurations 1.678 1.683 1.695 1.697 – – – – – – – – – 1.722 1.731 1.743 1.751 1.773 1.781 1.799 1.674 1.677 1.685 1.655 1.748 1.566 1.554 1.061 1.013 0.943 0.376 0.377 0.360 0.352 ?0.338 ?0.292 ?0.331 ?1.471 ?2.169 1.78 ?2.67 1.31 ?5.01 ?1.05 2.26 1.09 2.74 ?4.98 ?4.48 ?2.43 ?1.28 ?1.21 ?3.26 ?4.33 6.65 ?4.53
The optimized geometric parameters have been determined based on PW91-PAW calculations. The adsorption energies given in parentheses represent the RPBE values obtained at the PW91 optimized geometries. a FeI and FeII denote iron atoms in the ?rst and second layer, respectively. b This distance is measured perpendicular to the surface. c ? The relative distances are taken with respect to the bulk-optimized interlayer distance of 1.4326 A. d The adsorption energies Eads1, Eads2 are calculated with respect to energies of isolated H2 and H species
These results obtained using PW91-USPP calculations have also been con?rmed based on PW91-PAW calculations in connection with the seven layers slab and double side adsorption (model M3, see entries H2-1F in Table 3). In this case, we found similar weak binding energies for adsorbed H2 and the same decreasing trend of the binding energies with coverage. Particularly, the adsorption energy changes from 2.88 kcal/mol at u = 0.25 to 1.69 kcal/mol at u = 1.0. Both sets of calculations using M1 and M3 models indicate that a weak physisorbed state of molecular H2 can be identi?ed on Fe(1 0 0) surface. The strength of the binding energy for this state is dependent on surface coverage and decreases with the increase of coverage. As was shown in Section 3.4, this state could serve as a precursor to dissociative adsorption of H2. 3.5.2. H adsorption at low coverages (u = 1/4) The resulting geometric and energetic parameters for various binding con?gurations of H on the Fe(1 0 0) surface are given in Tables 2 and 3. These results are reported as function of speci?c surface sites and coverages. In the case of u = 0.25, the lowest coverage considered in this study, the adsorption con?gurations at 1F, 2F and 4F sites are represented in Fig. 2a–c. Based on the energetic data reported in Table 2, it follows that at low coverages the adsorptions at 2F and 4F sites are the most stable with
practically identical energies of about 60 kcal/mol. The adsorption energy of the on-top con?guration is smaller by 12 kcal/mol than the one of 4F and 2F sites. As will be shown in Section 3.6, both 4F and 2F con?gurations are stable minima while 1F con?guration is not a stable minimum but a second-order saddle point on the potential energy surface. Relative to isolated H2 molecular state, adsorption at 2F and 4F states is exothermic while adsorption at 1F state is endothermic as indicated by energies Eads1 in Tables 2 and 3. In the case of bare surface we have seen that opposite types of relaxations take place in the top two layers of the slab. The top ?rst layer relaxes inwards with an interlayer separation D12 = ?3.58% relative to the bulk value while the second layer relaxes outwards with a relative interlayer separation D23 = 2.09%. In the presence of adsorbed H, these values change function of the adsorption sites. For example, D12 vary from 1.98% in the case of adsorption at 2F site to ?4.94% for adsorption at 4F site. The energetic values determined for the case of the slab model with six layers and single side adsorption (model M1) were little changed relative to the case when calculations were done using the slab model with seven layers and double side adsorption (model M2). The corresponding results from model M2 are also indicated in Table 2. For the model M2, we notice, for example, that the adsorption at the 2F site
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Fig. 2. The adsorption con?gurations of H on Fe(1 0 0) surface: (a) one-fold (on-top) con?guration; (b) two-fold (bridged) con?guration; (c) four-fold hollow site con?guration. Panels (d) and (e) correspond to full coverage adsorption at two-fold and four-fold sites when all respective sites are occupied.
increases by only 0.14 kcal/mol to 60.17 kcal/mol while the energy at 4F site increases by 0.83 kcal/mol to 60.71 kcal/mol relative to energies determined for model M1. Consequently, in this case, the binding at 4F site becomes only slightly larger than the one at 2F site. These modi?cations in the binding energies are also re?ected by surface relaxations, which are slightly larger than those, obtained for model M1. For example, in the case of 4F site adsorption (model M2) the interlayer separations are D12 = ?6.12% and D23 = 4.81% versus ?4.94% and 4.28% obtained for M1 model. The lack of signi?cant changes when going from model M1 (with a sixlayer slab and single side adsorption) to model M2 (with a seven-layer slab and double side adsorption) indicates that at least in this regime of coverages the dipole–dipole interactions are well screened by the supercell models considered. In Table 3, we indicate the results obtained based on GGA-PAW calculations using PW91 (model M3) or RPBE (model M4) exchange-correlation functionals. A direct comparison of models M2 and M3 indicate that there are no major differences between the two sets of values. Both sets of models M2 and M3 show that 4F site is more stable than 2F site. Additionally, the PW91-PAW results are only slightly larger than those obtained using PW91-USPP calculations with the largest variation of 1.03 kcal/mol seen
for the case of adsorption at 4F site. The RPBE-PAW results given in parentheses in Table 3 indicate binding energies, which on average over the entire range of coverages are about 3.0 kcal/mol smaller than those determined using PW91-USPP (model M2). At u = 0.25 coverage however, the RPBE adsorption energies at 2F and 4F sites are 2.64 and 2.86 kcal/mol, respectively, smaller than the corresponding PW91-USPP results. These variations among the two functionals are typical to those obtained before for H adsorption on Fe(1 1 0) surface . The results presented in this section indicate that in the regime of low coverages H bonding at both bridge and hollow sites is possible. PW91-USPP method when applied to the six layer slab (model M1) predicts almost equal values for adsorption energies at 2F and 4F sites in agreement to previous ?ndings reported by Eder et al.  In the case of the seven layers slab model, PW91-USPP calculations indicate however a slightly higher stability for 4F site than 2F site. Similarly, the PW91-PAW method indicates that 4F site is more stable by 2.6 kcal/mol than 2F site. The slightly higher stability of H adsorbed at 4F site relative to 2F site was also found by Jiang and Carter based on PBE-PAW calculations in connection to a ?ve layers slab model . At u = 0.25 ML coverage, it was determined that the hollow
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
site is more stable than bridge site by 1.38 kcal/mol while at even lower coverages (u = 0.11 ML) this difference decreases to 0.92 kcal/mol. The calculated PW91-USPP (PW91-PAW) adsorption energies at low coverages found in this study, ranging from 59.87 to 60.71 kcal/mol (60.84– 61.74 kcal/mol), are in excellent agreement with the value of 59 kcal/mol estimated by Burke and Madix based on the analysis of desorption spectra of H2 . As a ?nal note, beside the 1F, 2F and 4F con?gurations described above we have also attempted to identify an asymmetric three-fold hollow site as proposed by Merrill and Madix . However, upon optimizations starting from such an asymmetric con?guration we found that it evolves towards the centered four-fold hollow con?guration. As a result, we were not able to con?rm the existence of such a three-fold hollow site adsorption con?guration. 3.5.3. Coverage dependence of H binding energies Another problem we analyzed in the current study was related to the dependence on surface coverage of the chemisorption properties of atomic H. This has been done using surface models M1–M4 introduced above for the coverage range 0.25 u 1. The corresponding results for various geometric parameters and adsorption energies are detailed in Tables 2 and 3 while representative con?gurations are depicted in panels (d) and (e) of Fig. 2. By inspecting the values of adsorption energies provided in these tables for speci?c types of sites it can be concluded that in the case of adsorption at 1F and 2F sites the corresponding adsorption energies decrease with the increase in coverage. The largest variations are seen to take place for the case of 1F sites. In this case, by increasing the coverage from 0.25 to 1 the corresponding binding energies decrease by 5.5– 5.6 kcal/mol depending on the particular model M1–M4 used in calculations. A similar decrease is obtained in the case of adsorption at 2F sites. However, for these sites the drop of adsorption energies is less signi?cant, with values of 2.21, 2.51, 2.56 and 3.18 kcal/mol for models M1–M4, respectively. This behavior suggests the existence of repulsive interactions between neighbor H atoms adsorbed at either 1F or 2F sites, respectively. In the case of adsorption at 4F hollow sites, we found a different behavior, namely the adsorption energy slightly increases with the coverage increase. The corresponding variation in the binding energy however is dependent on the level of calculations and the slab model used. Particularly, PW91-USPP calculations indicate increases of the binding energies of about 1.5 and 0.8 kcal/mol for the slab models with six and seven layers, respectively, while PW91-PAW and RPBE-PAW calculations show smaller variations of 0.27 and 0.16 kcal/mol, respectively. We notice that independent on the model used the increase of coverage at 4F sites leads to a decrease of the vertical height h(H-surf) of the H atom to the surface and of the bond distances r(FeI– H) between the H atom and the Fe atoms in the ?rst layer (see Tables 2 and 3). Concomitantly, the distance r(FeII–H)
between the H atom and the Fe atom beneath, in the second layer, increases suggesting opposite relaxations for the ?rst and second layers of Fe atoms. Overall, these results indicate that the most stable adsorption con?guration over the coverage range 0.25 u 1 is the 4F hollow site followed by the bridge 2F con?guration. These ?ndings are similar to those obtained by Jiang and Carter . In particular, they have determined an increase of the binding energy of about 0.46 kcal/mol for the hollow site with the increase of coverage from 0.25 to 1, while a decrease of about 1.84 kcal/ mol was found to take place for binding at bridge site over the same coverage range. Further extensions of these investigations were done by considering the case of atomic adsorption at mixed types of surface sites, such as 4F–2F, 4F–1F or 4F–2F–1F over the coverage range 0.25 u 3. Given the large number of calculations involved, we have analyzed these cases only based on PW91-USPP calculations and the slab model with six layers. The corresponding variations of the H adsorption energy are given in Fig. 3. Here the binding energies have been calculated with respect to isolated surface and isolated H2 molecule based on Eq. (1). In this ?gure, we have included also for completeness the results presented above corresponding to adsorption at single type of 1F, 2F and 4F sites. Additionally, we have included results corresponding to subsurface adsorption, which will be discussed in the next section. Representative con?gurations indicating adsorption at mixed sites are presented in Fig. 4. From Fig. 3, it can be clearly seen that in the range 0.25 u 1 the binding energy at 4F sites increases slightly with coverage while that at 2F or 1F sites decreases. Both energy values at 4F and 2F sites are positive indicating the existence of stable con?gurations (exothermic adsorption). However, the adsorption energies at 1F sites have negative
Fig. 3. Variation of the binding energy of H calculated based on Eq. (1) with respect to surface coverage. Notation 4F, 2F and 1F denote adsorption at respective single type of sites. Notations 4F–2F, 4F–1F and 4F–2F–1F indicate adsorption at mixed type of surface sites. The abbreviations 4F-H2 and S1 correspond, respectively, to adsorption of H2 molecule on a surface with preadsorbed H at all 4F sites and adsorption at subsurface S1 sites as described in text.
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Fig. 4. Representative adsorption con?gurations of H on Fe(1 0 0) for coverages larger than 1. The initial distribution of H atoms over different types of surface sites is as follows: (a) 4F(4)–1F(1); (b) 4F(4)–2F(1); (c) 4F(4)–2F(4); (d) 4F(4)–2F(1)–1F(1); (e) 4F(4)–2F(4)–1F(1); (f) 4F(4)–2F(4)–1F(4). Notation 4F(4) and 2F(4) denote adsorption at all four 4F and 2F sites, respectively.
values, which correspond to unstable con?gurations (endothermic adsorption) with respect to gas phase H2. In the coverage range, 0.5 u 1, the adsorption at mixed sites, such as 4F–2F is also possible and the corresponding binding energies are found to be slightly smaller than those obtained at single type 4F sites. In those cases when for a given coverage several adsorption con?gurations are possible, we have included here only the one with the most stable adsorption energy. As an example, in the case when two H atoms adsorb at 4F and 2F sites, the indicated value 4F–2F given in Fig. 3 of 7.92 kcal/mol corresponds to a large ? separation between the H atoms of about 3.28 A. In this con?guration, only one of the surface Fe atoms is bonded to both H atoms adsorbed at 2F and 4F sites, respectively. For the case of adsorption at neighbor 4F–2F sites, the ?nal ? H? ? ?H separation found is 2.46 A. In this case, there are two surface Fe atoms bonded to both H atoms. In this case, the corresponding binding energy is smaller than in previous case and has a value of 6.46 kcal/mol. Overall, the above presented analysis indicates that in the coverage range
0.25 u 1 adsorption at hollow sites leading to full occupation of 4F sites is energetically the most favorable process. These results are also consistent with experimental ?ndings obtained by Merrill and Madix . In their experimental studies Burke and Madix have determined a hydrogen saturation coverage on clean Fe(1 0 0) of 1.0 ? 0.1 ML  when the surface was exposed under UHV conditions to molecular H2. It presents interest however to determine if adsorption coverages higher than 1 ML can be achieved on Fe(1 0 0) surface, particularly if experiments would be done using atomic H under high pressure conditions. In order to study coverages u > 1 we considered that all 4F sites are occupied and analyzed further adsorptions at 1F or 2F sites. We illustrate such con?gurations in Fig. 4a and b for the case of u = 1.25 coverage. Here, all 4F sites are occupied and one additional H atom occupies either the 1F (4F(4)–1F(1)) or 2F (4F(4)–2F(1)) sites. By successive ?lling of the remaining 1F or 2F sites, we have determined the dependence of the binding energies in the range 1 u 2. As can be seen from Fig. 3, there is a
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
continuous decrease almost linear with coverage of the binding energies for mixed 4F–2F or 4F–1F cases, respectively. Among these con?gurations over the 1 u 2 coverage range, the most stable correspond to adsorption at 4F–2F sites while the mixed 4F–1F con?gurations are less favorable. A pictorial view of the 4F(4)– 2F(4) con?guration is given in Fig. 4c. We also note that the 4F(4)–1F(4) con?guration at u = 2 has a negative energy indicating an endothermic adsorption relative to gaseous H2. Beside this sequential scheme of ?lling of 4F sites followed by 2F or 1F sites a different scheme has also been considered. This consists in initial ?lling of 4F sites followed by a simultaneous adsorption at mixed 2F and 1F sites. We will denote this scheme as 4F–2F–1F. The smallest coverage where we analyzed such a combination corresponds to u = 1.5. In this case, there are six H atoms in the supercell distributed initially at four 4F sites, 4F(4), one 2F site, 2F(1) and one 1F site, 1F(1). We have analyzed two different cases of such mixed con?gurations with respect to positions of atoms adsorbed at 2F and 1F sites. In the ?rst case the H atoms are positioned at 2F and 1F sites, which are separated by one surface Fe atom. Upon optimization the H atoms remain at their respective sites with a ?nal separation ? between atoms at 2F and 1F sites of 3.26 A. The corresponding binding energy in this case was found to be 5.07 kcal/mol relative to isolated gas phase H2. In the second case analysed, we have considered that the H atoms are positioned at nearby 2F and 1F sites. This con?guration, denoted as 4F(4)–2F(1)–1F(1) is represented in Fig. 4d. As illustrated in this ?gure, it can be seen that upon optimization the H atoms from nearby 1F and 2F sites move towards each ? other leading to a ?nal separation H? ? ?H of 0.81 A. This state is similar to the on-top molecular state described earlier in this section and has a binding energy of 6.27 kcal/mol relative to isolated gas phase H2. The energy of this state of higher binding energy has been considered in Fig. 3 for 4F– 2F–1F states at u = 1.5. Similarly to previous case, for u = 1.75 coverage corresponding to a site distribution 4F(4)–2F(2)–1F(1), we have analyzed two different situations. In the ?rst one, the H atom at 1F site is separated by one surface Fe atom from the other two H atoms adsorbed at 2F sites. Upon optimization H atoms continue to occupy their original 4F, 2F, and 1F sites, respectively, and no recombination takes place. In this case, the distances between H atom adsorbed at 1F and those from ? 2F sites are 3.15 and 3.29 A, respectively. The binding energy of H atoms in this con?guration was found equal to 4.39 kcal/mol. The second case we analyzed corresponds to a con?guration where two H atoms are at nearby 2F and 1F sites while a second H atom is positioned at a 2F site separated by a surface Fe atom from the 1F site. In this case, we found that due to attractive interactions between H atoms from nearby 2F and 1F sites they recombine and as a result the corresponding binding energy increases to 5.47 kcal/ mol. Finally, for the same coverage u = 1.75 we have analyzed a con?guration having a site distribution 4F(4)–
3F(2). This con?guration is even less favorable than the other two 4F(4)–2F(2)–1F(1) con?gurations described above with a binding energy of 3.97 kcal/mol. The increased stability of mixed 4F–2F–1F con?gurations relative to 4F–2F con?gurations is also maintained at even higher coverages as can be seen from Fig. 3. For example, at u = 2, we ?nd that the mixed state 4F(4)–2F(2)–1F(2) has a binding energy of 4.86 kcal/mol while the binding of 4F(4)–2F(4) state is smaller with a value of 2.41 kcal/mol. Adsorption con?gurations for coverages above u = 2 have been studied by initial ?lling of all four 4F and 2F sites and sequential ?lling of the remaining 1F sites. Two such con?gurations are illustrated in panels (e) and (f) of Fig. 4. As in previous 4F–2F–1F cases, we observed that upon optimization H atoms from neighbor 1F and 2F sites tend to recombine leading to formation of adsorbed H2 molecules. Overall, as indicated in Fig. 3, as a result of coverage increase the binding energy per atom for 4F–2F–1F states decreases from 4.8 kcal/mol at u = 2 to about 2.1 kcal/mol at u = 3. When H recombination takes place the resulting structures can be considered as being composed from H atoms adsorbed at 4F sites and H2 molecules adsorbed at 1F sites. In order to test this possibility, we have analyzed a ?nal set of con?gurations denoted as 4F(4)–H2(m), with m = 1, 2 and 4, corresponding to a surface with all 4F sites occupied and with sequential occupation by H2 of 1F sites. The variation of the binding energies for 4F(4)–H2(m) structures is represented in Fig. 3. We notice that over the coverage range from 2 to 3 the variation of the binding energies for 4F–2F–1F and 4F-H2 curves is very close and that the curves are almost superimposed. This result supports the view that the ?nal states obtained from 4F–2F–1F con?gurations upon recombination of H from 2F to 1F sites are similar to those obtained by H2 molecular adsorption on a surface precovered with H at 4F sites. As a ?nal note we need to point out that the binding energies per atom indicated in Fig. 3 for the case when H2 recombination takes place include both the interaction between H atoms and surface atoms as well as the lateral interactions between H atoms. An alternative approach to get a measure of the interaction between the adsorbed H2 molecules and the surface with adsorbed H at 4F sites was obtained by recalculating the binding energies with respect to the energy of the surface with adsorbed H at all 4F sites plus the energy of isolated H2 molecule. In this case, the binding energies of H2 molecules were found to decrease to 1.78 kcal/mol for 4F(4)–H2(2) and 0.50 kcal/mol for 4F(4)– H2(4) con?gurations, respectively, indicating that H2 molecules can easily desorb from the surface. The general picture emerging from this analysis indicates that in the regime of low coverages occupation of 4F sites is energetically the most favorable. Once full coverage of 4F sites is achieved occupation of 2F sites followed by 1F sites is possible. In those cases when H atoms are positioned at neighbor 1F and 2F sites there is an attractive interaction between these atoms leading to H–H recombination and
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
formation of adsorbed H2 molecules. Such molecules however interact very weekly to the surface as re?ected by the adsorption energy of 0.5 kcal/mol at full 4F(4)–H2(4) coverage. As a result, desorption from these physisorbed states can easily take place. 3.5.4. Subsurface hydrogen absorption on Fe(1 0 0) surface In addition to H adsorption on Fe(1 0 0) surface, we have studied the absorption in subsurface layers. These investigations have been done using two sets of calculations, one using PW91-USPP method and a six-layers slab, and the second using PW91-PAW method in conjunction with a seven layers slab model. In both these cases the coordinates of the top three layers of the slab and of the H atom have been allowed to relax. A pictorial view of the resulting equilibrium con?gurations are given in Fig. 5 and the corresponding geometric and energetic parameters are detailed in Tables 2 and 3. As indicated in Fig. 5a, we have identi?ed a subsurface, tetrahedral like, equilibrium con?guration denoted as S1(1), positioned beneath a 2F surface con?guration. The results of PW91-USPP calculations show that in this state the H atom is positioned about ? 0.32 A beneath the top layer of the surface atoms. The bond
distances between H atom and the two neighbor Fe atoms in ? the ?rst layer are r(FeI–H) = 1.689 A while the distances to the other Fe atoms in the second layer are r(FeII– ? H) = 1.795 A. Relative to isolated Fe slab and isolated H atom the corresponding binding energy is 52.70 kcal/mol. By the increase of subsurface coverage (see Fig. 5b and c) the corresponding binding energy increases to 53.88 kcal/ mol. From the analysis of the absorption energies Eads1 given in Table 2 and calculated with respect to isolated H2 molecule, we observe that absorption at subsurface sites is still exothermic. Moreover, we notice that as a result of subsurface absorption and coverage increase some important deformations of the slab take place. For example, the relative interlayer distance between the ?rst and second layer D12 decreases by ?4.1% at u = 0.5 and by ?6.4% at u = 1.0, respectively, relative to bulk-optimized values, indicating an inward contraction of the ?rst layer towards the H atom. A reverse expansion is observed to take place for the D23 interlayer which increases by 6.96% and 13.0% at u = 0.5 and 1.0, respectively. The results obtained using PW91-PAW method and the seven layers slab are given for comparison in Table 3. They agree well with USPP data and only slight increases less than 1.05 kcal/mol (u = 0.25) of the binding energies are noticed.
Fig. 5. Subsurface absorption con?gurations of H on Fe(1 0 0) surface at coverages 0.25 (a), 0.5 (b) and 1.0 (c). Panels (d) and (e) correspond to tetrahedral subsurface absorption con?gurations in the second and below the second layers, respectively.
D.C. Sorescu / Catalysis Today 105 (2005) 44–65 Table 4 Calculated vibrational frequencies for H adsorbed at different sites at u = 0.25 and 1.0 coverages using PW91-USPP and PW91-PAW methods Con?guration 1F(1) 2F(1) 4F(1) 1F(4) 2F(4) 4F(4) Local minimum No Yes Yes No Yes Yes n (cm?1) (USPP) 1640, 1161, 1146, 1681, 1205, 1053, 283i, 283i 1095, 360 227, 209 331i, 330i 1089, 449 605, 290 n (cm?1) (PAW) 1670, 1183, 1193, 1738, 1221, 1127, 274i, 269i 1127, 368 296, 292 341i, 340i 1106, 448 226, 222
n (cm?1)  1953 1343 1151
The nature of each adsorption site as a local minimum on the potential energy surface is also indicated.
In addition to the subsurface state S1, we have also characterized other tetrahedral states denoted as S2 and S3 (see Fig. 5d and e) positioned in the second and below the second layer of surface atoms, respectively. As indicated in both Tables 2 and 3, these con?gurations have signi?cantly smaller binding energies of 45.06 (45.97) kcal/mol and 47.26 (47.41) kcal/mol, respectively, relative to S1 subsurface sites. Additionally, we notice the change in sign for Eads1 energies for these sites indicating that absorption becomes endothermal. The large decrease of the binding energies for S2 and S3 sites is also re?ected by some signi?cant vertical relaxations in subsurface layers. For example, in case of S2 state, in order to accommodate the adsorbed H atom both ?rst and second layers move outwards by 6.1% and 6.7%, respectively, while in the case of S3 state we see a signi?cant expansion taking place for the distance between the second and third layers by 8.2% (9.9%), and a contraction between ?rst and second layers by ?4.48% (?4.53%). 3.6. Vibration frequencies analysis An important aspect of previous experimental and theoretical studies [8,11] was the site assignment of adsorbed H based on the corresponding vibrational frequencies. Based on a 12-atom Fe cluster and MINDO/ SR calculations Blyholder et al.  have shown that there is an inverse relationship between the ligancy and the perpendicular vibrational frequency of H atom interacting with Fe surface. Particularly, it was found that vibrational frequencies vary as follows: 1953 cm?1 for 1F site, 1343 cm?1 for the 2F site and 1151 cm?1 for the 4F site . Based on EELS and TPD spectra, Merril and Madix  have assigned the high temperature b2 state to H adsorbed in a 4F hollow site. This state converts with increasing coverage to a low temperature b1 state that was suggested to have lower ligancy. The vibrational signature assigned to these states was 700 cm?1 for b2 state and 1000 cm?1 for b1 state . It is clear from these experimental and theoretical studies that the assignment of vibrational frequencies of H adsorbed at 4F site is different. Particularly, Blyholder et al.  have determined a vibrational frequency of 1151 cm?1 for H adsorbed at 4F site while Merril and Madix  suggest that this site should have a vibrational frequency of 700 cm?1.
We have determined the vibrational frequencies of adsorbed H at 1F, 2F and 4F sites using M1 (PW91-USPP) and M3 (PW91-PAW) models. These determinations were done for u = 0.25 and 1 coverages and the corresponding results are reported in Table 4. First, we note that not all the adsorption con?gurations correspond to local minima on the potential energy surface. In the case of H adsorption at 1F site, two imaginary frequencies were determined indicating that this site is a second-order saddle point. This characteristic is valid at both low and full coverages and has not been observed by Blyhoder et al. . In contradistinction, all vibrational frequencies of H adsorbed at 2F and 4F sites were positive, indicating that these sites are indeed local minima on the potential energy surface. These ?ndings related to the nature of the critical points on the potential energy surface are similar to those reported by Jiang and Carter based on PAW-PBE analysis . From the inspection of the vibrational frequencies at different sites we observe that these generally decrease with the increase of ligancy as determined previously by Blyholder et al. . Additionally, we observe that in the regime of low coverages the difference between vibrational frequencies at 2F and 4F sites is small. By increasing the coverage from u = 0.25 to 1, there is a more distinct difference among frequencies of these sites with values around 1205–1221 cm?1 for H adsorbed at 2F sites and between 1053 and 1127 cm?1 for H adsorbed at 4F sites. These values support the view expressed in Ref.  that in the regime of high coverages the vibrational band at 1000 cm?1 is characteristic to H adsorbed at the 4F hollow site. This is consistent also with the results obtained by Blyholder et al. . However, we have found no evidence in the case of low coverages of vibrational bands centered around 700 cm?1 as described by Merril and Madix . We have shown in previous section that for coverages above u = 1 it is possible to have adsorption con?gurations, such as 4F(4)–1F(1). For such states, we have also evaluated the corresponding vibrational frequencies and have determined that they are positive de?nite. In particular, a vibrational frequency of 1815 cm?1 was determined for H atom adsorbed at 1F site when the neighbor 4F sites are fully occupied. This result suggests that once all 4F sites are occupied further adsorption at 1F site becomes possible leading to con?gurations, such as 4F(4)–1F(1).
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Fig. 6. Iso-surfaces of the difference electron density for H adsorbed at (a) one-fold, (b) two-fold and (c) four-fold sites. The light regions enclosing the H atoms correspond to charge accumulation while darker regions correspond to charge depletion.
3.7. Charge distribution analysis An interesting ?nding of the adsorption calculations presented in Section 3.5 for coverages in the range 0.25 u 1 was that adsorption energies decrease rapidly for 1F sites and somewhat less signi?cant for 2F sites indicating the existence of repulsive interactions in these cases. In contradistinction, adsorption at 4F sites presents a small increase of the binding energy with coverage indicating the existence of a small attraction between H atoms. In order to gain more insight into the factors that determine this type of behavior, we have analyzed the difference electron densities, Dr ? rH?slab ? rslab ? rH (4)
where rH–slab, rslab and rH are the charge densities for the H/ Fe(1 0 0) slab, the bare slab and isolated H atom having the same position as the adsorbed H, respectively. The isosurfaces of the difference electron density corresponding to H adsorbed at 1F, 2F and 4F sites are given in Fig. 6. We found that in all these cases there is charge depletion from Fe surface atoms and charge accumulation on H adatoms. The amount of charge transferred however depends on the adsorption site. By integrating Dr around the H atom within ? an arbitrary sphere of radius 0.9 A, we found a charge transfer of about 0.06e for H adsorbed at 4F site and about 0.2e when H adsorbs at 1F site. These results indicate that electric charge ?ows from the surface to the adsorbed H atoms. As a result, a dipole moment will be generated and its amplitude increases with the separation of the H atom from the surface. By increasing the surface coverage, the electrostatic repulsion between H atoms increases. This is seen particularly for the case of adsorption at 1F sites where atoms are located far from the surface and somewhat less strongly for 2F sites. In the case of 4F sites, H atoms are close to the surface and consequently there is a much better screening of the resulted dipole moment by the surface electrons leading to negligible interactions with coverage
increase. A similar model to the one described here was shown to explain the modi?cations of adsorption energies of H on Ni(1 1 1), (1 0 0) and (1 1 0) surfaces . In particular, it was found that on Ni (1 0 0) surface there are signi?cant repulsive interactions among H atoms adsorbed at 1F sites while small attractive interactions were determined in the case of H atoms adsorbed at hollow sites. A further con?rmation of the charge redistribution effects described above can be obtained from the analysis of the work function. Particularly, it is expected that as a result of charge accumulation on H atoms the largest modi?cation of the work function will take place in the case of 1F sites while the smallest variation will be seen for 4F sites. We have determined the modi?cations of the work function for adsorption at 1F and 4F sites for the case of u = 0.25 and 1.0 coverages. The calculated results con?rm the above expectations, namely the change of the work function for 1F sites of 1.4 eV is signi?cantly larger than modi?cation for 4F sites of 0.02 eV. 3.8. Minimum energy paths for H diffusion 3.8.1. Surface diffusion We have shown in previous sections that the stable con?gurations of H on Fe(1 0 0) surface correspond to adsorption at bridge (2F) and hollow (4F) sites. Consequently, in order to characterize the surface diffusion we have considered the motion between these sites. The minimum energy paths corresponding to these processes are represented in Fig. 7. The potential energy surface describing the motion from 2F ! 4F con?gurations (see Fig. 7) has been investigated for the case of two different coverages. In the ?rst case corresponding to u = 0.25 there is a single H atom in the supercell which diffuses from 2F to 4F sites. The corresponding minimum energy potential path from Fig. 7 indicates the existence of a classical barrier of 1.93 kcal/mol. If zero-point energy corrections are also
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Fig. 7. Potential energy surface for diffusion of H from a two-fold site to the neighboring four-fold site at u = 0.25 and 0.75 coverages. In the case of u = 0.75 coverage, the 4F sites lateral to diffusion pathway are occupied.
included this barrier drops to 1.22 kcal/mol. In the region close to 4F site, corresponding to images 7–9, this potential does not present signi?cant variations indicating the existence of a shallow region. However, no additional minimum has been identi?ed in this region. Such a minimum denoted as displaced bridge-center was previously proposed by Blyholder et al.  and used by Merrill and Madix  to argue in favor of their proposed pseudo threefold site. The results presented here however do not support the existence of such a minimum. For the same diffusion pathway, we have analyzed the case when the adjacent 4F sites on both sides of the diffusion path are occupied. This case corresponds to u = 0.75 coverage. We found that in this case the corresponding diffusion barrier increases slightly to about 2.3 kcal/mol. From the results presented in Sections 3.5 and 3.6, it follows that the 1F site corresponds to a second-order saddle point on the potential energy surface. As a result, diffusion from 1F ! 2F and 1F ! 4F is a simple barrier-less downhill motion. Consequently, diffusion along 2F ! 1F ! 2F and 4F ! 1F ! 4F pathways require a signi?cantly larger amount of energy of about 12 kcal/mol than the one found for 2F $ 4F pathway, making the former pathways unfavorable energetic. The above results indicate that the diffusion barrier between 2F ! 4F sites is quite small with values of 1.93 and 1.77 kcal/mol for direct and reverse pathways. This barrier increases slightly with coverage to 2.30 kcal/mol (2.15 kcal/ mol) at u = 0.75 for direct (reverse) pathways. The calculated diffusion barriers discussed above correspond only to a classical hopping model of diffusion. In practical
cases, however, quantum tunneling effects should also be considered, particularly in the regime of low temperatures. 3.8.2. Surface to subsurface H diffusion We have shown in Section 3.5 that subsurface absorption of H at a tetrahedral site is possible on Fe(1 0 0) surface. In this section, we analyze the minimum energy pathway for diffusion from a surface adsorption site to a subsurface site. For this purpose, we have analyzed two different pathways. The ?rst one corresponds to diffusion from the surface 2F site to the S1 subsurface site while the second process corresponds to diffusion from the surface 4F hollow site to the subsurface S1 site. Both these pathways have been calculated using nudged elastic band method and PW91USPP. The surface has been modeled as a slab with 6 layers with the top three layers allowed to relax. The potential diagrams for the two indicated diffusion pathways are given in Fig. 8. We found that diffusion barriers for 2F ! S1 and 4F ! S1 processes are comparable having values of 9.55 and 7.52 kcal/mol, respectively. When comparing these values to those obtained for diffusion on the surface it is clear that H diffusion towards the subsurface is less favorable. However, the barriers to subsurface sites can be overcome with temperature increase making diffusion and adsorption to subsurface sites possible. Jiang and Carter have also analyzed hydrogen absorption in the subsurface layer of Fe(1 0 0) surface . They have determined a barrier of 8.76 kcal/mol for H to diffuse from surface to subsurface state, similar to values reported in this study. From the subsurface S1 state H atom can further diffuse into bulk. We have determined the classical diffusion
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Fig. 8. Potential energy surfaces for (a) diffusion of H from a two-fold to a subsurface S1 con?guration and (b) from a four-fold to a subsurface S1 con?guration.
barriers for a sequence of pathways involving tetrahedral adsorption sites S1 ! S2 ! S3 ! S4 positioned at successively deeper layers under the surface. The corresponding con?gurations, together with the calculated minimum energy pathways, are given in Fig. 9. From this analysis, we obtained that the largest diffusion barrier of 8.44 kcal/ mol is seen for the initial S1 ! S2 process. Once H atom reaches the second layer, further diffusion towards the bulk is possible with much lower activation energies. Namely, the barrier for S2 ! S3 step is about 1.89 kcal/mol while for S3 ! S4 process it is 2.1 kcal/mol. The overall picture emerging from these calculations is that the largest barriers for H diffusion from top of the surface towards the bulk are seen for surface to subsurface step as well as for diffusion among ?rst and second layers
underneath the surface. Deeper in the bulk the corresponding diffusion barriers decrease in magnitude to values close to those obtained in bulk, as will be shown in the next section. The calculated diffusion barriers presented above are expected to be only upper limits to the real diffusion barriers due to omission of quantum effects, such as tunneling. 3.9. Hydrogen absorption and diffusion in bulk Fe A ?nal problem considered in this study was related to absorption and diffusion of H atom in bulk bcc Fe. This problem has been studied based on PW91-USPP and PW91PAW methods with 2 ? 2 ? 2 and 3 ? 3 ? 3 supercells, containing 16 and 54 Fe atoms, respectively. The
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Fig. 9. Potential energy surface for H diffusion among different subsurface sites at increasingly larger depths below the surface. The atomic con?gurations corresponding to various local minima along the diffusing pathway are also indicated.
corresponding k-point grids of 6 ? 6 ? 6 and 4 ? 4 ? 4 have been used for the smaller and for the larger supercells, respectively. The ?rst objective of our calculations focused on description of the absorption properties of H at tetrahedral and octahedral sites, proposed as main absorption sites in previous experimental studies [37–39]. The atomic con?gurations of these sites are depicted in Fig. 10 while geometric and energetic parameters describing absorption properties are given in Table 5. In the case of 2 ? 2 ? 2 supercell, we found that the most stable absorption con?guration for H corresponds to tetrahedral site. At this site H atom has four nearest neighbors with whom it forms ? Fe–H bonds of about 1.667 A at PW91-USPP level or ? 1.654 A at PW91-PAW (see Table 5). For H adsorbed at octahedral site, there are two shorter Fe–H bonds of 1.575 ? ? (1.563) A and four other longer bonds of 2.018 (2.000) A as determined based on PW91-USPP (PW91-PAW) calculations. At both these theoretical levels we found that the binding energy at octahedral site is smaller than the one for tetrahedral site by 4.15 and 4.18 kcal/mol, respectively. Based on vibrational analysis of H adsorbed at tetrahedral and octahedral sites, we determined that tetrahedral site corresponds to a true minimum while octahedral site is a
second-order saddle point on the potential energy surface. These results are similar but not identical to those obtained before by Walch . He found an energetic difference of about 4.6 kcal/mol between the octahedral and tetrahedral sites but the octahedral site was characterized as a maximum on the potential energy surface instead of a saddle point as we found here. The dependence of the above-presented results on the decrease in H concentration was determined from calculations performed using a larger 3 ? 3 ? 3 supercell. The corresponding geometric and energetic parameters are given in Table 5. We found no modi?cation in the stability-order of tetrahedral and octahedral sites with the octahedral site being higher in energy by 3.39 (3.34) kcal/mol at USPP (PAW) levels than the tetrahedral site, respectively. When going from 2 ? 2 ? 2 to 3 ? 3 ? 3 supercells we found only small variations of the Fe–H bonds at tetrahedral site. These results indicate that at tetrahedral site H atom is well accommodated by crystalline lattice without introduction of signi?cant stresses. In contradistinction, for H adsorbed at octahedral site we determined larger relaxations. These are observed in the particular case of short Fe–H bonds which ? ? expand by 0.019 (0.02) A to 1.594 (1.583) A at USPP (PAW) levels while longer Fe–H bonds contract to about 1.999
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Fig. 10. Potential energy surface for diffusion of H in bulk bcc iron between tetrahedral sites. Case (a) corresponds to a pathway through a trigonal site (see the middle panel) while case (b) corresponds to a pathway through an octahedral site (see the middle panel).
? (1.989) A. These larger deformations indicate that H atom at octahedral site is in a more strained position than at tetrahedral site and as a result the former site is not stable. For the tetrahedral sites we evaluated the absorption energies using a relation similar to Eq. (1) where the energies of slab and slab–H systems are changed to those of the bulk and bulk–H systems, respectively. In the case of 2 ? 2 ? 2 supercell, we found that the corresponding absorption energies have values of ?1.87 kcal/mol (?4.51 kcal/mol) at PW91-USPP (PW91-PAW) levels.
The negative sign of these energies indicates in our sign convention that absorption of H in bulk bcc Fe is endothermic. In order to facilitate comparison with experimental data, we have also determined zero point energy corrections corresponding to H atom adsorbed at tetrahedral site and for isolated H2 molecule. When these corrections are included the heats of solution of ?4.19 kcal/ mol at PW91-USPP and ?7.05 kcal/mol at PW91-PAW levels are obtained. These values increase slightly to ?4.06 kcal/mol at PW91-USPP and ?6.75 kcal/mol at
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
Table 5 Calculated interatomic Fe–H distances (R(H–Fe)), relative energies (DE) and dissolution energies (DEdiss) for H adsorbed in bulk bcc Fe using PW91-USPP and PW91-PAW calculationsa Con?gurationb PW91-USPP R(H–Fe) (2 ? 2 ? 2) Structure Th Oh Tr (3 ? 3 ? 3) Structure Th Oh Tr Exp. 
a b c d
PW91-PAW DEc 0.00 4.15 2.45 0.0 3.39 2.16 DEdiss ?4.19 R(H–Fe) 1.654 (4) 1.563 (2), 2.000 (4) 1.593 (2), 1.649 (1) 1.661 (4) 1.583 (2), 1.989 (4) 1.603 (2), 1.668 (1) DE 0.00 4.17 2.49 0.0 3.34 2.15 DEdiss ?7.05
1.667 (4)d) 1.575 (2), 2.018 (4) 1.605 (2), 1.662 (1) 1.672 (4)d 1.594 (2), 1.992 (4) 1.615 (2), 1.678 (1)
Geometric parameters are given in Angstroms and energies in kcal/mol. Symbols Th, Oh and Tr denote the tetrahedral, octahedral and trigonal sites in the bulk. The relative energies are given with respect to the energy of H adsorbed at tetrahedral site. The number given in parentheses indicates the number of equivalent Fe–H bonds for each site.
PW91-PAW levels, in the case of the larger 3 ? 3 ? 3 supercell. The calculated heat of solution at PW91-PAW level compares very well to the experimental value of ?6.83 kcal/mol  while the value obtained based on PW91-USPP calculations is less accurate. Additionally, our calculated value of ?6.83 kcal/mol compare very well to the one of ?6.91 kcal/mol determined by Jiang and Carter based on PAW-PBE calculations . Overall, the results presented in this section and those related to bulk properties of iron described in Section 3.1, indicate that PW91-PAW method provides a more accurate description than PW91USPP calculations of various properties of bulk iron and for the interaction of H with the iron lattice. As a result, PAW calculations represent the method of choice among the methods analyzed here. Beside identi?cation of the absorption sites the values of the diffusion coef?cients of hydrogen in iron has been subject of signi?cant debate. Kiuchi and McLellan  have compiled a list of 62 sets of data for the diffusivity of H through bcc Fe. It was concluded that there exist signi?cant inconsistencies among various sets with values distributed anywhere between 1.16 and 16.5 kcal/mol. The main reason for the scatter in data has been attributed to the trapping of hydrogen by lattice defects [38,39] or to spurious surface effects [37,38]. Based on electrochemical methods Kuchi and McLellan have determined that in the regime of low temperatures ?408 to 808 H diffusion obey an Arrhenius dependence and has a low activation energy of 1.36 kcal/ mol. This energy corresponds to diffusion between adjacent tetrahedral sites. In the regime of higher temperatures 50– 5508, a slight increase in the activation energy to 1.6– 1.7 kcal/mol has been determined, effect associated to an increase of the octahedral site occupation. More recently, Hagi  has also determined the diffusion coef?cient of hydrogen in iron free of trapping by dislocations and impurities using electrochemical permeation measurements. He estimated an activation barrier of 1.07 ? 0.04 kcal/mol based on the trapping theory. These results indicate that
recent electrochemical measurements were able to provide more consistent values for the barrier of H diffusion in the range 1.07–1.36 kcal/mol. In this study, we have determined theoretically the diffusion barrier of H in bulk bcc Fe and compared this to results of electrochemical measurements. For this purpose, we have analyzed two different pathways for diffusion among neighbor tetrahedral sites. The ?rst one corresponds to a direct diffusion pathway among tetrahedral sites (see Fig. 10a) while the second corresponds to diffusion through an octahedral site (Fig. 10b). Nudged elastic band calculations have been done for the 2 ? 2 ? 2 supercell described above. In the case of the ?rst pathway we identi?ed a transition state corresponding to a trigonal site (see middle panel below Fig. 10a positioned at about 2.45 kcal/mol (2.49 kcal/mol) above the tetragonal site as determined by PW91-USPP (PW91-PAW) calculations. In the trigonal con?guration, H atom binds to three neighbor Fe atoms such that two of the Fe–H bonds are equal. The corresponding values of these bonds are given in Table 5. For the pathway through octahedral site, we found that the activation energy is equal to the difference between the energies of octahedral and tetrahedral sites, which is expected as the octahedral site is a saddle point on the potential energy surface. The calculated diffusion barriers depend also slightly on H concentration. For example, from optimizations of the tetrahedral and trigonal con?gurations in the 3 ? 3 ? 3 supercell (see Table 5) we found that diffusion barriers decrease to 2.16 kcal/mol (2.15 kcal/mol) at PW91-USPP (PW91-PAW) levels. When zero point energy corrections for the equilibrium (tetrahedral) site and for transition state (trigonal site) are included the corrected diffusion barriers are 1.16 kcal/mol (1.09 kcal/mol) at PW91-USPP (PW91PAW) levels. These calculated diffusion barriers agree quite well with previous experimental measurements [37,39] which found values in the range 1.1–1.36 kcal/mol. Additionally, our
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
calculations are also in good agreement with recent theoretical results determined by Jiang and Carter  who determined a zero-point-energy-corrected diffusion barrier of 0.97 kcal/mol. Moreover, in agreement with the results of Ref.  the most favorable diffusion pathway of H inside bulk Fe we found corresponds to hoping among neighbor tetrahedral sites through a path containing a trigonal site. The diffusion pathway through the octahedral site is energetically not favorable as this site corresponds to a second-order saddle point on the potential energy surface. Overall, the results presented in this section con?rm previous experimental ?ndings [37,39] that in iron free of dislocations and impurities there is a very low activation energy for H diffusion and that the most energetically channel for H diffusion among tetrahedral sites takes place through a trigonal site.
4. Conclusions We have performed plane-wave DFT calculations with PW91-USPP and PW91-PAW methods to calculate the adsorption and diffusion properties of hydrogen on Fe(1 0 0) surface and in bcc Fe bulk. The main conclusions of this study can be summarized as follows:
(f) (a) We have determined that dissociative chemisorption of H2 molecule on Fe surface can readily take place with a classical barrier of about 3.6 kcal/mol. A weakly physisorbed state above Fe surface atoms can serve as a precursor state to dissociative chemisorption. (b) Adsorption of H on Fe(1 0 0) at low coverages (u = 0.25) can take place at both 2F and 4F sites with binding energies in the range 59–61 kcal/mol. These values have been con?rmed using slab models with six and seven layers, where adsorption was done, respectively, on one and both sides of the respective slabs. No signi?cant differences between the adsorption energies calculated using PW91-USPP and PW91-PAW methods have been found. However, the use of RPBE functional leads to adsorption energies, which on average are about 3.0 kcal/mol smaller than those determined using PW91 functional. Relative to isolated H2 in gas phase, adsorption of H at 2F and 4F sites is exothermic while adsorption at 1F sites is endothermic. (c) A detailed picture of the variation of H binding energy on surface coverage in the range 0.25–3.0 has been determined. It was found that with the increase of coverage from 0.25 to 1 adsorption at 4F site is the most stable. Over the coverage range 1.0–2.0 in addition to fully occupied 4F sites adsorptions at 2F sites or mixed 2F–1F sites are energetically favorable. The increase of H coverage above 2.0 leads H2 recombination as a result of attractive interactions between H atoms from neighbor 1F and 2F sites. The resulted H2 species weakly adsorb at 1F sites in a con?guration parallel to
the surface. By the increase in coverage the binding energies of such con?guration decrease, fact that should facilitate H2 desorption. Beside adsorption at surface sites, atomic H can also absorb at subsurface sites in a tetrahedral con?guration. These sites are positioned just beneath the surface at positions below 2F surface sites. The corresponding binding energies were found to increase slightly with coverage reaching values of 53.88 kcal/mol when all subsurface sites are occupied. Other subsurface absorption sites in the second and third layers have been identi?ed to correspond also to tetrahedral sites but absorption at these sites has smaller binding energies than those obtained for subsurface S1 sites. Signi?cant relaxations of surface top layers have been observed as a result of hydrogen absorption. Relative to gas phase H2 absorption at S1 subsurface sites is exothermic while absorption at deeper sites S2 and S3 is endothermic. Diffusion of H atom on surface takes place with a small classical barrier of about 1.9 kcal/mol among the 2F and 4F sites. This barrier was found to increase slightly to 2.3 kcal/mol with the increase of coverage from 0.25 to 0.75. Diffusion among 2F sites or 4F sites passing through 1F sites has signi?cantly larger barriers of about 12 kcal/mol. From the surface sites H can diffuse to subsurface and further to the bulk sites. The largest barriers in these processes were obtained for the surface to subsurface step with diffusion barriers in the range 7.52–9.55 kcal/ mol, and for diffusion among ?rst and second layers underneath the surface with barriers of 8.44 kcal/mol. Beyond second layer the diffusion barriers decrease more signi?cantly to values close to those obtained in bulk Fe. Absorption and diffusion of H in bulk Fe were investigated and the most stable absorption con?guration was found to correspond to interstitial tetrahedral sites. The octahedral sites were found to correspond to a second-order saddle point on the potential energy surface. Relative to gas phase H2 absorption of H in bulk bcc Fe is endothermic. The PW91-PAW calculated heat of solution of ?6.75 kcal/mol for H adsorbed at tetrahedral sites was found to be in very good agreement to experimental value of ?6.83 kcal/mol while the PW91-USPP result is slightly less satisfactorily (?4.06 kcal/mol). Among tetrahedral sites H can diffuse with small activation energies of about 1.1 kcal/mol through a path through the trigonal site. Such low diffusion barriers support the view that mobility of H inside bulk bcc Fe is very high.
Acknowledgments We gratefully acknowledge the computational resources provided under a Challenge award received from the DOD
D.C. Sorescu / Catalysis Today 105 (2005) 44–65
High Performance Computing Modernization Of?ce, and by ARL MSRC supercomputer center.
¨  G. Ertl, H. Knozinger, J. Weltkamp (Eds.), Handbook of Heterogeneous Catalysis, VCH, Weinheim, 1997.  G.A. Somorjai, Introduction to Surface Chemistry and Catalysis, Wiley, New York, 1994.  (a) F. Fischer, H. Tropsch, Brennst. Chem. 4 (1923) 276; (b) R.B. Anderson, The Fischer–Tropsch Synthesis, Academic, Orlando, FL, 1984.  R. Gibala, R.F. Hehemann (Eds.), Hydrogen Embrittlement and Stress Corrosion Cracking: A Troiano Festschrift, ASM International, 1984.  F. Bozso, G. Ertl, M. Grunze, M. Weiss, Appl. Surf. Sci. 1 (1977) 103.  J. Benziger, R.J. Madix, Surf. Sci. 94 (1980) 119.  M.L. Burke, R.J. Madix, Surf. Sci. 237 (1990) 20.  P.B. Merrill, R.J. Madix, Surf. Sci. 347 (1996) 249.  W. Moritz, R. Imbihl, R.J. Behm, G. Ertl, T. Matsushima, J. Chem. Phys. 83 (1985) 1959. ?  A.M. Baro, W. Erley, Surf. Sci. 112 (1981) 759.  G. Blyholder, J. Head, F. Ruette, Surf. Sci. 131 (1983) 403.  T.J. Raeker, A.E. DePristo, Surf. Sci. 235 (1990) 84.  S.P. Walch, Surf. Sci. 143 (1984) 188.  M. Eder, K. Terakura, J. Hafner, Phys. Rev. B (2001) 115426.  D.E. Jiang, E.A. Carter, Phys. Rev. B 70 (2004) 064102.  G. Kresse, J. Hafner, Phys. Rev. B 48 (1993) 13115. ¨  G. Kresse, J. Furthmuller, Comput. Mater. Sci. 6 (1996) 15. ¨  G. Kresse, J. Furthmuller, Phys. Rev. B54 (1996) 11169.  D. Vanderbilt, Phys. Rev. B41 (1990) 7892.  G. Kresse, J. Hafner, J. Phys. Condens. Matter 6 (1994) 8245. ¨  E.G. Moroni, G. Kresse, J. Hafner, J. Furthmuller, Phys. Rev. B56 (1997) 15629.
     
    
    
¨ P.E. Blochl, Phys. Rev. B50 (1994) 17953. G. Kresse, D. Joubert, Phys. Rev. B59 (1999) 1758. H.J. Monkhorst, J.D. Pack, Phys. Rev. B13 (1976) 5188. M. Methfessel, A.T. Paxton, Phys. Rev. B40 (1989) 3616. G. Kresse, J. Hafner, Phys. Rev. B47 (1993) 558. ? (a) G. Mills, H. Jonsson, Phys. Rev. Lett. 72 (1994) 1124; ? (b) G. Mills, H. Jonsson, G.K. Schenter, Surf. Sci. 324 (1995) 305; ? (c) H. Jonsson, G. Mills, K.W. Jacobsen, Nudged elastic band method for ?nding minimum energy paths of transitions, in: B.J. Berne, G. Ciccotti, D.F. Coker (Eds.), Classical Quantum Dynamics in Condensed Phase Simulations, World Scienti?c, 1998, p. 385. (a) J.P. Perdew, Y. Wang, Phys. Rev. B45 (1992) 13244; (b) J.P. Perdew, J.A. Chevary, S.H. Vosko, K.A. Jackson, M.R. Pedersen, D.J. Singh, C. Frolhais, Phys. Rev. B46 (1992) 6671. S.H. Vosko, L. Wilk, M. Nusair, Can. J. Phys. 58 (1980) 1200. (a) J.P. Perdew, K. Burke, M. Ernzerhof, Phys. Rev. Lett. 77 (1996) 3865; (b) B. Hammer, L.B. Hansen, J.K. N?rskov, Phys. Rev. B59 (1999) 7413. D.C. Sorescu, D.L. Thompson, M.M. Hurley, C.F. Chabalowski, Phys. Rev. B 66 (2002) 035416. C. Kittel, Introduction to Solid State Physics, Wiley, New York, 1996. D.E. Jiang, E.A. Carter, Surf. Sci. 547 (2003) 85. Z.Q. Wang, Y.S. Li, F. Jona, P.M. Marcus, Solid State Commun. 61 (1987) 623. K.P. Huber, G. Herzberg, Molecular Spectra and Molecular Structure IV. Constants of Diatomic Molecules, Van Norstrand Rheinhold, 1979 . G. Kresse, J. Hafner, Surf. Sci. 459 (2000) 287. K. Kiuchi, R.B. McLellan, Acta Metall. 31 (1983) 961. ¨ J. Volkl, G. Alefeld, in: A.S. Nowick, J.J. Burton (Eds.), Diffusion in Solids, Academic Press, New York, 1975. H. Hagi, Mater. Trans. JIM 35 (1990) 112. J.P. Hirth, Metall. Trans. A 11 (1980) 861.