Applied Surface Science 256 (2010) 5376–5380
Contents lists available at ScienceDirect
Applied Surface Science journal homepage: www.elsevier.com/locate/apsusc
Monte Carlo model of CO adsorption on supported Pt nanoparticle A.V. Myshlyavtsev a,b, P.V. Stishenko a,* a b
Omsk State Technical University, Omsk, Russia Institute of Hydrocarbons Processing SB RAS, Omsk, Russia
A R T I C L E I N F O
A B S T R A C T
Article history: Available online 4 January 2010
For molecular simulations with thousands of atoms it is desirable to use a lattice gas model because it is fast and easy-to-use for computations. Unfortunately, simulation of adsorption on heterogeneous surfaces within this model is rather complicated due to a large variety of available adsorption site types. We propose the combined model with lattice representation of adsorbent atoms and arbitrary location of adsorbate atoms. Using this model simulation of CO adsorption on supported Pt nanoparticles has been performed. With the proposed approach the above-mentioned difficulties were successfully overcome. ß 2009 Elsevier B.V. All rights reserved.
Keywords: Adsorption Monte Carlo simulation Pt nanoparticle Carbon monoxide
1. Introduction In heterogeneous catalysis, adsorption and reaction processes usually occur on supported metal nanoparticles. To study such processes experimentally, a vast variety of surface science techniques came into the use in the last decades of the 20th century. Besides, considerable progress in manufacturing of model catalysts with well-defined properties like particle size, shape and separation has been achieved. Successful interpretation of the experimental results requires constructing mathematical models to provide accurate description of the systems under consideration. Application of the conventional mean-field models is rather limited here because of the peculiarities of the reaction performance on the nanometer scale, including the inherent heterogeneity of metal crystallites as well as spontaneous and adsorbateinduced changes of the shape and degree of dispersion of supported catalysts. Despite the potential power of molecular dynamics simulation method, the use of this technique for the analysis of the physical and chemical processes over the supported particles is also limited due to the short length and time scales typical of this approach. Under such circumstances, the use of stochastic simulations based on the Monte-Carlo technique, is almost inevitable [1–3]. Often as a model of the supported particle a top projection of the truncated pyramid has been considered, i.e., the pyramid is represented by a N N square lattice, where the central M M array of sites mimics the top facet (active surface of the particle)
* Corresponding author at: Omsk State Technical University, prospekt Mira, 11, g. OMSK 644050, Russia. Tel.: +7 3812 65 23 79; fax: +7 3812 65 23 79. E-mail address:
[email protected] (P.V. Stishenko). 0169-4332/$ – see front matter ß 2009 Elsevier B.V. All rights reserved. doi:10.1016/j.apsusc.2009.12.084
and the periphery corresponds to the support (or to the side facets with different catalytic properties). Despite the simplicity of such models, both predictable and unexpected results were obtained by modelling of different peculiarities the reaction performance over the supported particles. Such factors as reactant supply due to the diffusion over the support (spillover) [4], interplay of the reaction kinetics on different facets of the supported particle due to the adsorbed species diffusion between the facets [5], jump-wise reshaping of the active particle under the influence of the adlayer composition [6], oscillations and chaos in the catalytic reactions occurring on the supported catalysts [7], influence of geometric parameters of the supported particles on the reactivity and selectivity of the multi-route catalytic reactions [8], reciprocal effect of the surface morphology and adsorption-reaction processes over the catalytic particle [9], were examined in various studies. However, the theoretical models taking into account the dynamic change of the shape and surface morphology of the catalytic particles under the influence of adsorbed species and temperature action are practically lacking. One example of such a model is the one elaborated in [10–12]. It is based on the well known ‘‘solid-onsolid’’ (SOS) model [13]. As a model of the catalytic metal particle the SOS located on a neutral support was considered. Unfortunately, this model is too simple to consider a real supported nanocrystal. The goals of our study are (1) elaboration of the statistical lattice model of the supported catalyst particle using the realistic manybody interatomic potential, (2) elaboration of the non-lattice model of the adsorption taking into account the presence of different adsorption sites and (3) characterization of dependency between the chemical potential of adsorbate and the coverage of particle faces and adsorption sites.
A.V. Myshlyavtsev, P.V. Stishenko / Applied Surface Science 256 (2010) 5376–5380
2. The model Using the Monte-Carlo technique we have defined equilibrium state of support-nanoparticle-adsorbate system on different values of temperature and chemical potential for CO adsorption process on the Pt nanoparticle surface. This process is of a great interest for industry and surface science. Simulations were conducted with the Metropolis algorithm in a grand canonical ensemble. To represent location of nanoparticle and support atoms the lattice gas model was employed. This model is rather simple to use and implement. It is realistic enough for nanoparticles consisting of several thousands of atoms. This size is typical of nanoparticles used in industrial catalysts. Another advantage of this model is simple enough algorithmic accomplishment of almost any concept describing the occurrence of physical and chemical processes. In the case of modelling adsorption, it is important that adsorption energies obtained from experiments or more accurate simulations (i.e., ab initio or molecular dynamics) are usually linked to concrete types of adsorption sites on the definite faces. Therefore this data can be used only if the adsorbent crystal has enough strict structure to detect the face and to find a fitting adsorption site on it. For calculation of interaction energy between the nanoparticle atoms an empirical tight-binding potential has been employed. Potentials of this type have been used in many investigations by simulation of free and supported nanoparticles [14–17]. Atomic interactions potentials used in the model have defined a form of the system energy Hamiltonian:
E¼
0 1 sffiffiffiffiffiffiffiffiffiffiffiffiffi N X X X ffi @ Ri j Bi j A i¼1
j 6¼ i
(1)
j 6¼ i
Here N is a number of atoms into the particle, Rij and Bij—the repulsion and bound terms of interaction energy between i-th and j-th atoms respectively. j runs through all atoms distance to which from i-th atom does not exceed some specified value. Rij and Bij are calculated in different ways for different potentials. For the TB-SMA potential: Ri j ¼ Ai j e pi j ðri j =r0 1Þ 2
Bi j ¼ ji j e2qi j ðri j =r0 1Þ
(2) (3)
Here A, j, p, q and r0 – the coefficients defined by interacting atoms types, ri j – the distance between atoms. Values of these coefficients and the way they were obtained are specified below in the simulation details section.
5377
Table 1 TB-SMA parameters for Pt–Pt and Pt-support interactions. p and q are not specified for Pt-support because they are irrelevant for the first coordination sphere.
Pt–Pt Pt–support
A, eV
j, eV
p
q
r0, A˚
0.2975 0.2800
2695 2695
10.612 –
4.004 –
3.92 3.92
To decrease computational cost of energy calculation we took into account only the interactions up to the fifth coordination sphere. In units of distance it means 8.76 A˚. For calculation of interaction energy between the nanoparticle atoms and the support atoms the same kind of potential has been used, but interaction radius was limited by the first coordination sphere, or up to 3.92 A˚. Such approach to simulation of supported nanoparticles was successfully used in Ref. [18] with similar empirical potentials. The values of potential parameters for nanoparticlesupport interactions are also specified in the simulation details section. For realistic adsorption simulation on nanoparticle surface it is necessary to take into account many adsorption site types. It makes the use of lattice model for this purpose rather complicated. Although the number of all possible adsorption sites is finite, in practice it is hard to enumerate them all and create a complete lattice for adsorbent molecules. Every adsorption site type introduces in the lattice unit cell several points obtained by different spatial orientation of atoms which form an adsorption site. For example, for a bridge adsorption on (1 1 1) face there are 24 possible spatial orientations—3 orientations obtained by rotation around normal to each of 8 faces. Large unit cells with many points suitable only for the molecules adsorbed on some definite site types make simulation algorithms more complicated, slow and memory consuming. To avoid these difficulties instead of a lattice gas model we used a continuous space model to store coordinates of adsorbate molecules. Continuous model allows arbitrary location of particles. In the employed model each adsorbed molecule is represented as a single point in continuous space. This point designates location of adsorbate atom which takes part in adsorption bound and can be used for calculation of lateral interactions between adsorbate molecules. To describe interactions between the adsorbed molecules a simple hard-sphere potential was employed and the hard-sphere radius was taken as equal to the catalyst lattice period, for platinum it means 3.92 A˚. Of course, the real life interactions between the adsorbate molecules are much more complex, but quantitative calculation of exact data about
Fig. 1. Supported Pt nanoparticle (4033 atoms).
5378
A.V. Myshlyavtsev, P.V. Stishenko / Applied Surface Science 256 (2010) 5376–5380
Fig. 2. Adsorbed CO molecules on 4033 atom Pt nanoparticle (colour: green spheres designate CO molecules, white spheres: Pt atoms, orange spheres: Pt atoms employed in the occupied adsorption sites; black and white: light gray spheres designate CO molecules, white spheres: Pt atoms, dark gray spheres: Pt atoms employed in the occupied adsorption sites). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of the article.)
adsorption was not a goal of this work. So, for the sake of simplicity, we used the hard-sphere approximation. Coordinates of adsorption centers in the continuous space are being detected by searching suitable configurations on the nanoparticle surface in the lattice gas model. These configurations are given for each adsorption site type as a set of internal coordinates of each atom or empty lattice point which form an adsorption site. Computational resources required for searching given configurations in the lattice data structure are close to those for enumeration of atom neighbours.
Each adsorption site type corresponds to a known value of energy. This value is used to calculate contribution of the adsorbed molecule to the total system energy. This value is defined specifically for the site type bound energy and with the chemical potential of the adsorbed molecule. In this work we described and used in simulations the following types of adsorption sites (active centers): a-top, bridge, threefold hcp and fcc hollows on (1 1 1) face and a-top, bridge and fourfold hollow on (1 0 0) face. Atoms of these faces almost completely form surface of platinum nanoparticles investigated in this work.
Fig. 3. Equilibration of relative proportions of the occupied adsorption sites in the course of Monte Carlo simulation. Last stage of the simulated annealing procedure (m = 1.9 eV, T = 581 K).
A.V. Myshlyavtsev, P.V. Stishenko / Applied Surface Science 256 (2010) 5376–5380
5379
3. Simulation details For system energy calculation a second-moment approximation of tight-binding potential (TB-SMA) has been employed. The parameter values for the potential were taken from Cleri and Rosato work [19] where they were fitted to reproduce experimental values «of the cohesive energy, lattice parameters (by a constraint on the atomic volume), and independent elastic constants». The values of these parameters are shown in Table 1. Energy of interactions between the nanoparticle and the support atoms was calculated with the same potential, but only atoms in the first coordination sphere were taken into account. Potential parameters were adjusted to reproduce particle shape imposed by the TEM images [20]. Their values are also shown in Table 1. The values of p and q are not specified, because for the first coordination sphere the multiplier r/r0 1 near these parameters is always zero, so they do not affect the calculations results. It should be noticed that the lattice parameter for support was taken as equal to the nanoparticle one. We suppose that these simplifications do not significantly affect adsorption on particle surface but only shape of the particle matters. If more realistic support description is considered to be necessary and the information about the particle-support interaction energy will be available for the simulated support material, it can be introduced without modification of the proposed algorithm whereas the lattice model is still acceptable for the particle itself, which is certainly the case for big enough particles. The energy values for CO adsorption on each site type were taken from [21,22], where they were calculated with the molecular dynamics method in the density functional theory (DFT) framework. Simulation was started from searching for equilibrium shape for the supported nanoparticles. Initially a cubic particle was placed on the support and with Metropolis algorithm and simulated annealing method an equilibrium state was achieved. Simulations were performed for the particles consisting of 4033 atoms. The simulated annealing procedure consisted of three consequent steps: 6000 Monte Carlo steps (MCS) at 1977 K, 1000 MCS at 1163 K and 1000 MCS at 581 K. After that isotherms were built for the temperature values 465, 581 and 697 K in the range of the chemical potential (m) from 2.3 to 1.3 eV. The system was considered as the equilibrium after stabilization of relative proportions of occupied adsorption sites. 4. Results and discussion Fig. 1 shows supported 4033 atom Pt nanoparticle obtained in simulation with the Metropolis algorithm and the simulated annealing method. From this particle shape adsorption simulation was started. In the figure one can see that this shape is similar to the TEM images in [20]. In Fig. 2 one can see the same particle after equilibration with the adsorption process on different values of chemical potential. The system was considered to be equilibrium based upon the graph in Fig. 3. This graph shows a number of occupied adsorption sites on each step. Overall height of the curve represents the total number of adsorbed molecules. Each area under the curve corresponds to an adsorption site type. As one can see from the graph in Fig. 3, a number of performed simulation steps is sufficient to consider the system to be equilibrium. Fig. 4 shows the obtained adsorption isotherm with the parts of each adsorption site type for 4033 atom nanoparticle in the same manner as in Fig. 3. Similar distributions of occupied adsorption sites were obtained for 465, 581 and 697 K temperature values. In general, the obtained distribution of occupied adsorption sites conforms to the used values of adsorption energies–adsorption on the sites with lower energy is more probable. At low values of chemical potential adsorption occurs only on (1 0 0) faces, and
Fig. 4. Parts of each adsorption site type depending on the adsorbate chemical potential.
these faces become fully covered when m reaches 2.1 to 2.0 eV. Covering of (1 1 1) faces begins when m exceeds 2.0 eV and reaches full coverage at m = 1.7 eV. The character of the obtained isotherms is the same for all temperature values employed in the simulation. It was experimentally revealed in [23] that part of occupied bridge adsorption sites on (1 0 0) face is coverage dependent. We tried to determine if this effect can be reproduced in a framework of the proposed model. The graphs in Fig. 5 show the change of bridge adsorption contribution to the total coverage of (1 0 0) faces. One can see that at low coverage levels part of adsorption on the bridge sites is more significant than at high coverage and is about 50% for m = 2.3 eV. With the growth of chemical potential and therefore of coverage level, the part of bridge adsorption decreases down to 10–30%%. It does not conform to general
5380
A.V. Myshlyavtsev, P.V. Stishenko / Applied Surface Science 256 (2010) 5376–5380
that of lateral interactions between the adsorbed molecules on high coverage values. The last factor can be taken into account by introduction of more accurate potential for the adsorbed molecules in each particular case. We can conclude that combination of the lattice gas model for the adsorbent and the continuous space model for the adsorbate can be used to overcome difficulties caused by variety of different active center types on realistic nanoparticles. The data structures used for the lattice gas and continuous space models do not introduce any limitations on simulation algorithms and keep the model open to extension with new factors employing already known techniques and optimizations. References Fig. 5. Part of occupied bridge adsorption sites on (1 0 0) face depending on the chemical potential.
dependency of coverage on m, because energy for on the bridge site (2.27 eV) is lower than for the a-top sites (2.22 eV). In a framework of the employed model, it can be explained by higher requirements for configuration on the bridge adsorption site. It means that adsorption on the a-top site needs only one noncovered atom of Pt, whereas for the on-bridge adsorption two nearby Pt atoms are required. This difference becomes significant on high coverage in the presence of surface defects and near edges of the nanoparticle, because these defects induce appearance of domains with the a-top adsorptions. The reproduced dependency of bridge adsorption contribution is not quantitatively precise. We suppose that for more accurate results, more attention must be paid to mutual influence of the adsorbed molecules. A more precise potential instead of hardspheres one (e.g. Lennard–Jones) can be easily introduced and sufficiently improve the model accuracy. With the proposed model a new adsorption site type can be easily introduced if it is considered to be significant for the simulated process. Only description of this adsorption site in the internal coordinates and adsorption energy obtained from experiments or computations would be required. 5. Conclusion In the presented model lattice representation of the adsorbent nanoparticle and the continuous space model for location of adsorbate molecules were combined. It has allowed to simulate adsorption on a realistic Pt nanoparticle surface. The obtained results reproduced the expected coverage growth with the increase of chemical potential. Also known dependency of proportion of different occupied adsorption site types on the surface coverage was qualitatively reproduced. For different temperature values, similar results were obtained and character of the reproduced dependency was the same. It was shown, that coverage of different faces is determined by the energies of adsorption sites and, to some extent, by probability of their appearance on the face. Model is accurate enough for simulation of adsorption systems with a finite number of adsorption site types and with known adsorption energies in the thermodynamic equilibrium. Deviations of model predictions from the precise values can be caused by the influence of unaccounted adsorption sites, by the non-equilibrium state of the simulated system and by
[1] V.P. Zhdanov, B. Kasemo, Simulations of the reaction kinetics on nanometer supported catalyst particles, Surf. Sci. Rep. 39 (2000) 25–104. [2] V.P. Zhdanov, Impact of surface science on the understanding of kinetics of heterogeneous catalytic reactions, Surf. Sci. 500 (2000) 966–985. [3] V.I. Elokhin, A.V. Myshlyavtsev, Catalytic processes over supported nanoparticles: simulation, in: J.A. Schwarz, C.I. Contescu, K. Putyera (Eds.), Dekker Encyclopedia of Nanoscience and Nanotechnology, Marsel Dekker, New York, 2004, pp. 621– 632. [4] V.P. Zhdanov, B. Kasemo, Kinetics of rapid heterogeneous reactions on the nanometer scale, J. Catal. 170 (1997) 377–389. [5] V.P. Zhdanov, B. Kasemo, Monte Carlo simulation of the kinetics of rapid reactions on nanometer catalyst particles, Surf. Sci. 405 (1998) 27–37. [6] V.P. Zhdanov, B. Kasemo, Kinetics of rapid heterogeneous reactions accompanied by the reshaping of supported catalyst crystallites, Phys. Rev. Lett. 81 (1998) 2482–2485. [7] H. Persson, P. Thormahlen, V.P. Zhdanov, B. Kasemo, Simulations of the kinetics of rapid reactions on supported catalyst particles, Catal. Today 53 (1999) 273–288. [8] A.S. McLeod, L.F. Gladden, Relating metal particle geometry to the selectivity and activity of supported-metal catalysts: a Monte Carlo Study, J. Catal. 173 (1998) 43–52. [9] F. Gracia, E.E. Wolf, Monte Carlo simulations of the effect of crystallite size on the activity of a supported catalyst, Chem. Eng. J. 82 (2001) 291–301. [10] E.D. Resnyanskii, E.I. Latkin, A.V. Myshlyavtsev, V.I. Elokhin, Monomolecular adsorption on rough surfaces with dynamically changing morphology, Chem. Phys. Lett. 248 (1996) 136–140. [11] E.D. Resnyanskii, A.V. Myshlyavtsev, V.I. Elokhin, B.S. Bal’zhinimaev, Dissociative adsorption on rough surfaces with dynamically changing morphology: a Monte Carlo model, Chem. Phys. Lett. 264 (1997) 174–179. [12] E.V. Kovalyov, V.I. Elokhin, A.V. Myshlyavtsev, Stochastic simulation of physicochemical processes performance over supported metal nanoparticles, J. Comput. Chem. 29 (2008) 79–86. [13] J. Lapujoulade, The roughening of metal surfaces, Surf. Sci. Rep. 20 (1994) 195– 249. [14] A. Rapallo, G. Rossi, R. Ferrando, A. Fortunelli, B.C. Curley, L.D. Lloyd, G.M. Tarbuck, R.L. Johnston, Global optimization of bimetallic cluster structures. I. Size-mismatched Ag–Cu, Ag–Ni, and Au–Cu systems, J. Chem. Phys. 122 (2005) 194308. [15] F. Baletto, C. Mottet, R. Ferrando, Growth simulations of silver shells on copper and palladium nanoclusters, Phys. Rev. B 66 (2002) 155420. [16] B.B. Maranville, M. Schuerman, F. Hellman, Simulation of clustering and anisotropy due to Co step-edge segregation in vapor-deposited CoPt3, Phys. Rev. B. 73 (2006) 104435. [17] D. Cheng, W. Wang, S. Huang, Core-shell-structured bimetallic clusters and nanowires, J. Phys.: Condens. Mat. 19 (2007) 356217. [18] K.P. McKenna, P.V. Sushko, A.L. Shluger, Transient atomic configurations of supported gold nanocrystallites at finite temperature, J. Phys. Chem. C 111 (2007) 2823–2826. [19] F. Cleri, V. Rosato, Tight-binding potentials for transition metals and alloys, Phys. Rev. B 48 (1993) 22–33. [20] K.M. Bratlie, H. Lee, K. Komvopoulos, P. Yang, G.A. Somorjai, Platinum nanoparticle shape effects on benzene hydrogenation selectivity, Nano Lett. 10 (2007) 3097–3101. [21] H. Orita, N. Itoh, Y. Inada, All electron scalar relativistic calculations on adsorption of CO on Pt(1 1 1) with full-geometry optimization: a correct estimation for CO site-preference, Chem. Phys. Lett. 384 (2004) 271–276. [22] S. Yamagishi, T. Fujimoto, Y. Inada, H. Orita, Studies of CO adsorption on Pt(1 0 0), Pt(4 1 0), and Pt(1 1 0) surfaces using density functional theory, J. Phys. Chem. B 109 (2005) 8899–8908. [23] R. Martin, P. Gardner, A.M. Bradshaw, The adsorbate-induced removal of the Pt{1 0 0} surface reconstruction. Part II: CO, Surf. Sci. 20 (1995) 69–84.