Langmuir 1997, 13, 903-908
903
Monte Carlo Simulation Study for Adsorption of Dimers on Random Heterogeneous Surfaces† Tomoshige Nitta,* Hideki Kiriyama, and Takeshiro Shigeta Department of Chemical Engineering, Faculty of Engineering Science, Osaka University, Toyonaka, Osaka, 560, Japan Received October 30, 1995. In Final Form: February 27, 1996X Adsorption isotherms of dimers on random heterogeneous surfaces, based on the lattice gas model, are calculated by the grand canonical ensemble Monte Carlo simulation method. Random heterogeneous surfaces are presumed to consist of two active sites with different energies a and b which are randomly distributed on a simple square lattice, on which dimer molecules composed of the identical segments occupy two adjacent sites without lateral interactions. Calculated adsorption isotherms are used to test the applicability of two versions of the multisite occupancy model originally proposed by Nitta et al. (1984); one is characterized by the quasi-chemical approximation for segment-site pairs and another is by moleculebond pairs. The two model equations are found to be good only when the difference in rb ) exp{(b - a)/kT} is smaller than 10. For rb values greater than 50, some discrepancies are observed in such a way that the first model overestimates the energetically more stable states while the second model underestimates them.
1. Introduction Surface heterogeneity, which comes from different atoms, different slit width, and surface edges and cracks, is known to affect the adsorption equilibrium and catalysis on solid surfaces.1,2 The adsorption energy distribution function is an important property used to describe the surface heterogeneity; however, the surface topography is also an important factor used to characterize the surface heterogeneity. The patchwise heterogeneous surface3 and the random heterogeneous surface4 are two extreme models representing the surface topography, while the correlated energy surface proposed by Zgrablich et al.5-7 is an intermediate topography model. Although the size of adsorbate molecules is an important factor influencing the adsorption isotherms,1 there have been very few theories on multisite occupancy adsorption on heterogeneous surfaces.1,8-11 One significant problem for the adsorption of oligomers on heterogeneous surfaces is how to take into account the entropic defect due to the preferential adsorption of molecular segments on energetically more stable sites. The multisite occupancy model proposed by Nitta et al.9 used the quasi-chemical approximation (QCA) theory which was first introduced by Guggenheim12 for treating the local compositions occurring * To whom correspondence should be addressed: Tel and Fax: +81-6-850-6265. E-mail:
[email protected]. † Presented at the Second International Symposium on Effects of Surface Heterogeneity in Adsorption and Catalysis on Solids, held in Poland/Slovakia, September 4-10, 1995. X Abstract published in Advance ACS Abstracts, September 15, 1996. (1) Rudzinski, W.; Everett, D. H. Adsorption of Gases on Heterogeneous Surfaces; Academic Press: London, 1992. (2) Jaroniec, M.; Madey, R. Physical Adsorption on Heterogeneous Solids; Elsevier: Amsterdam, 1988. (3) Ross, S.; Olivier, J. P. On Physical Adsorption; Interscience: New York, 1964. (4) Hill, T. L. J. Chem. Phys. 1949, 17, 762. (5) Ripa, P.; Zgrablich, G. J. Phys. Chem. 1975, 79, 2118. (6) Riccardo, J. L.; Chade, M. A.; Pereyra, V.; Zgrablich, G. Langmuir 1992, 8, 1518. (7) Riccardo, J. L.; Pereyra, V.; Zgrablich, G.; Rojas, F.; Mayagoita, V.; Kornhauser, I. Langmuir 1993, 9, 2730. (8) Marczewski, A. W.; Derylo-Marczewska, M.; Jaroniec, M. J. Colloid Interface Sci. 1986, 109, 310. (9) Nitta, T.; Kuro-oka, M.; Katayama, T. J. Chem. Eng. Jpn. 1984, 17, 45. (10) Nitta, T.; Yamaguchi, A. J. Chem. Eng. Jpn. 1992, 25, 420. (11) Nitta, T.; Yamaguchi, A. Langmuir 1993, 9, 2618.
S0743-7463(95)00957-7 CCC: $14.00
in liquid mixtures due to the difference in intermolecular interactions. They applied the QCA theory to segmentsite pairs and derived an analytical expression for the isotherm equation. However, the validity or the applicability of the approximation has not been demonstrated so far. Very recently, Zgrablich and his co-workers13 reported the comparison of the multisite occupancy model with adsorption data generated by Monte Carlo (MC) simulations for a model system of dimers on the random heterogeneous surface. They concluded that the validity of the multisite occupancy model is limited in the case where the difference in the adsorption energies of active sites is small and claimed that the Fermi-Dirac (F-D) approach, which treats each energy level of the adsorbed dimers to be independent, is able to qualitatively represent the simulated isotherms of the dimers. However, their comparisons of the F-D theory are made primarily with the experimental data of gases on zeolites, instead of simulated data. The present work is a study to explore the applicability of the multisite occupancy model based on the lattice gas model by using the grand canonical ensemble MC simulation method for the model system of homogeneous dimers on the random heterogeneous surface. In section 2, we describe how we have made random heterogeneous surfaces and their characteristics in terms of bond correlations. Section 3 outlines the interaction model and the MC simulation method. In section 4, we will describe two theoretical models based on the quasi-chemical approximation, whose data will be compared with the simulation data. Comparisons of adsorption isotherms and internal energies are reported in section 5, while those of site coverages and bond coverages are given in section 6. Some concluding remarks will be given in the last section. 2. Random Heterogeneous Surfaces We use a surface of simple square lattices composed of site a and site b. The algorithm we used for making a random heterogeneous surface is briefly outlined, as follows. At first, all (12) Guggenheim, E. A. Mixtures; Clarendon Press: Oxford, 1952; Chapter 11. (13) Ramirez-Pastor, A. J.; Nazzarro, M. S.; Riccardo, J. L.; Zgrablich, G. Surf. Sci. 1995, 341, 249.
© 1997 American Chemical Society
904
Langmuir, Vol. 13, No. 5, 1997
Nitta et al.
Figure 1. Illustrative diagrams of two random surfaces: (a) fb ) 0.1; (b) fb ) 0.3. that have already been occupied. When a dimer occupies sites R and β, the interaction energy is given by -ab ) -a - b. The total energy of N molecules adsorbed is therefore expressed as
Table 1. Number of Bonds for the Random Surface theory
simulation
Maa Mab Mbb
fb ) 0.1 16200 3600 200
16190 3620 190
Maa Mab Mbb
fb ) 0.2 12800 6400 800
12790 6420 790
Maa Mab Mbb
fb ) 0.3 9800 8400 1800
9760 8480 1760
U ) -(Naaaa + Nabab + Nbbbb)
where Naa, Nab, and Nbb are the number of bonds occupied by dimer molecules. We used the standard Metropolis algorithm for three different steps: molecular creation, destruction, and rotation. The transition probability Pmn from state m to state n is given14,15 as
{
Mab ) Mfafb,
Mbb ) Mf 2b/2
Fn Fm
Pmn ) min 1,
of the lattices are set as site a and the number of site b lattices Mb is set at zero. Then we choose one site randomly; if it is the site a, then we change it to site b and add one to Mb. If it is not site a, then we again choose one site randomly and repeat the process until Mb g M fb, where M is the total number of sites and fb is the fraction of site b among M sites. Figure 1 shows two computer-generated random surfaces for fb ) 0.1 and 0.3 where black squares represent site b. The size of a surface is 100 × 100 lattices making M ) 10 000. When the surfaces are perfectly random, the theoretical numbers of bonds per site are given by the following equation:
Maa ) Mf 2a/2,
(2)
(1)
The numbers of bonds for the three surfaces for fb ) 0.1, 0.2, and 0.3 were calculated and divided by 4 (the coordination number) to get the MRβ value comparable to that in eq 1. Table 1 shows the numbers of bonds Maa, Mab, and Mbb for the three surfaces along with the numbers calculated from eq 1. They are not identical but almost the same as to consider the surfaces to be almost energetically random.
3. Interaction Model and Simulation Method The dimer molecule is composed of two identical segments which occupy two neighboring sites. The interaction energies between one segment and sites a and b are denoted by -a and -b, respectively, where a and b are positive. No lateral interaction energy is assumed; therefore, one molecule is like a hard dimer expelling other molecules to occupy the two sites
}
(3)
The ratio of the probability densities Fn/Fm of state n and state m is given for the three steps as follows.
(a) The molecular creation step
(
)
Fn ∆µs Un - Um M ) exp Fm (N + 1) kT kT
(4)
where M is the total number of sites on the surface, N is the number of dimer molecules existing on the lattice, ∆µs is the chemical potential difference on the surface, and Um - Um is the difference of internal energies between the trial state n and the present state m.
(b) The molecular destruction step
(
)
Fn N ∆µs Un - Um ) exp Fm M kT kT
(5)
(c) The molecular rotation pivoting at one end
(
)
Un - Um Fn ) exp Fm kT
(6)
(14) Nicholson, D.; Parsonagee, N. G. Computer Simulation and Statistical Mechanics of Adsorption; Academic Press: London, 1982. (15) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, 1987.
Monte Carlo Study for Adsorption of Dimers
Langmuir, Vol. 13, No. 5, 1997 905
Details of the simulation are as follows: After one trial of rotation, 10 trials of creation or destruction were executed. For equilibration of the system, 440 MMCSs (mega Monte Carlo steps) were discarded and 1100 MMCSs were normally used for ensemble averages. The minimum value of success we set for creation and destruction trials was 20 000; all runs exceeded this criterion. The quantities calculated in the present MC simulations are summarized as follows
The new variables Na and Nb are solved by using the free energy minimization principle subject to the constraint eq 14. The expression of eq 13 for g with eq 14 is equivalent to an intuitive expression for site coverages θa and θb, the Langmuir type equation given by
(
θb/(1 - θb)
) exp
θa/(1 - θa)
)
b - a ≡ rb kT
(15)
surface coverage θ ) 2〈N〉/M
(7)
site coverages θa ) 〈Na〉/Ma;
θb ) 〈Nb〉/Mb
The parameter rb represents the relative strength of adsorption on site b compared to that on site a. From the above relations, the adsorption isotherm equation for a random heterogeneous surface is derived as
(8)
( )
ln Kaf ) - ln 2θ + 2 ln
bond coverages θRβ ) 〈NRβ〉/Mab
internal energy (10)
where Na and Nb are the number of sites a and b, respectively, occupied by dimers.
where f is the fugacity of a dimer molecule in gas phase and Ka is the adsorption equilibrium coefficient for a hypothetical homogeneous surface composed of only site a. For the derivation of eq 16, the reader should refer to the previous paper of Nitta et al.9 or the textbook of Rudzinski and Everett.1 The relation between Kaf and ∆µs is given as Kaf ) exp(∆µs/kT). The site coverage θa is calculated by solving eq 17 with respect to Y,
4. Theoretical Models
fa
The canonical partition function Q composed of N molecules distributing over M sites is expressed by
Q ) jsNg exp(-E/kT)
M! 1 N!(M - 2N)! MN
∑R
Yθ ≡ θa/(1 - θa)
ln g ) ln g0 -
N*ij!(Mij - N*ij)!
∑i ∑j ln
(19)
Nij!(Mij - Nij)!
subject to
Naa + Nab + Nbb ) N
(20)
The above assumptions yield almost the same expressions as the first version of QCA; they are summarized as follows:
the Langmuir type assumption θRβ/(1 - θRβ)
(13)
(18)
Equation 17 is a variant of constraint eq 14. Another version of quasi-chemical approximation is considered. Instead of the number of adsorbed sites, the number of bonds Nij adsorbed on site i and adjacent site j could be a variable for taking into account the entropic effect of preferential adsorption; that is
(12)
N*R!(MR - N*R)! ln
(17)
where Y is defined as
When the surface is heterogeneous, more active sites are preferentially adsorbed for the system to be energetically more stable, which results in the decrease in the number of configurations for distributing molecules on lattice sites. This entropic effect due to the preferential adsorption on more active sites has not been solved exactly; however, one approach for this problem is the quasi-chemical approximation (QCA) originally proposed by Guggenheim12 for treating lateral interactions between unlike molecules and later applied to segment-site interactions by Nitta et al.9 The expression for g by using QCA is given as
ln g ) ln g0 -
rbY Y + fb )1 1 + θY 1 + rbθY
(11)
where jsN is the molecular partition function for a molecule vibrating in three directions on lattice sites, g is the combinatory factor representing the number of ways to distribute N dimer molecules over M heterogeneous sites, and E is the total energy. In the case where the surface is homogeneous, the Flory-Huggins equation for dimers is used for the combinatory factor, which is denoted by g0.
g0 )
(16)
(9)
where Rβ ) aa, ab, bb, and
U ) -〈Naa〉aa - 〈Nab〉ab - 〈Nbb〉bb
θa 1 - θa
θaa/(1 - θaa)
(
) exp
)
Rβ - aa ) rRrβ kT
(21)
NR!(MR - NR)! the adsorption isotherm
where NR denotes the number of occupied sites R and the asterisk represents the random distribution for the homogeneous surface. This term has been introduced for normalzing the g factor. The constraint equation for NR is
Na + Nb ) 2N
(14)
(
ln Kaf ) -ln 2(1 - θ) + ln
)
θaa 1 - θaa
(22)
the definition of a variable Y Yθ ≡ θaa/(1 - θaa)
(23)
906
Langmuir, Vol. 13, No. 5, 1997
Nitta et al.
and the determining equation for Y rabY rbbY Y f 2a + 2fafb + f 2b )1 1 + θY 1 + rabθY 1 + rbbθY
(24)
Equation 24 is solved with respect to Y; then bond coverage θaa is calculated through eq 21. Therefore, the right-hand side of eq 22 is evaluated from surface coverage θ, giving the relation between the fugacity f and θ. It is noted that the set of equations from 21 to 24 resembles the isotherm equation of the F-D approach proposed by Ramirez-Pastor et al.13 5. Adsorption Isotherms and Internal Energies Figure 2 shows adsorption isotherms for three random heterogeneous surfaces with fb ) 0.1, 0.2, and 0.3 The ordinate is the surface coverage, and the abscissa is the logarithm of Kaf. The keys represent the simulation data points; the dotted and the solid lines represent the theoretical isotherms of the first and the second versions of the QCA theory, respectively. The parameter rb, defined by eq 15, represents the relative strength of adsorption on site b compared to that of site a. When rb equals unity, which implies that the surface is homogeneous, the FloryHuggins equation is a good expression for representing the dimer adsorption on the homogeneous surface. When the parameter rb is smaller than approximately 10, the two theories give almost the same curves and are good approximations for representing the simulation data of surface coverage against the fugacity. However, when rb is greater than 100, the first version of the theory (shown by the dotted lines) predicts a lager θ than the simulation data in the low surface coverage region and a smaller θ in the middle region. The crossing point (in θ) is approximately the surface fraction fb. On the other hand, the second version of the theory (solid lines) predicts a smaller θ in the low surface coverage region and a larger θ in the middle region. It is interesting to see that the simulation points would be well represented if we mixed the two theoretical isotherms. Figure 3 shows the negative internal energy of the system against the surface coverage θ for the surface of fb ) 0.3. When parameter rb is large, the θ values obtained for the first version of QCA theory (dotted lines) are larger than the simulation data points, while the θ values obtained for model 2 (solid lines) are usually smaller. Figure 4 shows the results for the surface of fb ) 0.1. The predictions of the two approximations deviate from the simulations in the region of lower surface coverage, and the discrepancy seems to be magnified compared to the previous case, fb ) 0.3.
Figure 2. Adsorption isotherms of dimer on three random heterogeneous surfaces: (...) model 1; (s) model 2.
6. Site Coverages and Bond Coverages Figure 5 shows the comparisons of the site coverages of model 1 (ordinate) with those of the simulations (abscissa) for fb ) 0.3. The site coverages θa and θb in the ordinate are calculated from eqs 15, 17, and 18 by specifying the surface coverage θ, which is obtained from a simulation run, and the corresponding site coverages θa and θb are given for the abscissa; these points should be on the straight line if the theory is satisfactory. In the case of large rb, θb values calculated by model 1 are always larger than the simulation values, while θa values calculated by model 1 are smaller. This is attributable to the fact that model 1 treats each site independently; therefore, a larger number of the active site b is favorable for each segment even when two sites of b are not bonded. In the case where fb ) 0.1, the
Figure 3. Internal energy vs surface coverage for fb ) 0.3: (...) model 1; (s) model 2.
discrepancies for θb between the simulation and model 1 become larger than in the case where fb ) 0.3. Figure 6 shows similar comparisons of bond coverages of model 2 (ordinate) with those of simulations (abscissa) for fb ) 0.3. The ordinate bond coverages θaa, θab, and θbb are calculated from eqs 21, 23, and 24 by specifying the surface coverage θ. At first glance, we see that the values
Monte Carlo Study for Adsorption of Dimers
Langmuir, Vol. 13, No. 5, 1997 907
Figure 4. Internal energy vs surface coverage for fb ) 0.1: (...) model 1; (s) model 2.
Figure 6. Comparisons of the bond coverages θaa, θab, and θbb between simulations and model 2 for fb ) 0.3.
Figure 7. Two types of active bonds b-b: (O) site a, (b) site b; (a) one effective bond b-b, (b) two effective bonds b-b.
Figure 5. Comparisons of the site coverages θa and θb between simulations and model 1 for fb ) 0.3.
of θbb, bond coverages on more active sites, exceed unity in the abscissa and approach unity as the surface coverage θ approaches unity. Since the values of θbb obtained from the simulation are evaluated by 〈Nbb〉/Mbb, the value of Mbb should be responsible for the strange result of the values of θbb exceeding unity. The theoretical value of Mbb defined by eq 1 is the number of bonds b-b per active site since the number of bonds b-b evaluated on the model surface was divided by 4 (the coordination number) as stated in section 2.
However, the effective number of bonds b-b available for dimer adsorption on adjacent active sites b-b is not as simple as presumed in eq 1. Figure 7 shows two types of bonds b-b which have two bonds b-b but differ in the effective number for dimer adsorption; that is, the effective number of bonds b-b is 1 for type a while it is 2 for type b. This is because the occupation of one site inevitably excludes the other three bonds issued from the site from being occupied. The value of Mbb defined by eq 1 has presumed that all bonds b-b were type a even though the simulated random surface had many bonds of type b. Therefore, we will need to have some measure for characterizing all bonds on the surface, which may be called the bond correlation or the lattice connectivity. At present it is open to question. It is noted here that the theoretical number of bonds, Maa, Mab, and Mbb defined by eq 1, is the limiting value, which is appropriate when θ approaches unity. This is because the state of the maximum entropy is realized at full surface coverage since the total energy of the system does not change with changes in the configurations of the dimers when they cover the whole surface. 7. Concluding Remarks Monte Carlo simulations were performed to obtain adsorption isotherms of dimers on random heterogeneous
908
Langmuir, Vol. 13, No. 5, 1997
surfaces consisting of two different energy sites. When the surface is homogeneous, the Flory-Huggins equation is found to be a good expression for representing the adsorption isotherms of dimer molecules. The data obtained with two versions of the quasi-chemical approximation (QCA) approaches to represent the combinatory factor for distributing dimer molecules on heterogeneous lattices are compared to the simulation data. The first version of QCA, which treats each active site
Nitta et al.
independently, is found to overestimate the occupancy of more energetically stable sites, while the second version of QCA, which treates each bond independently, underestimates them. It is noteworthy that the simulated adsorption isotherms might be well represented if we could mix the isotherms predicted from the two versions of QCA. More theoretical work is needed to represent the adsorption isotherms calculated by computer simulations. LA950957T
