Langmuir 1997, 13, 1307-1316
1307
A Novel Two-Parameter Biexponential Distribution To Map Surface Heterogeneity: Genesis of Common Distributions and Isotherms† S. D. Prasad‡ Physical Chemistry Division, National Chemical Laboratory, Pune 411 008, India Received October 31, 1995. In Final Form: January 20, 1997X A generalized two-parameter biexponential distribution (BIED) is suggested to map surface heterogeneity which in turn can be visualized as a linear combination of three popular site energy distributions, viz., the negative exponential, positive exponential, and constant. Taking advantage of the linearity property of the Mellin transform operator, closed form expressions for isotherms and their pressure derivatives are derived. Depending on the values of R and β parameters, the relative weights given to each of the component distributions are different, and a wide spectrum of isotherm and isosteric heat behavior is possible. When R and β parameters are equal, the BIED is symmetrical and the isosteric heat vs surface coverage shows a linear fall resembling the constant distribution. A method based on the product of the surface coverage and pressure derivative is used to discriminate between the BIED and constant distributions. The symmetrical BIED also can approximate the Gaussian distribution and when R and β are unequal, even the skewed Gaussian distribution. This is true only as far as isotherms are concerned, but their isosteric enthalpy of adsorption vs surface coverage show markedly different behavior. A limited applicability of the BIED to the experimental system/Ar/Rutile and the data of Drain and Morrison is shown. It appears that the data can be adequately fitted by the Skewed-Gauss distribution, as has been originally done by Sokolowski et al. or by the BIED. The local isotherm employed is the Langmuir model.
I. Introduction The quantitative description of surface heterogeneity is a problem, which has not been totally solved until recently.1-4 Attempts have been made in the past to postulate a priori theoretical forms for the site-energy density distribution, such as the negative exponential,5,6 Gaussian,7 skewed Gaussian,8 constant,9,10 Gamma,11 positive exponential,12 sinusoidal,13 Maxwell-Boltzmann,14 log normal distribution,15 etc. The agreement of the experimental isotherm or calorimetric differential heat with that computed from the chosen prior distribution is thought to validate the distribution.1-4 Alternatively, methods involving inverting the adsorption integral † 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. ‡ E-mail:
[email protected]. X Abstract published in Advance ACS Abstracts, February 15, 1997.
(1) Steele, W. A. The Interaction of Gases with Solid Surfaces; Pergamon Press: Oxford, 1974. (2) Rudzinski, W.; Everett, D. H. Adsorption of Gases on Heterogeneous Surfaces; Academic Press: New York, 1992. (3) Clark, A. The Theory of Adsorption and Catalysis; AcademicPress: New York, 1970. (4) Jaroniec, M.; Madey, R. Physical Adsorption on Heterogenous Solids; Elsevier: New York, 1988. (5) Halsey, G. D.; Taylor, H. S. J. Chem. Phys. 1976, 64, 1762. (6) Prasad, S. D.; Doraiswamy, L. K. Phys. Lett. 1983, 94A, 219. (7) Ross, S.; Olivier, J. P. On Physical Adsorption, Wiley (Interscience): New York, 1964. Sanford, C.; Ross, S. J. Phys. Chem. 1954, 58, 288. (8) Rudzinski, W.; Jaroniec, M.; Sokolowski, S. Phys. Lett. 1974, 48A, 171. Dubinin, M. M.; Radushkevich, L. V. Dokl. Akad. Nauk SSSR, 1947, 55, 331. Stoeckli, H. F. J. Colloid Interface Sci. 1977, 59, 184. (9) Temkin, M.; Levich, V. Zh. Fiz. Khim. 1941, 14, 296. Temkin, M.; Levich, V. Zh. Fiz. Khim. 1946, 20, 1441. (10) Prasad, S. D.; Doraiswamy, L. K. Chem. Phys. Lett. 1983, 97, 31. (11) Smutek, M. Surf. Sci. 1975, 52, 465. (12) Prasad, S. D.; Doraiswamy, L. K. Chem. Phys. Lett. 1983, 99, 129. Temkin, M. I. Kinet. Catal. USSR 1967, 5, 1005. (13) Hill, T. L. J. Chem. Phys. 1949, 17, 762. (14) Kindl, B.; Wojciechowski, B. W. J. Chem. Soc., Faraday Trans. 1 1973, 69, 1926. Hsu, C. C.; Rudzinski, W.; Wojciechowski, B. W. J. Chem. Soc., Faraday Trans. 1 1976, 72, 453. (15) Hoory, S. E.; Prausnitz, J. M. Surface Sci. 1967, 6, 377.
S0743-7463(95)00544-0 CCC: $14.00
equation for the overall isotherm through integral transform methods,16-19 have also been suggested, after postulating a suitable local isotherm (most often of the Langmuir type).1-4 A main difficulty which is commonly encountered is that many diverse functional forms for the energy-distribution integrate out to give roughly the same overall isotherm, within the experimental error. Clearly one should seek something more than mere numerical agreement with the experimental data. However, a few distribution functions are more popular than the others.1-4 Considerable experimental evidence exists (in terms of explaining the differential heat vs surface coverage and isotherms) for the validity of the three-site energy distributions, viz., (a) the positive exponential12,20-23 distribution, (b) the negative exponential distribution,4,5, 20-23 and (c) the constant distribution.9,21 These give rise to (a) the Temkin negative power isotherm,12 (b) the Freundlich isotherm4,5 (also known as the power-law isotherm), and (c) the Temkin-Pyzhev isotherm.4,5 In addition, the Dubinin-Radushkevich isotherm, which has remarkable success in analyzing microporous adsorbents, corresponds to a skewed-Gaussian distribution.8,21,23 There is a school of thought that the Gaussian7 distribution is perhaps the most fundamental (as a result of the central limit theorem which states that a number of independent distributions added on will give rise to a Gaussian).24 However, many of the experimental siteenergy distributions clearly point to an asymmetry (16) Prasad, S. D.; Doraiswamy, L. K. Phys. Lett. 1977, 60A, 11. (17) Prasad, S. D.; Doraiswamy, L. K. Chem. Phys. Lett. 1984, 108 (1), 83. (18) Landman, U.; Montroll, E. W. J. Chem. Phys. 1976, 64, 1762. Jagiello, J. Langmuir 1994, 10 (8), 2778. (19) Sips, R. J. Chem. Phys. 1950, 18, 1024. (20) Rudnitsky, L. A.; Alexeyev, A. M. J. Catal. 1975, 37, 232. (21) Bhat, Y. S.; Prasad, S. D.; Doraiswamy, L. K. J. Catal. 1984, 87 (1), 10. (22) Bhat, Y. S.; Prasad, S. D.; Doraiswamy, L. K. AIChE J. 1985, 31, 1682. (23) Bhat, Y. S.; Prasad, S. D.; Doraiswamy, L. K. Stud. Surf. Sci. Catal. 1984, 19, 419. (24) Ash, R. A. Basic Probability Theory; John Wiley & Sons: New York, 1970; p 171.
© 1997 American Chemical Society
1308 Langmuir, Vol. 13, No. 5, 1997
(skewed) when plotted against the heat of adsorption, implying the validity of the skewed-Gaussian distribution (Dubinin-Radushkevich isotherm).8 The foregoing discussion clearly shows the need for finding a generalized distribution which will give rise to the genesis of the most known experimental isotherms and site-energy distributions. In addition it should be mathematically tractable; i.e., it is desirable to have closed form expressions for the adsorption isotherms, pressure derivatives, differential heat, etc. We propose a generalized two-parameter (a, β) biexponential distribution to achieve this objective. This exponential distribution will be shown to be a linear combination of three distributions, viz., (a) the negative exponential distribution, (b) the positive exponential distribution, and (c) the constant distribution. The generalized isotherm in turn will be a linear combination of three popular isotherms, i.e., (a) the Freundlich isotherm, (b) the Temkin negative power law isotherm, and (c) the Temkin-Pyzhev isotherm. In addition, when the a and β parameters of the distribution are chosen to be equal (symmetric case), itwill be shown to approximate the Gaussian distribution, provided that the full width at half maximum (fwhm) is chosen appropriately. Lastly, when a and β are chosen to be not equal, then the distribution can fairly approximate the behavior of the skewed-Gaussian distribution. II. Formulation The most interesting scenario in heterogeneous catalysis corresponds to the case when the limits (upper and lower) of the heats of adsorption are finite, in contrast to the (0,∞) limits chosen often for mathematical convenience.5,19 This mainly arises from the Balandin Volcanoe principle,3 which clearly shows that there is an optimal distribution of heats for the catalytically most active sites.3 With infinite upper limit, the integral equation becomes singular and the Stieltjes transform method can be used for inverting the integral equation.19 The case of finite limits is more difficult but also has been solved completely.17 In our present study the lower limit value is 33.47 kJ mol-1 (8 kcal‚mol-1). The upper limit of the heat of adsorption is 133.88 kJ‚mol-1 (32 kcal‚mol-1), which is a representative value in many chemisorption systems.3,21-23,25,30 In the previous studies10,12,22 the limits Q1 ) 17 kcal‚mol-1 (71 kJ‚mol-1) and Q2 ) 31 kcal‚mol-1 (129.6 kJ‚mol-1) were employed, which were representative of the CO/Ni system.25-27 Even for this well-studied system, wide variability in the isosteric heat vs coverage has been observed experimentally, illustrating the difficulties involved in the choice of Q1 and Q227 in any theoretical computation of the isotherms. The analysis employs the random patch model (RPM) of the heterogeneous surface, which has enjoyed great popularity, because of its ready applicability and which at the same time is mathematically tractable. As is customary, each of these patches constitutes a separate thermodynamic entity, and adsorption equilibrium prevails between the gas and surface phases on an individual patch. The other model, which is equally successful, is (25) Broker, F. J.; Wedler, G. Discuss. Faraday Soc. 1966, 41, 95. Wedler, G.; Schroll, G. Z. Phys. Chem. (Frankfurt) 1973, 85, 216. (26) Rideal, E.; Sweet, F. Proc. Roy. Soc. London 1960, A257, 291. (27) Wedler, G. Chemisorption: An experimental Approach; Butterworths: London, 1976; p 39. Brenman, D.; Graham, M. J. Discuss.Faraday Soc. 1966, 41, 95. (28) Beeck, O. Adv. Catal. 1951, 2, 151. Beeck, O. Discuss. Faraday Soc. 1950, 8, 118. (29) Castright, R. D.; Dumesic, J. A. J. Catal. 1994, 148, 771. (30) Somorjai, G. A. Chemistry in Two Dimensions: Surface; Cornell University Press: Ithaca, NY, 1981.
Prasad
the random-site model developed by Hill,13 Steele,31 Pieorotti,32 Rudzinski,33 and co-workers which only considers total random distribution of adsorption sites ignoring all spatial correlations. Interaction effects between adsorbed molecules is neglected, and the Langmuir isotherm (described often as a local isotherm) describes the adsorption on an individual patch. It is further assumed that with the variation of temperature, site-energy distribution is basically unaffected. This is certainly not true, if we keep in mind that thermal roughening can be significant in many model systems. With these assumptions, the overall adsorption isotherm is given by4-6
θ (T,p) ) θli(T,p) )
∫QQ
2
1
θli(T,p) δ(Q) dQ
1 1 + (b0/p) exp(-Q/RT)
(1a) (1b)
where θli(T,p) denotes the Langmuir local isotherm, δ(Q) is the site-energy density distribution, Q is the heat of adsorption, and Q1 and Q2 are the limits of the heat of adsorption. b0 is the entropy factor corresponding to the standard state of the adsorbed phase at half coverage. Equation 1 is the most fundamental one in the RPM description of adsorption isotherms. The biexponential distribution δ(Q) is defined as
δ(Q) ) CN [exp(RQ2) - exp(RQ)][exp(-βQ1) exp(-βQ)] (2) where CN is the normalization constant determined by
∫QQ
CN
2
1
δ(Q) dQ ) 1
(3)
with the following definitions: exp(RQ2) ) a2; exp(RQ1) ) a1; exp(-βQ1) ) b1; exp(-βQ2) ) b2 The distribution defined by eq 2 can be simplified and becomes
δ(Q) ) CN[a2b1 - a2 exp(-βQ) - b1 exp(RQ) + exp(R - β)Q] (4) and the normalization constant CN equals
CN ) 1/∆
(5a)
where
∆ ) a2b1(Q2 - Q1) - (a2/β)(b1 - b2) (b1/R)(a2 - a1) + [(1/(R - β))](a2b2 - a2b1) (5b) R*β Equation 4 can be thought of as a linear combination of the three distributions. The first, second, and third terms in parentheses denote the constant, the negative exponential, and the positive exponential distributions, respectively. The last (fourth) term becomes a positive exponential distribution or negative distribution accordingly as R > β or R < β in that order. The parameters R and β appearing in the BIED have a physical significance. They are the pressure exponents (31) Steele, W. A. J. Phys. Chem. 1963, 67, 2016. (32) Pierotti, R. A.; Thomas, H. E. Trans. Faraday Soc. 1973, 70, 125. (33) Rudzinski, W.; Lajtar, L. J. Chem. Soc., Faraday Trans. 2 1981, 77, 153. Rudzinski, W.; Lajtar, L.; Patrykiejew, A. Surf. Sci. 1977, 67, 195. Rudzinski, W.; Patrykiejew, A.; Lajtar, L. Surf. Sci. 1978, 77, L655.
Mapping Surface Heterogeneity
Langmuir, Vol. 13, No. 5, 1997 1309
experimentally observed, when the generalized isotherm corresponding to the BIED reduces to the Temkin negative power isotherm and Freundlich isotherm in the appropriate range of pressures. As will be shown later, the weights given to each of these components distributions, viz., positive-exponential or negative exponential, depends on the relative values of R and β. The pressure exponents are seldom greater than unity, i.e., 0 < R < 1, 0 < β < 1. For reduction of the number of variables by 1, many simulations have been carried out with R + β ) 1. For this case, the BIED is essentially a one-parameter distribution (psuedo two parameter). A considerable spectrum of isotherm behavior is predicted for this case. Since the probability density distribution δ(Q) has to integrate to unity, the relative weights of terms, corresponding to each of these distributions can be found very simply. But before doing that, we define the absolute total weight as
|∆| ) |a2b1(Q2 - Q1)| + |a2(b1 - b2)/β| + |b1(a2 - a1)/R| + |(a2b2 - a1b1)/(R - β)| (6a) and the relative weights
wpos - wneg ) |exp{(R - β)Q2}/β - exp(R - β)Q1/R + exp(RQ2 - βQ1)(1/R - 1/β)|/∆| (6h) In conclusion, for R and β, the weight of the positive exponential distribution is higher than that of the negative exponential distribution. It is desirable to have the right range of values of |R| < 1 , |β| < 1, so that the experimentally observed pressure exponents of the isotherm agree, when the BIED is defined. Alternatively, experimentally observed pressure exponents β and R can be used to define a unique BIED. In addition, we should know the isosteric heat variation with surface coverage to parametrize the heat. It will be shown that this in turn can reveal more information about R and β values. It is also to be brought in mind that wmix gives a weight for another positive exponential distribution when R > β
wpos ) |(a2b2 - a1b1)/(R - β)|/|∆|
(6i)
and when R < β
wcons ) |a2b1(Q2 - Q1)|/|∆|
(6b)
wneg ) |(a2b2 - a1b1)/(R - β)|/|∆|
wneg ) |a2(b1 - b2)/β|/|∆|
(6c)
wpos ) |b1(a2 - a1)/R|/|∆|
(6d)
At first glance it may seem that the BIED is analogous to the exponential distribution with a maximum,20,34 as proposed by Rudnitsky and Alexeyev, viz.20,34
Dividing eqs 6c and 6d by 6b, we get the following relations for the ratio of the relative weights
wneg/wcons ) [1 - exp{β(Q1 - Q2)}]/β(Q1 - Q2) (6e) wpos/wcons ) [1 - exp{R(Q1 - Q2)}]/R(Q2 - Q1) (6f) when R f 1, β f 1, the exponential term in both these equations become vanishingly small, and hence wcons > wpos, wneg. On the other hand when R f 0, β f 0, expanding the exponentials in eqs 6c and 6f and retaining only the first-order terms, we get
wpos = wneg = wcons i.e., all the component distributions have equal weights in the BIED and are important. Similarly the ratio of wmix to wcons can be written as
wmix/wcons ) {exp(-β(Q2 - Q1)) exp(-R(Q2 - Q1))}/[(R - β)(Q2 - Q1)] (6g) when R f 1, and β f 1, it is obvious that this mixed term is not of much significance. When R ) β, the mix term is not defined, and this simply adds another term to the weight of the constant distribution
when
Finally some comments are in order as to the relative importance of wpos and wneg. By starting with eqs 6c and 6d, we can easily deduce the following result.
R f 0, β f 0
by expanding the exponents in eq 6g and retaining only first-order terms, it can be shown that
wmax/wcons ) 1, R * β Here again all the terms become equally important. In conclusion, for most of the R and β values of real significance, the weight of the constant distribution term wcons predominates.
δ(Q) ) exp(-r1Q) δ(Q) ) exp(r2Q)
(6j)
Q < Qc Q > Qc
r1 and r2 are the parameters defining the distribution and Qc denotes characteristic parameter of the distribution. At one time, either the positive exponential or negative exponential type distribution dominates. In contrast to these, in BIED, both these constituent distributions contribute appreciably at all times. Besides, the constant distribution is not given any importance at all in their work. When R ) β, the BIED can resemble a Gaussian distribution, especially for low values of R and β, which is not possible with their exponential distribution.20,34 In parts a and b of Figure 1 the normalized distributions are again grouped into two sets, depending on the relative values of R and β. In the former set R < β, whereas in the latter set R > β. Notice skewness of the plots and the locations of their means (first absoute moments). When R ) β exactly, the distribution becomes symmetric and eqs 5 and 6 need modification. In Figure 1c the symmetrical case R = β is considered. The limiting regions of behavior for R ) β can be expressed as two types.
case 1: R ) β f 0 After expanding the exponents in eq 2 and keeping only the first powers in Q, the distribution becomes parabolic
δ(Q) ) CNRβ(Q2 - Q) (Q - Q1)
(7)
The distribution is again symmetric with respect to Qmax. The curves corresponding to the quadratic distribution, namely, 1 and 3 in Figure 1c can approximate the Gaussian distribution provided that the fwhm of the biexponential distribution equals that of the Gaussian. The fwhm of (34) Reference 2, pp 121-124.
1310 Langmuir, Vol. 13, No. 5, 1997
Prasad
a
b Figure 2. Plot illustrating the fitting of the Gaussian trial function for the symetric BIED with R ) β.
c
Figure 3. Plot of the weights of the constituent constant, positive and negative exponential distributions comprising the biexponential distribution plotted against R; β ) 1 - R. Notice the insignificant weight given to the mix distribution.
Figure 1. Plot of the normalized site-energy density function BIED (biexponential distribution function) vs the heat of adsorption Q; (a) R < β; (b) R > β; (c) R ) β.
the Gaussian is given by the well-known expression 0.5
fwhm ) 2(2 ln 2) σ
(8)
the probability density function becomes rectangular. This is illustrated in curve 5 of Figure 1c. By differentiating δ(Q) (eq 2) with respect to Q and setting it equal to zero, and keeping in mind that R ) β, we easily obtain
Qmax ) (Q1 + Q2)/2
(10)
Keeping in mind the definition of the biexponential distribution, eq 2 can be rewritten as
δ(Q) ) CN[λ - exp(-R)] [λ - exp(R)]
(11)
where σ denotes the variance. A reference to Figure 2 indeed shows that the Gaussian distribution can also be approximated to a satisfactory degree by the biexponential. But the inflection point, which is present in the former is totally absent in the latter. The other possibility is
where, is a dummy variable to show mathematically the symmetry of the δ(Q) distribution
case 2: R ) β f 1
and
After a few elementary manipulations, we transform δ(Q) to
δ(Q) ) CN exp[R(Q2 - Q1)][1 exp[R(Q - Q2)]][1 - exp(R(Q1 - Q))] (9) If e-5 ) 0, for Q1 + 5 < Q < Q2 - 5, the terms inside the large braces can be put equal to unity. For these Q values,
λ ) exp(Q2 - Q1)/2
(12)
Q ) Qm +
(13a)
δ(Qm + ) ) δ(Qm - )
(13b)
Since eq 13b is obvious, by virtue of eq 11, the symmetry of the site energy distribution with respect to Qm is proved. In Figure 3 the weights for each of the constituent distributions of the biexponential distribution defined by
Mapping Surface Heterogeneity
Langmuir, Vol. 13, No. 5, 1997 1311
Table 1. Table of Isotherms for the Biexponential Distribution Expressed as a Linear Combination of the Positive and Negative Exponentials δ(Q) ) ω1 exp(R′Q/RT) + ω2 exp(-β′Q/RT); ω1 ) -b1; ω2 ) -a2
θ ) Fpos + Fnega
Fpos
Fneg p < p2 [RTω2/(β′ - 1)] (p/b0)β′ F(1,1-β′;2-β′;-x1-1)x1β′-1 x2β′-1F(1,1-β′;2-β′;-x2-1)
[RTω1/(R′ + 1)](b0/p)R′ [x2-(1+R′)F(1,R′+1;R′+2;-x2-1) x1-(1+R′)F(1,R′+1;R′+2;-x1-1)] [RTω1/R′](b0/p)R′ [x2-R′F(1,-R′;1-R′;-x2) - R′π/sin(R′π) x1-(1+R′)(R′/1 + R′)F(1,R′+1;R′+2;-x1-1)]
p2 < p < p1 (RTω2/β′) (p/b0)β′ [β′π/sin(β′π) + β′/(β′ - 1)x1β′-1F(1,1-β′;2-β′;-x1-1)-x2β′F(1,β′;β′+1;-x2)] p > p1 (p/b0)β′(RTω2/β′)F(1,β′;β′+1;-x1)x1β′ - F(1,β′;β′+1;-x2)x2β′]
[RTω1/R′](b0/p)R′ [F(1,-R′;1-R′;-x2)x2-R′ - F(1,-R′;1-R′;-x1)x1-R′] θli ) 1/(1 + b0/p exp(-Q/RT); a
x1 ) b0 exp(-Q1/RT)/p;
x2 ) b0 exp(-Q2/RT);
R′ ) RRT;
p1 ) b0 exp(-Q1/RT); p2 ) b0 exp(-Q2/RT)
The constant and mix experssions(Fcons, Fmix) are not shown. F(a,b;c;z) is a hypergeometric function of z.
eqs 6a-f are plotted. As is obvious, the constant distribution has the largest weight. The plot is made as a function of the a parameter and the β parameter is simply β ) 1 - R. As R increases, the constant distribution reaches its maximum weight at R ) 0.5. On the other hand, the weights of the negative and the positive exponential distribution components increase and decrease monotonically. At all times the contribution of the constant distribution component remains predominant. Later it will be shown from the analysis of the isosteric heat vs surface coverage plots that a good percentage of the isotherms generated computationally also show linear fall with the surface coverage, like in the case of “pure” constant distribution. But this rule is not “universally” valid as will be shown later, implying that the weights of the other “components” positive exponential to negative exponential are also important. The weight of the mix distribution is, however, minimal, as is easily seen from the plot. III. Derivation of the Generalized Isotherm The Mellin transform method for the exact derivation of adsorption isotherms for the positive exponential distribution and negative exponential distribution with finite limits of heats of adsorption has been employed extensively by the author and co-workers in a number of previous publications.10,22-23 Therefore, we present only the minimum of details, which are given in the Appendix. Since the Mellin transform is a linear operator, the observation that the biexponential distribution can be split into a linear combination of positive exponential and negative exponential forms in Q (see eq 4) is very useful. The Mellin transform can be applied to each of these functions in turn.22-23 The isotherm expressions derived for the positive exponential and the negative exponential distributions along with the relevant variables are presented in Table 1. It is a trivial exercise to compute the isotherm expression for the constant distribution, so we omit the details. The various pressure ranges mapped in Table 1 are necessary because the hypergeometric series converges only when the absolute value of the argument is less than unity, and for larger values, techniques of analytical continuation are necessary (see Appendix).35,10,22-23 If Fpos, Fneg, and Fconst denote the corresponding functions for the positive exponential, negative exponential, and constant distributions, then the overall isotherm can (35) By pure negative and positive exponential distributions, we mean general functional forms of the type δ(Q) ) CN exp(-βQ) and δ(Q) ) CP exp(RQ), respectively, CN and CP being the normalization constants.
be written as
θ(T,p) ) CN[a2b1Fcons - a2Fneg - b1Fpos + Fmix] (14) The last term, Fmix, will be a positive exponential or negative exponential depending on the relative values of R and β; i.e., R > β or R < β give the possible outcomes. The exponent of this last fourth term, Fmix, will of course be different from that of Fpos and Fneg. IV. The Limiting Forms of the Generalized Isotherm If we recall the definitions of p1 and p2 in Table 1, the most interesting range corresponds to p2 < p < p1. From our computations p1 ∼ 108 exp(-8) ∼ 3.35 × 104 Torr ∼4.472 MPa at RT ) 1 kcal‚mol-1, p2 ∼ 108 exp(-32) Torr ∼ 1.266 × 10-6 ) 0.1688 mPa. It is precisely in this range that the preponderance of power law isotherms occur. The isotherm corresponding to the negative exponential distribution can be very well approximated by10,22-23
θ(T,p) ∼ A1 + A2pβ′
(15)
whereas the positive exponential distribution gives10,22-23
θ(T,p) ) B1 - B2p-R′
(16)
the R′ and β′ exponents depend on the values of R and β chosen to define the distribution. Each of these isotherms has been extensively employed in the literature. The former is the well-known Fruendlich isotherm and the latter is the Temkin negative power law isotherm. The logarithmic form given by the constant distribution eq 15. (Temkin-Pyzhev isotherm), i.e.
(
θ ) RT/(Q2 - Q1) ln
)
b0/p + exp(Q2/RT)
b0/p + exp(Q1/RT)
(17)
can also be reduced to a power-law isotherm in a very narrow range of pressure. The overall isotherm will be the sum of power law relations in p, as well as the logarithmic terms given by eq 17. But in actual practice, the power law form is still a very good approximation as is evident from Figure 4 where in three sets of isotherms are plotted with various values of R and β. We have considered the exact (in a computational sense) isotherms of the biexponential distribution and also considered the power law approximation. The power law expression provides such an excellent fit, the curves are identical and are not shown separately. This is rather
1312 Langmuir, Vol. 13, No. 5, 1997
Prasad
Figure 4. Biexponential distribution function (BIED) isotherm plots for three valves of R; β ) 1 - R. The power-law approximation fits are so good that they are indistinguishable.
Figure 5. Plot of the isosteric heat vs surface coverage θ for the biexponential distribution for various R.
surprising as the power-law form given by eq 13 is true only for a “pure” negative exponential distribution.35 The “pure” positive exponential distribution should give rise to a power law isotherm with a negative pressure exponent (see eq 14). The weights of each of these component distributions are nonzero, in the actual biexponential distribution. But it appears that only the positive pressure exponent of the power law approximation dominates. It may be argued that the “fit” may be valid only from empirical considerations; hence reading too much into the sign of the pressure exponents can be misleading. However, the foregoing discussion may be considered as illustrating the ambiguities involved in identifying a unique site-energy distribution from the experimental isotherm data. V. The Isosteric Heat Behavior An expression for the overall isosteric heat can be derived from eq 1 by using the Clasius-Clapeyron equation1-3,36-38
(∂ ln p/∂T)θ ) Qiso/RT2
(18)
Hence we maintain constancy of the overall mean surface coverage θ and not the local surface coverage θli while differentiating. Using eq 18 and also the following relation for each patch
(∂ ln p/∂T)θli ) -1/p(∂θli/∂T)p/(∂θli/∂p)T
(19)
We get the following expression, after averaging over δ(Q)
Qiso )
∫QQ θli(1 - θli) Q δ(Q) dQ/ ∫QQ θli(1 - θli) δ(Q) dQ 2
1
2
1
(20)
The Qiso refers to the overall isosteric heat of the heterogeneous surface and Q in the numerator of eq 20 refers to the isosteric heat of the local ith patch; θli is as before the local surface coverage (isotherm). In Figure 5 we have plotted the isosteric heat behavior of the biexponential distribution as a function of the surface coverage θ. The three curves described correspond to the cases (1) R > β; (2) R ) β, and (3) R < β, respectively. In (36) Everett, D. H. Proc. Chem. Soc. 1957, 38. (37) Everett, D. H. Trans. Faraday Soc. 1950, 46, 453, 942, 957. (38) Hill, T. L. J. Chem. Phys. 1950, 18, 246.
Figure 6. Isosteric heat vs surface coverage for the “pure” normalized positive exponential distribution and the negative exponentail distribution with R ) β ) 0.5. Notice that the fit of the isosteric heat to eq 22 is indistinguishable for the latter.
the case of (1) the curve is concave, while for (3) it is convex, case (2) displays a near linear plot corresponding to almost zero curvature. This is very similar to the behavior of the constant distribution, which almost for all ranges of surface coverage predict a linear fall of the isosteric heat with the surface coverage.
Qiso ) Q° - bθ
(21)
where Q° denotes the initial heat at zero coverage and θ the gross (overall) surface coverage. Reference to Figure 5 shows that this is indeed true for the biexponential distribution; also discrimination between a constant distribution and the biexponential distribution (when R ) β) is really hard, illustrating the difficulties in identifying a unique site-energy distribution. Curves 1 and 3 however show appreciable curvature. In the former R > β, so the constituent negative exponential distribution dominates, while for the latter, the positive exponential part predominates. In Figure 6 the isosteric heat vs surface coverage is plotted for a pure “positive exponential distribution” with the same R value and a negative exponential distribution with the same β value.25 These pure distributions are separately normalized.35 The biexponential distribution has both R and β parameters and corresponding weights for each of these pure distributions. Comparison of the positive exponential distribution with the biexponential
Mapping Surface Heterogeneity
Langmuir, Vol. 13, No. 5, 1997 1313
Figure 7. Plot of the isosteric heat vs surface coverage for the BIED for a wide range of R, β, and θ.
distribution (Figure 5) shows that the isosteric heat behavior of the former is appreciably different. A reference to Figure 6 also shows that the isosteric heat behavior of the pure negative exponential can be neatly fitted with an equation of the following form, as is well-known
Qiso ) -a log(bθ)
Figure 8. Plot of the pressure derivative product vs surface coverage for the range of small surface coverage R ) 0.1-0.4 for the biexponential distribution. The product is on a log scale.
(22)
For the data presented in Figure 6, a ) 4.567 and b ) 2.902 × 10.-2 Unfortunately for the biexponential distribution, no such simple parametrization of the isosteric heat surface coverage curve is possible. In Figure 7 isosteric heat vs θ plot is made for an extended range of pressure. The plot is qualitatively similar to that presented in Figure 5 (curves 1 and 2). A simple parametrization of the isosteric heat θ curves does not seem to be possible. In sharp contrast to this, when R ) β ) 0.5 the isosteric heat falls as a neat linear function of surface coverage, as illustrated in Figure 7 (cf. the midcurve 2 in Figure 5) even over a very extended pressure range. The constant distribution is the only other well-known distribution, which shows this behavior. Thus it is difficult to distinguish between the constant and biexponential distributions on the basis of the differential heat behavior alone. In the next section, we investigate the pressure derivative method of discriminating between rival distributions.
Figure 9. Same as Figure 8 but for the range R ) 0.5-0.9. The product is on a log scale.
VI. Pressure Derivative Method The pressure derivative method has been suggested earlier to discriminate between rival distributions. Since the constant distribution and the biexponential distribution both predict a linear fall of the differential heat with the surface coverage, an alternative method involving the pressure derivative product,39-41 which has been employed earlier,39-41 is suggested. The pressure derivative product is defined as
d 2 (θ /2) ) θ dθ/dp dp
(23)
where θ denotes the mean surface coverage θ. Since the θ increases and dθ/dp decreases with pressure, one would expect the product to display a maximum, when plotted as a function of pressure or surface coverage. The plot of (39) Prasad, S. D.; Doraiswamy, L. K. Chem. Phys. Lett. 1984, 104 (4), 315. (40) Datar, A. S.; Prasad, S. D.; Doraiswamy, L. K. Chem. Phys. Lett. 1989, 159, 337. (41) Datar, A. S.; Prasad, S. D. Langmuir 1991, 7, 1310.
Figure 10. Comparison of the pressure derivative product for the constant and the biexponential distributions: curve 1, R ) β; curve 3, R ) β ) 0.001.
dθ/dp against surface coverage is, however, more informative. The location of the maximum, θmax in the coverage axis is shown to lead to useful information.39-41 In Figures 8 and 9 θ dθ/dp are plotted against θ for two sets of R and β for the biexponential distribution. In Figure 10 the pressure derivative products of the constant and the biexponential distribution are plotted against the
1314 Langmuir, Vol. 13, No. 5, 1997
Prasad
Figure 11. Plot of the pressure derivative product for the Gaussian distribution.
Since the Gaussian is a popular distribution, it is very interesting to see how far this distribution can approximate the behavior of the biexponential distribution. The fact that the numerical fit is good has been already shown (see Figure 2). It is also of interest to study the behavior of the pressure derivative product, which is shown in Figure 11. Since the mean of the symmetrical biexponential distribution (BIED) coincides with that of the Gaussian, the only parameter of interest is the σ. The plots are shown for two values of σ, viz., σ ) 4, σ ) 5. For smaller values of σ the highest isosteric heat (found from extrapolation to θ f 0) is about Q = 24 kcal‚mol-1, which is unphysical. Notice the sharp dependence of the pressure derivative product on the value of σ, the variance of Gauss distribution, in Figure 11. Thus the pressure derivative product of the Gaussian is indeed different. Similar conclusions that apply to the skewed-Gaussian approximation also apply to the unsymmetrical BIED (which occurs when R * β) and hence are not shown separately. VIII. Application of Bied to Experimental Data
surface coverage. The pressure derivative product of the BIED is smaller than that of the constant distribution. For best comparison purposes the case R ) β is considered. When R ) β ) 1, the weight of the constituent constant distribution (wcons) is largest in BIED. Even for such a case, the pressure derivative product of the BIED is several-fold smaller than that of the constant distribution as is evident from Figure 10. When R ) β ) 0.001, the BIED reduces to the parabolic distribution. But even for such a case, the pressure derivative product is much smaller than that of the constant distribution (see Figure 10). Another fact which enables us to distinguish between the BIED and the constant distribution is that for the constant distribution, the location of θmax is given by
An attempt is made to fit the BIED to physical adsorption data of Drain and Morrison,43 for the system Ar/rutile at 85 K. Even though, for these data, the lower and upper limits are 1.54 and 4.2 kcal‚mol-1 (as employed by Sokolowski et al.42), it is doubtful whether the Langmuir isotherm is very rigorous. This is because a cursory examination of Drain and Morrison’s data reveals that the data double backs on itself, and hence a BET type local isotherm is perhaps more appropriate, rather than the concave Langmuir isotherm. Neverthless Sokolowski et al.42 has fitted their data satisfactorily with a Langmuir type local isotherm, along with a skewed-Gaussian siteenergy distribution. It may look very strange that the hypothetical limits employed in the work, Q1 ) 8 kcal‚mol-1 and Q2 ) 32 kcal‚mol-1 are far removed from those mentioned above. But our focus is on chemisorption, and the temperature of 503 K (RT ) 1 kcal‚mol-1) is much higher than the 85 K employed in physical adsorption. If we compute Q/RT, this works out to be 9.76 and 24.7 for Drain’s data in contrast to 8 and 32 for our work, so these hypothetical limits are not entirely unrealistic, especially when scaled by the temperature of the experiment. In the literature, fine grained chemisorption data are not available for heterogeneous surfaces, as commomly as in the case of physisorption. In Figure 12, Drain and Morrison’s data43 are fitted to the skewed-Gauss distribution and the BIED. By making use of the condensation approximation, Sokolowski et al.42 had already shown that the overall isotherm will have a Dubinin-Radhuskevich form with a corresponding skewed-Gaussian distribution as the site-energy distribution
θmax ) 1/(Q2 - Q1)
θ(T,p) ) exp[-B(RT ln(p/p0))2]
(25)
δ(Q) ) 2B(Q - Q1) exp[-B(Q - Q1)2]
(26)
p0 ) b0 exp(-Q1/RT)
(27)
Figure 12. Plot of Drain and Morrison’s data fitted to the skewed-Gaussian and BIED. The experimental data are for the system Ar/rutile at 85 K.
(24)
i.e., on the spread of the heats of adsorption. For the BIED, this θmax depends in addition on the(R and β) values. In conclusion the pressure-derivative product constitutes a useful procedure to distinguish between the two distributions. VII. A Gaussian Approximation to the Biexponential Distribution When R ) β the biexponential distribution is a symmetrical one when plotted against heat (see Figure 1).
Here B is the characteristic parameter of the distribution with Q1 the lower limit of the heat of adsorption, and b0 (42) Sokolowski, S.; Jaroniec, M.; Cerofolini, G. F. Surface. Sci. 1975, 47, 429. (43) Drain, L. E.; Morrison, J. A. Trans. Faraday Soc. 1952, 48, 840.
Mapping Surface Heterogeneity
Langmuir, Vol. 13, No. 5, 1997 1315
the Hobson constant, and p0 the saturated vapor pressure. The b0 can also be computed using Hobson’s relation
b0 ) (1.76 × 104)(MT)0.5
(28)
b0 ) 1.02 × 106 in our studies as well as in the previous work; M is the molecular weight and T is the temperature.42 Our best values for fitting the skewed-Gauss distribution are B ) 1.569 × 10-6 cal-2 mol2, Q1 ) 1.54 kcal‚mol-1, Q2 ) 3.6 kcal‚mol-1, saturation capacity Va ) 719.4 cm3, p0 ) 112 Torr. The parameters for the BIED were R ) 10-2 and β ) 4.594 with identical limits of the heats of adsorption. The isotherms for the BIED were computed using the analytical formulas given in the Table 1, along with an appropriate normalization constant. The resulting plots are displayed in Figure 12. The fit of the skewed-Gauss distribution is best seen in the lower regions of pressure, but the BIED also seems to fit the data equally well. Some concluding remarks are in order. The fitting of an isotherm alone does not guarantee that the site-energy distribution is correct; further investigations of differential heats and specific heats are necessary to validate it further. Interaction effects are neglected;41 the Langmuir local isotherm may fail to decribe adsorption in many real situations. These are fairly obvious limitations. In spite of all this, it is desirable to look for an isotherm which has only maximal two parameters and which can also probably explain the genesis of many popular isotherms such as the Temkin isotherm, Freundlich isotherm, etc. IX. Conclusions A novel two-parameter (R,β) biexponential distribution (BIED) is proposed to quantitatively map surface heterogeneity. This BIED when expanded becomes a linear combination of three distributions: (a) negative exponential distribution; (b) positive exponential distribution; (c) constant. Since these distributions give rise to common isotherms and also explain the observed differential heat surface coverage behavior, it is interesting to study the BIED’s properties. Depending on the choice of the R and β parameters, the relative weights given to each of these distributions and their component isotherms are different. When the R and β parameters are equal in BIED, the BIED is symmetrical with respect to the heat of adsorption. The BIED can even provide an excellent approximation to the Gaussian distribution by choosing the full width at half maximum (fwhm) of the former to be a suitable multiple of the variance of the latter. In addition, when R and β are not equal, the distribution is skewed and can resemble the skewed-Gaussian distribution (the DubininRadhuskevich isotherm). The analysis neglects interaction effects among adsorbed molecules and uses the Langmuir local isotherm along with the random-patch model, with realistic heats of adsorption. For the derivation of the adsorption isotherms, the well-established Mellin transform technique employed in previous publications is used. The component distributions (a) and (b) of the BIED give rise to hypergeometric series in pressure and other adsorption heat limits, etc. The simplification of these functions to power-law expressions in various pressure ranges is illustrated. The differential heat behavior in various surface coverage regions is also studied. Depending on the relative values of R and β parameters, the differential heat vs surface coverage shows a wide spectrum of behavior ranging from that of the positive exponential to the negative exponential. When R and β are chosen equal,
the differential heat shows linear fall with the surface coverage like the constant distribution, making the discrimination between the constant and BIED distribution difficult. A method based on the product of surface coverage and the pressure derivative is suggested to discriminate between the BIED, constant, and Gaussian distributions, which show qualitatively different behavior in each case, making discrimination between the rival distributions possible. An attempt is made to fit the BIED to Drain and Morrison’s physical adsorption data 43 corresponding to the Ar/rutile system at 85 K. The fit is adequate, as is true for the skewed-Gaussian distribution, originally employed bySokolowski43 et al., along with a Langmuir local isotherm. X. Appendix In this Appendix, we give a short derivation of the adsorption isotherms corresponding to BIED. We present the Mellin transform method,44,45 which has been used earlier.10,22-23 Hence we give only the essential details below: We first break the BIED into two components (see eq 4), omitting the constant part, which is trivial. The Fmix boils down to either the positive or negative exponential distribution depending on the R and β values. The BIED can be written as a linear combination of the positive or negative exponentials:
δ(Q) ) -b1 exp(RQ) - a2 exp(-βQ)
(A1)
The gross isotherm is
θ(T,p) )
∫QQ θ(T,p)li δ(Q) dQ 2
(A2)
1
θli(T,p) ) 1/[1 + (b0/p) exp (-Q/RT)]
(A3)
We now give the Mellin transform method for deriving the adsorption isotherms corresponding to eq A2 for the negative and positive exponential distribution. We propose the following substitutions:
(b0/p) exp(-Q/RT) ) x; R′ ) RRT; β′ )βRT p1 ) b0 exp(-Q1/RT); p2 ) b0 exp(-Q2/RT); x1 ) p1/p; x2 ) p2/p With these substitutions, the negative and positive exponenetial contributions to the overall isotherm can be transformed as
∫xx xβ′-1/(1 + x) dx
(A4a)
∫xx x-R′ -1/(1 + x) dx
(A4b)
Fneg ) (p/b0)β′RT Fpos ) (b0/p)R′RT
1
2
1
2
we recall that,44 for a δ′ (x) defined as
δ′ (x) ) (1 + ax)-n )0
0 < x < x1 x > x1
(A5a)
The Mellin transform Mµ is defined as (44) Tables of Integral Transforms; McGraw-Hill: New York, 1954; Vol. 1, pp 308 and 310. (45) Higher Transcedental Functions; McGraw-Hill: New York, 1953.
1316 Langmuir, Vol. 13, No. 5, 1997
Mµ)
∫0∞xµ-1δ′(x) dx
Prasad
(A5b)
The following identities also come in handy in simplifying the final result
For eq A4a this works out as
Mµ ) µ-1x1µ F(n,µ;µ+1; -ax1)
(A6)
Γ(v) Γ(1-v) ) Π/sin(Πv)
(A10)
Γ(1+v) ) vΓ(v)
(A11)
Similar results can be derived for the positive exponential distribution eq A4. In the high-pressure region, eq A6 can be readily used, at both the limits x1 and x2 as the argument of the hypergeometric series is less than unity, Thus
Finally for the most useful working range of pressure p2 < p < p1 we have
Fpos ) (RT/µ) (b0/p)µ [F(1,-µ;1-µ;-x2)x2-µ -
Fneg ) (RTω2 /β′) (p/b0)β′ [β′π/sin(β′π) +
F(1,-µ;1-µ;-x1)x1µ (A7) Whenever the argument of the hypergeometric series exceeds unity in magnitude, the techniques of analytical continuation are used to extend the range of validity,45 i.e., of the hypergeometric function
F(a,b;c;z) ) B1(-z)-a F(a,1-c+a;1-b+a;z-1) + B2(-z)-b F(b,1-c+b;1-a+b;z-1) (A8) Since on the right-hand side the argument is less than unity, the usual expansion of the hypergeometric series is valid. All that remains to be shown isto evaluate B1 and B2 which are given by45,10
B1 ) Γ(c) Γ(b-a)/Γ(b) Γ(c-a); B2 ) Γ(c) Γ(a-b)/Γ(a) Γ(c-b) (A9)
β′/(β′ - 1)x1β′-1 F(1,1-β′;2-β′;-x1-1) x2β′ F(1,β′;β+1;-x2)] (A12) We keep in mind that ω1 ) -b1 and ω2 ) -a2.
Fpos ) [RTω1/R′] (b0/p)R′ [x2-R′ F(1,-R′;1-R′;-x2) R′π/sin(R′π) - x1-(1+R′) (R′/1+R′) × F(1,R′+1;R′+2;-x1-1)] (A13) We recall that depending on the values of Q1, Q2, x1, x2, p1, and p2 may drastically change and each of the three pressure ranges presented in Table 1 may be applicable along with the requisite analytical formulas. LA950544S
