Order−Disorder Phase Transitions in Adsorbed Films. I. Monolayer

Oct 25, 2007 - In the case of films ordering into the c(2 × 2) commensurate phase at low temperatures, ... In particular, the monolayer films undergo...
0 downloads 0 Views 460KB Size
15664

J. Phys. Chem. C 2007, 111, 15664-15676

Order-Disorder Phase Transitions in Adsorbed Films. I. Monolayer and Bilayer Films of Square Symmetry† A. Patrykiejew* and S. Sokołowski Department for the Modelling of Physico-Chemical Processes, Faculty of Chemistry, MCS UniVersity, 20031 Lublin, Poland ReceiVed: April 13, 2007; In Final Form: June 29, 2007

Canonical as well as grand-canonical ensemble Monte Carlo simulation methods are used to study the disordering of thin, monolayer, and bilayer films formed by Lennard-Jones fluids on the (100) plane of the face-centered-cubic (fcc) crystal. Two classes of systems, ordering into the (1 × 1) and c(2 × 2) commensurate phases at low temperatures, are considered. It is demonstrated that the mechanism of the disordering of adsorbed films ordering into the (1 × 1) commensurate phase at low temperatures occurs gradually and in a layerby-layer mode. Only in the case of large misfit between the adsorbate atoms and the surface lattice do the commensurate films undergo a sharp first-order transition, which occurs at different temperatures in different layers. In the case of films ordering into the c(2 × 2) commensurate phase at low temperatures, the mechanism of disordering strongly depends on the film thickness as well as on the misfit between the adsorbate and substrate lattice. In particular, the monolayer films undergo a continuous order-disorder transition belonging to the universality class of the two-dimensional Ising model, whereas bilayer films may disorder via a single first-order or continuous transition, depending on the relative size of adsorbed atoms and the surface lattice unit cell.

I. Introduction Phase transitions and ordering phenomena in monolayer as well as multilayer films adsorbed on well-defined crystal surfaces have been attracting a great deal of interest over the last four decades.1-7 Strong experimental evidence exists stemming from adsorption studies of simple molecules on metals,8-13 the basal plane of graphite,14-19 lamellar metal dihalides,20-22 and on other substrates (e.g., MgO)24-26 that the structure and thermodynamic properties of adsorbed monolayers are determined primarily by the competing effects due to the adsorbate-adsorbate interaction and the substrate corrugation potential. In general, the adsorbate-adsorbate interaction tends to establish an adsorbed monolayer structure of hexagonal symmetry and the natural interatomic distance as determined by the properties of the interaction potential. Alternatively, the substrate corrugation potential favors the formation of commensurate structures in which the adatoms occupy adsorption sites that result from the periodic variation of the potential field generated above the crystal surface. Of course, not only the lateral modulations of the substrate potential field but also the total strength of the adsorbate-substrate interaction influences the behavior of adsorbed films. Weakly adsorbed films are likely to lose any form of ordering because of desorption that sets already at low temperatures. In particular, when the adsorbate does not wet the solid surface, the stability of adsorbed layers is very limited. When, however, the substrate potential is strong enough to support the formation of a stable monolayer as well as multilayer adsorbed films, the corrugation potential comes into play. In such systems, the structure of possible commensurate adsorbed layers depends on the symmetry properties †

Part of the “Keith E. Gubbins Festschrift”. * Corresponding author. E-mail: [email protected].

of the surface lattice as well as on the relative size of adsorbed atoms and the surface lattice unit cell, which can be conveniently measured by the so-called dimensional incompatibility parameter,21 defined as

I ) (a - r1)/r1

(1)

where a is the surface lattice constant and r1 is the distance between adsorbate atoms in the surface phase. Of course, in the case of a rectangular symmetry of the surface lattice, one can define two different incompatibility parameters I1 and I2, by taking different lattice constants as a. The role of dimensional incompatibility on the properties of adsorbed monolayers was demonstrated experimentally in the case of adsorption of simple adsorbates (e.g., rare gases, methane) on lamellar metal dihalides.23 The surface of lamellar metal dihalides of the general formula MeX2 forms a honeycomb lattice of the lattice constant depending on the sizes of metal and halide atoms. It is thus possible to obtain substrates of different surface unit lattice cells and of the same symmetry. It was demonstrated that the critical as well as the triple point temperatures in monolayer films reach maximum values for vanishing misfit between the adsorbate and the adsorbent. This behavior is characteristic of the commensurate structures, and it was also observed in lattice gas models of adsorbed monolayer films.27 The formation of commensurate phases is usually possible when the misfit between adsorbate and adsorbent is sufficiently small and when the surface corrugation potential is strong enough to localize adsorbate atoms over adsorption sites, as defined by the positions of the corrugation potential minima. When these conditions are not fulfilled, various incommensurate phases may appear.7,28-30 One expects that the formation of the second layer on top of the commensurate monolayer film should also lead to the

10.1021/jp0728813 CCC: $37.00 © 2007 American Chemical Society Published on Web 10/25/2007

Order-Disorder Phase Transitions

J. Phys. Chem. C, Vol. 111, No. 43, 2007 15665

formation of a commensurate phase, at sufficiently low temperatures at least, because the monolayer becomes a source of the corrugation potential of the symmetry properties determined by the structure of the surface lattice. This prediction is not always correct, as shown by X-ray diffraction study of bilayer oxygen film on graphite,31 which demonstrated that the top layer forms a lattice of rectangular rather than the expected hexagonal symmetry. When the monolayer undergoes the commensurate-incommensurate (C-IC) or the order-disorder (O-D) phase transition before the adsorption in the second layer sets in, the formation of the second layer may considerably affect the structure of the first layer and different scenarios are possible. For example, Monte Carlo simulation study of adsorption on the (100) plane of a model fcc crystal demonstrated32 that when the misfit between adsorbate and adsorbent is sufficiently large the (1 × 1) monolayer ordered structure may undergo the C-IC transition leading to the formation of a haxagonally ordered floating incommensurate phase. The second layer also forms a hexagonal lattice and stabilizes hexagonal ordering in the first layer. It is also possible that the presence of the second layer may lead to the recovery of the commensurate structure in the first layer. Such situations can be expected when the surface potential felt by the atoms adsorbed in the second layer is still strong. The main goal of this work is to elucidate the effects of the misfit between adsorbate and adsorbent on the properties of the monolayer as well as bilayer films formed on a square lattice of sites and to investigate the mechanism of disordering in the films ordering into the (1 × 1) and c(2 × 2) structures. The paper is organized as follows. In Section 2, we present the model and Monte Carlo methods used in our study. Sections 3 and 4 present the results obtained for monolayer and bilayer films, respectively. The paper concludes in Section 5, which summarizes our findings. II. Model and Monte Carlo Methods Here, we use a simple fluid model, assuming that it consists of atoms interacting via the Lennard-Jones (12,6) potential

u(r) )

{

4[(σ/r)12 - (σ/r)6] r e rmax r > rmax 0

(2)

which is cut at the distance rmax ) 3.0σ and we do not use any long-range corrections.33 Throughout this paper, we use the parameter  as a unit of energy. The fluid is assumed to be in contact with a wall, being a perfect (100) plane of an fcc crystal, forming a square lattice characterized by the unit vectors a1 and a2 of the same length a, taken as a unit of length. The potential field due to the solid substrate felt by the fluid atom is a periodic function of the two-dimensional vector r* )(x*,y*) (x* ) x/a, y* ) y/a), which represents the position of an adatom over the surface plane, and hence can be written in the form of a Fourier series34

V*(z*, τ*) ) *gs[Vo(z*) +

∑ Vg(z*) fg(τ*))

(3)

g*0

In the above, the sum runs over all nonzero reciprocal lattice vectors of the length g, z* ) z/a is the distance from the surface, the parameter /gs ) gs/ measures the strength of the interaction between a fluid atom and a single atom of the solid and a vast majority of results presented here has been obtained for /gs ) 1.0. In the case when the interaction between an adsorbate atom and an individual atom of the solid is also

TABLE 1: Dimensional Incompatibilities and the Heights of Potential Barrier for Diffusion for the Systems of Different σ* values σ*

ordered structure

I

V/D

0.82 0.84 0.86 0.88 0.90 0.92 0.93 0.94 1.20 1.22 1.24 1.26 1.28 1.30 1.31 1.32 1.33 1.34

(1 × 1) (1 × 1) (1 × 1) (1 × 1) (1 × 1) (1 × 1) (1 × 1) (1 × 1) c(2 × 2) c(2 × 2) c(2 × 2) c(2 × 2) c(2 × 2) c(2 × 2) c(2 × 2) c(2 × 2) c(2 × 2) c(2 × 2)

0.11229 0.08580 0.06055 0.03645 0.01342 -0.00861 -0.01927 -0.02971 0.07489 0.05727 0.04022 0.02371 0.00771 -0.00779 -0.01537 -0.02283 -0.03017 -0.03741

1.613 1.596 1.578 1.559 1.541 1.523 1.513 1.504 1.266 1.248 1.230 1.212 1.194 1.177 1.168 1.159 1.150 1.142

represented by the Lennard-Jones (12,6) potential, one has analytical expressions for the Fourier components Vo(z) and Vg(z) as well as for the functions fg(τ).34 A quantitative measure of the surface potential corrugation is the height of potential barrier between adjacent sites, calculated as

/gs ) |V*(z/min, τ/min) - V*(z/sp, τ/sp)|

(4)

where τ/min ) and τ/sp are the locations of the potential minima over the site and the saddle point, respectively, and z/min and z/sp are the corresponding distances from the surface at which the fluid-surface potential reaches a minimum. Of course, the magnitude of VD depends on σ* ) σ/a and changes linearly with /gs. Table 1 summarizes the values of VD for the systems considered in this work. This Table also gives the values of the misfit between the adsorbate and the adsorbent for the systems ordering into the (1 × 1) as well as into the c(2 × 2) commensurate phases at low temperatures. In the calculations of the misfit, we have taken r1 as equal to the length of the lattice constant of the fcc crystal formed by an atom interacting via the Lennard-Jones potential at zero temperature 1.0964σ*35, and the lattice constant of the ordered state was assumed to be a and x2a in the case of the (1 × 1) and c(2 × 2) structures, respectively. Similar to earlier works,36-38 we use Monte Carlo methods in the canonical as well as grand-canonical ensembles.39,40 All calculations have been performed using square simulation cells of the size in x and y directions equal to L* × L*, with L* ranging between 10 and 40, with standard periodic boundary conditions applied in x and y directions. In the case of simulations in the canonical ensemble, the parallel tempering method has been used41,42 allowing us to study a series of thermodynamic states, that is, at different temperatures, in a single run. The method is very attractive because it considerably reduces the number of Monte Carlo steps needed to reach equilibrium as well as to obtain reliable estimations of averages, by a factor of 10-2. A single Monte Carlo step (MCS) consisted of N (N being the number of particles in the system) randomly chosen trials to change the system state by attempted displacement of a randomly chosen particle, and also by a randomly chosen three-dimensional vector from the cube of the side dmax (dmax was adjusted periodically

15666 J. Phys. Chem. C, Vol. 111, No. 43, 2007

Patrykiejew and Sokołowski

to keep the acceptance ratio at about 40%). The number of MCS’s performed during a single run varied between 105 and 106, and a similar number of MCS’s was used to equilibrate the system. When the Monte Carlo method in the grand-canonical ensemble was used, we performed isothermal runs starting from a dilute, gas-like and then increased the chemical potential up to the values when two adsorbed layers were formed. In all cases, the adsorption runs have been followed by the desorption runs at the same temperature. One should note that whenever the system undergoes a first-order phase transition the adsorption-desorption isotherm exhibits a hysteresis loop, resulting from the metastability effects. The transition points have been located assuming that the effects due to metastability are more or less the same in both phases. Of course, we are aware that this method is very crude and does not provide reliable results when the hysteresis loops are large, but our main aim has been to obtain a qualitative rather than quantitative picture of the phase behavior of the systems studied. A very tedious calculations of the system free energies are necessary in order to determine the locations of the phase transition points accurately.40 The recorded quantities included the average number of particles (in the case of simulations carried out in the grandcanonical ensemble), the average adsorbate-adsorbate, 〈egg〉, and the adsorbate-solid, 〈egs〉, interaction energies (per particle), the density profiles of the adsorbate, averaged over the entire surface wall, n(z*), and the in-plane pseudo two-dimensional radial distribution functions, g(r), calculated separately for each atomic layer, using the same procedure as applied in refs 43 and 44. In the case of Monte Carlo simulations in the canonical ensemble, we have also recorded the heat capacity. The internal structures of the adsorbed films have also been probed by the calculation of the following bond-orientational order parameters38,45,46

Ψk,l )

〈|

1

Nb,l



∑ ∑ exp[ikφij]| ,, k ) 4,6 i∈l j∈l

(5)

that are suitable to detect the formation of layers of square and hexagonal symmetry. In the above equation, φij is the angle between the “bond” joining the nearest neighboring atoms i and j (both located in the lth layer) and the fixed reference axis, Nb,l is the number of such “bonds” in the lth layer for each atomic layer, and 〈...〉 means the averaging over all configurations generated during a single Monte Carlo run. The first nearest neighbors of each molecule have been determined using the in-plane cutoff, defined as the location of the first minimum at the radial distribution function.43 Moreover, we have also monitored the behavior of the Fourier transform of the local density in each layer

|Fq(l)| )

|

1

∑ exp(-iq‚ri)| Nl

Nl i)1

(6)

which allows us to determine the long-range order corresponding to the presence of the lattice structure of the assumed commensurate phase characterized by the reciprocal lattice vector q.47,49 In the case of the (1) phase we have taken q ) (2π/ 0.5L)(1, 1) and q ) (π/0.5L)(1, 1) in the case of the c(2 × 2) phase. In the above equation, Nl is the number of atoms in the layer l and ri ) (xi, yi) is the two-dimensional vector defining the location of the ith atom with respect to the surface.

To perform finite size scaling analysis, we have performed calculations using simulation cells of different size L × L, with L ranging between 10 and 40 and calculated the susceptibilities as well as the fourth-order Binder cumulants40 of the abovedefined order parameters Ψk,l and |Fq(l)| using the definitions

χop ) NlT*[〈op2〉 - 〈op〉2]

(7)

and

Uop(L) ) 1 -

〈op4〉 3〈op2〉2

(8)

where op marks either Ψk,l or |Fq(l)|. It should be stressed that the simulation method used in this work has a serious drawback because of the assumption that the simulation cell size is fixed in the directions parallel to the solid surface. In reality, solid-like phases undergo thermal expansion upon the increase of temperature. By keeping the simulation cell size fixed we cannot reduce the planar stress, which considerably influences the stability of solid-like adsorbed layers, just the same as it influences the melting of bulk materials.48 One should note, however, that the in-plane simulation cell size variation is hard to implement in our model because of the assumed properties of the fluid-solid potential and the necessity of using periodic boundary conditions. This difficulty could be (at least partially) overcome by taking into account the thermal expansion of the substrate and using temperaturedependent fluid-solid potential. Alternatively, one can also try to use molecular dynamics simulation methods in which the substrate atom’s vibrations are explicitly taken into account. III. Monolayers In this section, we present and discuss the results obtained for monolayer films, which at sufficiently low temperatures form the commensurate (1 × 1) or c(2 × 2) phases. To elucidate the effects of misfit between the adsobate atoms and the surface lattice, we have performed a canonical ensemble MC simulation for a series of systems characterized by different sizes of adsorbate atoms, σ*, and at the surface densities corresponding to completely filled commensurate phases, equal either to Fs ) N/M ) 1 (N is the number of adsorbed atoms and M ) L × L is the number of adsorption sites), in the case of the (1 × 1) structure, or to 0.5, in the case of the c(2 × 2) structure. Figure 1 presents the plots of temperature changes of the order parameters |Fq(1)| (parts a and c) and Ψ4(1) (parts b and d) obtained for the systems ordering into the (1 × 1) structure and characterized by different sizes of adsorbate atoms, ranging between σ* ) 0.82 and 0.93. As expected, the disordering of completely filled monolayers ordered into the (1 × 1) structure occurs gradually. The same situation takes place in lattice gas models when the disordering of the (1 × 1) structure under the condition of a fully filled lattice is considered.50 The inspection of density profiles, radial distribution functions, and snapshots of configurations has demonstrated that the disordering is accompanied by the promotion of adsorbed atoms to the second layer and that even at high temperatures, where the order parameters already assume rather low values of the order of 0.3, the atoms in the monolayer are highly localized. In particular, the radial distribution functions recorded at high temperatures do not exhibit the behavior characteristic to fluid-like phases but rather show the structure characteristic of partially ordered systems. This is illustrated by the results obtained for the system with σ* ) 0.86 and given in Figure 2,

Order-Disorder Phase Transitions

Figure 1. Temperature changes of the order parameters |Fq(1)| (parts a and c) and Ψ4(1) (parts b and d) obtained for the systems characterized by different σ* (shown in the figure) and the simulation cell size 16 × 16 × 10.

J. Phys. Chem. C, Vol. 111, No. 43, 2007 15667

Figure 3. Plots of the susceptibility of the order parameter |Fq(1)| vs temperature for monolayer films formed by atoms of different sizes and for different simulation cell sizes, L × L × 10 (the values of L are shown in the figure).

Figure 4. Plots of temperature dependencies of the order parameters |Fq(1)| and Ψ4(1) for the system with σ* ) 0.94, obtained for two different sizes of the simulation cell L × L × 10 (the values of L are shown in the figure).

Figure 2. Radial distribution functions at different temperatures (shown in the figure) for the system characterized by σ* ) 0.86 and the simulation cell size 16 × 16 × 10. Note that the subsequent curves have been shifted by 1.0 for clarity. The inset shows the density profiles obtained for the same system at two different temperatures (shown in the figure).

which present the radial distribution functions and the density profiles at different temperatures. Quite-similar results have been obtained for other systems with adsorbate atoms of σ* between 0.82 and 0.93. The observed behavior of the radial distribution functions is a direct consequence of rather-high values of the potential barrier between adjacent minima of the surface potential (cf. Table 1). The adsorbed atoms are forced to occupy the positions close to the locations of the surface potential minima (adsorption sites) even at the temperatures at which the film is already disordered. It should be noted that in the systems characterized by σ* ∈ [0.82, 093] the total potential energy of the system changes smoothly with temperature and hence we do not observe any anomalies in the behavior of the heat capacity. In fact, the heat capacity versus temperature curves exhibit only broad maxima and do not show any finite size effects, as expected for the systems that do not undergo a phase transition. Similarly, the

susceptibilities of the order parameters |Fq(1)| and Ψ4(1) also exhibit finite maxima without any appreciable effects resulting from the finite size of the simulation cell. Figure 3 presents the temperature changes of χFq(1) obtained for the systems of different σ* and quite-similar results have been obtained for χΨ4(1). In the case of larger adsorbate atoms, with σ* ) 0.94, we have found a different behavior. Namely, the system exhibits a sharp first-order transition at T* ) 0.183 ( 0.019, rather than a gradual disordering as observed in systems of smaller adsorbate atoms discussed previously. Figure 4 presents the temperature changes of the order parameters |Fq(1)| and Ψ4(1), which demonstrate the presence of well-developed hysteresis loops. The displayed plots demonstrate that the low-temperature phase is the commensurate (1 × 1) structure and that the hightemperature phase still shows a considerable square ordering. The inspection of snapshots has shown that the high-temperature phase is the incommensurate phase consisting of large commensurate domains separated by walls in which the adsorbed atoms are hexagonally ordered (see Figure 5). Also, the radial distribution functions recorded for that phase have indicated the presence of small, but quite-well seen, peaks corresponding to hexagonal ordering. The incommensurate phase disorders gradually when the temperature increases and the disordered phase exhibits a partial localization of adatoms, just the same as that observed in systems discussed previously. The stability of the commensurate (1 × 1) structure depends on the misfit between the adsorbate atoms and the substrate surface lattice, and one expects that it should be the highest

15668 J. Phys. Chem. C, Vol. 111, No. 43, 2007

Patrykiejew and Sokołowski

Figure 7. Temperature dependencies of the order parameter |Fq(1)| (part a) and |Ψ4(1)| (part b) for two systems of different σ* ordering into the c(2 × 2) commensurate phase obtained for different sizes of the simulation cell L × L × 10 (the values of L are shown in the figure). Figure 5. x-y projection of the configuration that shows the structure of the incommensurate phase obtained for the system with σ* ) 0.94 at T* ) 0.17 and the simulation cell of the size 20 × 20 × 10.

Figure 6. Plots of temperatures at which the susceptibilities of the order parameters |Fq(1)| and Ψ4(1) as well as the heat capacity reach a maximum vs the misfit between the adsorbate atom size and the substrate surface lattice size for monolayer films ordering into the commensurate (1 × 1) structure at low temperatures. Note that the results for the system with σ* ) 0.94 correspond to the location of the first-order phase transition. The numbers in the figure show the values of σ* for all of the systems considered.

when the misfit is the smallest. Figure 6 shows the plots of the temperatures at which the heat capacity, as well as the susceptibilities of the order parameters |Fq(1)| and Ψ4(1) reach maximum values, plotted against the natural misfit between the adsorbate and the adsorbent, as defined by eq 1. One readily notes that in the case of positive misfit exceeding about 0.04, the locations of these maxima are only very-weakly affected by the strain induced by the periodic substrate potential. This, again, results from the fact that the corrugation potential “pins” the small atom into the registry positions more efficiently than in the case of larger adsorbate atoms. One should also note that the maxima at the curves shown in Figure 6 do not occur at the zero misfit but are shifted to about 0.04. Quite-similar effects were found for the experimentally determined23 dependences

between critical and triple-point temperatures and the misfit for monolayer films formed by noble gases and simple molecules (CH4) on lammelar dihalides. The fact that the curves show a much-steeper fall for negative misfit is a direct consequence of the form of the fluid-fluid interaction potential. Even a small compression of the film leads to a large increase of the fluidfluid interaction energy and hence to a lowering of the stability of the commensurate structure. The second series of calculations has been carried out for the systems characterized by σ* ∈ [1.2,1.34], which order into the c(2 × 2) commensurate structure at sufficiently low temperatures. This structure has the surface density F ) 0.5 and the atatoms occupy one of two equivalent and interpenetrating sublattices of the surface lattice. Thus, the first nearestneighbor distance between the adatoms is equal to x2a. Within the formalism of the two-dimensional lattice Ising model, this commensurate structure corresponds to the antiferromagnetic phase, which is known to undergo a continuous order-disorder transition.51 One can also expect that the monolayer films, formed on surfaces characterized by sufficiently strong and a highly corrugated surface potential field, should also belong to the universality class of the two-dimensional Ising model. To verify the above hypothesis, we have performed a series of calculations for the systems of adatoms of different σ* and of different simulation cell sizes. The behavior of the order parameter |Fq(1)| obtained for different sizes of the simulation cell have demonstrated that the systems with σ* between 1.20 and 1.33 all exhibit a continuous order-disorder transition. Parts a and b of Figure 7 show two examples of temperature dependencies of the order parameters |Fq(1)| and Ψ4(1), respectively, obtained for the systems characterized by different values of σ* eqaul to 1.20 and 1.30, and for the simulation cells of different sizes. It is seen that the behavior of the order parameter |Fq(1)| shows the finite size effects characteristic of the continuous order-disorder phase transition, whereas the order parameter Ψ4(1) behaves differently and does not exhibit any appreciable finite size effects though demonstrates that the square ordering is severely reduced when the system undergoes the order-disorder transition. The data presented in Figure 7 indicate that the bond-orientational order parameter, Ψ4, is not an appropriate order parameter for the

Order-Disorder Phase Transitions

J. Phys. Chem. C, Vol. 111, No. 43, 2007 15669

Figure 8. Plots of the susceptibility of the order parameter |Fq(1)| vs temperature for the system characterized by σ* ) 1.20 and different sizes of the simulation cell (shown in the figure). The inset gives the log-log plots of the maximum values of χFq(1) against the size of the simulation cell for the system with different σ*. The dashed lines have been plotted assuming the value of the critical exponent’s ratio γ/ν ) 1.75, which corresponds to the exact value of the two-dimensional Ising model.

phase transition observed here. The reason is that the phase transition occurring in the systems considered is an orderdisorder transition of the same type as that known from twodimensional lattice models on square lattices rather than the melting transition in two dimensions.7,52 Rather-convincing evidence that the order-disorder transition belongs to the universality class of the two-dimensional Ising model provides the finite size scaling analysis of the susceptibility of the order parameter |Fq(1)|. Figure 8 shows the plots of χFq(1) versus temperature obtained for the system with σ* ) 1.20 and different simulation cell sizes. Quite-similar results have been obtained for other system with σ* ∈ [1.22,1.32]. In the inset of Figure 8, we present the log-log plots of the maximum of χFq(1) (χFq(1),max(L)) versus the simulation cell size, which demonstrates that in all cases the results are consistent with that predicted for the two-dimensional Ising model value of the ratio of the critical exponents γ/ν ) 1.7551 (with γ and ν being the susceptibility and the correlation length critical exponents, respectively) and resulting from the finite size scaling theory relation40

χmax ∝ Lγ/ν

(9)

It should be emphasized that all of the lines given in the inset of Figure 8 have been drawn assuming the exact value of the critical exponents ratio. The best fit to the data have given the values of γ/ν ) 1.75 ( 0.02. Also, the plots of the fourthorder cumulants of the order parameter |Fq(1)| for different sizes of the simulation cell have fixed points are very close to the values U* ) 0.61 ( 0.01, consistent with the value for the twodimensional Ising model.53 Of course, to study the orderdisorder phase transition of the c(2 × 2) phase one could also use the order parameter defined as the difference in the densities of adatoms on different sublattices, just as it is done in studies of lattice gas models.7 On the other hand, the susceptibility of the bond-orientational order parameter Ψ4(1) does not scale with the simulation-cell size in a way characteristic of the second-order phase transition (see parts a and b of Figure 9). The observed behavior is a consequence of the fact that the liquid-like phase possesses a partial square order and the adsorbed atoms remain partially

Figure 9. Parts a and b present the plots of the susceptibility of the bond-orientational order parameter |Ψ4(1)| vs temperature for the systems characterized by different σ* valuesand different sizes of the simulation cell (shown in the figure). Parts c and d show the plots of the order parameter |Fq(1)| (part c) and the bond-orientational order parameter |Ψ4(1)| (part d) vs temperature, obtained for different sizes of the simulation cell (shown in the figure) for the system with σ* ) 1.34. The inset of part c shows, in enlarged scale, the region in which the first-order phase transition takes place.

localized on adsorption sites of a square lattice. Therefore, the parameter Ψ4(1) is not a proper order parameter to describe the phase behavior of the systems considered here. A partial localization of adatoms at the temperatures above the phase transition point results from the rather-high corrugation of the surface potential and has been confirmed by the calculations of radial distribution functions. The system of σ* ) 1.34 undergoes two, rather than one, phase transitions. The first transition, which takes place at the temperature of T* ≈ 0.16, transforms the commensurate c(2 × 2) monolayer into the incommensurate solid-like phase, whereas the second phase transition, occurring at a considerably higher temperature of T* ≈ 0.71, leads to the formation of a disordered liquid-like phase. Part c of Figure 9 shows the temperature changes of the order parameter |Fq(1)| obtained for different sizes of the simulation cell. The inset shows, in enlarged scale, the region of the first transition. Part d of Figure 9 shows the behavior of the bond-orientational order parameter, Ψ4, over the temperature region in which the first-order phase transition occurs. It is seen that in small systems with L ) 12, 16 as well as 20 the order parameters change continuously, but in a way characteristic of the first-order transition,50,54,55 and the curves cross each other. The lack of hysteresis loops, often appearing in computer simulations of systems undergoing a first-order phase transition, results directly from the parallel tempering method applied here. Each simulation run is carried out over a certain temperature range, from below to above the transition point, and the configurations generated at each temperature are exchanged randomly. In the case of a sufficiently large simulation cell (L ) 40), we do observe the existence of metastable states and a well-developed hysteresis loop. One should also note that the scaling plots of the order-parameter susceptibilities (see Figure 10) are quite consistent with the scaling relation40

χmax ∝ Ld

(10)

15670 J. Phys. Chem. C, Vol. 111, No. 43, 2007

Figure 10. Plots of the susceptibilities of the order parameter |Fq(1)| (part a) and the bond-orientational order parameter |Ψ4(1)| (part b) vs temperature for the system characterized by σ* ) 1.34 and different sizes of the simulation cell (shown in the figure) in the region of the first-order phase transition. The insets give the log-log plots of the maximum values of χFq(1) (part a) and χΨ4(1) (part b) against the size of the simulation cell.

where d ) 2 is the system dimension, predicted to hold in the case of first-order phase transitions. Parts c and d of Figure 9 also show that the transition leads to a large drop in the magnitude of the bond-orientational order parameter, Ψ4(1), whereas the order parameter |Fq(1)| drops only from about to 0.97 to about to 0.85. This suggests that the adsorbed atoms in the high-temperature phase do not exhibit large displacements from the positions in the commensurate phase, whereas the orientational ordering changes significantly. The structure of that phase has been determined by a direct inspection of snapshots of configurations. Figure 11 presents two examples of configurations, both recorded at the temperature above the transition point and different sizes of the simulation cell. Part a of Figure 11 shows the configuration obtained for L ) 20, and part b of Figure 11 shows the configuration recorded for L ) 40. In both cases, one finds the same structure of decagonal order in which each adsorbed atom has five nearest neighbors. Such aperiodic structures are known to be characteristic of quasicrystals formed by two- and three-component

Patrykiejew and Sokołowski

Figure 11. Parts a and b show the snapshots of configurations showing the structure of the incommensurate monolayer phase for the system characterized by σ* ) 1.34 recorded for two different sizes of the simulation cell with L ) 20 (part a) and 40 (part b). The dashed lines in part a show the substrate surface lattice. For a better visualization of the structure, we have put the lines joining the nearest-neighbor atoms in parts of the depicted configurations.

alloys.56 Such forms of order have never been observed in onecomponent systems. A more-detailed study of the conditions under which such two-dimensional quasicrystals can exist will be published elsewhere. Upon the increase of temperature, we observe the continuous, presumably Ising-like, phase transition, that occurs at T* ≈ 0.71, leading to the formation of a disordered liquid-like phase. Figure 12 presents the plots of the susceptibility of the order parameter |Fq(1)| against the temperature obtained for different simulation cell sizes, and the inset shows the scaling plot, which supports the conclusion that the phase transition belongs to the universality class of the two-dimensional Ising model. Of course, the susceptibility of the bond-orientational order parameter, Ψ4(1), does not exhibit the finite size effects characteristic of the

Order-Disorder Phase Transitions

Figure 12. Plots of the susceptibilities of the order parameter |Fq(1)| vs temperature for the system characterized by σ* ) 1.34 and different sizes of the simulation cell (shown in the figure) in the region of the second-order phase transition. The inset gives the log-log plot of the maximum values of χFq(1) against the size of the simulation cell.

second-order phase transition, the same as the previously discussed systems of lower σ*. Our claim that the high-temperature phase transition is Isinglike is based on the observation that the structure of the incommensurate structure disorders at the temperatures below the Ising-like transition point. The resulting phase retains a considerable degree of square ordering because of the corrugation potential, and hence it is expected to undergo an Ising-like transition. Of course, the O-D Ising-like transition could be also studied using the order parameter defined in terms of the sublattice densities. It seems, however, that the order parameter (|Fq(1)|) used here is well-suited to determine the location and nature of this phase transition. IV. Bilayers The bilayer films formed by small atoms of σ* ∈ [0.82,0.93] ordering into the (1 × 1) structure in each layer behave in a way similar to monolayers, and each layer disorders gradually when the temperature increases. Figure 13 presents the temperature changes of the order parameter |Fq(n)| for both layers, n ) 1 and 2 (part a), of the conjugated susceptibility, χFq(n) (part b), calculated separately for each layer, and of the heat capacity (part c), obtained for two systems characterized by σ* ) 0.86 and 0.92. The results for the systems with other σ* values between 0.82 and 0.93 are qualitatively the same and hence not shown here. As expected, no appreciable finite size effects are observed, the same as in the monolayer films discussed previously. In all cases, the stability of the first layer, adjacent to the solid surface, appears to be enhanced by the presence of the second layer, whereas the disordering of the second layer occurs at temperatures lower than the disordering of monolayer films. Part d of Figure 13 presents the plots of the temperature at which the susceptibility of the order parameter |Fq(n)| reaches a maximum in layers 1 and 2 versus the misfit between the adsorbate and the adsorbent, and together with the results given in Figure 6 supports the above statement. The results obtained for the bilayer film formed by the atoms of σ* ) 0.94 have demonstrated that the second layer undergoes a sharp first-order phase transition at the temperature of about 0.14, that is, below the first-order transition temperature found

J. Phys. Chem. C, Vol. 111, No. 43, 2007 15671

Figure 13. Plots of the order parameter |Fq(n)| (part a) of the susceptibility χFq(n) (part b) and the heat capacity (part c) against temperature, calculated separately for each layer, n ) 1 and 2, and obtained for bilayers formed by atoms of different values of σ* and of different sizes of the simulation cell (shown in the figure). Part d gives the plots of temperatures at which the susceptibilities of the order parameters |Fq(n) (n ) 1,2) reach a maximum vs the misfit between the adsorbate and the adsorbent for bilayer films ordering into the commensurate (1 × 1) structure at low temperatures. Note that the result for the system with σ* ) 0.94 and n ) 2 corresponds to the location of the first-order phase transition.

in the case of the monolayer film, whereas the first layer disorders at considerably higher temperature of about 0.42, and the disordering is a continuous process. This shows that the presence of the second layer considerably affects the properties of the first layer and stabilizes the commensurate phase, even when the second layer is already disordered. The inspection of configurations recorded below and above the temperature at which the second layer undergoes the phase transition has shown that the transition in the second layer is of the same nature as the commensurate-incommensurate phase transition observed in the monolayer film. In particular, the structure of the incommensurate phase is like that depicted in Figure 5 and consists of commensurate domains separated by walls in which the adsorbate atoms exhibit hexagonal ordering. Alternatively, the first layer does not undergo the commensurate-incommensurate phase transition, but rather disorders gradually as the temperature increases. The above discussion is illustrated by the plots of the order parameter |Fq(n)| (n ) 1, 2) given in Figure 14. The mechanism of disordering of bilayer films formed by atoms of σ* ∈ [1.20,1.34] is different and also changes with σ*. We begin the presentation of results with those obtained for the system of adatoms characterized by σ* ) 1.20. A perfectly ordered bilayer commensurate film has been prepared as a starting configuration. Then the annealing experiment has been carried out by performing Monte Carlo simulation over a wide range of temperatures, starting from T* ) 0.1 up to T* ) 0.5. Next, the freezing run has been carried out by decreasing the temperature from 0.5 back to 0.1. It occurred that during the annealing run both layers disorder at the same temperature of about 0.41 and a liquid-likely bilayer is formed. The freezing run, however, has not led to the restoration of the commensurate bilayer structure, but to the incommensurate phase in which the density of the first layer is considerably higher than the density of the second layer. Figure 15 presents the temperature changes

15672 J. Phys. Chem. C, Vol. 111, No. 43, 2007

Figure 14. Temperature changes of the order parameter |Fq(n)| (n ) 1and 2) for the bilayer film formed by atoms of σ* ) 0.94. The results prsented here were obtained for the simulation cell of L ) 16. The dashed vertical line marks the first-order phase transition that takes place in the second layer.

Figure 15. Temperature changes of the order parameter |Fq(n)| (n ) 1and 2) during the annealing and freezing runs for the bilayer formed by atoms of σ* ) 1.20. The starting configuration has been a perfect commensurate bilayer of the surface density Fs ) 1.0. The inset shows the configuration recorded at T* ) 0.2 during the freezing run. The atoms located in the first and second layer are shown as open and filled circles.

of the order parameter |Fq(n)| (n ) 1, 2) obtained during the annealing and freezing runs and shows that the bilayer liquid freezes into incommensurate structure. It should be noted that this incommensurate phase exhibits neither square nor hexagonal ordering. The inset of Figure 15 shows a snapshot recorded at T* ) 0.2 and shows that the second layer is only partially filled. The above results suggest that the commensurate phase is not the only stable state for the monolayer film and that the bilayer film does not order into the commensurate c(2 × 2) phase. To verify this prediction, we have performed the grandcanonical ensemble Monte Carlo simulation at two different temperatures of T* ) 0.4 and 0.6. Figure 16 shows the adsorption-desorption isotherms obtained and demonstrates that the monolayer film undergoes a first-order commensurateincommensurate phase transition at T* ) 0.4, whereas the isotherm recorded at the higher temperature of T* ) 0.6 does not show any trace of that transition. The density profile shown

Patrykiejew and Sokołowski

Figure 16. Adsorption-desorption isotherms recorded at two different temperatures (shown in the figure) for the system characterized by σ* ) 1.20. The insets show the density profiles recorded at µ* ) -6.0 and at two different temperatures.

in the inset of Figure 16, both recorded at the same value of the chemical potential, equal to µ* ) - 6.0, demonstrates that at T* ) 0.4 the density is equal to 1.083 and practically all particles are located at the first layer, whereas at T* ) 0.6 the density is lower and equal to 1.002, but some of the adsorbed atoms are already promoted to the second layer. This result also supports our observation of the commensurate-incommensurate transition, which occurs in the first layer. An increase of σ* to 1.22 leads to quite-different behavior because the commensurate bilayer film becomes a stable state at low temperatures. In particular, the annealing and freezing runs show that the disordering of the commensurate bilayer occurs via a single phase transition. The transition temperature, which is equal to about T* ) 0.53, is, however, considerably lower than the temperature of the order-disorder transition observed in the monolayer film, which takes place at T* ≈ 0.6. Moreover, the order-disorder transition occurring in the monolayer film has been demonstrated to be a continuous phase transition belonging to the universality class of the twodimensional Ising model, whereas the bilayer undergoes a firstorder phase transition (see Figure 17). In particular, we have found that the susceptibility of the order parameter |Fq(n)| scales with the system size in agreement with the law given by eq 10, as demonstrated by the results given in inset a to Figure 17. Inset b of Figure 17 shows the changes in the densities of both layers (evaluated by the integration of density profiles), recorded for the simulation cell of the size 30 × 30 and demonstrates that the phase transition is accompanied by an increase (decrease) of the density in the first (second) layer. Thus, the transition leads to the transfer of adatoms from the second to the first layer. The reader should note that the density of the first layer reaches a value of about 0.54 at the temperature above the transition point. A further increase of temperature above 0.6 leads to a gradual decrease of the first layer density because of desorption. Quite-similar results have been obtained for the system with σ* ) 1.24 and 1.26. The scaling plots for the order parameter |Fq(n)| susceptibility have demonstrated that eq 10 is satisfied quite well. In particular, for σ* ) 1.24 the exponents for the first and the second layer have been found to be equal to 2.06 and 2.08, respectively, and for σ* ) 1.26 we have obtained the values equal to 1.95 and 2.03. Although the bilayer systems are not strictly two-dimensional, the development of fluctuations in the direction perpendicular to the surface are suppressed and

Order-Disorder Phase Transitions

Figure 17. Temperature changes of the order parameter |Fq(n)| for the first (n ) 1) and second (n ) 2) layer for a bilayer film formed by the adsorbate atoms of σ* ) 1.22 obtained for different sizes of the simulation cell. Inset a shows the scaling plot of the susceptibility χFq(n), and inset b presents the layer densities as functions of the temperature.

Figure 18. Scaling plots of the susceptibility χFq(n) for bilayers formed by adsorbate atoms of different sizes (shown in the figure). The numbers at the data give the slopes of the fit to the power law.

the correlation length cannot exceed the interlayer distance. Therefore, these systems are bound to behave as effectively twodimensional. In the case of still-larger adsorbate atoms of σ* ) 1.28 as well as 1.30, the mechanism of the disordering changes and instead of a first-order transition we have found a sequence of two continuous transitions. When σ* ) 1.28, these transitions occur at the temperatures of about T/2 ) 0.78 in the second layer and at T/1 ) 0.795 in the first layer, whereas in the case of σ* ) 1.30 the transitions occur at T/2 ) 0.80 and at T/1 ) 0.81. In both systems, the order-disorder transitions have been found to belong to the universality class of the two-dimensional Ising model. In particular, the susceptibilities of the order parameters |Fq(1)| and |Fq(2)| scale with the system size as predicted by the eq 9. Figure 18 presents a comparison of the appropriate scaling plots obtained for the systems with σ* ) 1.26, 1.28, and 1.30 and demonstrates that the behavior of the systems with σ* ) 1.28 and 1.30 is quite consistent with the assumption that γ/ν ) 1.75. One has to ask the question why the bilayers formed by atoms of σ* between 1.22 and 1.26 disorder via a single discontinuous (first-order) phase transition while the films formed by larger atoms of σ* ) 1.28 and 1.30 disorder in a layer-by-layer mode via two continuous phase transitions? A plausible answer to this question can be proposed by considering the changes in the properties of adsorption isotherms and the structure of such

J. Phys. Chem. C, Vol. 111, No. 43, 2007 15673

Figure 19. Adsorption-desorption isotherms obtained for the system characterized by σ* ) 1.22 at different temperatures (shown in the figure). The numbers show the values of the order parameter |Fq(1)| in different phases. One should note that in the case of bilayer films the order parameter |Fq(2)| assumes values very close to the corresponding values of |Fq(1)|.

bilayer films at different temperatures. From the adsorptiondesorption isotherms calculated using the grand-canonical ensemble Monte Carlo simulation, for the system with σ* ) 1.22 (see Figure 19), it follows that at low temperatures, below 0.4, the condensation of the second layer on top of the commensurate monolayer leads to the formation of the commensurate bilayer structure, of the total number surface density equal to about 1.0. At the temperature T* ) 0.4, the formation of the second layer on top of the commensurate monolayer is accompanied by the transformation of the first layer into a denser, incommensurate structure. Also, the second layer assumes a similar incommensurate structure, and both layers reach a density of about 0.54 at the onset of the third-layer condensation. At still higher temperatures of 0.5 and 0.55, the commensurate-incommensurate transition has been found to occur already in the first adsorbed layer, prior to the condensation of the second layer. The above-described changes in the behavior of adsorption isotherms are confirmed by the observed changes of the order parameter |Fq(1)| in the first (n ) 1) and second (n ) 2) layer at the phase transition points, which are also shown in Figure 19. In particular, the isotherm at T* ) 0.4 demonstrates that the condensation of the second layer leads to a sudden decrease of the order parameter |Fq(1)| from about 0.92 to about 0.42, and the order parameter |Fq(2)| also takes upon the value close to 0.42. This clearly indicates the presence of the C-IC phase transition induced by the formation of the second layer. The canonical ensemble simulation, carried out at the density corresponding to a completely filled commensurate bilayer, shows that the film disorders at a temperature of about 0.53, whereas the simulation in the grand-canonical ensemble demonstrates that already at T* ) 0.4 the bilayer forms an incommensurate phase of higher density. This observation can be explained as follows. The canonical ensemble simulation artificially stabilizes the commensurate structure. The transformation to the incommensurate phase would require the formation of incompletely filled incommensurate bilayer phase, the same as that observed in the case of the system with σ* ) 1.20. However, in the case of larger adsorbate atoms of σ ) 1.22, as well as 1.24 and 1.26, the stability of the commensurate phase is higher and hence the formation of an incompletely filled incommensurate structure is more difficult. The mechanism of the disordering of the second layer in bilayer films formed by atoms of σ* ) 1.32, 1.33, and 1.34 is,

15674 J. Phys. Chem. C, Vol. 111, No. 43, 2007

Patrykiejew and Sokołowski

Figure 20. Part a shows a snapshot of the top-layer configuration for the system characterized by σ* ) 1.34 at T* ) 0.1. The atoms within walls are shown as black circles. Part b is the local density profile recorded at the same temperature.

at first sight, quite similar to the behavior of the second layer in the already discussed monolayer film formed by atoms of σ* ) 1.34 (cf. parts c and d of Figure 9). Thus, during the annealing runs, in which the starting configuration was assumed to be a perfectly ordered commensurate bilayer, we have observed discontinuous jumps of the order parameters, |Fq(2)| and Ψ4(2). These discontinuities occur at rather-low temperatures equal to about 0.19 when σ* ) 1.32, about 0.14 when σ* ) 1.33, and about 0.05 when σ* ) 1.34. The jumps of the order parameters suggest the presence of a first-order phase transition. The inspection of snapshots has demonstrated that the hightemperature phase in the second layer has a quite-different structure from the incommensurate monolayer, which exhibits decagonal ordering (cf. Figure 11). In the present case, the second layer assumes the structure in which stripes of the commensurate domains are separated by highly localized walls running along the direction diagonal to the surface lattice symmetry axes (see part a of Figure 20). The atoms in the walls are located directly above the saddle points between a pair of atoms from the layer adjacent to the surface and hence at considerably larger distance from the surface than the atoms forming the commensurate domains, located over the sites formed by four atoms from the first layer (see Figure 20b). The atoms within stripes are also displaced from the centers of surface potential minima by a displacement vector τk ) (∆xk, ∆yk) with ∆xk ) ∆yk, depending on the distance from the wall k (see Figure 21 a). The length of this displacement vector decreases with distance from the wall separating commensurate domains (see Figure 21b). Freezing runs have shown that the commensurate structure is not recovered even at very-low temperatures. Moreover, the ground-state energies of bilayer films of stripe-like structure have been found to be lower than those of the commensurate bilayer films. For example, when σ* ) 1.34 these two energies are equal to -9.56 and -9.54, respectively, and though the difference between them is not large it is nevertheless sufficient to make the film with the stripelike second layer more stable. The results obtained for different values of σ* (equal to 1.32, 1.33, and 1.34) have shown that the distance between adjacent walls decreases when σ* becomes larger. It should be emphasized that the top-layer density remains equal to the density of the perfectly ordered commensurate layer. All of these facts allow us to conclude that the appearance of the stripe-like structure in the second layer is induced by the strain resulting from the a large misfit. The commensurate c(2 × 2) structure in the second layer is not a stable phase even at the ground state. Therefore, the discontinuous jumps of the order parameters do not represent a true phase transition but are due

Figure 21. Part a shows the positions of atoms in the second layer with respect to the atoms in the first layer. Part b shows the length of the displacement vector τk plotted against the distance from the wall (k) obtained for the system with σ* ) 1.33 at T* ) 0.05 and the simulation box size L ) 40.

to the relaxation from the metastable commensurate structure, taken as a starting configuration, to the stable stripe-like structure. One should note that the presence of the second layer causes the first layer to retain a commensurate c(2 × 2) structure up to rather-high temperatures. As a consequence, the first layer does not undergo the C-IC transition at all. Instead, it exhibits a continuous Ising-like order-disorder transition to the disordered phase, as proven by the finite size scaling of the susceptibility of the order parameter |Fq(1)|. Namely, we have found that the scaling relation (9) is satisfied and γ/ν ) 1.75. V. Summary and Concluding Remarks We have performed extensive Monte Carlo simulation study of disordering of the commensurate monolayer and bilayer films formed by atomic Lennard-Jones fluids on the (100) face of an fcc crystal. In particular, we have been interested in the changes in the mechanism of the order-disorder transitions taking place in the films characterized by different misfit between the size of the adsorbed atoms and the size of the surface lattice unit cell. The results reported have demonstrated that in the case of systems that order into the (1 × 1) phase at low temperatures, the disordering of the ordered structure occurs gradually and that the stability of the commensurate film depends primarily on the misfit between the adsorbate and the substrate. In particular, the commensurate films exhibit the highest stability for small positive natural misfit. In the case of a sufficiently large misfit (σ* ) 0.94 and I ) -0.02971), the commensurate monolayers have been found to undergo a discontinuous (first-order) phase transition to the incommensurate phase of the domain-wall structure in which rectangular commensurate domains are separated by walls. Within the walls, adatoms are hexagonally packed. Of course, the formation of such incommensurate structures is possible only when the height of potential barrier between adjacent surface potential minima is low enough. We have performed Monte Carlo calculations for different values of the parameter /gs, lower as well as higher than 1.0,58 and found that an increase of /gs to 1.1 stabilizes the commensurate structure and excludes the possibility of the formation of incommensurate structure for the monolayer film formed by atoms of σ* ) 0.94.

Order-Disorder Phase Transitions Alternatively, when /gs is lowered to 0.9 the stability of the commensurate film is reduced significantly and already at the temperature of T* ) 0.06 the incommensurate phase has been observed to appear. In the case of bilayer films ordering into the (1 × 1) structure in each layer, the disordering occurs in a layer-by-layer mode. In general, the presence of the second layer considerably stabilizes the commensurate structure in the first layer. In particular, the bilayer film formed by atoms of σ* ) 0.94 does not exhibit the commensurate-incommensurate phase transition in the first layer at all, whereas a similar transition takes place in the second layer (cf. Figure 14). When the adsorbed atoms are large enough to exclude mutual occupation of the adjacent potential minima and the adsorbed monolayers as well as bilayers form the c(2 × 2) commensurate structure, the mechanism of the disordering changes and depends on the misfit between the adsorbate and the substrate lattice. In general, when the misfit is small, the commensurate monolayer films exhibit a continuous order-disorder transition belonging to the universality class of the twodimensional Ising model, as demonstrated by the finite size scaling analysis for the simulation data. In the case of a sufficiently large negative misfit, the commensurate film undergoes a first-order phase transition to an incommensurate phase of decagonal order in which each adsorbed atom has five first-nearest neighbors (cf. Figure 11). Of course, the stability of the commensurate structure depends on the misfit between the adsorbate and the adsorbent and, similar to the case of monolayers ordering into the (1 × 1) structure. The orderdisorder transition temperature reaches a maximum value when the misfit equals about -0.02. The behavior of bilayer films ordering into the c(2 × 2) structure at low temperatures shows some peculiarities associated with the appearance of the C-IC phase transitions observed in the systems with sufficiently large positive misfits, exceeding about 0.02. Namely, instead of a sequence of two Ising-like phase transitions occurring in each layer we observe a single first-order phase transition involving both layers. It has been shown by grand-canonical ensemble MC simulation that the bilayer film does undergo the C-IC phase transition, which involves both layers. Alternatively, when the misfit is sufficiently small (I ∈ (-0.02, 0.02)) the bilayer films undergo two second-order phase transitions, both belonging to the universality class of the two-dimensional Ising model. This can be attributed to a particularly high stability of commensurate structures. Of course, the corrugation potential felt by the second layer is lower than that experienced by the atoms located in the first layer, and hence the second layer disorders at a lower temperature than the first layer. When the negative misfit exceeded about -0.02, the second layer was found to exhibit a stripe-like structure (cf. Figure 20a). The stripe-like phase shows the presence of axial modulations along the diagonal of the surface lattice symmetry but retains the density of the commensurate c(2 × 2) phase. One can expect that in the case of thicker films similar stripe-like structures may appear in the top layers, whereas the remaining layers should retain the commensurate structure. The above hypothesis is supported by our results, which have demonstrated that the presence of the second layer stabilizes the commensurate structure in the first layer, even at the temperatures at which the second layer is already disordered. Of course, the properties of the first layer are much-more affected by the surface corrugation potential than the higher layers because the surface potential, given by eq 3, decreases with the distance from the

J. Phys. Chem. C, Vol. 111, No. 43, 2007 15675 surface as z*-4, whereas the Fourier coefficients Vg(z*) for g * 0 decay exponentially with z*. The disordering and melting of thick adsorbed films, with the number of layers up to 12, will be discussed in our future work. Acknowledgment. This work has been partially supported by the European Community under grant no. MTDK-CT-2004509249. References and Notes (1) Dash, J. G. Films on Solid Surfaces; Academic Press: New York, 1975. (2) Ordering in Two Dimensions; Sinha, S. K, Ed.; North-Holland: Amsterdam, 1980. (3) Phase Transitions in Surface Films; Dash, J. G., Ruvalds, R., Eds.; Plenum: New York, 1980. (4) Binder, K.; Landau, D. P. In AdV. Chem. Phys.; Lawley, K. P., Ed.; Wiley: Chichester, 1989; p 91. (5) Phase Transitions in Surface Films 2; Taub, H., Torzo, G., Lauter, H. J., Fain, S. C., Jr., Eds.; Plenum: New York, 1991. (6) Bruch, L. W.; Cole, M. W.; Zaremba, E. Physical Adsorption: Forces and Phenomena; Clarendon Press: Oxford, 1997. (7) Patrykiejew, A.; Sokołowski, S.; Binder, K. Surf. Sci. Rep. 2000, 37, 207. (8) Cohen, P. I.; Unguris, J.; Webb, M. B. Surf. Sci. 1976, 58, 429. (9) Unguris, J.; Bruch, L. W.; Moog, E. R.; Webb, M. B. Surf. Sci, 1979, 87, 415. (10) Grunze, M.; Kleban, P. H.; Unertl, W. N.; Rys, F. S. Phys. ReV. Lett.1983, 51 582. (11) Glachant, A.; Jaubert, M.; Bienfait, M.; Boato, G. Surf. Sci. 1981, 115, 219. (12) Moog, E. R.; Webb, M. B. Surf. Sci. 198, 148, 338. (13) Rosmeyer, C.; Girardet, C.; Zeppenfeld, P.; George, K.; Bu¨chel, M.; Comsa, G. Surf. Sci. 1994, 313, 251. (14) Thomy, A.; Duval, X.; Regnier, J. Surf. Sci. Rep. 1981, 1, 1. (15) Venables, J. A.; Sequin, J. L.; Suzanne, J.; Bienfait, M. Surf. Sci. 1984, 145, 345. (16) Marx, R.; Wassermann, E. F. Surf. Sci. 1982, 117, 267. (17) Zhang, Q. M.; Kim, H. K.; Chan, M. H. W. Phys. ReV. B 1986, 33, 5149. (18) Gangwar, R.; Colella, N. J.; Suter, R. M. Phys. ReV. B 1989, 39, 2459. (19) D’Amico, K. L.; Bohr, J.; Moncton, D. E.; Gibbs, D. Phys. ReV. B 1990, 41, 4368. (20) Enault, A.; Larher, Y. Surf. Sci. 1977, 62, 233. (21) Nardon, Y.; Larher, Y. Surf. Sci. 1974, 42 299. (22) Matecki, M.; Thomy, A.; Duval, X. J. Chim. Phys. 1974, 71, 1484. (23) Larher, Y. In Surface Properties of Layered Structures; Benedek, G., Ed.; Kluwer: Dortrecht, 1992; p 261. (24) Coulomb, J. P.; Sullivan, T. S.; Vilches, O. E. Phys. ReV. B 1084, 30, 4753 (25) Jordan, J. L.; McTague, J. P.; Hastings, J. B.; Passell, L. Surf. Sci. Lett. 1985, 150, L82. (26) Sidoumon, M.; Angot, T.; Suzanne, J. Surf. Sci. 1992, 272, 347. (27) Patrykiejew, A.; Borowski, P. Phys. ReV. B. 1990, 42, 4670. (28) Bak, P. Rep. Prog. Phys. 1982, 45, 587. (29) den Nijs, N. In Phase Transitions and Critical Phenomena; Domb, C., Lebowitz, J. L., Eds.; Academic Press: London, 1988; Vol. 12, p 219. (30) Patrykiejew, A.; Sokołowski, S.; Binder, K. Surf. Sci. 2002, 512, 1. (31) Morishige, K J. Chem. Phys. 1991, 95, 2667. (32) Patrykiejew, A.; Sokołowski, S.; Zientarski, T.; Binder, K. Surf. Sci. 1999, 421, 308. (33) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications; Academic Press: San Diego, CA, 1996. (34) Steele, W. A. Surf. Sci. 1973, 36, 317. (35) Broughton, J. Q.; Gilmer, G. H. J. Chem. Phys. 1983, 79, 5095. (36) Salamacha, L.; Patrykiejew, A.; Sokołowski, S.; Binder, K. J. Chem. Phys. 2004, 120, 1017. (37) Salamacha, L.; Patrykiejew, A.; Sokołowski, S.; Eur. Phys. J. E 2006, 18, 425. (38) Patrykiejew, A.; Sokołowski, S. J. Chem. Phys. 2006, 124, 194705. (39) Allen M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, 1987. (40) Landau, D.P.; Binder, K. A Guide to Monte Carlo Simulation in Statistical Physics; Cambridge University Press: Cambridge, 2000. (41) Tesi, M. C.; Janse van Rensburg, E. J.; Orlandini, E.; Whittington, S. G. J. Stat. Phys. 1996, 82, 155.

15676 J. Phys. Chem. C, Vol. 111, No. 43, 2007 (42) Hansmann, U. H. E. Chem. Phys. Lett. 1997, 281, 140. (43) Vishnyakov, A.; Neimark, A. V. J. Chem. Phys. 2003, 118, 7585. (44) Salamacha, L.; Patrykiejew, A.; Sokołowski, S.; Binder, K. Eur. Phys. J. E 2004, 13, 261. (45) Strandburg, K. J. ReV. Mod. Phys. 1988, 60, 161. (46) Patrykiejew, A.; Zientarski, T.; Binder, K. Acta Phys. Polon. 1996, 89, 735. (47) Rappaport, N. J. The Art of Molecular Dynamics Simulation; Cambridge Univ. Press: Cambridge, 1995. (48) Tartaglino, U.; Tosatti, E. Surf. Sci. 2003, 532-535, 623. (49) Binder, K. In Cohesion and Structure of Surfaces; De Boer, F. R., Pettifor, D. G., Eds.; Elsevier: Amsterdam, 1995; Vol. 4. (50) Rzysko, W.; Patrykiejew, A.; Binder, K. Phys. ReV. B 2005, 72, 165416. (51) Binder, K.; Landau, D. P. Phys. ReV. B 1980, 21, 1941.

Patrykiejew and Sokołowski (52) Nelson, D. R. In Phase Transitions and Critical Phenomena; Domb, C., Lebowitz, J. L., Eds.; Academic Press: London, 1988; Vol. 7, p1 (53) Kamieniarz, G.; Blo¨te, H. W. J. J. Phys. 1993, 26, 201. (54) Ferrenberg, A. M.; Landau, D. P.; Binder, K. Phys. ReV. E 1998, 58, 3353. (55) Vollmayr, K.; Reger, J. D.; Scheucher, M.; Binder, K. Z. Phys. B: Condens. Matter 1993, 91, 113. (56) Quasicrystals: An Introduction to Structure, Physical Properties and Applications; Suck, J.-B., Schreiber, M., Ha¨usller, P., Eds.; Springer: Berlin, 2002. (57) Patrykiejew, A.; Sokołowski, S. To be submitted for publication. (58) The potential barrier between adjacent sites, V/D, is proportional to / gs.