2188
Ind. Eng. Chem. Res. 1998, 37, 2188-2192
Monte Carlo and Cellular Automata Modeling of CO Oxidation on a Catalytic Surface Including the Eley-Rideal Step and CO Diffusion P. Sreekumar, V. K. Jayaraman, and B. D. Kulkarni* Chemical Engineering Division, National Chemical Laboratory, Pune 411008, India
We present the results of cellular automaton and Monte Carlo simulations for the case of CO oxidation on a catalytic surface with simultaneous adsorption, reaction, and diffusion, including the Eley-Rideal step in the reaction mechanism. The results show that the carbon monoxide diffusion on the catalytic surface does not alter the qualitative nature of phase transitions. Also for any given adsorption to diffusion ratio the cellular automaton predicts slightly lower values of the phase transition point. The automaton has considerable speed advantage and corroborates the Monte Carlo simulation results. 1. Introduction Reactions occurring on a catalyst surface exhibit a rich variety of behavioral patterns. As a consequence, experimental and theoretical efforts to characterize surface catalyzed reactions have been the subject matter of many publications. Models of differing levels of rigor have been proposed to predict the qualitative and quantitative behavior of surface catalyzed reactions (Beusch et al., 1972; Hugo and Jakubith, 1972; McCarthy et al., 1975; Sheintuch and Schmitz, 1977; Lagos et al., 1979). Interest in microscopic modeling of surface reactions gathered momentum with the simple model for oxidation of CO over a platinum surface proposed by Ziff et al. (1986), who considered the catalytic surface as a twodimensional square lattice with periodic boundary conditions. They assumed that a CO molecule from the gas phase adsorbs onto a vacant site, while oxygen gets dissociatively adsorbed onto two neighboring vacant sites. The adsorbed CO and oxygen react and instantaneously desorb, vacating the two sites. The Monte Carlo simulations, based on this simplified picture of the surface, showed two kinetic phase transitions. For low values of gas phase partial pressure of CO, the surface is completely covered with oxygen (yCO < y1), while for higher values of the partial pressure, yCO, the surface was found to be completely covered by carbon monoxide. Thus, a reactive phase exists only in a narrow window of CO partial pressures. The phase transition at y1 was first order, whereas the second transition at y2 was second order. These interesting results opened up a series of detailed microscopic studies. This includes studies on monomer-dimer, dimer-dimer, and trimer-trimer models and various simulation procedures (Fichthorn et al., 1987; Meakin and Scalapino, 1987; Chopard and Droz, 1988; Kohler and ben-Avraham, 1991; Ziff et al., 1991; Ziff and Brosilow, 1992; Pan et al., 1995). An excellent review of application of the Monte Carlo simulation method for the study of reaction processes has appeared recently (Albano, 1996). In the past, in addition to mean field analysis, detailed Monte Carlo simulations were performed to get deeper insight into the underlying mechanisms of the * To whom correspondence should be addressed.
surface catalyzed reactions. While the Monte Carlo approach is more rigorous than the mean field approach, it deals with the sites one at a time and the updating is done sequentially. The recently developed cellular automaton (CA) modeling approach provides an interesting and computationally cost and time saving alternative to Monte Carlo simulations. CA was originally introduced by Von-Neumann (1966), as a possible idealization of biological systems with the specific purpose of modeling biological reproduction. Since then, it has been applied and reintroduced for a wide variety of purposes (Hartmann and Tamayo, 1990; Sahimi and Stauffer, 1991; Zygourakis et al., 1991a,b). Cellular automaton is discrete in space, time, and state. The state of any particular site in the lattice is determined by the state of its neighborhood. The automaton evolves in discrete time steps subject to a set of transition rules. These rules have to be formulated by careful analysis of the system dynamics. The updating is synchronous, and this distinct feature of the cellular automaton allows parallel processing. This results in substantial reduction in the computational costs. Due to the fast and efficient parallel environment in CA we can afford to conduct the simulations in larger matrices. Also the large number of unsuccessful attempts inherent in the Monte Carlo approach are not there in CA simulations. This distinct feature enhances the computational speed of CA in any computing environment considerably. A number of studies using mean field, Monte Carlo, and cellular automata have been made in the recent past to understand the surface catalyzed reactions using oxidation of CO on platinum as an example (Carberry, 1976; Hlavacek and van Rompay, 1981; Engel and Ertl, 1982; Ertl, 1983; Moller et al., 1986; Chopard and Droz, 1989; Levermore and Boghosian, 1989; Mai et al., 1990; Mai and von Niessen, 1991, 1992, 1993; Tambe et al., 1994). A comprehensive treatment taking account of all three steps, viz., adsorption, reaction, and diffusion, together and incorporating the Eley-Rideal step, however, has not been attempted so far. The reaction mechanism of CO oxidation has been investigated thoroughly, and various reaction schemes have been proposed (Carberry, 1976: Hlavacek and van Rompay, 1981). The Eiley-Rideal step has in the past been included by various authors, and its importance
S0888-5885(97)00526-5 CCC: $15.00 © 1998 American Chemical Society Published on Web 04/03/1998
Ind. Eng. Chem. Res., Vol. 37, No. 6, 1998 2189
Figure 1. Algorithm for Monte Carlo simulation.
has also been discussed (Meakin, 1990; Tambe et al., 1994). The incorporation of this step also has the advantage of generalizing the overall reaction scheme. The ZGB reaction scheme with the Eley-Rideal step (will henceforth be called the ERZGB model) can be written as follows:
CO + S f COS
(1)
O2 + 2S f 2OS
(2)
COS + OS f CO2v + 2S
(3)
CO + OS f CO2v + S
(4)
The fourth step, viz., the Eley-Rideal step, has been included in the earlier Monte Carlo simulations by Meakin (1990) and CA simulations by Tambe et al. (1994). In this paper we have incorporated the diffusion step to the ERZGB model (the ZGB model with the inclusion of the Eley-Rideal step) with the aim of studying its effect on the qualitative and quantitative behavior of the system. Simulations based on the Monte Carlo approach as well as cellular automaton have been carried out with the aim of comparing the results qualitatively and quantitatively. 2. Simulation Procedure 2.1. Monte Carlo Simulation. We first simulated the ERZGB model as a function of the adsorption to
diffusion ratio by Monte Carlo approach. All simulations were done with a two-dimensional lattice (50 × 50). The simulations begin by picking a random number that decides whether an adsorption or diffusion step occurs (this depends on the adsorption to the diffusion ratio). If adsorption is the chosen step, a site is randomly chosen next. If the site is vacant, we choose CO or O2 with another random number. This depends on the mole fraction of CO or O2 in the gas phase reservoir above the lattice. If CO was chosen, it adsorbs onto the site. All four nearest sites are checked for oxygen atoms, and if present, the reaction occurs with a randomly chosen oxygen atom. If oxygen was the chosen molecule, a second site is randomly picked, and if empty, the oxygen molecule adsorbs dissociatively to the two sites. All six neighbors are checked for the presence of a CO, and if present, one CO-O pair is reacted randomly. If the second site is occupied by a CO, the trial ends. If the randomly picked site is occupied by an oxygen atom and the chosen molecule is CO, the Eley-Rideal reaction occurs and the site is vacated. On the other hand, if the chosen mechanism in the beginning of the trial is diffusion, a CO in the lattice is randomly chosen. If its randomly chosen neighbor is occupied by an O or CO, the trial ends. If it is empty, the CO diffuses to that site. The simulation as described above was conducted until 104 MC steps. This procedure has been depicted in Figure 1. 2.2. CA Simulations. In the CA simulations, the lattice is divided into Margolous blocks consisting of 2
2190 Ind. Eng. Chem. Res., Vol. 37, No. 6, 1998
Figure 2. Cellular automata transition rules for CO oxidation on catalytic surface with the Eley-Rideal step.
× 2 cells. Such division correctly accounts for the reaction stoichiometry (Mai and von Niessen, 1991; Tambe et al., 1994). Any cell inside a Margolous block is a neighbor of three other cells. The adsorption and reaction steps are carried out with the help of CA transition rules. The probability for adsorption and reaction transitions depends on the partial pressure of CO in the gas phase, the classical statistical weights of the individual configurations, and the present state of a Margolous block. The transition rules along with the respective probabilities are shown in Figure 2. As an illustration, we consider a Margolous block which is empty. The various possible transitions are as follows: (1) One, two, three, or four CO molecules can be adsorbed with probabilities 4yCO, 12yCO2, 24yCO3, 24yCO4, respectively; or (2) one or two oxygen molecules can get dissociatively adsorbed with probabilities 12yO2 or 24yO22, respectively. Simultaneous adsorption, reaction coupled with the diffusion step are carried out similar to the strategy followed by Mai and von Niessen (1993). An example of the temporal evolution of the system is illustrated in Figure 3. As can be seen from the figure, a driving velocity pointing to one of the lattice directions is randomly attributed to each CO molecule. All configurations of CO molecules pointing to a vacant site
are randomly reduced to a pair consisting of a vacant site and a CO molecule pointing toward it, and the location of the vacant site and the CO molecule is exchanged. This diffusion step is carried out for all such pairs. The reaction step is then carried out synchronously within Margolous blocks. After one sweep through the lattice, the blocks are shifted by one cell to the right followed by a shift to a cell down, then left, and up to obtain all possible reactive configurations in a Margolous block. Following the reaction step, a rotation step consisting of a random rotation of CO molecules clockwise or anticlockwise is performed. Adsorption and reaction including the Eley-Rideal step are carried out in Margolous blocks subject to the transition rules described in Figure 2. The evolution of the system can then be described as follows:
|t + 1〉 ) (Rot. Θ R Θ D)m Θ (R Θ A)n|t〉 In the above equation the randomly occurring order of events m and n depend on the adsorption to diffusion ratio. In the above CA simulations we have used a lattice dimension of 1000 × 1000. The results were verified by expanding the lattice space to 2000 × 2000. For any given value of yCO, several runs were performed
Ind. Eng. Chem. Res., Vol. 37, No. 6, 1998 2191
Figure 3. Temporal evolution of lattice space. Table 1. Phase Transition (Critical yCO) Values for Different A/D Ratios A/D model
no diffusion
A/D ) 1
A/D ) 1/10
A/D ) 1/40
MC CA
0.4971 0.4315
0.5624 0.5040
0.62862 0.5219
0.6474 0.5438
with different initial random number seeds. Every such run consisted of several thousand lattice updates which were increased several-fold near the phase transition point. The simulation results are shown in Table 1. Both MC and CA results indicate that with a decrease in the adsorption to diffusion ratio the phase transition point is shifted to the right, thereby increasing the reactive regime in the lattice. This is due to the fact that, with an increase in CO mobility, there is an increase in the reactive events and the steady state coverage of O and CO decrease. In comparison with the earlier ZGB model (Ziff et al., 1986), the phase transition value for any particular adsorption to diffusion ratio (ADR) is lower for the ERZGB model. These results also conform to the earlier results of the ERZGB model without the diffusion step by Meakin (1990) and Tambe et al. (1994). A comparative look at the MC and CA diffusion results indicates that, for any given particular ADR ratio, the CA predicts lower values for the phase transition point than the corresponding MC simulations. This is due to the different ways in which the neighborhoods are defined in CA and MC. Both in MC and in CA simulations we found incorporation of diffusion has not altered the qualitative nature of the phase transition. In other words, the phase transition remains first order irrespective of whether or not the diffusion step is incorporated. This conforms to the earlier results obtained by Mai and von Niessen (1993) for the ZGB model. 3. Conclusions In summary, the present work comprehensively analyzes the effects of adsorption, diffusion, and reaction using MC and CA simulations. Both simulation results indicate that oxygen does not poison the surface, if the Eley-Rideal step is included, irrespective of whether
CO molecules diffuse or not. Thus, we can conclude that diffusion of CO does not alter the qualitative nature of the results. Both MC and CA simulations predict lower values for the phase transition point for the ERZGB model when compared with the original ZGB model. In general there is qualitative and near quantitative agreement between both approaches for any particular value of the adsorption to diffusion ratio. We found CA simulations to be severalfolds faster than MC simulations. Thus, the advantages of fast CA calculations can, in the future, be applied to characterize various other surface catalyzed reaction mechanisms. Acknowledgment Prof. L. K. Doraiswamy during his tenure at the National Chemical Laboratory, Pune, India, emphasized not only experimental research but also theoretical analysis and simulations. It is in this spirit that we wish to dedicate this work to him. We also gratefully acknowledge financial support from the Department of Science and Technology, India. Literature Cited Albano, E. V. The Monte Carlo simulation method: A powerful tool for the study of reaction process. Heterog. Chem. Rev. 1996, 3, 389. Beusch, H.; Fieguth, P.; Wicke, E. Thermisch und Kinetisch Verursachte Instabilitaten Im Reaktionsverhalten Einzelner Katalysatorkorner. Chem.-Ing.-Technol. 1972, 44, 445. Carberry, J. J. Chemical and Catalytic Reaction Engineering; McGraw-Hill: New York, 1976. Chopard, B.; Droz, M. Cellular Automata Approach to Nonequilibrium Phase Transitions in a Surface Reaction Model: Static and dynamic properties. J. Phys. A: Math. Gen. 1988, 21, 205. Chopard, B.; Droz, M.; Kolb, M. Cellular Automata Approach to Nonequilibrium Diffusion and Gradient Percolation. J. Phys. A: Math. Gen. 1989, 22, 1609. Engel, T.; Ertl, G. In The chemical physics of solid surfaces and heterogeneous catalysis; King, D. A., Woodruff, D. P., Eds.; Elsevier: Amsterdam, 1982; Vol. 4. Ertl, G. In Catalysis: Science and technology; Anderson, J. R., Boudart, M., Eds.; Springer: Berlin, 1983; Chapter 3. Fichthorn, K. A.; Ziff, R. M.; Gulari, E. In Catalysis; Ward, J. M., Ed.; Elsevier: Amsterdam, 1987; Chapter 3. Hartman, H.; Tamayo, P. Reversible Cellular Automata and Chemical Turbulence. Physica D 1990, 45, 293.
2192 Ind. Eng. Chem. Res., Vol. 37, No. 6, 1998 Hlavacek, V.; van Rompay, P. Current Problems of Multiplicity, Stability and Sensitivity of States in Chemically Reacting Systems. Chem. Eng. Sci. 1981, 36, 1587. Hugo, P.; Jakubith, M. Dynamisches Verhalten und Kinetik der Kohlenmonoxid-Oxidation am Platin- Katalysator. Chem.-Ing.Technol. 1972, 44, 383. Kohler, J.; ben-Avraham, D. The Dimer-Trimer Model for Heterogeneous Catalysis: A New Universality Class. J. Phys. A: Math. Gen. 1991, 24, L621. Lagos, R. E.; Sales, B. C.; Suhl, H. Theory of Oscillatory Oxidation of Carbon Monoxide over Platinum. Surf. Sci. 1979, 82, 525. Levermore, C. D.; Boghosian, B. M. In Cellular Automata and Modelling of Complex Physical Systems; Manneville, P., Boccara, N., Vichniac, G. Y., Bideaux, R., Eds.; Springer Proceedings in Physics; Springer: Heidelberg, 1989. Mai, J.; von Niessen, W. Cellular Automaton Approach to a Surface Reaction. Phys. Rev. A 1991, 44, R6165. Mai, J.; von Niessen, W. A Cellular Automaton model with Diffusion for a Surface Reaction. Chem. Phys. 1992, 165, 57. Mai, J.; von Niessen, W. Diffusion and Reaction in Multicomponent Systems via Cellular Automaton Modelling: A + B2. J. Chem. Phys. 1993, 98, 2032. Mai, J.; von Niessen, W.; Blumen, A. The CO + O2 Reaction on Metal Surfaces. Simulation and Mean-Field Theory: The influence of diffusion. J. Chem. Phys. 1990, 93, 3685. McCarthy, E.; Zahradnik, J.; Kuczynski, G. C.; Carberry, J. J. Some unique aspects of CO oxidation on supported Pt. J. Catal. 1975, 39, 29. Meakin, P. Simple Models for Heterogeneous Catalysis with a Poisoning Transition. J. Chem. Phys. 1990, 93, 2903. Meakin, P.; Scalapino, D. J. Simple Models for Heterogeneous Catalysis: Phase Transition-like Behavior in Non-Equlibrium Systems. J. Chem. Phys. 1987, 87, 731. Moller, P.; Wetzl, K.; Eiswirth, M.; Ertl, G. Kinetic Oscillations in the Catalytic CO Oxidation on Pt(100): Computer Simulations. J. Chem. Phys. 1986, 85, 5328.
Pan, H. Y.; Wang H. J.; Zhao Z. S. Dimer-dimer surface reaction of Albano’s type: A cellular automaton approach. J. Phys. A: Math. Gen. 1995, 28, 4279. Sahimi, M.; Stauffer, D. Efficient Simulation of Flow and Transport in Porous Media. Chem. Eng. Sci. 1991, 46, 2225. Sheintuch, M.; Schmitz, R. A. Oscillations in Catalytic Reactions. Catal. Rev.sSci. Eng. 1977, 15, 107. Tambe, S. S.; Jayaraman, V. K.; Kulkarni, B. D. Cellular Automaton Modelling of a Surface Catalysed Reaction with Eley-Rideal Step: The Case of CO oxidation. Chem. Phys. Lett. 1994, 225, 303. Von-Neumann. In The theory of self-reproducing automata; Burks, M. W., Ed.; University of Illinois: Urbana, IL, 1966. Ziff, R. M.; Brosilow, B. J. Investigation of the First-Order Phase Transitions in the A-B2 Reaction Model using a ConstantCoverage Kinetic Ensemble. Phys. Rev. A 1992, 46, 4630. Ziff, R. M.; Gulari, E.; Barshad, Y. Kinetic Phase Transitions in an Irreversible Surface Reaction Model. Phys. Rev. Lett. 1986, 56, 2553. Ziff, R. M.; Fichthorn, K.; Gulari, E. Cellular automaton version of the AB2 Reaction Model Obeying Proper Stoichiometry. J. Phys. A: Math. Gen. 1991, 24, 3727. Zygourakis, K.; Bizios, R.; Markenscoff, P. Proliferation of Anchorage-Dependent Contact Inhibited Cells: I. Development of Theoretical Models Based on Cellular Automata. Biotechnol. Bioeng. 1991a, 38, 459. Zygourakis, K.; Markenscoff, P.; Bizios, R. Proliferation of Anchorage-Dependent Contact Inhibited Cells: II. Experimental Results and Validation of Theoretical Models. Biotechnol. Bioeng. 1991b, 38, 471.
Received for review July 23, 1997 Revised manuscript received November 25, 1997 Accepted January 5, 1998 IE9705261
