Density Functional Theory and Kinetic Studies of Methanation on Iron

Jun 29, 2007 - John M. H. Lo, and Tom Ziegler*. Department of Chemistry, University of ..... Julien Couble , Daniel Bianchi. Applied Catalysis A: Gene...
2 downloads 0 Views 357KB Size
11012

J. Phys. Chem. C 2007, 111, 11012-11025

Density Functional Theory and Kinetic Studies of Methanation on Iron Surface John M. H. Lo and Tom Ziegler* Department of Chemistry, UniVersity of Calgary, Alberta, T2N 1N4 Canada ReceiVed: March 20, 2007; In Final Form: May 14, 2007

The production of methane from CO and H2 on a clean Fe(100) surface was studied using periodic DFT calculations in conjunction with kinetic modeling. It was found that the adsorption and subsequent dissociation of CO and H2 at a finite temperature are significantly influenced by the presence of coadsorbed CO and H. The surface coverage of CO is controlled by entropy, and a maximum of 50% coverage can be attained. The kinetic modeling adopting the computed activation barriers for the elementary steps involved in methanation revealed that CH is a major C1 species on the surface. The production of CH4 was found to be more favored at a high reaction temperature and H2 partial pressure but suppressed by a high pressure of CO.

1. Introduction Fischer-Tropsch (FT) synthesis is an environmentally friendly and economically promising process that produces high-quality chemicals and liquid fuels from synthesis gas (CO and H2).1-3 This process has been widely used since its first discovery in the 1920s,4 especially in Germany and Japan during the second World War and also in the United States afterward. Currently, the FT technology is utilized in Malaysia and South Africa to generate various forms of valuable diesels and petroleum products with low sulfur content. It is anticipated that the FT synthesis will play a more important role in the future petroleum industry due to the increasing global demands of liquid fuels and the presence of large reserves of natural gas. Consequently, many experimental5-11 and theoretical12-16 studies have recently been devoted to unraveling the complex mechanisms of this process through which a large array of linear, branched, and oxygenated hydrocarbons are formed. One of the obstacles for the mechanistic studies of FT synthesis stems from the differences in product selectivity of various FT active catalysts. The commercial FT synthesis is often performed on the surface of iron, cobalt, or nickel catalysts supported on either aluminum or titanium oxide;1,3 frequently, an alkali promoter such as potassium is added to enhance the selectivity.2 It has been observed that while iron and cobalt favor the production of long-chain hydrocarbons, nickel shows a high selectivity toward methane;17 this diversity of functionalities has been attributed to the different surface and electronic structures of the active sites on these catalysts.18,19 In order to account for the abundant kinetic and spectroscopic data from the experiments, various reaction mechanisms and kinetic models comprising different elementary reaction steps have been developed.20,21 Among these are three principal types of mechanisms, the carbene mechanism,4 the hydroxylcarbene mechanism,22 and the carbonyl insertion mechanism.23 These mechanisms, despite the different proposed steps of * To whom correspondence should be addressed. E-mail: ziegler@ ucalgary.ca.

hydrogenation and CC bond formation, share the same general elementary reaction steps

CO + 3H2 f CH4 + H2O

(1)

1 CO + 2H2 f (CnH2n) + H2O n

(2)

CO + H2O a CO2 + H2

(3)

The carbene mechanism is the currently most accepted scheme due to the support by the current experimental studies,8,24,25 although there is evidence that several mechanisms may operate simultaneously during the FT process.7 The carbene mechanism consists of six elementary steps, (1) the adsorption and dissociation of CO and H2, (2) formation of H2O and CO2 from the adsorbed C, O, and H, (3) the hydrogenation of surface C to yield various alkyl groups CHn, (4) polymerization of surface CHn, (5) the formation of alkanes via the terminal hydrogenation, and (6) β-elimination of the surface alkyl and subsequent hydrogenation to produce olefins. In the present study, the methane production from CO and H2, which constitutes the initiation step of the FT process, on an iron surface was investigated. Iron is among the most common metals used as FT catalysts because of its advantages from the economical perspective and its comparable activities to more expensive cobalt and ruthenium catalysts. However, while the process of CH4 formation on the latter catalysts and the influence due to the structures of catalysts have been thoroughly studied,14,26 it is only recently that similar investigations have been carried out on iron.15,16 The present study affords the first extensive DFT investigation concerning the kinetics of methane production from CO and H2 catalyzed by iron. The present work explored the adsorption and dissociation of CO and H2 on Fe(100) and the subsequent hydrogenation reactions that lead to the formation of H2O, CHn, and CO2 (Steps 1-6). Fe(100) possesses a similar surface energy to the thermodynamically more stable Fe(110) surface,27 but it is believed that the surface relaxation and reconstruction during the surface cleavage account for an extra stability of open

10.1021/jp0722206 CCC: $37.00 © 2007 American Chemical Society Published on Web 06/29/2007

Methanation on an Iron Surface

J. Phys. Chem. C, Vol. 111, No. 29, 2007 11013

Fe(100). In all calculations performed in this work, a clean Fe(100) was assumed. The current study was limited to C1 units so as to simplify the calculations. In addition to discovering the stable configurations of adsorbed CHn species and computing the potential energy surfaces corresponding to their formations, the lateral interaction due to the presence of coadsorbed species was also studied. Kinetic modeling was performed to examine the dynamics and temperature dependence of methanation on Fe(100) by utilizing the computed data. 2. Computational Details 2.1. Calculations of Activation Energies of Elementary Steps. The electronic structure calculations performed in this work were done using the periodic density functional theory implemented in the Vienna ab initio simulation package (VASP).28-30 The fully nonlocal Vanderbilt ultrasoft pseudopotentials (usPPs),31,32 which describe the electron-ion interaction, and plane-wave basis sets with a kinetic energy cutoff of 400 eV were employed in the spin-polarized calculations, in which the exchange and correlation were described by the proposed functional of Perdew and Wang,33 referred to as PW91, in the generalized gradient approximation (GGA).34 A Methfessel-Paxton smearing function of 0.2 eV was utilized in all calculations.35 The lattice parameter of body-centered Fe(100) was determined by optimizing the smallest asymmetric unit with various sets of Monkhorst-Pack k-point grids36 in the reciprocal space. The computed lattice constant of 2.8553 Å, the bulk modulus of 156 GPa, and the associated magnetic moment of 2.30 µ0, based on the 12 × 12 × 12 k-point mesh, are all in very good agreement with the experimental values of 2.8665 Å37, 170 GPa,38 and 2.22 µ0,39 respectively. Individual spin-polarized calculations on the atoms C, H, and O and molecular systems CO, H2O, H2, and CHn were performed in a 10 Å × 10 Å × 10 Å box to determine their structures and energetics in the gas phase. Reasonable agreement with the experimental geometries was achieved in these calculations; for instance, the computed equilibrium CO bond distance of 1.144 Å is only 1.4% longer than the experimental value, 1.128 Å.40 The Fe(100) surface was modeled with a five-layer slab of Fe with a 2 × 2 supercell structure. Only the adsorption of molecules on one side of the slab was considered; therefore, a vacuum of 10 Å was introduced between two repeated slabs in order to minimize the induced dipole interaction. During a geometry optimization, the adsorbates and first two layers of Fe atoms were allowed to relax using the RMM-DIIS algorithm,41 while the remaining Fe layers were frozen at their optimized lattice positions. Geometry optimizations were completed when all of the forces were smaller than 0.01 eV/Å. Binding energies of adsorbates were calculated according to the following equation

Eads )

Eslab + NEfrag - Efrag+slab N

(4)

where Efrag is the energy of isolated fragment, Eslab is the energy of the metal slab, and Efrag+slab is the total energy of the adsorbate-slab system. The parameter N refers to the number of adsorbate molecules on the surface. Notice that in the calculations concerning H2, Efrag was defined as (1/2)EH2. The Brillouin zones were sampled by a smaller 7 × 7 × 1 k-point grid in all supercell calculations in order to enhance the computational efficiency. Testing calculations revealed that

this k-point mesh yields the energy, which differs only by 6 × 10-4 eV from that obtained by using a larger 12 × 12 × 1 mesh. It is, however, that a slightly larger lattice constant of 2.8730 Å (about 0.6% longer) resulted, which may have led to the undesirable lattice strain when this k-point mesh was employed in conjunction with the lattice constant (i.e., 2.8553 Å) determined with the 12 × 12 × 12 mesh. This adverse effect would be significant in the calculations where the metal slabs are frozen but could be mostly compensated by relaxing the surface during the geometry optimization. Since surface relaxation was allowed in all calculations of the present work, it is believed that the influence on the computed geometries and energies of adsorbates should be minimal. Vibrational frequency calculations were also performed within the harmonic approximation by computing the Hessian matrix using numerical double force differentiation with atomic displacements of (0.02 Å. A partial frozen phonon approximation was applied in which one degree of freedom, along z-direction, for each Fe atom on the surface was included in the Hessian calculations. Comparisons with the cases where the full frozen phonon approximation was applied showed that artificial imaginary modes may be introduced when the whole slab is assumed frozen. To determine the activation energies of the elementary steps occurring on the Fe(100) surface, transition-state structures were located by means of the nudged elastic band (NEB) method implemented in VASP.42 This method simultaneously minimizes the intermediate states along the reaction path with respect to the total energy, with the restriction of atomic motions to the hyperplane perpendicular to the reaction path. The initial guesses of intermediate images were generated by linear interpolations between the start and end configurations. 2.2. Kinetic Model of FT Synthesis. Utilizing the activation energies supplied by the DFT calculations for various adsorbed species, the kinetics of CO hydrogenation on the Fe(100) surface was simulated. It has been shown recently by Jansen and coworkers in a surface study of NO on Rh(100) and Rh(111) that lateral interactions play a crucial role in the surface coverage, rates of reactions, and thus the reaction mechanism.43-45 In order to account for the influence of neighboring adsorbates on the barriers of reactions, the current model adopted a p(2 × 2) lattice as one unit of the reaction surface of Fe(100), each of which had four available active sites. The number of such units was determined by the total surface area of the catalyst. Every reaction unit was considered independent of each other in the calculations. At most, four species might have simultaneously adsorbed on a reaction unit, and the number of species determined the local coverage. It was assumed that reactions, that is, dissociation and bond formation and breaking, could only happen within a reaction unit. To simplify the modeling, no diffusion of surface species was allowed. Activation barriers of reactions were computed in accordance with the number of other coadsorbates present in the reaction unit; that is, the activation energy for a reaction was calculated for different surface coverages of the reaction unit. In the kinetic model, the rate constant of an elementary reaction was computed by the Eyring equation46

k)

(

)

kBT ∆Gact exp h RT

(5)

where kB, h, and R are Boltzmann constant, Planck’s constant, and the gas constant, respectively. T is the absolute simulation temperature in Kelvin. ∆Gact ) ∆Eact - T∆Sact is the change in Gibb’s free energy in the activation process. For the adsorption

11014 J. Phys. Chem. C, Vol. 111, No. 29, 2007

Lo and Ziegler

or desorption of a gaseous molecule, such as CO, ∆Sact was computed and included in order to account for the loss of the translational and parts of the rotational and vibrational entropies of the gaseous molecules when they approached to the surface.47 On the other hand, ∆Gact was set equal to the activation energy ∆Eact for reactions between adsorbed molecules. This assumption was based on the fact that the vibrational entropy is small (for instance, Svib ) 0.003 kcal/mol for CO at 400 K, as deduced in the present work), and compared to the translational and rotational counterparts (which are, respectively, 14.9 and 4.8 kcal/mol at 400 K), constitutes only a small fraction of the total entropy of a molecule. Neglecting its contribution did not cause a significant change in the computed activation barriers of the surface reactions. In addition, it was not possible to accurately evaluate the rotational entropy of an adsorbed molecule; therefore, it is assumed that all of the rotational entropy was lost when the molecule was adsorbed. In the case of reactions of adsorbed species, it is assumed that they possess a similar mobility and, in turn, similar rotational entropies. Consequently, the contribution to ∆Sact due to the rotational motion was ignored. The reaction rate constants of all of the elementary steps included in the process of CO hydrogenation were used to formulate the system of rate equations that describes the kinetic behavior of different species on Fe(100). Reaction units with various combinations of surface species were all considered distinct and taken into account in the modeling by individual differential equations. That is, for example, the units containing one adsorbed CO and two adsorbed CO molecules are two different species and have different reaction rates of additional CO adsorption and dissociation. Through this approach, the lateral effects of coadsorbed molecules on reaction rates could be explicitly dealt with. Three categories of reactions were considered in the present model. The first category concerns the adsorption and desorption of molecules such as H2, CO, CO2, H2O, and CH4

H2(g) + [4*] a [2*H + 2*] CO(g) + [4*] a [*CO + 3*] CO2(g) + [4*] a [*CO2 + 3*] H2O(g) + [4*] a [*H2O + 3*] CH4(g) + [4*] a [*CH4 + 3*]

(6)

The second category includes the dissociation processes associated with the molecularly adsorbed species

[*CO + 3*] a [*C + *O + 2*] [*CO2 + 3*] a [*CO + *O + 2*] [*H2O + 3*] a [*OH + *H + 2*]

(7)

The final category contains the hydrogenation reactions of surface carbon-containing species

[*C + *H + 2*] a [*CH + 3*] [*CH + *H + 2*] a [*CH2 + 3*] [*CH2 + *H + 2*] a [*CH3 + 3*] [*CH3 + *H + 2*] a [*CH4 + 3*]

(8)

In the reactions mentioned above, [*] and [*M] represent a vacant site and an adsorbed species M, respectively. To

incorporate the steric interaction induced by the neighboring adsorbed molecules, subsequent additions of CO and H2 were permitted in all possible reaction units, provided that one vacant site (for CO) or two vacant sites (for H2) were present. The system of rate equations derived according to these elementary reaction steps on Fe(100) were solved numerically for various temperatures and pressures by means of the Gear’s method48 implemented in the program suite of MATLAB.49,50 3. Results and Discussions 3.1. Electronic Structure Calculations of Surface Species. 3.1.1. Adsorption of H and CO. Due to its significance in the metallurgy and steel industry, hydrogen adsorption on a metal surface, especially iron, has been the topic of extensive investigations using various experimental techniques such as low-energy diffraction (LEED),51 X-ray photoelectron spectroscopy (XPS),52 and high-resolution electron energy loss spectroscopy (HREELS).53 From the theoretical perspective, several calculations have been performed to analyze the relative stability and preference of hydrogen adsorption sites on various low-index Fe surfaces.54-56 Sorescu has recently studied H adsorption on Fe(100) using GGA-PW91/usPP and GGA-PW91/ PAW methods57 and identified the most stable configuration of surface H as the hollow position with a binding energy of 7.2 kcal/mol. Jiang and Carter have also predicted the hollow site as most stable, with an adsorption energy of 8.1 kcal/mol.58 The present work confirms their finding that the four-fold site is the most favorable for surface H. The higher binding energy of 8.6 kcal/mol is attributed to the smaller calculated H-surface distance (3.44 versus 3.82 Å57). The adsorption of two and four H atoms on a p(2 × 2) Fe unit, which corresponds to, respectively, 0.50 and 1.00 ML coverage, was also studied. The results show no significant coverage dependence of the binding energy. The full coverage of Fe(100) only marginally increased the binding energy of H by 0.05 kcal/mol, in agreement with the prediction of Sorescu that H is more strongly bound by 0.3-1.5 kcal/mol when the surface coverage increases.57 It has been identified from the temperature-programmed desorption (TDP) and XPS spectra that there are four chemisorption modes of the CO molecule on the Fe(100) surface, (1) on-top, (2) bridge, (3) hollow, and (4) tilted configurations.59-61 The most stable configuration was determined by near-edge X-ray absorption fine structure (NEXRFS) to be a tilted CO molecule at the hollow position, with an angle of 45 ( 10° from the surface normal.62 Sorescu and co-workers, using GGAPW91/usPP and GGA-PBE/usPP, demonstrated that CO prefers the tilted configuration at the hollow position, and the corresponding binding energy is between 43.8 and 46.7 kcal/mol.63 A later study by Bromfield et al. afforded the same conclusions, although a much higher binding energy (≈59 kcal/mol) was reported.64 Very similar results were obtained in the present work. All of the four possible adsorption geometries were considered. The most stable configuration at low coverage was found to be the four-fold state with CO tilted by 49° from the normal. Its adsorption energy of 46.9 kcal/mol and vibrational frequency of 1165 cm-1 agree very well with both experiments and the results from Sorescu et al.63 The other three configurations were found, corresponding to the transition states. It was also found that at low coverage, the trend of relative stability was E(twofold) ≈ E(one-fold) < E(four-fold), where the one-fold state was only 1 kcal/mol more stable than the two-fold state. Compared to the case of surface H, there exists a more pronounced lateral interaction between coadsorbed CO mol-

Methanation on an Iron Surface ecules in the p(2 × 2) unit cell. The binding energy of the second CO slightly decreased to 46.2 kcal/mol, and the molecule adopted the tilted configuration on the hollow site diagonal to the first CO. The small change in energy is attributed to the competition for surface density from the Fe atoms between the two adjacent CO molecules. The subsequent CO additions, however, become less likely although still thermodynamically favorable. The adsorption energies of CO at 0.75 and 1.00 ML of coverage are 34.9 and 30.9 kcal/mol, respectively. Bromfield et al. have investigated the CO adsorption on Fe(100) at saturation (i.e., 1.00 ML) and obtained a consistent binding energy of 1.5 eV.64 As pointed out in their work, the third and fourth additions of CO alter the surface structure of CO-adsorbed Fe(100); CO molecules have to adopt the less favorable upright configuration when the coverage is higher than 0.5 ML because of the sterics, thus eliminating the Fe-O interaction and weakening the Fe-C bond. 3.1.2. The Influence of Lateral Interactions on the Dissociation and Adsorption of CO and H2. The presence of coadsorbed molecules on the surface influences both the thermodynamics and kinetics of the adsorption and dissociation of CO and H2. CO molecules readily underwent molecular adsorption on clean Fe(100), with merely a small barrier of 0.14 kcal/mol, which is in line with the experimental observation that CO adsorption takes place at 150 K59 at low coverage. However, the adsorption barrier increased with the surface coverage; at 0.50 ML, the barrier was 4.1 kcal/mol, while it became 13.0 kcal/mol as full surface coverage was achieved. The higher barrier is due to the fact that the adsorbed CO molecules induce a significant electrostatic repulsion that hinders the approach of additional CO. In addition to the larger adsorption barrier, the presence of coadsorbed CO altered the adsorption mode of CO molecules; at high coverage (>0.50 ML), CO molecules in the tilted configuration at the hollow sites were converted to the upright configuration so as to minimize the mutual repulsion and the competition of electron density from the surface Fe atoms. The neighboring molecules on the surface affected not only the molecular adsorption but also the dissociation of the adsorbed CO molecules. The calculated activation energy for CO dissociation at 0.25 ML was 25.6 kcal/mol, consistent with the experimental value of 26.2 kcal/mol62 reported by Moon et al. The corresponding activation energies at higher coverage were also computed, and an increasing trend was noticed. At 0.50 ML, the activation barrier was slightly lowered to 25.3 kcal/mol, and the reaction was favored thermodynamically by about 2 kcal/mol. These results are not in agreement with the conclusion of Bromfield et al.64 that this reaction has a higher activation energy and is endothermic. Nevertheless, the dissociation at even higher coverage (0.75 ML) turned out to be marginally endothermic by 0.1 kcal/mol, yet the activation energy increased rapidly to 32.6 kcal/mol. To complete the picture, the dissociation of CO in the presence of adsorbed C and O was studied. It was discovered that this reaction possesses a huge barrier of activation (35.6 kcal/mol) and is largely endothermic (by 10.4 kcal/mol). It can hence be concluded that high surface coverage stabilizes the molecularly adsorbed CO molecules and inhibits their dissociation. In contrast to the case of CO, the adsorption process of H2 is not influenced much by the surface coverage. It has been pointed out by Boszo et al.51 that H2 interacts with the Fe(100) surface via dissociative adsorption. For a clean Fe(100) surface, the addition of H2 to two adjacent hollow sites has the estimated barrier of 3.5 kcal/mol.57 The present work yielded a consistent value of 3.2 kcal/mol for this process. The addition of a second

J. Phys. Chem. C, Vol. 111, No. 29, 2007 11015 H2 was found to be exothermic by 15.3 kcal/mol, but the activation energy was slightly decreased to 2.1 kcal/mol. Replacing surface H with CO did not change the barrier significantly; in the presence of adsorbed CO, the activation energy for the dissociative addition of a second H2 was increased to 2.3 kcal/mol. Indeed, the rather coverage-independent activation energy of H2 is not surprising; the small size of the hydrogen molecule makes it insensitive to the steric repulsion from surface adsorbates. 3.1.3. Adsorption of Other Intermediate Species. In order to understand the details of the mechanism involved in the FT process, the structures of various proposed surface intermediates such as CHn (n ) 0-4), O, CO2, and H2O were investigated at 0.25 ML of coverage. C and O were the dissociation products of CO, and subsequent hydrogenation of C and O produced many different kinds of CHn species and water.1 When CO recombined with surface O, CO2 was generated. Three possible configurations of surface carbide, namely, onefold, two-fold, and four-fold states, were investigated. Among these states, the four-fold state was the most stable configuration, with an adsorption energy of 183.0 kcal/mol. The two-fold state was less stable than the four-fold state by 42.3 kcal/mol, and the on-top state was the least stable, possessing the binding energy of 118.9 kcal/mol. These results illustrate that carbon, compared to CO, was strongly bound to the Fe(100) surface, as observed by previous experiments.65 The carbide surface of iron is believed to be an active phase that catalyzes the synthesis of hydrocarbons in the FT process.66 Similar to carbon, the strongest binding of O appeared at the hollow position at 0.25 ML, with a computed adsorption energy of 79.2 kcal/mol. This value is in good agreement with the value recently deduced by Hafner and co-workers67 using the DFT method for low-surface coverage. The adsorption at on-top and bridging sites were less stable; their corresponding binding energies were 35.1 and 15.1 kcal/mol higher in energy than that of the four-fold state, respectively, and they correspond to transition states, as revealed by the vibrational frequency calculations. For both C and O, subsurface adsorption is also available, but the corresponding binding energies are much smaller compared to those for surface adsorption, as shown by Sorescu.16 Therefore, they were excluded in the kinetic modeling described in the following sections. Methylidyne CH was adsorbed on the Fe(100) surface in a linear fashion; CH attained an upright configuration and formed a zero angle with respect to the surface normal. Adsorption states of CH on one-fold, two-fold, and four-fold positions were analyzed. The results reveal that the hollow site is most favorable, with the binding energy of 163.5 kcal/mol. The previous studies of Sorescu16 also yielded the same conclusions, but the cluster calculations by Anderson predicted the most stable configuration at the bridging position.68 In the present work, one-fold and two-fold states were found to be less stable than the four-fold state, and they both correspond to transition states, contrary to the finding by Sorescu that the on-top state is the second local minimum configuration of adsorbed CH. In the case of CH2, again the adsorption on three different adsorption sites was considered. Similar to C and CH, the most stable configuration of adsorbed CH2 was at the hollow position. The calculated adsorption energy was 100.4 kcal/mol, which was much smaller than the CH and C binding energies because of the weaker interaction between sp2-hybridized C and the surface Fe atoms. The vibrational frequency analysis showed that this state is a local minimum on the potential energy surface.

11016 J. Phys. Chem. C, Vol. 111, No. 29, 2007

Lo and Ziegler

TABLE 1: Calculated Structures and Adsorption Energies (in Kilocalories per Mole) for Various Surface Intermediates at the Most Stable Sites on Fe(100) coverage

r(Fe-X)a

C (four-fold) CH (four-fold) CH2 (four-fold) CH3 (two-fold)

0.25 0.50 1.00 0.25 0.25 0.25 0.25

2.048 2.048 2.048 2.048 2.095 2.144 2.167

CH4 (one-fold)

0.25

3.021

CO (tilted)

0.25

1.970 2.250 1.970 2.250 1.972 2.250 2.317 2.144 2.186 2.101 2.017 2.315

configuration H (four-fold)

0.50 0.75

CO2 (four-fold) O (four-fold) OH (four-fold) H2O (two-fold)

1.00 0.25 0.25 0.25 0.25

r(X-H)b

1.107 1.139 1.096 1.112 1.094 1.100

Θc

110.75

Eads 8.56 8.39 8.57 183.00 163.47 100.43 41.83

101.24 102.70 108.43 110.02 110.94

120.63 0.980 0.982

r(C-O)

1.38 1.321

46.87

1.320

46.19

1.320

34.92

1.238 1.223 1.342

30.98 23.22 79.16 95.45 8.48

a

X stands for C in the cases of C, CO, CHn, and CO2, H in the case of H, and O in the cases of O, OH, and H2O. b X stands for C in the cases of C, CO, CHn, and CO2 and O in the cases of OH and H2O. c The angle of HCH except for CO2 and H2O.

In addition to the four-fold state, the two-fold state, of which the CH2 plane is perpendicular to the Fe-Fe bridge, was also determined to be a minimum that lies 14.3 kcal/mol higher in energy than the four-fold state. It is anticipated that this twofold state is favorable because C adopts the four-coordinated tetrahedral geometry. On the other hand, the one-fold state is a transition state according to the vibrational frequency calculations and has a binding energy of 72.7 kcal/mol. The adsorption of methyl on the Fe(100) surface was studied for the on-top, bridged, and hollow positions. Unlike the other CHn species, CH3 at the hollow site was unstable. The onefold and two-fold states were similar in stability. The two-fold state has the adsorption energy of 41.8 kcal/mol and is only 1.6 kcal/mol more stable than the one-fold state. The vibrational frequency calculations prove that the two-fold state is the local minimum. Due to the three-fold symmetry of the methyl radical, there are a number of possible rotational isomers of adsorbed CH3. Calculations indicate that these structures are stable and have essentially the same adsorption energies. Methane activation has been extensively investigated theoretically on the surfaces of nickel, platinum, and some heavy transition metals19,69-73 because of its relevance to the steam reforming process. It has been shown from these studies that methane is mainly physisorbed on the surface via the weak van der Waals interaction due to the strong C-H bond. The strength of adsorption of CH4 on metal surfaces does not exceed 5 kcal/ mol.71 The present work investigated the one-fold, two-fold, and four-fold states of adsorbed CH4. Since CH4 may coordinate to the metal surface via one, two, or three hydrogen atoms, there were, in total, nine possible configurations for CH4 adsorption. Among them, only six configurations were stable with respect to the free CH4 molecule, but the corresponding adsorption energies were all less than 2 kcal/mol. It was found that CH4 prefers to adsorb at the on-top position where two hydrogen atoms point to the surface, denoted by one-fold/two-fold state, and the associated binding energy was only 1.4 kcal/mol. It has to be emphasized that, due to the limitations of the DFT methods implemented in the VASP program, such a weak van der Waals

interaction cannot be accurately described. Therefore, the present value of the binding energy should be considered with great caution. The hydroxyl radical OH is an intermediate for the production of H2O and methanol during the FT synthesis and is formed from the hydrogenation of surface O atoms resulting from the CO dissociation. The calculations performed in the current study indicate that the configuration of OH adsorbed on the Fe-Fe bridge is the most stable, with the OH bond tilted by 69° relative to the surface normal. This state has a high binding energy of 95.5 kcal/mol, which agrees well with the value deduced by Hafner and co-workers.74 The strong adsorption of OH on the Fe(100) surface may result from the bonding interaction of O lone pairs and the vacant 3d orbitals of bridging Fe atoms. Contrary to the strong binding of OH on the Fe(100) surface, a H2O molecule was only weakly adsorbed at the two-fold site via the oxygen atom. The estimated adsorption energy for H2O on the Fe-Fe bridge was only 8.5 kcal/mol, comparable to the binding energy of the H atom. Stable configurations were also found at the one-fold and four-fold sites, but they both had similar adsorption energies to that of the two-fold states. Therefore, it is expected that a H2O molecule is mobile and diffuses rapidly on the Fe(100) surface. CO2 is a product of the water-gas shift reaction that occurs alongside of the CO hydrogenation in the FT synthesis (eqn 3). This side reaction is particularly catalyzed by Fe but not other FT-active catalysts, and it serves to manipulate the amount of H2 in the complicated FT reactions by converting H2O into H2. Choe, Park, and Huh have performed cluster calculations of the CO2 dissociation on the Fe(111) surface and noticed that the bent structure at the short-bridge position is most favorable.75 In the present work, however, the most stable configuration was found to be at the hollow site, where the OCO angle is 120.6° and the adsorption energy is 23.2 kcal/mol. Both the one-fold and two-fold states are unstable with respect to the free CO2 molecule. All of the structural parameters and energetics concerning the surface intermediates at their most preferred configurations are listed in Table 1.

Methanation on an Iron Surface

J. Phys. Chem. C, Vol. 111, No. 29, 2007 11017

Figure 1. Energy profiles of hydrogenation steps (a) C + H a CH, (b) CH + H a CH2, (c) CH2 + H a CH3. and (d) CH3 + H a CH4. Relative energies are given in kcal/mol.

3.1.4. Hydrogenation of Surface Intermediates. Note that CO and H2, after being adsorbed on the Fe(100) surface, undergo dissociation to yield C, O, and H, which generates different surface intermediates in the recombination reactions of these species. Having the most stable configurations of the surface intermediates determined, the paths of their formation were computed. Starting from surface carbide, four hydrogenation processes were investigated, CH, CH2, CH3, and CH4 formation, while for surface O atoms, the hydrogenation leading to, respectively, OH and H2O were studied. In addition, the oxidation of unreacted CO to give CO2 was also examined as it constitutes an important part of the FT synthesis. Figure 1a-d shows the calculated minimum energy paths for the consecutive hydrogenations of C that lead to CH4. In the first hydrogenation, C + H f CH, both C and H are assumed to be at adjacent hollow sites. During the hydrogenation, H moves over the two-fold site and forms CH, with the upright configuration at the four-fold site. The minimum energy path for this exothermic process has the activation barrier of 14.4 kcal/mol, which agrees with the barrier of 14.9 kcal/mol determined by Sorescu.16 The second hydrogenation (Figure 1b) corresponds to the addition of hydrogen to methylidyne to yield CH2. Before the C-H bond formation, both CH and H are at the adjacent fourfold sites. Along the reaction profile, the H atom migrates from a hollow site and forms a partial bond with CH, which is slightly tilted in the transition state. The H shift continues over the bridging site until CH2 is formed at the neighboring hollow site. This reaction possesses a barrier of 14.6 kcal/mol, which is very similar to the first hydrogenation. However, this process is endothermic, and the barrier for backward reaction is only 5.1 kcal/mol. The third hydrogenation step, as the second one, is an endothermic process, although to a lesser extent. This reaction requires CH2 and H at two adjacent four-fold sites; however, the HCH plane has to be oriented so that it is facing the H at the neighboring hollow position. The computed minimum energy profile shows that H first shifts toward CH2, forming a CH2-H complex at the transition state, followed by the shifting

of CH3 to the most stable two-fold state, as illustrated in Figure 1c. The higher activation barrier, 18.6 kcal/mol compared to the first two hydrogenation steps, is not unexpected because of the transition state involving the shift of both CH2 and H. In the final step of methanation, the methyl group at the bridging position interacts with H at the diagonal hollow site to form CH4 that lies on the one-fold site. The reaction path considered in the present work assumes that H is not at the hollow site next to CH3 because of the steric repulsion between the methyl hydrogen and surface hydride. The most stable configuration of CH3 and H coadsorbed on a p(2 × 2) unit is that CH3 sits at the two-fold site which is farthest from H, as depicted by Figure 1d. During the course of hydrogenation, CH3 shifts toward the on-top site and rotates to a proper orientation to form a pseudotetrehedral CH4 transition state which, afterward, leads to CH4 on the one-fold site when H moves closer to CH3. This reaction profile, unlike the second and third hydrogenation steps, is largely exothermic and has a barrier of 10.8 kcal/mol. The current study also investigated the hydrogenation of surface O atoms that resulted from the dissociation of CO, which leads to the formation of hydroxyl and H2O. Both O and H preferred the adsorption sites at the four-fold position, while OH favors the bridging site. Therefore, the minimum energy profile that connects these species was calculated. The computed activation barrier for the hydrogenation of O was 17.5 kcal/ mol, whereas the backward reaction possessed a similar barrier of 18.6 kcal/mol; the comparable activation barriers suggest that there may exist an equilibrium between O, H, and OH on the Fe(100) surface. This agrees with the observations from Bernasek that an ordered layer of OH on Fe(100) undergoes disproportionation to produce H2 and H2O at about 300K.76 The surface hydroxyl may be hydrogenated to generate H2O that is chemisorbed at the bridging site. In this process, H at a hollow position shifts toward OH, which lies on the adjacent Fe-Fe bridge and forms an O-H bond. It is assumed that the hydrogen of OH is pointing away from the surface hydride so as to minimize the repulsion. Employing these configurations, the barrier for OH + H f H was determined to be 18.0 kcal/

11018 J. Phys. Chem. C, Vol. 111, No. 29, 2007 mol, which is slightly lower than the first step of the hydrogenation of O. Furthermore, the second hydrogenation step was endothermic; the reverse process had a small barrier of 1.8 kcal/ mol. Consequently, the formation of H2O is essentially not favorable; high reaction temperatures are necessary to initiate this reaction. Finally, the oxidation of surface CO to yield CO2 was also examined. The initial state consisted of a tilted CO and O, which were simultaneously at the most stable four-fold configuration. The reaction proceeded via the migration of O over the twofold site and the formation of a CdO bond with the neighboring CO. This process had a very high barrier for both forward and backward reactions; the calculated barrier for CO2 formation was 26.3 kcal/mol, while that for the CO2 decomposition was 20.5 kcal/mol. Hence, as the formation of H2O, the production of CO2 is significant only at high reaction temperatures. 3.1.5. Magnetic Properties of Fe(100) with Adsorbed Species. It has been observed that the formation of the surface leads to the dramatic increase in the average magnetic moment of Fe.16 In the present calculations involving the five-slab model, the Fe(100) surface is found to be ferromagnetic, as in the case of bulk Fe, having the enhanced average magnetic moment of 2.71 µ0 compared to the bulk value of 2.30 µ0. Detailed analysis shows that the largest increase occurs on the top layer whose magnetic moment is 3.052 µ0, while for the two successive sublayers, the increments are much smaller (2.457 and 2.620 µ0, respectively). The adsorption of CHn species causes the demagnetization of the Fe(100) surface, which is possibly due to the charge transfer from the d band of surface Fe atoms to the adsorbates. For example, the adsorption of a C atom at the hollow site of the Fe(100) surface results in the decline of the average magnetic moment for the surface Fe atoms from 3.052 to 2.835 µ0. Interestingly, this demagnetization effect on the Fe atoms at both the top and first sublayer is strongly dependent upon the types and adsorption geometries of the adsorbates. For instance, it is noticed that the adsorption of CH2 at the four-fold site causes a greater demagnetization of the surface Fe atoms (2.807 µ0) than the C atom adsorbed at the same position does. On the other hand, the reduction of the magnetic moment is more pronounced on the Fe atoms that are directly bonded to the adsorbates. In the case of CH adsorption at the hollow site, the surface Fe atoms and the subsurface Fe atom beneath the hollow site undergo a substantial quenching of magnetization (-0.207 and -0.144 µ0, respectively). The other Fe atoms at the first sublayer, however, only experience a small demagnetization (-0.04 µ0). It has also been observed that the extent of demagnetization of the surface Fe atoms due to the adsorption is influenced by the surface coverage. For a free Fe(100) surface, the magnetic moment of the Fe atoms at the top layer is 3.052 µ0. The adsorption of a single CO molecule at a hollow site of a p(2 × 2) supercell, corresponding to the 0.25 ML of coverage, reduces the magnetization to 2.916 µ0. The adsorption of an additional CO molecule further lowers the magnetic moment to 2.734 µ0. At the 1.0 ML of coverage where all hollow sites are occupied by CO, the resulting magnetic moment is 2.353 µ0, which is close to the corresponding value for the bulk Fe. This negative dependence found in the present case of CO adsorption on the Fe(100) surface is consistent with the results deduced from the work of Hafner et al.,77 although an opposite trend has been reported for the oxygen adsorption on the same surface.67 These results imply that whether a metal surface is magnetized or demagnetized in the course of adsorption is intimately correlated

Lo and Ziegler to the nature of the adsorbate as well as the electronic structure of the surface. 3.2. The Kinetics of Methanation. 3.2.1. The Phenomenological Kinetic Model. The conventional rate equation kinetic model is based on the Langmuir-Hinshelwood approach, which assumes that all sites (in this case, the 2 × 2 unit cells of hollow sites) are energetically homogeneous and the surface species are randomly distributed. This phenomenological model, served as a simplification of the real surface of the catalyst, which may contain surface irregularities and defects, has been applied with success in several studies on the kinetics of hydrogenation reactions of CO and the polymerization on iron and cobaltbased-catalyzed Fischer-Tropsch synthesis.5,10 Nevertheless, recent studies have shown that islands of adsorbates may be formed on the surface of the catalyst during the course of the reactions of surface species; this apparent site preference would greatly alter the rates of reactions and thus the reaction mechanism.78 These calculations also found the fluctuating surface populations of different species with time, which is indicative of the competitive adsorption of gaseous reactants. The kinetic model developed in the present work is different from the aforementioned model, in a way that a p(2 × 2) supercell instead of a single hollow site is treated as the minimum unit. This approach affords an advantage that the lateral interaction of coadsorbed atoms and molecules could be explicitly taken into account since more than one species may be adsorbed on a “unit” at a time. In this case, the adsorption of a species is ultimately dependent upon the local environment around the adsorption site. The activation energy associated with such adsorption is then changed when different coadsorbed species are present. Accordingly, selective adsorption is introduced, and the assumption of random surface occupation can be partially removed. However, the macroscopic randomness of these (2 × 2) units on the Fe(100) surface is still retained unless the complete stochastic approach that describes every step of surface reactions (adsorption, dissociation, recombination, and desorption) on the Fe(100) surface by statistics is employed. The details of such a method can be found in the reference of Battaile and Srolovitz.79 The current model distinguishes the (2 × 2) units with different compositions and calculates the variation of their concentrations in time. The total population of a particular species, such as CO, is then given by summing over the concentrations of different types of units that contain *CO. The resulting concentrations correspond to the so-called “mean-field” concentrations derived from solving the phenomenological rate equations. While the mean-field concentration only provides the information concerning the average surface population of the species at various instances, the concentrations of individual units allow for the analysis of their time evolution on the surface and provide insight into the fluctuations of surface species. For example, the comparison of the concentrations of units containing *CO and (*CO + *H), respectively, helps acquire the details of the adsorption preference of CO in the presence of H and the ratio of CO populations with and without the neighboring H, thus yielding the proper description of the competitive adsorption of CO and H at the four-fold site. It is worth noting that, whereas the kinetic model employed in this study includes partially the site selectivity of adsorbates and the surface fluctuations due to the lateral interactions, which are missing in the conventional rate equation model, it is still an approximation that models the Fe surface as separated, noninteracting (2 × 2) pieces. The most accurate description of the dynamics of the methanation on the Fe(100) surface can

Methanation on an Iron Surface

J. Phys. Chem. C, Vol. 111, No. 29, 2007 11019

TABLE 2: Calculated Forward and Backward Activation Energies (in Kilocalories per Mole) of Elementary Steps in CH4 Formation on an Fe(100) Surfacea reaction step

Ef

Eb

COg a COs COg + COs a 2COs COg + 2COs a 3COs COg + 3COs a 4COs H2g a 2Hs H2g + 2Hs a 4Hs COg + 2Hs a COs + 2Hs H2g + COs a COs + 2Hs COs a Cs + Os 2COs a Cs + Os + COs 3COs a Cs + Os + 2COs Cs + Os + COs a 2Cs + 2Os COg + Cs + Os a Cs + Os + COs COg + Cs + Os + COs a Cs + Os + 2COs Cs + Hs a CHs CHs + Hs a CH2s CH2s + Hs a CH3s CH3s + Hs a CH4s Os + Hs a OHs OHs + Hs a H2Os COs + Os a CO2s H2Os a H2Og CO2s a CO2g

0.14 4.05 13.83 13.03 3.17 2.10 2.18 2.28 25.60 25.25 32.64 35.56 0.23 14.36 14.39 14.64 18.61 10.82 17.45 21.36 26.31 13.19 27.28

47.01 49.55 26.23 31.80 18.42 20.03 50.40 20.26 33.90 27.53 32.58 25.26 38.09 24.42 18.78 5.07 14.51 24.40 18.55 5.14 20.49 4.71 4.06

a Subscript s stands for surface species, and g stands for gaseous species.

be achieved only when the whole Fe surface is considered as a single unit and all possible surface reactions including neighboring effects are simultaneously taken into account; however, this approach is practically impossible (unless the statistical treatment is used). On the basis of the conditions that only the hollow sites be present in the (2 × 2) units and no adsorption of the gaseous species on the on-top or bridge sites be allowed, the present model should be able to provide reliable data regarding the kinetics of CH4 formation on the Fe(100) surface. 3.2.2. The Adsorption and Dissociation of CO on Fe(100). The computed activation barriers for the elementary steps involved in the methanation during the FT synthesis are listed in Table 2. By applying these values and eq 5, rate constants for the reaction steps in Table 2 can be obtained for various reaction temperatures. The rate equations describing the elementary reactions are then formulated and numerically solved utilizing the computed rate constants. First of all, the surface coverage of CO and the temperature dependence of its dissociation were investigated. The simulations included the reaction steps where CO participates, that is, the reactions 1-4 and 9-14 in Table 2. It was noticed that adsorption of CO is thermodynamically favorable regardless of the surface coverage. It was experimentally observed, however, only at low temperatures and high CO exposure.59 At temperatures higher than 220 K, CO started to desorb from the surface. A further increase in temperature triggered the dissociation of CO to produce surface carbide and oxide. To understand the temperature effects on the adsorption mechanism of CO, the entropy of CO calculated using GGABP86/TZP method80,81 in the ADF program82-84 was included in ∆Gact in eq 5 and calculations were performed at several temperatures between 150 and 473 K. The partial pressure of 0.05 atm for CO was assumed. As seen in Figure 2, the adsorption of CO on the Fe(100) surface demonstrated a significant temperature dependence. At very low temperatures (Figure 2a), CO was predominantly adsorbed molecularly. Due to the low CO exposure, the surface was only half covered

(corresponding to 100% coverage by [2*CO]). No significant dissociation of CO was noticed, as observed experimentally by Bernasek and co-workers,59 because of the high dissociation barrier and insufficient thermal excitation. Surface coverage increased with temperatures; at 300 K (Figure 2b), the percentages of surface coverage due to [3*CO] and [4*CO], which correspond to 75 and 100% coverage in a p(2 × 2) unit cell, respectively, rose rapidly, while that of [2*CO] dropped exponentially. Meanwhile, the dissociation of CO became competitive with CO absorption, as reflected by the dominance of [*CO + *C + *O] at the end of simulation. The preference of CO dissociation was more manifested at higher temperatures. Figure 2d reveals that [*CO + *C + *O] prevailed on the surface, with a small coverage by [2*CO] at 473 K. Higher reaction temperatures activated the CO dissociation, which possesses a high reaction barrier, and favored the formation of *C and *O. Simultaneously, the adsorption of CO was suppressed at high temperatures due to the large entropic contribution of CO that made the adsorption unfavorable on the free-energy surface, resulting in the diminishing CO coverage on the Fe(100) surface, as shown in Figure 2a-c. It is worth noting that the equilibrium of CO adsorption and dissociation was more rapidly attained at higher temperatures; the maximum coverage by [*CO + *C + *O] was achieved within 2.8 h at 473 K, while it took 100 times longer when temperature was dropped to 300 K. In addition to temperature, the reactions of CO adsorption and dissociation on Fe(100) also displayed pressure dependence. Higher surface coverage by CO was expected when CO gas of higher pressure was utilized because of the higher frequency of collisions. Similar results were observed in the present study. Figure 3a-d shows the variation of surface coverage by different species as a function of the partial pressure of CO. The simulations were carried out for 106s (∼278 h) at 300 K. Two major changes are noticed from Figure 3, the increase in [3*CO] and [4*CO] and the drop of [*CO + *C + *O]. Note that these three species indicate the total surface coverage of Fe(100) surface. A sharp growth of [4*CO] with increasing CO pressure, accompanied by the decline of [*CO + *C + *O], illustrates that a larger area of Fe(100) is fully covered by adsorbed CO at high CO pressure. The high CO pressure also leads to the higher coverage by [2* CO + *C + *O], which results from the adsorption of CO to [*CO + *C + *O]. According to Figure 3a,d, it is deduced that the surface coverage changes from 72.5% for pCO ) 0.05 atm to 87.5% for pCO ) 5.00 atm. The calculations of the compositions of surface species reveal the dramatic increase of CO coverage from 27.5 to 60.5% for the CO pressures of 0.05 and 5.00 atm, respectively. At high CO pressure, it was also found that CO dissociation is suppressed; the coverage by C and O on the surface is reduced from 22.5 to 13.5%. This phenomenon has also been inferred by the kinetic studies of iron-catalyzed FT synthesis5,85 and is possibly due to the competition with CO adsorption for the vacant sites necessary for the CO dissociation. 3.2.3. The Formation of Carbon Filaments on Fe(100) in the Absence of H2. It has been well-known that Fe is an active catalyst for the water-gas shift reaction (eq 3), which enables the use of low H2/CO syngas in the FT synthesis by converting CO and H2O produced in the course of the reaction into CO2 and H2. In the absence of H2, Fe instead activates the Boudouard reaction in which CO is disproportionated into CO2 and C, resulting in the growth of a carbon filament86,87 on the Fe surface. Nevertheless, this reaction is extremely slow at low temperatures and low CO concentration,88 and it is thus

11020 J. Phys. Chem. C, Vol. 111, No. 29, 2007

Lo and Ziegler

Figure 2. Surface coverage of CO on Fe(100) at different temperatures. (a) 150; (b) 300; (c) 400; and (d) 473 K. Only configurations with a significant temperature dependence are shown here, 2CO, 3CO, 4CO, CO + C + O.

considered not significant in the FT process operating at low temperatures. Experiments concerning CO adsorption also provide evidence that only a modest amount of CO dissociation occurs below 300 K.59 To unveil the importance of the Boudouard reaction in the FT synthesis, reaction steps 21 and 23 in Table 2 were added into the reaction scheme of CO adsorption and dissociation. Calculations were performed at several temperatures ranging from 150 to 500 K, which is at the higher end of the operating temperatures of industrial FT reactors. The results are plotted in Figures 4 and 5, respectively. As shown in Figure 4, the coverage of Fe(100) by carbide is extremely low for temperatures below 275 K, regardless of the applied partial pressure of CO. The low carbide content is due partly to the suppressed CO dissociation and partly to the slow Boudouard reaction. In the low-temperature regime, the most abundant surface species is adsorbed CO, and only a negligible amount of CO2 is produced. As the temperature rises, the surface

coverage by carbide increases because of the enhanced CO dissociation and CO2 formation, both of which have rather high activation energies. For a temperature below 400 K, the extent of carbide deposition is apparently independent of the partial pressure of CO, as illustrated by Figure 5. For instance, only a small increase in the amount of surface carbide is observed at 400 K when the CO partial pressure is increased from 0.05 to 5.00 atm. Only above 400 K does a significant dependence of carbide formation on CO pressure become apparent; at these temperatures, the dissociation of CO, which leads to the deposition of C on the surface, turns to be more favored compared to the competing CO desorption, which is suppressed by the increased CO pressure. About 50% of available sites are occupied by C when 5.00 atm of CO is utilized at 450 K, and the coverage increases to 60% at 500 K. The participation of the Boudouard reaction in carbon deposition can be evaluated by monitoring the surface species

Methanation on an Iron Surface

J. Phys. Chem. C, Vol. 111, No. 29, 2007 11021

Figure 3. Surface coverage of CO on Fe(100) at 300 K with different partial pressures. (a) 0.05; (b) 0.5; (c) 1.0; and (d) 5.0 atm. Only four configurations are shown here, 3CO, 4CO, CO + C + O, and 2CO + C + O.

such as [*C], [2*C], [*CO + *C], [2*CO + *C], and [3*CO + *C], which are formed only during the recombination of *CO and *O

[*CO + *C + *O] f [*CO2 + *C] f [*C]

(9)

[*CO + 2*C + *O] f [*CO2 + 2*C] f [2*C] (10) It is noticed that [*CO + *C + *O] remains the dominant species at low temperatures (