Adsorption of Spherical Molecules in Probing the ... - ACS Publications

Kazakh-American University, Satpaev Str., 18 a, Almaty 480013, Kazakhstan,. Laboratoire Environnement et Mine´ralurgie, INPL and CNRS UMR 7569, BP 40...
0 downloads 0 Views 225KB Size
Langmuir 2002, 18, 2075-2088

2075

Adsorption of Spherical Molecules in Probing the Surface Topography. 1. Patchwise Heterogeneity Model Vadim Sh. Mamleev,† Fre´de´ric Villie´ras,*,‡ and Jean-Maurice Cases‡ Kazakh-American University, Satpaev Str., 18 a, Almaty 480013, Kazakhstan, Laboratoire Environnement et Mine´ ralurgie, INPL and CNRS UMR 7569, BP 40, 54 501 Vandoeuvre les Nancy Cedex, France Received February 12, 2001. In Final Form: October 16, 2001 The Fowler-Guggenheim local isotherm corresponding to the lattice model of localized adsorption on a patchwise surface was used to evaluate the statistical distribution of adsorption sites with respect to their adsorption energies from an experimental adsorption isotherm. An accurate low-pressure argon adsorption isotherm at 77 K on muscovite was used as initial information. The calculations were carried out in two ways that imply different surface constructions. The first one is based on the classical hypothesis of Langmuir, who assumed a finite manifold of adsorption sites and wrote an overall isotherm as a series summing the additive contributions of the different sites. The second is derived on the assumption of an infinite manifold of sites, when the sum turns into an integral. Both representations are equivalent only with the exact adsorption model. However, as shown, in a number of cases they can give similar results even with the approximate model. The influence of the ω parameter expressing the reduced energy of lateral interactions was tested for different ω values ranging between 0 and 4 in units of RT. Muscovite possesses crystallinity. Therefore, one may expect that its surface consists of recurring adsorption sites. Peaks on the energy distribution could mean at least a short range ordering (SRO) on this surface. A derivative of a local isotherm for sites of each kind is known to be a bell-shaped curve. The greater the energy of lateral attraction is, the narrower is the bell; that is, the Langmuir local isotherm (the case of ω ) 0) corresponds to the widest bell. Hence, if the hypothesis concerning SRO holds true, an energy distribution obtained with ω ) 0 must be discrete, inasmuch as a broadening of the calculated partial distributions can lead only to a deteriorating of the approximation accuracy of local derivatives and, correspondingly, of the overall derivative and the overall isotherm. Indeed, the distributions calculated with a small ω represent sums of δ-like peaks. With a small ω, both algorithms show the same results; therefore, any numerical artifacts are excluded. However, SRO on the surface is not proven yet, because with a large ω the algorithms lead to continuous distributions rather than to discrete ones. The relative approximation error of the overall isotherm with ω ) 0 is not too large (=2-3%). However, the large approximation error of the overall derivative (=15%) unambiguously indicates the unacceptability of the assumption ω ) 0 (local Langmuir isotherm). Thus, to make final conclusions about the structure of the surface and validate such an analysis, it is necessary to improve the patchwise model, to test the method by using an independent molecular simulation, and to determine realistic values of the energies of lateral interactions being different and not equal for adsorption sites of each kind.

Introduction The problem of probing the topography of heterogeneous surfaces by analyzing physical adsorption of gases was considered repeatedly.1-6 However, the accuracy of the information derived from such an analysis is still open to question. The history of the problem leads us to the work of Irving Langmuir,7 who was the first in generalizing his own wellknown isotherm describing localized adsorption on a homogeneous surface to the case of a heterogeneous * To whom correspondence should be addressed. E-mail: [email protected]. † Kazakh-American University. ‡ Laboratoire Environnement et Mine ´ ralurgie, INPL and CNRS UMR 7569. (1) Roginskii, S. Z. Adsorption and Catalysis on Heterogeneous Surfaces; Academy of Sciences of the USSR: Moscow, 1948 (in Russian). (2) House, W. A. Colloid Sci. 1983, 4, 1. (3) Jaroniec, M.; Bra¨uer, P. Surf. Sci. Rep. 1986, 6, 65. (4) Sircar, S.; Myers, A. L. Surf. Sci. 1988, 205, 353. (5) Jaroniec, M.; Madey, R. Physical Adsorption on Heterogeneous Solids; Elsevier: Amsterdam, 1988. (6) Rudzin˜ski, W.; Everett, D. H. Adsorption of Gases on Heterogeneous Surfaces; Academic Press: London, 1991. (7) Langmuir, I. J. Am. Chem. Soc. 1918, 40, 1361. (8) Zeldowich, Ya. B. Acta Physicochim. URSS 1935, 1, 961. (9) Roginskii, S. Z.; Todes, O. M. Acta Physicochim. URSS 1946, 21, 519.

surface made up of a set of adsorption sites with different energies. In the absence of lateral interactions, the overall isotherm of localized adsorption on a heterogeneous surface has the form7 N

FiP/[P + A0 exp(-i/RT)] ) n(P) ∑ i)1

na

(1)

where n is the number of adsorbed molecules per unit of surface area, na is a limiting amount of admolecules per unit of surface area, P is the gas pressure in a bulk phase, N is the number of distinguishable adsorption sites on the surface, Fi is the fraction of sites with an adsorption energy i (1 e i e N), A0 is the Langmuir constant, R is the gas constant, and T is temperature. It is obvious that eq 1 implies a finite diversity of adsorption sites on the surface. Later on, Roginskii, Zeldowich, and Todes1,8,9 modified eq 1, assuming the variety of adsorption sites tends to infinity. Their equation was of the integral form

na

∫

max

min

F()P/[P + A0 exp(-/RT)] d ) n(P)

(2)

where min and max are the minimum and maximum adsorption energies.

10.1021/la010233j CCC: $22.00 © 2002 American Chemical Society Published on Web 02/16/2002

2076

Langmuir, Vol. 18, No. 6, 2002

Mamleev et al.

To generalize the problem, it is rational to write eqs 1 and 2 in another form. In the case of a discrete distribution, the main equation and the natural constraints are expressed as follows: N

Fiθ(P,ηi) ) naΘ(P) ) n(P) ∑ i)1

na

(3)

N

Fi ) 1 ∑ i)1

(4)

Fi g 0 (1 e i e N)

(5)

Having replaced the sums by integrals, for a continuous distribution, we obtain

∫ηη

na

max

min

F(η) θ(P,η) dη ) naΘ(P) ) n(P)

∫ηη

max

min

F(η) dη ) 1

F(η) g 0 (ηmin e η e ηmax)

(6) (7) (8)

where θ(P,η) (or θ(P,ηi)) is a local isotherm corresponding to sites of a defined kind, Θ(P) ) n(P)/na is an overall isotherm describing a total degree of the monolayer occupancy, and η (or ηi) is a reduced adsorption energy related to  by the equality η )  - RT ln A0. To recalculate the reduced energy η into the absolute scale of , the Langmuir constant A0 should be determined. For this purpose some information is necessary in addition to an overall isotherm, for example, a number of isotherms measured at different temperatures or data about heat of adsorption. However, to see relative mutual disposition of different fragments of a distribution, it is enough to carry out the calculations on the scale of η. The function n(P) is determined experimentally and therefore includes errors. For this reason, neither eq 3 nor eq 6 can be exact, even if the physical model, on the basis of which these are written, is correct. Despite the conventional use of the equality symbol, eq 3 or eq 6 can be solved only approximately by the best fit of their righthand and left-hand sides. In the early publications of Russian schools,1,8,9 the change from eq 3 to eq 6 has not been substantiated at a physical level10 and apparently is dictated for mathematical convenience. Handling by integrals is somewhat easier than by sums, since the large number of methods of the functional analysis allows eq 6 to be solved. The mathematical techniques for the analysis of the series like eq 3 are much less developed. At the same time, the statistical derivation of isotherms of localized adsorption implies a finite manifold of adsorption sites.11 Even if the number of different types of sites is very large, say N ) 6.02 × 1023 per kilogram of an adsorbent, the solution of eq 6 from the standpoint of the lattice model should be written10,11 in the form of the sum N

F() )

Fiδ( - i) ∑ i)1

(9)

(10) Lopatkin, A. A. Theoretical Fundamentals of Physical Adsorption; Moscow State University: Moscow, 1983 (in Russian). (11) Tovbin, Yu. K. Theory of Physical Chemistry Processes at GasSolid Interface; Nauka: Moscow, 1990 (in Russian).

where δ(x) is the Dirac δ function. Equation 9 corresponds to the series in eq 3. Because of experimental errors and the inaccuracy of the model, an inversion of eq 3 with a large N loses physical sense. It is hardly possible to obtain realistic information about the energy distribution with N > 10. That is why the analysis of heterogeneity presupposes either an averaging of adsorption energies for different sites by eq 3 or a possibility of representation of the energy distribution in the form of a smooth histogram by eq 6. Langmuir7 intuitively ascribed to the averaging in eq 3 the chemical sense. Indeed, one can suppose that similar surface clusters with definite geometry (steady combinations of superficial atoms) form typical adsorption sites. They correspond to approximately equal adsorption energies that can be attributed to one peak on a real energy distribution F(). Thus, determination of the sets {Fi} and {i} means the chemical identification of the clusters. The approach of Roginskii, Zeldowich, and Todes1,8,9 ignores the chemical identification. In other words, it implies that clusters having different conformations and chemical natures can correspond to the same interval of energies ranging from  to  + d. When developing the numerical methods for evaluation of F(), most authors1-6,12-18 adhered to this approach. If the function F() exhibits legible peaks, their assignment to certain clusters is still possible. However, in the situation where the energy distribution resembles “white noise”, the analysis of surface cluster organization obviously becomes unobtainable. If a model for local isotherms describing adsorption on certain sites would be exact, one should admit that solving the integral equation is preferable in comparison to inverting the series. However, the exact model is absent. Because of its inaccuracy, the peaks on a distribution can be unresolved, even if typical adsorption centers really exist in nature and could be detected with an exact model. Moreover, for reproducibility of a solution of a Fredholm integral equation of the first kind like eq 6, regularization (smoothing) is necessary.12-18 The regularization is able to eliminate the peaks totally. Note that regularization is physically lawful19,20 only in the cases when a searched for physical function, according to a priori information, is really a smooth one. Since such information is inaccessible, the use of eq 6 instead of eq 3 actually has no advantages. Actually, a preference for use of either eq 3 or eq 6 depends on a priori information about the structure of a surface. For surfaces carrying recurring adsorption sites, the use of eq 3 is preferable. For example, it is obvious to assume that energy distributions should be discrete in the case of surfaces with expressed crystalline structures.21 The detection of peaks of F() for such surfaces would mean at least the existence of short range ordering (SRO) in clusters. (12) House, W. A.; Jaycock, M. J. Colloid Polymer Sci. 1978, 256, 52. (13) Merz, P. J. Comput. Phys. 1980, 38, 64. (14) LamWan, J. A.; White, L. R. J. Chem. Soc., Faraday Trans. 1991, 87, 3051. (15) Von Szombathely, M.; Bra¨uer, P.; Jaroniec, M. J. Comput. Chem. 1992, 13, 17. (16) Koopal, L. K.; Vos, C. H. W. Langmuir 1993, 9, 2593. (17) Papenhuijzen, J.; Koopal, L. K. In Adsorption from Solutions; Ottewill, R. H., Rochester, C. H., Smith, A. L., Eds.; Academic Press: London, 1983. (18) Jagiello, J. Langmuir 1994, 10, 2778. (19) Tikhonov, A. N.; Arsenin, V. Ya. The Methods of Solution of Incorrect Problems; Nauka: Moscow, 1979 (in Russian). (20) Tikhonov, A. N.; Goncharskii, A. V.; Stepanov, V. V.; Yagola, A. G. Regularization Algorithms and a priori Information; Nauka: Moscow, 1983 (in Russian). (21) Villie´ras, F.; Michot, L. J.; Bardot, F.; Cases, J. M.; Franc¸ ois, M.; Rudzin˜ski, W. Langmuir 1997, 13, 1104.

Probing Surface Topography

Langmuir, Vol. 18, No. 6, 2002 2077

In general, the necessary condition for surface crystallinity is crystallinity of the solid itself. Unfortunately, the correctness of eq 3 is disputable even in the case of crystalline adsorbents. For instance, the Monte Carlo simulation22-24 of experimental argon adsorption data25,26 for the surface of crystalline TiO2 agrees better with the hypothesis of an amorphous surface organization than with the hypothesis of a crystalline organization. The calculated energy distribution shows the absence not only of long range ordering (LRO) natural for crystalline surfaces but even of SRO. Numerical experiments22-24 are not yet accurate enough to make unequivocal conclusions. Besides, the Monte Carlo method allows the investigation of only a finite number of hypothetical variants of superficial structures. It is possible to conclude that some of these variants correspond better to the experiment than other ones, but it is impossible to accept them as ones strictly corresponding to reality. One should bear in mind that eq 6 belongs to the class of ill-posed (or incorrect19,20) mathematical problems. The small deviations in the right-hand side of this equation can induce large deviations of the function F(). For this reason, one can find a multitude of possible energy distributions due to deviations between simulated and experimental isotherms (of 5% or more).22-24 These distributions can be both discrete and continuous ones. We note that, as opposed to the Monte Carlo simulation, on the basis of analysis of adsorption isotherms, a number of authors16,17,27,28 even for amorphous surfaces found the distributions F() representing five to seven narrow peaks. These results are still disputable. For example, the strict analysis of errors15 for these evaluations shows that the number of significant peaks is apparently less. In any case, if SRO really exists on the surfaces, both approaches, that is eqs 3 and 6, must give comparable results. The application of eqs 3 and 6 becomes problematic in systems with significant lateral interactions between admolecules. These equations can include the simplest corrections for the lateral effects. For systems with patchwise heterogeneity, the isotherm of Fowler-Guggenheim,29 which takes into account lateral interactions within the framework of the lattice model in the meanfield approximation, is often used. This isotherm is written as follows:

θ(P,η) ) P/[P + exp(-η - Zwθ/RT)]

(10)

where w is the nearest-neighbor interaction energy defined so that w is positive for any attractive interaction, and Z is the coordination number of the two-dimensional surface lattice. If w ) 0, eq 10 leads to eqs 1 and 2 with the local Langmuir isotherm as the kernel. Equation 10 is derived under the assumption that a heterogeneous surface consists of macroscopic uniform patches. The existence of LRO means a homothetic structure of a surface or crystallinity. However, eq 10 was often used2,3 even for surfaces of amorphous adsorbents. In any case, it is a very crude approximation, since the (22) Bakaev, V. A. Surf. Sci. 1988, 198, 571. (23) Bakaev, V. A. Langmuir 1992, 8, 1372. (24) Bakaev, V. A. Langmuir 1992, 8, 1379. (25) Morrison, J. A.; Los, J. M.; Drain, L. E. Trans. Faraday Soc. 1951, 47, 1023. (26) Drain, L. E.; Morrison, J. A. Trans. Faraday Soc. 1952, 48, 316. Drain, L. E.; Morrison, J. A. Trans. Faraday Soc. 1952, 48, 840. (27) Bra¨uer, P.; House, W. A.; Jaroniec, M. Thin Solid Films 1982, 97, 369. (28) Sacher, R. S.; Morrison, I. D. J. Colloid Interface Sci. 1979, 70, 153. (29) Fowler, R. H.; Guggenheim, E. A. Statistical Thermodynamics; University Press: Cambridge, 1949.

parameter w has some averaged (effective) value for a surface. In the next article30 we shall improve this approximation. Correction for w * 0 is justified only for systems in which distances between admolecules’ positions exceed their collision diameter. It is difficult, at least within the framework of the localized adsorption model, to take into account correlation between admolecules, since it is obvious that simultaneous localization of molecules in two adjacent centers is impossible, if the distance between these centers is less than a collision diameter. Finally, as implied in the lattice model, the admolecules’ location relative to the surface must be determined. One can imagine a real heterogeneous adsorbent as a rough surface consisting of a set of “dimples” (or “holes”). The holes are characterized by various binding energies. If the difference max - min much exceeds RT, in the holes with the binding energy max the adsorbate will be in a localized state. As regards the localization of the adsorbate in the holes with the energy min, it is defined by the oscillation amplitude of a surface potential in the proximity of the holes. If this amplitude is commensurable with the magnitude of RT, the adsorbate can hardly be considered as being localized but mobile. For wide energy distributions (max - min = 10RT), one can expect a situation with localized adsorption on max sites and mobile adsorption in the proximity of min sites. In such a case, the calculated F() function will have a correct shape only for high adsorption energies. We have enumerated only some of the difficulties restricting the research concerning the heterogeneity problem when applying eqs 3 and 6 even for elementary systems. Nevertheless, such research, in our opinion, is not a hopeless task. In the past decade precise techniques for measurement of adsorption isotherms at 77 K over a very wide range of pressures (=10-6 to 103 Torr) have advanced considerably.21,31,32 The density of experimental data points is so large that it is possible to construct not only isotherms but also their first derivatives with respect to the logarithm of pressure. The derivatives of isotherms are more sensitive to the effects of heterogeneity in comparison to isotherms. In particular, in some cases the derivatives show several peaks. Without any mathematical consideration, one may conclude that these peaks correspond to different adsorption sites (or surface clusters). Obviously, these experimental successes necessitate efforts in the development of new mathematical methods for the analysis of energy distributions. Initially it is necessary to accept a few hypotheses that will be verified below. First, under adsorption of spherical molecules such as Ar on crystalline adsorbents, such as clay minerals, the energy distribution should be almost discrete. The ribbed surface of crystals such as clay minerals is frequently observed using electronic microscopy. Flat faces of monocrystals are unlikely to be amorphous to such an extent that these lose their cluster structure in common with predictions24 of the Monte Carlo method for TiO2. Second, the number of surface potential minimums, which are distinguishable for their different depths, can (30) Villie´ras, F.; Mamleev, V. Sh.; Nicholson, D.; Cases, J. M. Submitted to Langmuir. (31) Michot, L. J.; Franc¸ ois, M.; Cases, J. M. Langmuir 1990, 6, 677. (32) Villie´ras, F.; Michot, L. J.; Cases, J. M.; Berend, I.; Bardot, F.; Franc¸ ois, M.; Ge´rard, G.; Yvon, J. In Equilibria and Dynamics of Gas Adsorption on Heterogeneous Solid Surfaces; Rudzin˜ski, W., Steele, W. A., Zgrablich, G., Eds.; Studies in Surface Science and Catalysis, 104; Elsevier: Amsterdam, 1997; p 573.

2078

Langmuir, Vol. 18, No. 6, 2002

much exceed the number of peaks that can be observed in a distribution resulting from applying mathematical models such as eqs 3 and 6. However, the models can predict the number of some of the strongest sites. The last assumption is necessary for the reason that, with the exception of rare cases, the two-dimensional lattices of an adsorbate and an adsorbent do not coincide. Thus, in a close packed monolayer, adsorbed molecules cannot be localized in all adsorption sites. Furthermore, having occupied a part of adsorption sites, a localized adsorbate by itself can generate new surface potential minimums. At the same time, one should keep in mind that, with increasing pressure, the occupancy of adsorption sites occurs sequentially: strong sites are filled first followed by weak ones. That is why the weak sites do not influence adsorption on strong ones at low pressure and, consequently, an energy distribution in the high-energy region can be evaluated with good accuracy. At high pressure, a part of the weak sites are overlapped by admolecules adsorbed by strong sites. It is legitimate to assume that the weakest sites, remaining at the final stage of monolayer filling, are indiscernible. Eventual peaks in the lowenergy region can be ascribed to such weak sites. Third, although the localization of an adsorbate is possible only in strong sites, the lattice model conserves its phenomenological correctness and gives plausible estimations for the adsorption energy of weak sites. Actually, adsorption isotherms in elementary models of localized and mobile adsorption, relative to their shape, differ insignificantly, and by variation of effective parameters, these can be easily matched. The discrete distribution model for crystalline adsorbents seems to be more realistic for physical reasons and can be accepted a priori for analysis of the surface heterogeneity of monocrystals. The purpose of the present paper is to justify the use of the series like eq 3 instead of the integral like eq 6 for systems with SRO. The series can include more exact corrections for lateral effects21,30 and allows an interpretation of chemical properties of adsorption sites. To show the equivalency of both approaches, by solving eq 6 it is necessary to establish the existence of SRO [distinct peaks on the distribution F()] and then to confront the results of calculation for continuous and discrete distributions. Argon Adsorption Isotherm on Muscovite. Muscovite is a mica phyllosilicate of ideal half unit cell Si3AlO10Al2(OH)2K. The mica material originates from Aitik orebody (Sweden) and was obtained from Boliden Mineral Ab. Pure muscovite was handpicked from these materials while using a microscope, ground in a pebble mill, and wet-sieved to minus 75 µm fraction. A fraction minus 5 µm was then prepared using microsieves in an ultrasonic bath.33 An experimental argon adsorption isotherm was recorded at 77 K using the quasi-equilibrium device described in refs 21, 31, and 32. Argon N56 (purity >99.9996%) was provided by Alphagaz (France). Prior to the experiment, powder was outgassed during 18 h at 150 or 250 °C under a residual pressure of 10- 3 Torr. The Ar adsorption isotherm at 77.3 K on muscovite is shown in Figures 1 and 2 on conventional and logarithmic pressure scales. According to the definition of Brunauer et al.,34 the present isotherm belongs to the typical second type of isotherms. (33) Hanumantha Rao, K.; Cases, J. M.; Barre`s, O.; Forssberg, K. S. E. In Mineral Processing, Recent Advances and Future Trends; Mehrotra, S. P., Shekhar, R., Eds.; Allied Publishers Limited: New Delhi, 1995; p 29.

Mamleev et al.

Figure 1. Original isotherm of Ar adsorption on muscovite at 77.3 K and the isotherm of monolayer adsorption according to the BET theory on the scale of reduced pressure. P0 is the pressure of the saturated vapors of Ar. P0 ) 197.45 Torr at the level of Nancy (France).

Figure 2. Data of Figure 1 plotted on the scale of the pressure logarithm (pressure in Torr). The enlarged fragment illustrates the density and accuracy of the experimental points in the isotherm.

The setup used for the quasi-equilibrium volumetric experiments has been described in detail earlier.31,32 The reproducibility of the experiment has been examined repeatedly. Within the pressure range 10-6 < P < 100 Torr, its relative error with respect to measured adsorbed quantities did not exceed =0.2%. In comparison with reproducibility, a constant bias is more complicated for evaluation, because the latter, apart from experimental magnitudes, includes calculated ones. The absolute error in the measured adsorption is minimal at low pressures. However, the relative systematic error is determined mainly by the inaccuracy of the initial fragment of the isotherm. Hence, an objective judgment about the relative error can be obtained by investigating the initial isotherm fragment (Figure 3). A spurious discontinuous adsorption jump at the initiation pressure in the region 3 × 10-7 < P < 10-6 Torr is observed. However, this jump does not exceed =0.1% of the measured maximum adsorbed quantity. In the region 3 × 10-7 < P < 10-3 Torr, the measurement error on pressure is still rather large (=10%) but decreases quickly with increasing pressure. In the region 2 × 10-3 < P < 10 Torr, it apparently does not differ considerably from reproducibility. In the interval 3 × 10-7 < P < 2.15-3 Torr, the isotherm is the concave curve (see Figure 3). At 2.15-3 < P it becomes the convex curve. If lateral interactions at low pressures (34) Brunauer, S.; Deming, L. S.; Deming, W. S.; Teller, E. J. Am. Chem. Soc. 1940, 62, 1723.

Probing Surface Topography

Langmuir, Vol. 18, No. 6, 2002 2079

Figure 4. Original and corrected derivative of the isotherm.

Figure 3. Initial fragment of the isotherm. Illustration of its liner extrapolation toward zero pressure (a). The points used in the calculations with the corrected and uncorrected isotherm (b). The short nonlinear interval corresponding to equilibration of the ionization detector was excluded from the calculations.

are negligible, both for a homogeneous and for a heterogeneous surface an overall isotherm must be convex.1 This classical assumption seems logical. Therefore, most probably, the concavity of the isotherm is explained on account of the inaccuracy of the experiment. To estimate the systematic error (see eq 15), we have found the tangency point (=2.2-3 Torr) of the isotherm with the straight line corresponding to Henry’s law. The difference between the measured isotherm and that from Henry’s law within the pressure interval 5.86 × 10-4e P e2.43 × 10-4 was used for the evaluation of the systematic bias (see Figure 3b). In fact, it is the error of uncertainty in the measurement of the isotherm within this interval. It was established (see below) that this error very weakly influences the energy distribution calculated. The pressure interval Pmin ) 5.86 × 10-4 e P e Pmax ) 7.82 Torr incorporating the Henry’s law region up to the beginning of multilayer adsorption (Figures 2 and 3) was chosen for further analysis. A feature of the quasi-equilibrium experiment is that it includes (Figure 2) a very large number of experimental data points (=2000). Therefore, the problem of isotherm numerical differentiation does not arise. The differentiation was carried out by calculation of central differences in intervals incorporating 10 experimental data points. After such an operation over the pressure interval Pmin e P e Pmax, the relative error in the derivative did not exceed 2%. For calculations using the isotherm, 150 data points were selected, diminishing the amount of experimental data points by a factor of 10 (Figure 4). The decrease of the number of points does not influence a calculated energy distribution but considerably decreases computing time. We note also that not the isotherm but its derivative was used to check for the correctness of the calculated distributions described below. Because of its sensitivity,

the derivative is more convenient for visual control of discrepancies between the theoretical and experimental data. Numerical Procedure. At the present time, a theory that allows one to separate specific contributions of monolayer and multilayer adsorption, strictly speaking, is missing. However, the Monte Carlo numerical simulations22-24 of adsorption on oxides such as TiO2 facilitate this problem. First, at least in the case of adsorption on oxides, the monolayer is practically completely saturated23,24 (say, =95%) at rather low relative pressures: P/P0 = 0.1, where P0 is the saturation pressure of an adsorbate. Second, curiously enough, the simple Brunauer-Emmet-Teller model (BET), which is categorically unacceptable for the description of adsorption on homogeneous surfaces, gives rather exact asymptotic values24 in the case of heterogeneous surfaces, namely

n(P) = nB(P)(1 - P/P0)

(11)

where nB(P) is an experimental isotherm of multilayer adsorption. For calculation of adsorption energy distributions, we can formulate our own approach as follows. First, when solving the inverse problems, one can limit, in principle, an isotherm by the region of monolayer adsorption. However, to magnify information concerning superficial effects, one should increase the considered pressure interval up to appreciable multilayer adsorption. At the same time, because of the lack of a strict theory, the multilayer adsorption region (say, nB/n > 1.1-1.2) should be disregarded. Second, while solving the inverse problems, the value of na, which expresses a maximal number of “sitting places” in a monolayer, is much more important in comparison to the asymptotic adsorption at P f P0. Fortunately, the BET theory predicts the value of na with good accuracy. Moreover, at the present time, the BET approximation is the only well-tested method for na determination. A compromise in the choice of the upper pressure limit Pmax can hardly be solved rigorously. We were guided solely by an intuition, having restricted Pmax in order to include the maximum corresponding to monolayer saturation on the derivative and to eliminate the minimum related to multilayer adsorption (see Figure 4). Within the interval Pmin e P e Pmax, the correction of an isotherm was implemented by means of eq 11. As shown by numerical experiments, with such a correction, the calculated distributions are only slightly altered by normalization conditions 4 and 7, whereas, in the absence

2080

Langmuir, Vol. 18, No. 6, 2002

Mamleev et al.

of the correction, these are considerably perturbed in the low-energy region.35 The normalization, in turn, can be valid only with known na. As shown in Figure 2, the corrected isotherm has a distinct plateau, namely, na ) 1.83 cm3/g at standard temperature and pressure (STP). Note that, in the absence of information about na, one cannot hope35 to obtain a realistic shape for the calculated adsorption energy distribution in a low-energy region. In our opinion,35 the “penalty” method is the simplest and most convenient for the normalization. It consists of solving the following optimization problems: M

N

{1 - [na/n(Pj)]∑Fiθ(Pj ,ηi)}2 + ∑ j)1 i)1

∆R1({Fi},{ηi}) ) M-1

N

R1(1 M

{1 - [na/n(Pj)]∫η ∑ j)1

Fi)2 w min ∑ i)1

(12)

ηmax

∆R1[F] ) M-1

∫ηη

max

min

Ω[F] )

∫ηη

max

min

[qF2 + (dF/dη)2] dη

(16)

where q is a weighting coefficient. Instead of the original variational problem (eq 13), a perturbed problem is solved, namely

F(η) θ(Pj ,η) dη}2 +

∆R,R1[F] ) ∆R1[F] + RΩ[F] w min

F(η) dη)2 w min (13)

where R is the regularization parameter. The standard method for solving problem 17 with constraint 8 consists of discretization of the integrals and application of the NNLS algorithm (non-negative least squares) of Lawson and Hanson.37 The NNLS algorithm is a time-consuming method. For computational speedup, the IRA (improved regularization algorithm) was offered,35 in which the nonnegativity constraint was ensured by the barrier function method. We note that this method35,36 can be applied in solving both problems (12 and 13). However, the IRA possesses a number of shortcomings. In particular, it requires a good initial approximation for the function F(η) and does not guarantee a convergence with a bad initial approximation. Besides, it allows one to find a minimum of the discrete analogue of functional 17 only approximately, while the NNLS algorithm ensures reaching an exact minimum, although after a large number of iterations. The difference between the approximate minimum obtained by the IRA and the exact minimum obtained by the NNLS algorithm usually lies in a fifth decimal digit. Nevertheless, to reach the exact minimum as closely as possible, the classical NNLS method, slightly improved in the standard package of Tikhonov et al.,20 has been preferred. The calculations were carried out with M ) 150. When computing continuous distributions, 101 points on the energy grid were distributed evenly within the interval -7.5 e η/RT e 6.5. Zero weighting factor (q ) 0) was used in Tikhonov’s stabilizer. With such a value, the stabilizer minimally perturbs a solution in the neighborhood of the boundaries ηmin and ηmax. When calculating continuous distributions, an initial approximation for the function F(η) can be preset arbitrarily; in particular, F(η) ) 0 was used in the computations. It is worthwhile to note that, in discrete distribution calculations, the initial values of the components {Fi} should not be zeros.36 Initial values for {ηi} were used at equidistant points over the range -7.5 e ηi/RT e 6.5, and initial values for {Fi} were set as Fi ) 1/N. With such an initial choice, the algorithms for solving problems 12 and 17 showed a convergence in all cases.

min

R1(1 -

A solution for similar problems was considered earlier.35,36 In particular, an algorithm for evaluation of discrete distributions with constraint 5, but in the absence of normalization 4, that is with R1 ) 0, was described in detail.36 Its extension to the case R1 > 0 does not give rise to special complications, and computing details will not be discussed here. The calculation algorithm for continuous distributions requires some explanations. To attain stability, in solving ill-posed problems, the Tikhonov’s stabilizer is often applied.19,20 In particular, the first-order stabilizer has the following form:

where n(Pj) are the points on an overall isotherm given at M values of pressure Pj (1 e j e M), P1 ) Pmin, PM ) Pmax, and R1 is the penalty parameter. The quantity, the conditional minimization of which is envisaged according to eqs 12 and 13, is the average square of the relative deviation between theoretical and experimental isotherms and not the absolute deviation; that is, under constraints 4 and 5 (or 7 and 8), the objective function is of the form M

[1 - ntheor(Pj)/nexp(Pj)]2 ∑ j)1

∆ ) M-1

(14)

where ntheor(Pj) and nexp(Pj) are the values of the theoretical and experimental isotherms given for the M points. Such a choice for the objective function is justified, as it provides for an almost constant chance relative error over the whole pressure interval Pmin e P e Pmax. Both terms on the right-hand sides of eqs 12 and 13 are of the quadratic form. It can be proven19,20 that, at the minimum of ∆R1, these terms are strictly monotone functions of R1. If R1 f ∞, the normalization constraints (eqs 4 and 7) are fulfilled exactly. Thus, R1 can be determined by the data of accuracy in the normalization. Numerical experiments show that, with the parameter R1 ) 10, errors in the normalization do not exceed 0.1% both for discrete and for continuous distributions. As mentioned above, an initial fragment of the isotherm limits the systematic experimental error. Therefore, this error can be roughly estimated by means of linear extrapolation toward zero pressure of the initial isotherm fragment (Figure 3). By using fourteen initial points (see Figure 3b), we have 14

[1 - nexp(Pj)/nHenry(Pj)]2/150}1/2 ) ∑ j)1

1/2 ={ ∆systematic

0.0187 (15) where nHenry(Pj) are the isotherm values of the extrapolated line. (35) Mamleev, V. Sh.; Bekturov, E. A. Langmuir 1996, 12, 441.

(17)

(36) Mamleev, V. Sh.; Bekturov, E. A. Langmuir 1996, 12, 3630. (37) Lawson, C. L.; Hanson, R. J. Solving Least Squares Problems; Prentice Hall: Englewood Cliffs, NJ, 1974.

Probing Surface Topography

Langmuir, Vol. 18, No. 6, 2002 2081

Figure 5. (a) Continuous energy distributions calculated with (R ) 10-4) and without (R ) 0) regularization. (b) Cumulative distributions calculated for the models of continuous (CD) and discrete (DD) distributions with different N values. The parameter of lateral interactions ω is equal to 0.0 (the Langmuir model).

Figure 6. Same as in Figure 5, but with ω ) 0.5.

Results and Discussion In Figures 5-12 and Table 1, the calculated data obtained with the kernel expressed by isotherm 10 are represented for various values of the parameter ω ) Zw/ RT. With ω > 4, the model predicts the phase transition of the first kind in an adsorbed film. In the supercritical region (ω > 4), the model, strictly speaking, loses a physical sense. For this reason, the calculations were carried out under the condition ω < 4. Since F(η) in the form of a δ function cannot be displayed graphically, it is convenient to compare discrete and continuous distribution models by using the cumulative distribution (integral distribution) that has the following form:

F)

∫ηη

min

F(η) dη

(18)

For better visualization of the amount of the strongest adsorption sites, it is rational to plot (1 - F) versus -η; this increasing function vanishes at the limit ηmax. The continuous distributions were calculated with and without regularization. When solving problem 17, the main difficulty is the choice of the suitable regularization parameter R. Let D(R) be the value of ∆1/2 at the minimum of functional 17. The choice of R, through the discrepancy D(R), commonly consists of solving the following transcendental equation:19,20 1/2 D(R) ) ∆random

(19)

1/2 where ∆random is the average random experimental error. However, the selection of R from eq 19 is justified only in

Figure 7. Same as in Figure 5, but with ω ) 1.0.

the absence of constraints 7 and 8, when D(0) ) 0. The constrained minimum of ∆1/2 with R ) 0 is not equal to zero. The overall error of the approximation consists of two addends. The first addend is equal to the random 1/2 . The second one is the errors of the experiment, ∆random systematic deviation between the theoretical and experimental isotherms that is caused both by the systematic

2082

Langmuir, Vol. 18, No. 6, 2002

Mamleev et al.

Figure 10. Regularized continuous distributions calculated for the patchwise model with ω ) 3.0, 3.5, and 3.9.

Figure 8. Same as in Figure 5, but with ω ) 1.5.

Figure 11. Dependence of the relative mean-square error of the approximation of the isotherm upon the ω parameter used in the calculations.

Figure 9. Regularized (a) and nonregularized (b) continuous distributions calculated for the patchwise model with ω ) 2.0 and 2.5.

experimental error and by the imperfections of the model. To an approximation, these kinds of errors can compensate each other, but generally being averaged in many experiments, the systematic deviation D(0) is equal

Figure 12. Overall derivative corresponding to the models of discrete (a) and continuous (b) distributions with different ω values.

Probing Surface Topography

Langmuir, Vol. 18, No. 6, 2002 2083

Table 1. Results of the Calculation of Discrete Distributions with Different ω and N Valuesa ω)0

ω ) 0.5

ω ) 1.0

ω ) 1.5

) 0.032 255 for the continuous distribution with R)0

) 0.024 253 for the continuous distribution with R)0

) 0.015 859 for the continuous distribution with R)0

) 0.007 988 0 for the continuous distribution with R)0

∆1/2

∆1/2

∆1/2

N)2

N)4

N)5

N)7

∆1/2 ) 0.032 169

∆1/2 ) 0.024 178

∆1/2 ) 0.015 828

∆1/2 ) 0.007 971 9

-ηi/RT -3.9426 0.39488

Fi

-ηi/RT

0.16283 1.6672 na ) 1.8300

-4.1604 -1.4885 0.73503 8.8642

Fi

-ηi/RT -4.3937 -2.4829 -1.2097 0.85415 2.6879

0.11985 0.17365 1.5240 0.012384 na ) 1.8299

Fi

-ηi/RT -4.7016 -3.3726 -2.1767 -0.93245 0.10746 1.0506 2.5961

0.083272 0.068454 0.18121 1.3026 0.19559 na ) 1.8311

Fi

0.046149 0.053914 0.063809 0.18349 0.11176 1.0080 0.36162 na ) 1.8287

N)8

N)8

N)8

N)8

∆1/2 ) 0.032 170

∆1/2 ) 0.024 178

∆1/2 ) 0.015 828

∆1/2 ) 0.007 971 7

-ηi/RT

Fi

-3.9427 0.162 81 -2.4054 0.136 50 × 10-4 -1.3871 0.121 37 × 10-4 -1.3751 0.819 88 × 10-5 -0.730 37 0.305 72 × 10-4 0.002 723 4 0.427 03 × 10-7 0.394 93a 0.694 26a 0.394 93a 0.972 87a na ) 1.8300 a

∆1/2

-ηi/RT -4.1604 -1.4886 -1.4107 1.9199 0.735 04a 0.735 03a 4.4021 9.0578

Fi 0.119 85 0.173 51 0.139 72 × 10-3 0.298 88 × 10-4 1.2695a 0.254 43a 0.940 50 × 10-4 0.012 411 na ) 1.8300

-ηi/RT

Fi

-4.3937 0.083 272 -2.4828 0.068 464 -1.2100 0.179 35 -1.1710 0.001 850 9 0.854 17a 0.332 20a 0.854 17a 0.970 36a 2.6880a 0.032 711a 2.6880a 0.162 87a na ) 1.8311

-ηi/RT

Fi

-4.7002 0.046 324 -3.3532 0.055 760 -2.0931 0.069 819 -0.826 46 0.205 25 0.788 82 0.581 73 1.2762 0.540 18 2.6581a 0.101 55a 2.6577a 0.228 97a na ) 1.8296

Corresponding to the degeneracy.

to D h (0) > 0. If only one measurement is carried out, then one may use the evaluation D h (0) ≈ D(0). On the basis of additivity of the errors, it is reasonable to select R from the condition 1/2 D(R) = ∆random + D(0)

(20)

A solution resulting from such a regularization should be reproducible for isotherms measured in various experiments. Condition 20 is fulfilled with R ) 10-4, and the corresponding regularized distributions are shown in Figures 5-10. As observed with low ω values (ω ) 0.0-1.5), the obtained solution represents pronounced peaks. It is noteworthy that with R ) 0 both minimization variants (eqs 12 and 13) give almost equivalent results. The calculated minimum in eq 12 is less than the minimum in eq 13 (Table 1). This means that an exact solution of any of the variational problems (eqs 12 and 13) is a discrete distribution. This fact is supported not only by the minimum values of ∆1/2 but also by the degeneracy of discrete distributions that is observed in calculations with different values of N (Figures 5-8 and Table 1). In the case of continuous distributions, the values of ∆1/2 are higher because it is impossible to calculate precisely the integral in eq 6 with a finite density of the energy grid. However, with increasing numbers of knots on the energy grid, the minimums of ∆1/2 for the calculated discrete and continuous distributions must approach each other. The degeneration property of energy distributions has already been demonstrated in the case of the artificial example.36 The numerically generated isotherm, used as an example, was physically incorrect, and the conditional minimum of ∆ corresponded to the very large error. Upon rather exact approximation of the real isotherm, the detected distribution degeneracy is quite unexpected. It contradicts the usual viewpoint that the best solution for eq 6 is always a smooth function.

One can rather simply explain a δ-like solution in the case of surfaces having only one kind of adsorption site. It is known that the derivative of a local isotherm represents a bell-shaped curve.21 Let us estimate some parameters of a derivative of an isotherm. One can readily find the expression for the derivative of isotherm 10 with respect to the pressure logarithm as follows:

dθ/d ln P ) θ(1 - θ)/[1 - ωθ(1 - θ)]

(21)

The maximum of function 21 is reached for θ ) 1/2 and ln P ) -η/RT - ω/2. Thus, the height of the derivative is equal to

H ) 1/(4 - ω)

(22)

In the case of the Langmuir isotherm (ω ) 0) at the level of the half-height of the bell (Hh ) H/2 ) 1/8), the extents of occupancy are expressed as θ1,2 ) 1/2 ( x1/4-1/8, that is θ1 ) 0.854 and θ2 ) 0.146. In terms of the chemical potential, the energy width of the derivative is expressed by the simple formula

∆µ/RT ) ln P1 - ln P2 ) ln[(θ1/θ2)(1 - θ2)/(1- θ1)] - ω(θ1 - θ2) (23) where ∆µ is the difference of the chemical potentials of the adsorbate at the points θ1 and θ2. With the Langmuir local isotherm (ω ) 0), the halfwidth of the derivative is maximal: ∆µ/RT ) 3.525. In the presence of lateral attraction, the bell becomes narrower. For the same values of θ1 and θ2, but with ω ) 4.0, we have ∆µ/RT ) 0.697. If, because of lateral effects, the experimental derivative is narrower than that predicted by the Langmuir isotherm, the best solution of problem 13, with the Langmuir equation as a kernel, will be a δ

2084

Langmuir, Vol. 18, No. 6, 2002

Mamleev et al.

function. Indeed, any wide energy distribution, such as a Gaussian one, can lead only to an extension of the bell of the derivative and, consequently, can only increase the deviation between the theoretical and experimental derivatives. The appearance of several δ-like peaks in a distribution, for example, two peaks (Figure 5), means that a variety of adsorption sites is present on the surface with corresponding bells of derivatives narrower than those predicted by the Langmuir local isotherm (see Figure 12a). The relative approximation error of the overall isotherm with the Langmuir kernel may seem not too large (=3%, see Figure 11). However, the large approximation error of the overall derivative (=15%, see Figure 12a) unambiguously indicates physical unacceptability of the local Langmuir isotherm. In formal mathematical terms, the discrete distributions are the natural result of exact minimization 12 under constraint 5 with a finite M number. It is known that, for occurrence of the ∆ minimum, the fulfillment of Kuhn-Tucker conditions is necessary.38 These conditions have the following form:

∂∆/∂Fi - λi ) 0 (1 e i e N)

(24)

∂∆/∂ηi ) 0 (1 e i e N)

(25)

λiFi ) 0 (1 e i e N)

(26)

where λi is a non-negative Lagrange multiplier. Let us imagine first a variant of unconditional minimization 12, that is without relationships 5 and 26 and with λi ) 0 (1 e i e N) in eq 24. To interpolate M points of an experimental isotherm by a theoretical one, it is enough to define N ) M/2. Then equations 24 and 25 form a system of N ) M/2 algebraic equations to determine M/2 values of Fi and M/2 values of ηi. If, by chance, all M/2 components Fi obtained from this system are positive, the solution satisfies the physical sense and we have the distribution of N ) M/2 δ-like peaks (in our case N ) M/2 ) 75). Obviously, it is the maximum of the possible number of peaks, because, after adding one more term in eq 3, the system of equations 24 and 25 becomes an undetermined one and cannot be solved without constraints 5 and 26. At the same time, it is known13 that, in the absence of constraint 5, small high-frequency random errors in the input data cause solution oscillations that contain large positive and negative Fi components.13 As the unconditional minimum of ∆ contains negative Fi, a theoretical isotherm determined at the conditional minimum (under constraint 5) cannot interpolate an experimental one; it can only approximate it. In other words, one can hardly find a solution containing a large number of δ-like peaks, for example, 75; in practice, their number will be less. With large N, the exact conditional minimization of ∆ will lead to a degenerated solution, meaning that some Fi components are either equal to zero or correspond to identical energy values (Figures 5-8 and Table 1). In such a situation, the N number can be taken to be arbitrarily large, even N > M; all the same, the solution will exhibit a constant number of the same peaks. One should pay attention to the fact that the smooth fragments in continuous distributions (Figure 9) calculated without regularization (R ) 0) are connected to the use of the finite energy grid. Increasing the grid density will inevitably lead to a set of δ functions. (38) Gill, P. E.; Murray, W.; Wright, M. H. Practical Optimization; Academic Press: London, 1981.

The more systematic deviations between theoretical and experimental isotherms there are, the lower the number of peaks in a distribution will be. In particular, if a model contradicts a given integral isotherm at the qualitative level, only one δ-like peak will remain in the distribution, for example, as for the artificial isotherm considered earlier.36 In practice, one can express a solution as a sum of δ functions only when their number is low. It follows from Figures 5-8 and Table 1 that with ω ) 0.0-1.0 it is possible to find the exact minimum of ∆. With ω ) 1.5, computer round-off errors make it difficult to find. Figure 8 and Table 1 show that, despite the proximity of the two values of ∆1/2, the solutions calculated using ω ) 1.5 with N ) 7 and N ) 8 differ appreciably. The two values of ∆1/2 are so close that the final results become dependent upon the initial approximation for {Fi} and {ηi}. Comparing distributions calculated with different ω values,one may see that the number of peaks in the distributions grows when increasing ω. With ω ) 1.5 the value of N ) 6 is explicitly not sufficient to reach maximum accuracy for isotherm approximation (Figure 11), and the continuous distribution describes the isotherm more precisely in comparison to the discrete one. However, with N ) 7 or 8, the discrete distributions obtained provide better accuracy than the continuous ones. The diminution of ∆ with increasing ω can be explained by considering the contributions of separate distribution fragments to an overall isotherm. The form of an isotherm, in the neighborhood of any arbitrary point with coordinate ln Px, is defined by the form of an energy distribution within a limited energy interval that hardly exceeds the span (∆µ counted from the point ηx ) -RT(ln Px + ω/2). Correspondingly, a fragment of the isotherm within the interval ln Px - ∆µ/RT < ln P < ln Px + ∆µ/RT should be defined with good accuracy by a fragment of the distribution within the interval of the energies ηx - 2∆µ < η < ηx + 2∆µ. The known condensation approximation (CA)39,40 implies that ∆µ f 0; then

(1/na) dn/d ln P|Px)exp(-ηx/RT-ω/2) ) F(ηx) By differentiating both parts of eq 6 and by restricting the limits of integration in accordance with our assumptions, we have

(1/na) dn/d ln P|Px)exp(-ηx/RT-ω/2) ≈ +∆µ F(η) (dθ(P,η)/d ln P)|P )exp(-η /RT-ω/2) dη ∫ηη-∆µ x

x

x

x

Under the assumption that the function F(η) within the region of the integration is a constant, with the Langmuir isotherm (ω ) 0), one can obtain

(1/naRT) dn/d ln P|Px)exp(-ηx/RT) ≈ F(ηx)th[∆µ/(2RT)] ) 0.943F(ηx) Even if the function F(η) is a constant within the interval =7RT, the CA gives the error =6%. Thus, this method is hardly suited even for semiquantitative analysis of energy distributions. The assumption of a lack of correlation of points on an isotherm, which underlies this method, is too rough. If this assumption were valid, a theoretical isotherm would describe an experimental isotherm with exactitude. (39) Cerofolini, G. F. Surf. Sci. 1971, 24, 391. (40) Harris, L. B. Surf. Sci. 1971, 26, 667.

Probing Surface Topography

Figure 13. Discrete cumulative distributions calculated with and without correction of the overall isotherm by Henry’s law with ω ) 0.0, 0.5, and 1.0. It is seen that the correction slightly perturbs the distributions only in the high-energy region.

As a matter of fact, within the pressure logarithm interval =2∆µ/RT, points on an isotherm are considerably correlated; therefore, the theoretical isotherm is not an arbitrarily flexible function. The reason for the correlation is that fragments on an isotherm, which have the interval =2∆µ/RT on a logarithmic pressure scale, correspond to defined fragments of the interval =4∆µ in an energy distribution. The number of independent (strictly speaking, weakly dependent) fragments in a distribution (the number of peaks) with decreasing ∆µ (magnification of ω) grows. Simultaneously, with increasing ω, the accuracy of matching of an experimental isotherm with a theoretical one is improved. Figure 11 confirms this conclusion, except for ω ) 3.9, which is very close to the critical point ω ) 4.0. At this point, the magnification of approximation error is due to the loss of accuracy in the numerical integration of eq 6 caused by sharp oscillations of the kernel when changing values of P. This effect is also exhibited in oscillations of the calculated derivative of the isotherm (Figure 12b). It is of interest to what extent the constant bias affects the results of the evaluations. It is mentioned above that, for maximal elimination of the constant bias, we have used the linear extrapolation of the experimental isotherm at very low pressure (Figure 3). The calculation results for the corrected isotherm are represented in Figures 1315. The same approximation error with the extrapolated isotherm is reached at the lower values of ω (compare Figure 11 with Figure 15). With a large ω it is possible to approximate with high accuracy both the original and corrected (extrapolated) isotherms. The correction with Henry’s law causes only small perturbations of the cumulative distributions in the highenergy region, but, as a whole, a rather large constant

Langmuir, Vol. 18, No. 6, 2002 2085

Figure 14. (a) Discrete (DD) and continuous (CD) distributions and (b) the cumulative distribution calculated with and without correction of the overall isotherm by Henry’s law with ω ) 1.5. It is seen that the perturbations are small enough.

Figure 15. Dependence of the relative mean-square error of the approximation of the isotherm upon the ω parameter after recalculation of the isotherm by using Henry’s law.

bias does not significantly influence the calculated distributions and does not violate their quantitative features. Thus, neither the form of the distribution nor the error of the approximation is a criterion for selecting ω. Only the assumption that at low pressures lateral effects are small allows one to give preference to extrapolation of the isotherm. The extrapolation leads to the decrease of the ω value necessary for approximating the isotherm within the chance error of the experiment. Although solely discrete distributions always correspond to an exact minimum of ∆, the number of peaks in distributions is reproducible only for ω ) 0. From Figures 6-9 one can see that some peaks at ω > 0 merge to form common peaks, as a result of regularization. Thus, to control the reproducibility of a solution, it is necessary to use regularization even with small values of ω. When calculating continuous distributions, a main principle of regularization consists of smoothing the solution. Very small perturbations imported into func-

2086

Langmuir, Vol. 18, No. 6, 2002

Mamleev et al.

Figure 16. Cumulative distributions calculated for regularized continuous (CD) and discrete (DD) distributions with different N and ω values.

Figure 17. (a) Approximation of the full-length isotherm by using the isotherm truncated by a factor of ≈1.8. (b) Regularized continuous distributions. (c) Steplike cumulative distributions resulting from the exact minimization of ∆.

tional 13 due to regularization cancel the degeneracy; namely, the δ functions in distributions turn into bellshaped functions similar to the Gaussians. This indicates that within the framework of experimental accuracy, it is impossible to determine widths of distributions for different sites. One should be satisfied with information concerning only positions of peaks and peak areas. Therefore, the evaluation algorithms of discrete distributions,21,32,36 at least in a number of cases, can be used with the same success as that for continuous algorithms.12-18 Regularization in computing discrete distributions can be organized by suitable choice of the N value. In Figure 16 the calculation results of continuous and discrete distributions are shown with different N values that equal the numbers of bells in regularized continuous distributions. With ω ) 0.5 and N ) 3, the similarity of the obtained results does not raise doubts. With ω ) 1.0 and 1.5, the value of N ) 4 is obviously insufficient. The largest peaks in Figures 7 and 8 are too wide. Their approximation by using only one δ function causes large errors. However, with N ) 5 (Figure 16), the discrete and continuous cumulative distributions are in close agreement. With ω ) 0.0 and 0.5, the average calculated deviation exceeds the constant experimental bias (Figure 11). Therefore, with such parameters the model cannot be considered as being satisfactory. The calculation variants with ω > 1.5 also raise doubts. According to the wellknown opinion of the fundamental study of Drain and Morrison,26 up to Θ = 0.6, lateral interactions can be disregarded at all. To check this hypothesis, we have restricted the overall isotherm within the interval 0 < Θ < 0.5 and have repeated the calculations with ω ) 0.0, 0.5, 1.0, and 1.5 for the truncated isotherm (Figures 15, 17, and 18).

The results shown in Figure 17 for the Langmuir isotherm represent indisputable proof of the manifestation of lateral attraction at low monolayer coverage. As before, the best solution is the sum of two δ functions, and the cumulative distributions calculated for the full-length and truncated isotherms coincide with each other with excellent accuracy. However, as mentioned, the distribution obtained with the local Langmuir isotherm does not provide a satisfactory approximation of the overall isotherm and should be rejected. Using proof by contradiction, one may conclude that evaluation of a realistic distribution should allow for lateral attraction even at Θ < 0.5. One of the reviewers noted ref 41 where adsorption is considered with negative values of ω (repulsion). Probably the repulsion can lead to stepwise adsorption isotherms.41 The steps can be falsely interpreted as different kinds of sites. However, we exclude such an effect with respect to the experiment under consideration. Indeed, the bells of local derivatives with negative ω should be wider than those with ω ) 0, but even with ω ) 0, it is impossible to approximate the isotherm with satisfactory accuracy. In other words, the monotone increase of the error (see Figures 11 and 15) contradicts the hypothesis concerning the negative ω. In our opinion, the strongest adsorption sites are formed by interparticle contacts.30 For these sites lateral interactions can be significant. However, according to both our calculations and the common viewpoint, their fraction is small enough (=2%). The next sites in descending order of adsorption energy are expected to be point vacancy defects on cleavage planes. (41) Ramirez-Pastor, A. J.; Bulnes F. M.; Riccardo J. L. Surf. Sci. 1999, 426, 48.

Probing Surface Topography

Langmuir, Vol. 18, No. 6, 2002 2087

sponding to the high-energy peaks deviates more and more from the Langmuir isotherm. Actually, one should keep in mind that ω can vary for different adsorption sites. However, in the model under consideration, the ω parameter is accepted as being the same for all sites. The model used is very imperfect. Its basic assumption is an identity of adsorption sites within each patch. Surfaces of natural monocrystals seldom conform to such an assumption. In a much more general situation, one could use a model that takes into account the real spatial disposition of close adsorption sites with respect to each other.30

Figure 18. Two regularized distributions corresponding to the smallest values of ω with which the satisfactory approximation of the overall isotherm is reached. The considerable truncation of the isotherm (Figure 17a) practically does not change the calculated distributions in the high-energy region.

If so, these are separated from each other at distances considerably exceeding the collision diameter of Ar. Hence, for these sites, one may accept ω ) 0. The regularized solutions with ω ) 1.0 and 1.5 (Figure 18) can be selected as some compromise. On one hand, these distributions provide a satisfactory approximation of the overall isotherm. On the other hand, these correspond to the minimum value of ω, which, in turn, is required according to the hypothesis about the independence of strongest adsorption sites. Note that, after the truncation of the isotherm, the form of the distributions in the high-energy region is completely conserved. Thus, the correction of the isotherm by eq 11 does not influence the high-energy fragments. The factor, which really is of great importance for such an analysis, is the energy of lateral interactions. Its value predetermines the calculated distribution. One may conclude that a realistic distribution contains at least four peaks (Figure 18). Three of them may correspond to adsorption sites of an individual kind, while the largest peak certainly incorporates sites of different kinds, because, under adsorption on the weak sites, admolecules undergo strong lateral interactions. At the same time, to approximate the overall isotherm when using a large ω, one should assume within the low-energy region either a continuous distribution or a distribution consisting of many δ-like peaks. The patchwise model in itself neither confirms nor refutes SRO on the surface. However, since the surface of muscovite is explicitly crystalline, a discrete distribution seems more realistic. Magnification of ω leads to the splitting of the δ-like peaks into a series of peaks, and simultaneously it improves the approximation accuracy. However, with increasing ω, the local isotherm corre-

Conclusions Two observations made in this work contradict the foundational concepts accepted in the theory of physical adsorption. First, contrary to the conventional opinion that a distribution of adsorption sites with respect to their energies should be a continuous and smooth function, in a number of cases the algorithms of precise minimization of the relative deviation between theoretical and experimental isotherms demonstrate a discrete distribution representing a sum of δ functions. Second, contrary to the well-known assertion that lateral interactions can be disregarded until the extent of the monolayer occupancy, Θ, reaches =0.6, when restricting the overall isotherm by the interval 0 < Θ < 0.5, the Langmuir local isotherm predicts an unrealistic energy distribution. Since the restriction of Θ in the calculations does not influence the form of the distribution, it is possible to conclude that bad resolution of the distribution is really related to the neglect of lateral interactions. The hypotheses put forward at the beginning of the article are worthy of additional discussion. First, the discreteness of the energy distributions detected with low and moderate values of the parameter of lateral interactions, ω, really testifies in favor of SRO of the surface under study. If the surface consists of a finite variety of adsorption sites, and the derivative for a model local isotherm is wider than derivatives for true local isotherms for different sites, the calculated distribution should be discrete. Second, the peaks in a distribution really have distinct positions in the high-energy region, while the form of the distribution within the low-energy region cannot be resolved with good accuracy. The peak corresponding to the lowest energy has the largest width and the largest area; therefore, it by no means can be represented as a δ function. This means that it incorporates sites of different kinds. A realistic local isotherm for these sites certainly should correspond to a large ω. However, with a large ω, a theoretical integral isotherm becomes a flexible function and the evaluation of an energy distribution becomes very sensitive to any errors. Third, as shown, a good approximation of the experimental isotherm can be reached despite the absence of a strict model for the local isotherm. It is legitimate to assume that a distribution is continuous within the region of the peak corresponding to the lowest energy. Note that, in the algorithm destined for the calculation of discrete distributions, a few addends should be assigned for this peak. However, this is necessary only for reaching a good formal approximation of the overall isotherm. Under physical interpretation, the set of these δ functions should be considered as one peak. The purpose of the paper is not reached to the end. Indeed, the hypothesis concerning SRO of the surface does not contradict the original integral equation that gives

2088

Langmuir, Vol. 18, No. 6, 2002

the discrete distributions with small ω. If a solution with the Langmuir local isotherm would be a continuous function, we could assert that the surface, for sure, is amorphous. At the same time the discreteness of the distributions obtained with some ω values does not prove SRO, inasmuch as values of ω acceptable for evaluating a realistic distribution still require additional proofs. Since the absence of information concerning correct energies of lateral interactions is seen to be the fatal obstacle in such an analysis, the next step of the theory should consist of extraction of their values from an experiment. Molecular simulation of adsorption on a crystalline surface by the methods of Monte Carlo or molecular

Mamleev et al.

dynamics would be useful in order to make sure of the correctness of solving of the inverse problems. For increasing the resolution of the energy distributions in the computations, one should allow for information reflecting real crystallography of surfaces. This will be the subject of the second article.30 Acknowledgment. This work was implemented according to Project INTAS-2000, Ref.No: 00-505. V.Sh.M. is grateful to the Ministry of National Education, Research, and Technology of France for the grant making his visit to J.-M.C.’s laboratory possible. LA010233J