J. Phys. Chem. C 2009, 113, 7723–7727
7723
Study of Carbon Monoxide Oxidation on CeO2(111) Using Ultra Accelerated Quantum Chemical Molecular Dynamics Md. Khorshed Alam,† Farouq Ahmed,† Katsuyoshi Nakamura,† Ai Suzuki,‡ Riadh Sahnoun,§ Hideyuki Tsuboi,† Michihisa Koyama,§ Nozomu Hatakeyama,† Akira Endou,† Hiromitsu Takaba,§ Carlos A. Del Carpio,† Momoji Kubo,| and Akira Miyamoto*,†,‡,§ Departments of Applied Chemistry and of Chemical Engineering and Fracture and Reliability Research Institute, Graduate School of Engineering, Tohoku UniVersity, 6-6-11-1302 Aoba, and New Industry Creation Hatchery Center, Tohoku UniVersity, 6-6-10 Aoba, Sendai 980-8579, Japan ReceiVed: October 8, 2008; ReVised Manuscript ReceiVed: March 6, 2009
Ceria has attracted intensive interest in the past decade because of its vital role in emerging technologies for exhaust gas purification in automobiles. In this study, we have investigated the process of conversion of CO to CO2 via the creation of an oxygen vacancy on the ceria surface. The study was conducted using a new ultra accelerated quantum chemical molecular dynamics method. Through this simulation, we have demonstrated that a high-energy colliding CO molecule is adsorbed on the ceria, pulling up an O atom from the ceria surface to form a CO2 molecule. This molecular dynamics simulation of CO oxidation has been carried out for the first time using quantum chemical methods. Introduction The dramatic increase of fossil fueled vehicles around the world and the consequent adoption of stringent standards on gas emissions have led to extensive research activities in this area. Automobile companies are striving to purify automobile exhaust emissions. In this regard, particularly the improvement of catalytic performance and the design of new three-way catalysts (TWCs) are currently topics of intense research. Ceria (CeO2) has emerged as one of the most important promoters of heterogeneous catalytic reactions, and it is used as a component of modern emissions-control catalysts in automobiles. Moreover, it represents the most important application of rare earth oxides to the catalytic process.1 In automotive control systems, ceria is a critical component of the so-called TWCs that simultaneously convert hydrocarbons, CO, and NOx present in exhaust gases into H2O, CO2, and N2. One of the important roles that ceria plays in TWCs is to act as an oxygen buffer by absorbing and releasing oxygen, maintaining the partial pressure of oxygen for optimum catalytic activity. Another important role that ceria plays in TWCs is to promote the conversion of CO to CO2 through oxidation involving a lattice oxygen. The oxidation of CO as well as reduction of NO2 can be catalyzed with ceria, without the need for precious metals.2 However, there is an intrinsic difficultly as to how to determine the geometric and electronic changes of atoms in the system because of its complexity. Yang et al. used density functional theory (DFT) calculations to determine the CO adsorption properties on ceria (111) and (110) surfaces.3 They showed that there is only a weak adsorption of CO on the CeO2(111) surface, while in the case of the (110) surface weak and few strong adsorption sites were observed. Huang et al. investigated the CO adsorption and oxidation on the most stable (111) and (110) * To whom correspondence should be addressed. Fax: +81-22-795-7235. E-mail:
[email protected]. † Department of Applied Chemistry. ‡ New Industry Creation Hatchery Center. § Department of Chemical Engineering. | Fracture and Reliability Research Institute.
ceria surfaces by means of DFT+U calculations.4 They demonstrated that CO undergoes physisorption on the (111) surface and chemisorption on the (110) surface of ceria. Di Monte et al. found that ceria was reduced in three-way catalysts by CO that is contained in automotive exhaust gases.5 They observed that CO oxidation proceeds through a water-gas shift reaction (WGSR). Sayle et al. investigated the role of oxygen vacancy formation of carbon monoxide.6 The formation of lattice oxygen vacancy in ceria, which plays an important part in the oxidation of CO,7,8 is associated with the reduction of Ce(IV) to Ce(III).9 Nolan et al. published a DFT+U study on the interaction between CO and the (111), (110), and (100) surfaces of ceria that is the first ab initio study in which the stoichiometric and reduced ceria surfaces involved in this important catalytic reaction are consistently described.10 Widmann et al. studied the dynamics of the CO oxidation reaction on Au/CeO2 catalyst by applying a temporal analysis of products (TAP) technique.11 Within this framework in the present study we applied our in-house code for ultra accelerated quantum chemical molecular dynamics (UA-QCMD) to investigate the oxidation mechanism of CO to CO2, creating an oxygen vacancy in the ceria surface, and also to visualize the electronic and structural change undergone in this system. Computational Methods Our novel UA-QCMD simulator consists of two parts. The first part is an original tight-binding (TB) simulator.12-14 In this simulator, an electronic structure calculation is performed by solving the Schro¨dinger equation HC ) εSC (H, C, ε, and S refer to the Hamiltonian matrix, eigenvector, eigenvalue, and overlap integral matrix, respectively) with the diagonalization condition CTSC ) I (I refers to the unit matrix). To determine the off-diagonal elements of H, Hrs, the following corrected distance-dependent Wolfsberg-Helmholz formula15 was used:
K Hrs ) Srs(Hrr + Hss) 2
(1)
To solve the Schro¨dinger equation in this simulator, parameters for H are used. In electronic structure calculations using the
10.1021/jp8088963 CCC: $40.75 2009 American Chemical Society Published on Web 04/08/2009
7724
J. Phys. Chem. C, Vol. 113, No. 18, 2009
Alam et al.
TABLE 1: Determined Coefficients of a Single ζ Parameter in a Slater-Type Atomic Orbital element
AO
a0
a1
a2
a3
a4
a5
s p s p s p d f
2.6000 2.6009 2.1450 1.9664 2.7165 2.6592 3.1519 4.5582
0.1045 0.1440 0.0802 0.5181 0.0332 0.0010 0.3368 0.3363
-0.0056 -0.0610 0.0185 1.0817 0.0208 -0.0127 0.0000 0.0000
-0.0056 0.0322 0.0007 1.1475 0.0041 0.0093 0.0000 0.0000
-0.0012 -0.0012 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
0.0000 -0.0005 0.0000 0.0000 0.0000 0.0000 0.0000 0.0000
C O Ce
TABLE 2: Determined Coefficients of the Hrr Parameter element
AO
b0
b1
b2
b3
b4
b5
s p s p s p d f
-10.6984 -8.81270 -28.6066 -14.2851 -5.53030 -0.06120 -8.43560 -7.19290
-11.8180 -11.7230 -15.0328 -14.8159 -5.07280 -4.62470 -4.86090 -8.16630
-1.4167 -1.6433 -1.4510 -1.6244 3.4414 3.8922 7.5711 3.1559
0.0161 0.0032 -0.4949 0.0154 -4.3070 -4.4481 -9.4497 -1.5706
0.0477 0.0885 0.0419 0.0434 1.5907 1.6406 3.5251 0.0885
-0.0042 0.0110 0.0000 0.0000 -0.1994 -0.2082 -0.4442 0.0000
C O Ce
TABLE 3: Calculated Binding Energy of CeO2 (Bulk Model) binding energy (eV) system CeO2
DFT
TB
-179.7
-172.2
TB simulator, the total energy of a system is obtained using the following equation: occ
E)
∑
N N ZiZje2 + Eijrepul(rij) r ij i)1 j)i+1 i)1 j)i+1 N
nkεk +
k)1
N
∑∑
∑∑
(2)
where the first, second, and third terms on the right-hand side refer to the molecular orbital (MO) energy, Coulombic energy, and exchange-repulsion energy, respectively. The first term on the right-hand side of eq 2 is rewritten as follows:
functions of atomic charges. ζr and Hrr are calculated by the polynomial functions of atomic charges described by eqs 5 and 6, respectively. In eqs 5 and 6, Zi corresponds to the atomic charge on atom i. The parameters regarding ζ, i.e., a0, a1, a2, a3, a4, and a5 in eq 5, and regarding Hrr, i.e., b0, b1, b2, b3, b4, and b5 in eq 6, were adjusted to reproduce the binding energies and electronic structures of each reactant obtained by firstprinciples calculations; these are summarized in Tables 1 and 2, respectively. 5
occ
occ
occ
∑ nkεk ) ∑ ∑ nk(Ckr)2Hrr + ∑ ∑ ∑ nkCkrCksHrs
k)1
k)1 r
k)1 r
(3)
ζr ) a0 +
∑ ak(Zi)k
s
5
where the first and second terms on the right-hand side refer to the monoatomic contribution to the binding energy and the diatomic contribution to the binding energy, respectively (nk is the number of electrons occupied in the kth molecular orbital). The binding energy calculated from the second term of eq 3 is used for the determination of the DAB parameter in the Morsetype two-body interatomic potential function described in the following equation:
Hrr ) b0 +
∑ bk(Zi)k
(6)
k)1
The second part of our methodology involves the classical MD simulator. Details of this method are described elsewhere.17,18 To renew the Morse-type function and the atomic charges, TB calculations are performed at the beginning and at successive intervals during the MD simulation process. Results and Discussion
EAB ) DAB{exp[-2βAB(RAB - R*AB)] 2 exp[-βAB(RAB - R*AB)]}
(5)
k)1
(4)
where EAB, DAB, βAB, RAB, and R*AB are the interatomic potential energy between atoms A and B, the binding energy between atoms A and B, the factor for the potential curve, the interatomic distance between atoms A and B, and the equilibrium interatomic distance between atoms A and B, respectively. To set H and S in our TB simulator, the exponent of a Slatertype atomic orbital (AO),16 denoted as ζr, and valence-state ionization potentials (VSIPs) for the valence shell of an AO of Ce, C, and O atoms are necessary. The former parameters are used to calculate the S matrix and Hrs in eq 1. The latter ones are used for the diagonal element of H (Hrr or Hss in eq 1). In our TB simulator, these are represented by the polynomial
Before undertaking the study of the interaction of CO on CeO2(111), we performed a preliminary calculation on the CeO2 bulk model to validate our methodology. To do that, we calculated the binding energy of the bulk model using the DFT module in the Dmol3 program as well as our TB method. The energy calculations by DFT were performed at the density approximation level using the VWN function19 and basis function of double numerical polarization. The calculated binding energies in the bulk ceria by UA-QCMD and DFT are shown in Table 3. The accuracy of the used parameters can be further evidenced by comparing the bond populations of Ce-O by UA-QCMD and by DFT, as shown in Table 4. Comparison of the results from the two different methods validates our UAQCMD methodology.
Carbon Monoxide Oxidation on CeO2(111)
J. Phys. Chem. C, Vol. 113, No. 18, 2009 7725
Figure 2. Change in the bond energies and bond lengths between the C atom of CO and O(a) and O(b) in ceria observed during the simulation.
Figure 1. (I) Initial structure of the CeO2-CO system whose lattice parameters are a ) 6.6271 Å, b ) 7.6523 Å, and c ) 20.00 Å. Here the black spheres represent oxygen atoms, the white spheres represent Ce, and the gray sphere represents C. The bottom layers are fixed, while the upper layers are relaxed. (II) Snapshot at 370 fs where in the dynamics simulation CO comes close to the ceria surface. Due to thermal dynamics the positions of the upper layer atoms are shifted. (III) Snapshot at 3140 fs where CO interacts with the surface. The topmost O-Ce-O trilayer in this snapshot is slightly shifted as in Figure 1(II). (IV) Final structure (4270 fs). Here CO removes an oxygen atom and forms CO2, creating an oxygen vacancy in the ceria surface.
Figure 3. Time profile of atomic charge analysis of each cerium atom of the first layer in the surface, of the C atom in CO during the simulation, and of each oxygen atom of the first layer in the surface.
TABLE 4: Calculated Bond Population of CeO2 (Bulk Model) bond population Ce-O
DFT
TB
0.22
0.19 Figure 4. Change of total atomic charges of Ce atoms (first layer) in the surface and the O(a) atom during the simulation.
The (111) ceria surface was modeled with supercell slabs consisting of six layers.20 The theoretical equilibrium lattice parameter used for the thickness of the vacuum separating the slabs is 11.82 Å. As the initial structure, a single CO molecule was placed on the free surface above the Ce and O atoms, where the carbon atom in the CO faced the surface as shown in Figure 1(I). The initial distance between CO and the surface was 6.17 Å. In the initial configuration, lower layers are fixed and the other layers are relaxed. Simulations were performed under a temperature of 873 K for a model constituted of 26 atoms in the unit cell. To study the dynamics of the CO molecule on the CeO2(111) surface, we performed a single-point TB calculation, and the atomic charges and Morse potentials were updated at intervals during the MD calculation. The MD calculations were carried out for 4270 steps with a 1.0 fs time increment. The initial velocity of CO was set to 1000 m/s, and the movement direction was toward the ceria surface. After dynamics calcula-
tion we again transferred this model for a single-point TB calculation and resumed the MD simulations. Four snapshots taken from the simulation at different time steps are shown in Figure 1. In this calculation, we focused on three specific atoms, O(a), O(b), and C as shown in Figure 1 to discuss the bond-breaking, bond formation, electron transfer, and structural change at the CeO2-CO interface during the simulation. During the simulation bond lengths, bond energies, and structural change were observed. The mean velocity of CO at 873 K is about 1000 m/s. According to the MD simulation at first, the CO molecule comes close to the surface as shown in Figure 1(II), and at 3140 fs, the CO is adsorbed by the O(b) atom on the surface as shown in Figure 1(III). During the simulation, the Ce-O(b) bond in ceria was completely broken in 4270 steps as shown in Figure 1(IV). On the other hand, a new bond was formed between the liberated O and the C atom
7726
J. Phys. Chem. C, Vol. 113, No. 18, 2009
Alam et al.
TABLE 5: Calculated Bond Lengths and Bond Energies of the C Atom to Two Oxygen Atoms [O(a) and O(b)] in CO2 C-O(a) time (fs) 0 1370 1870 2940 4270
C-O(b)
bond length (Å)
bond energy (eV)
bond length (Å)
bond energy (eV)
1.14 1.16 1.14 1.17 1.18
-8.53 -7.92 -8.25 -7.82 -7.81
6.17 1.62 1.41 1.19 1.17
0 -1.66 -3.75 -7.37 -7.92
of CO as shown in Figure 1(IV). The formation of the CO2 molecule created an oxygen vacancy in the ceria surface and occurred at the 4270th step. The calculated bond lengths and bond energies for CO2 are 1.17 and 1.18 Å and -7.91 and -7.82 eV, respectively. Figure 2 shows the change of the bond distance between C-O(a) in CO and C-O(b) with the simulation time. The initial bond distance of CO was 1.14 Å, and the C-O(b) distance in ceria was 6.17 Å. During the simulation, the C-O(a) bond distance remained almost constant around 1.1 Å and the C-O(b) bond distance decreased from 6.17 to 1.17 Å. The C-O(a) bond distance undergoes a slight change while C-O(b) decreases with time. Figure 2 also shows the change of the bond energies between C and O(a) and O(b). It was found that the bond energy between C and O(a) fluctuated from -8.55 to -7.80 eV, while the bond energy between C and O(b) was zero at 0 fs; the energy increased with simulation. The bond energy between C and O(a) changes only marginally but that between C and O(b) increased with time; after 4270 fs, both bond energies were about -7.80 eV. Figure 3 shows that the charge of the C atom of CO gradually changes with the simulation time. During the dynamics simulation when CO starts to interact with the surface oxygen, the charge of the C atom also starts to increase. At the transition point (2940 fs), the highest electropositivity value was observed for C. At the same time, the electronegativity of O(b) in the ceria decreased drastically compared to that of other oxygen atoms in the ceria surface. In Figure 3, we also analyzed the change in the atomic charges with time. This figure shows the time profile of the atomic charges of the individual atoms of Ce and O of the surface layer. The figure indicates that the charges of Ce(a), Ce(b), and Ce(d) suddenly decreased at 2940 fs. On the other hand, the atomic charge on the C atom and O(b) atom gradually increased from the beginning of the simulation. Although the charge difference is not significant, it shows that the Ce atoms were reduced from Ce4+ to Ce3+. The reduction of Ce is apparent when the temporal change of the charges of Ce and O atoms on the surface is plotted, as shown in Figure 4. It can be observed that the total charge of Ce atoms in the first layer gradually decreases with time. To discuss the electron transfer between atoms during the CO oxidation process, we analyzed the change in the atomic charges of the ceria with time. Figure 4 shows the time profile of atomic charges on Ce (total charge of Ce atoms in the first layer) and O(a) during the oxidation process. This figure indicates that the atomic charge on O(a) increased gradually from the initial stage of the simulation. On the other hand, the atomic charge on the Ce atoms gradually decreased. This electron transfer was found when the Ce-O(b) bond broke by reacting with CO. This implies that some Ce atoms were reduced from Ce4+ to Ce3+ and thus the electron transfer took place from O(b) to Ce. Further detailed analysis showed that the electrons were transferred from the p-orbital of the O atom to the f-orbital
of the Ce atom. This reduction reaction is related to the specific characteristics of the Ce element, which has two oxidation states, Ce3+ and Ce4+.21 Oxygen vacancies play a central role in the defect chemistry of CeO2. We have therefore estimated the oxygen vacancy formation energy on the surface of ceria. The calculated oxygen vacancy formation energy is -3.87 eV, which is in good agreement with experimental results.6,22 During the molecular dynamics simulation, CO comes close to the surface and is adsorbed on the surface as shown in Figure 1(II,III). In these two snapshots the topmost O-Ce-O trilayer seems to have been shifted during the simulation; however, this is an instantaneous effect stemming from the vibrations and change of position of the atoms during the MD simulation. We also analyzed the change in the internal bond energies and bond lengths of CO and C-O(b) with time. Table 5 shows variations of the bond lengths and bond energies. The intermolecular bond length increased from 1.14 Å in CO in a vacuum to 1.18 Å,3 and the variation also happened in the case of the bond energy from -8.53 to -7.81 eV. On the other hand, the C-O(b) bond length significantly changed from 6.17 to 1.19 Å when CO was strongly interacting with the surface at 3140 fs as shown in Figure 1(III). Finally, the adsorbed CO molecule pulled up one oxygen atom from the surface, leading to the formation of a CO2 molecule. The bond lengths and the bond energies between C and O(a) and O(b) also indicate the formation of the CO2 molecule. Conclusion We have analyzed the oxidation of CO on the ceria (111) surface and unveiled the factors involved in this chemical transformation. We have found that the reaction follows the “Eleay-Rideal (ER)” mechanisms and observed a change of the electronic state as well as a structural change of the system. Finally, we have confirmed that our new tight-binding quantum chemical molecular dynamics method is very effective in investigating the fundamental aspects of TWCs. References and Notes (1) Kasˇper, J.; Fornasiero, P.; Graziani, M. Catal. Today 1999, 50, 285–298. (2) Rodriguez, J. A.; Jirsak, T.; Sambasivan, S.; Fischer, D.; Maiti, A. J. Chem. Phys. 2000, 112, 9929–9939. (3) Yang, Z.; Woo, T. K; Hermansson, K. Chem. Phys. Lett. 2004, 396, 384–392. (4) Nolan, M.; Watson, G. W. J. Phys. Chem. B 2008, 110, 16600– 16606. (5) Di Monte, R.; Kasˇpar, J. Top. Catal. 2004, 28, 47–57. (6) Sayle, T. X. T.; Parker, S. C.; Catlow, C. R. A. Surf. Sci. 1994, 316, 329–336. (7) Jin, T.; Okuhara, T.; Mains, G. J.; White, J. M. J. Phys. Chem. 1987, 91, 3310–3315. (8) Sanchez, M. G.; Gazquez, J. L. J. Catal. 1987, 104, 120. (9) Yao, H. C.; Yu Yao, Y. F. J. Catal. 1984, 86, 254–265. (10) Nolan, M.; Watson, G. W. J. Phys. Chem. B 2006, 110, 16600– 16606. (11) Widmann, D.; Leppelt, R.; Behm, R. J. J. Catal. 2007, 251, 437– 442.
Carbon Monoxide Oxidation on CeO2(111) (12) Elanany, M.; Sasata, K.; Yokosuka, T.; Takami, S.; Kubo, M.; Miyamoto, A.; Imamura, A. Stud. Surf. Sci. Catal. 2002, 142, 1867–1876. (13) Sasata, K.; Yokosuka, T.; Kurokawa, H.; Takami, S.; Kubo, M.; Imamura, A.; Shinmura, T.; Kanoh, M.; Selvam, P.; Miyamoto, A. Jpn. J. Appl. Phys. 2003, 43, 1859–1864. (14) Elanany, M.; Selvam, P.; Yokosuka, T.; Takami, S.; Kubo, M.; Imamura, A.; Miyamoto, A. J. Phys. Chem. B 2003, 107, 1518–1524. (15) Carzaferri, G.; Forss, L.; Kamber, I. J. Phys. Chem. 1989, 93, 5366– 5371. (16) Kotani, M.; Amemiya, A; Ishiguro, E.; Kimura, T. Tables of Molecular Orbitals; Maruzen: Tokyo, 1963. (17) Takaba, H.; Hayashi, S.; Zhong, H.; Malani, H.; Suzuki, A.; Sahnoun, R.; Koyama, M.; Tsuboi, H.; Hatakeyama, N.; Endou, A.; Kubo, M.; Del Carpio, C. A.; Miyamoto, A. Appl. Surf. Sci. 2007, 254, 7955– 7958.
J. Phys. Chem. C, Vol. 113, No. 18, 2009 7727 (18) Malani, H.; Hayashi, S.; Zhong, H.; Sahnoun, R.; Tsuboi, H.; Koyama, M.; Hatakeyama, N.; Endou, A.; Kubo, M.; Del Carpio, C. A.; Miyamoto, A. Appl. Surf. Sci. 2007, 254, 7608–7611. (19) Vosko, S. H.; Wilk, L.; Nusair, M. Can. J. Phys. 1980, 58, 1200– 1211. (20) Huang, M.; Fabris, S. J. Phys. Chem. C 2008, 112, 86438648. (21) Rajendran, A.; Takahashi, Y.; Koyama, M.; M.; Kubo, M.; Miyamoto, A. Appl. Surf. Sci. 2005, 244, 24–38. (22) Chen, Y.; Hu, P; Lee, M. H.; Wang, H. Surf. Sci. 2008, 602, 1736–1741.
JP8088963
