Dynamics of a monolayer of nitrogen physisorbed on graphite

Dynamics of a monolayer of nitrogen physisorbed on graphite

231 Surface Science 154 (1985) 231-253 North-Holland, Amsterdam DYNAMICS GRAPHITE Gianni CARDINI Department Received OF A MONOLAYER and Seamus ...

1MB Sizes 0 Downloads 56 Views

231

Surface Science 154 (1985) 231-253 North-Holland, Amsterdam

DYNAMICS GRAPHITE Gianni

CARDINI

Department Received

OF A MONOLAYER

and Seamus

of Chemistv, 10 September

OF NITROGEN

ON

F. O’SHEA *

University of Lethbridge, 1984; accepted

PHYSISORBED

Lethbridge,

for publication

Alberta, Canada TlK

30 November

3M4

1984

The dynamics of a crystalline monolayer of molecular nitrogen physisorbed on the basal plane of graphite, which was assumed to be rigid, has been studied by means of lattice and molecular dynamics. The study focusses mainly on the registered solid in which the centres of mass adopt the fi xfiR30” structure, but some results are also reported for denser monlayers, including one where the centre-of-mass lattice is non-triangular. Harmonic dispersion curves are reported for two models of the N,-N, interaction combined with Steele’s representation of the molecule-surface interaction. The one-phonon density of states is given for one of the models. Molecular dynamics results for the same model at two temperatures, 5 and 17 K, are analysed in terms of the dynamical structure factor, and one-phonon approximation to it, and the out-of-plane motions are treated in a manner similar to the one-phonon approximation. In-plane and out-of-plane motions are largely decoupled in all cases studied, and the in-plane motions show evidence of marked anharmonicity. The out-of-plane motions appear to exhibit Fermi-Pasta-Ulam type recurrences, although this is not analysed in detail.

1. In~oduction The study of physical adsorption of gases on solids [1,2] has grown rapidly in recent years, partly due to improvements in equipment, partly due to technological interest, but also because adsorbed systems are now recognized as exhibiting a very rich phenomenology [3]. The use of well characterized, relatively homogeneous substrates, such as exfoliated graphites [4], combined with simple adsorbates has made it possible to characterize the phases of adsorbates, and to clarify important concepts arising as a result. For the most part these studies have been based on the adsorption of the rare gases on graphites [3]. Recently reported experiments [5,6] have shown that it is now possible, by inelastic scattering of beams of light atoms, to determine vibrational dispersion curves for monolayers and multilayers physically adsorbed on surfaces.

* To whom correspondence

should

be addressed.

0039-6028/85/$03.30 @ Elsevier Science Publishers (North-Holland Physics Publishing Division)

B,V.

232

G. Cardmi, S.F. O’Shea / Nitrogen physlsorbed on graphite

The use of simple molecules, such as nitrogen, as adsorbates leads to a further increase in the number of phases which can be detected. Solid-solid phase transitions involving a change in rotational ordering were first detected for nitrogen on graphite [7-lo], and it has since been established that these order-disorder transitions can occur both for commensurate and incommensurate solids [ll]. Theoretical work [10,12-171 has focussed on determining the nature of the ordered solids, and on the nature of the rotational transition. Only one calculation [12] is reported on the lattice dynamics of such an overlayer, and a recent simulation study [10,16] has examined the rotational dynamics in the system in the region of the solid-solid transition. The possibility that scattering experiments of the type [5,6] mentioned above could be used to determine the vibrational dynamics of nitrogen overlayers makes it worthwhile to consider the vibrational behaviour over a range of temperatures. The work reported here deals mainly with very low temperatures, less than about half the transition temperature, and complements the simulation work [16] which considered the dynamics mainly at higher temperatures. Two important problems arise in studying the vibrations of the overlayer, one dealing with the interaction potentials to be used, and one with the dynamical methods. The problem of potentials is twofold, since the molecule-molecule interactions and moleculeesurface interactions need to be addressed separately. As the rotational dynamics are central to our concern the molecule-surface potential needs to be one which represents the variation in energy with molecular orientation over the surface. The only potential of this kind reported in the literature is due to Steele [18], and uses an atom-atom representation picturing the molecule as a diatomic and the solid as an array of atoms. The atom-atom interactions are of the Lennard-Jones 12-6 type with u = 3.36 A and E/K = 31. The MD work of Talbot et al. [16] used this potential and a model for the nitrogen-nitrogen interaction due to Chung and Powles (CP) [19], which also involves an atom-atom LJ model with u = 3.32 A and c/K = 35, combined with a molecular charge distribution approximated by three point charges. The potential used here is the Sq model of Murthy et al. [20], combining atom-atom LJ (u = 3.318 A, c/K = 35.6) core interactions with a five-charge electrostatic model. Negative charges of 4.0469 electron units are located on the molecular axis 0.6527 A on either side of the centre, positive charges of 5.2366 reside on each nucleus, and a neutralizing negative charge of 2.3794 is located on the centre itself. Detailed quantum mechanical calculations of the molecular charge distribution were used to parameterize the point charge approximation, and the LJ terms were adapted from a closely related empirical potential devised [21] to model the low temperature solid. The extra detail in the electrostatics leads to better agreement between the calculated and experimental librational frequencies for the low temperature, low pressure phase of solid nitrogen ((Y-N2) while maintaining the already satisfactory translational frequencies [20]. The CP model was derived to fit the

G. Cardini, SF. O’Shea / Niirogen physisorbed on graphite

233

properties of the fluid, and has been used successfully for this purpose. In our lattice dynamical calculations for solid (Y-N~we find that the CP model gives a good account of the translational modes, but gives consistently high librational frequencies, a defect it shares with many other simple models of the nitrogen-nitrogen interaction ]20]. We chose to use the 5q model in our simulations, but calculated the lattice dynamics of overlayers using both the CP and 5q models. When the layer is held in registry with the graphite substrate the difference between the models should be smaller for the overlayer than for the three dimensional solid [22]. Steele’s molecule-surface potential was used throughout, in its explicit form for the corrugated surface, but also as the basis for a Fourier expansion of the potential 118,231, the leading (4 = 0) term of which gives the energy of interaction with an equivalent flat surface. Throughout this work we use the rigid substrate approximation which treats the graphite as an array of atoms fixed in their equilibrium positions. Although this is not by any means perfect [24,22], it has been used in all previous studies of this system, and the alternative, which involves simultaneously treating the overlayer and substrate dynamics, is vastly more complicated [24]. With the present potentials and experimental data we do not feel that the increased effort is justified, and have elected to retain this approximation. The predictions of the lattice dynamical calculations are expected to be most accurate in the limit of zero temperature, and to become less accurate with increasing temperature [25]. The molecular dynamics method is classical but fully anharmonic, and should show clearly the increasing anharmonicity as the temperature increases [26]. The method of lattice dynamics has advantages which partly offset the restriction to harmonic behaviour. These include the explicit representation of the vibrations in terms of displacement coordinates, and the ease with which the natural symmetry of the dynamical system can be used to advantage during the calculations [25]. In the present case the symmetry proved to be of unusual interest. When studying the layer supported by a flat surface we found that the results exhibit an unexpected symmetry which simplifies the analysis and also provides useful insights into the corrugated case. The vibrations fall into two symmetry types, those which are symmetric in the plane of the overlayer and those which are antisymmetric. The origin of this symmetry is dynamical rather than structural, and results from the use of the harmonic approximation. Unlike the true surface potential, which distinguishes up from down, the dynamical surface potential is quadratic in the displacement from the equilibrium height above the surface. This artifical symmetry, which holds for all q values, is lost when the corrugation is added because the corrugation introduces terms which couple lateral and vertical motions through the mixed second derivatives of the potential with respect to displacements. The MD results are more difficult to interpret. The quantity most frequently reported in these circumstances is the dynamical structure factor S(q, w), which is the power spectrum of the wavevector dependent

density &)

operator

[26-281

= Ce-@r,(‘),

q is the wavevector, r,(t) is the position of atom j at time t, and the sum is over all atoms in the system. S(q, o) depends on the composite motions in-plane because the two-dimensional character of q filters out contributions normal to the plane. Although there is some coupling of motions in-plane and out-of-plane it is expected to be weak enough that it will not give rise to measurable intensity in S(q, a), and this expectation seems to be fulfilled. S(q, w) provides useful information regarding inelastic scattering from the system, and hence is of considerable interest with a view to testing the calculations against experiment. Some motions contribute strongly to S(q, w), while others can give rise to negligible intensity. From this point of view S(q, w) is not necessarily as useful for comparisons with lattice dynamical calculations, especially since the assignment of the peaks in S(q, w) is often difficult. They can be assigned by comparing the frequencies with those of the although this is unreliable when the corresponding harmonic modes, anharmonic shifts are comparable with the separation between harmonic frequencies. The one-phonon approximation to S(q, w) [29] gives another route to assignments. Using a sampling function am, a combination of coordinates representative of the mode of interest, one can construct the operator $ =CW”([) n

e-@(K),

where (R,) is the mean position of the centre of molecule n during the calculation. This operator enhances the intensity of the phonon represented by w,,(t) relative to those of other phonons, and in the presence of weak to moderate anharmoni~ity this is often sufficient to correlate the peaks in S(q, w) with the harmonic phonons. For in-plane translations, for example, one can use w,(t)

= q*%(t),

for longitudinal

modes, and

w, = (z x q)*u,(t), for transverse modes. z is the unit normal displacement of the centre of molecule n approach can be extended [30] to display example, using w, = Z’U,.

to the surface, from its mean the out-of-plane

and in, is the position. This motions, for

One can use as librational coordinate the difference between the cosine of the appropriate angle at its instantaneous value and at its reference value. For motions in-plane the angle is that between the projection of the molecular axis

G. Cardrni, S. F. O’Shea / Nitrogen physisorbed on graphite

235

(a) Y

Fig. 1. The unit cell in real (a) and reciprocal (b) space. b =&a for the triangular centre-of-mass arrangement used here almost exclusively. The labelling of the special points in the reciprocal lattice follows that used in orthorhombic systems. The Cartesian coordinate axes are parallel to the symmetry glide lines.

on the plane and the x-axis (fig. l), and for out-of-plane motions the angle is that between the molecular axis and the normal to the surface; these are equivalent to the usual angular coordinates (0, $) of the spherical polar coordinate system.

2. Calculations The structure (fig. 1) was assumed to have the centre of mass on a triangular lattice, fi x fiR30°, with the molecules parallel to the basal plane, and the molecules on the two sublattices were assumed to be related by glide lines, as was found experimentally [ll] and in earlier calculations [12-15,171. The molecular axis directions were obtained by minimizing the energy subject to these constraints. This was done for a corrugated surface obtained by summing explicitly over the carbon atoms of the surface and interior of the graphite (z = 3.3238 A), and repeated for a flat surface (z = 3.3402 A), obtained using the q = 0 term in Steele’s expansion of the surface potential [18]. In the latter case we maintained the triangular condition (b = fia), and the calculations were done at the commensurate spacing (a = 4.26 A) and the lattice constant was also optimized. Results were obtained both for the 5g and CP potentials, each combined with the Steele [18] model for the molecule-surface interaction. Dispersion curves were constructed along the special directions r --, X -+ S -+ Y + r (fig. lb), and the vibrational density of states was calculated for the 5q model with the corrugated and flat (a = 4.26 p\) surfaces, using a uniform dense mesh of about 24000 points within the reduced Brillouin zone. For the flat surface, where the in-plane and out-of-plane modes belong to distinct symmetry types, we have separated their contributions to the density of states.

236

G. Cardini, S. F. O’Shea / Nitrogen physisorbed on graphite

For the simulations, a system of 126 molecules (9 x 7 unit cells each containing 2 molecules) with rectangular periodic boundary conditions was used, and data obtained at nominal temperatures of 5 and 15 K. The equations of motion were integrated in the microcanonical ( E, V, N) ensemble using the Verlet algorithm with constraints [31], and energy conservation was obtained to at least 5 significant figures. The system was first thermalized at 15 K for 3500 timesteps of 5 X lo-l5 s, and then 8200 timesteps were used ta accumulate data for subsequent analysis. The temperature average for the run was 17.3 f 0.3 K. The starting configuration for 5 K was obtained from the final configuration at 15 K by cooling in stages (12, 10 and 8 K) of 2000 steps, and then a further 3500 steps were run at 5 K. Two successive runs of 8200 timesteps gave an average temperature of 5.2 k 0.2 K for the data accumulation stage. During the generation of the power spectra a Gaussian filter having a FWHM of 0.8 cm-’ was used to smooth the final results.

3. Results and discussions The optimum angle 0 between the molecular axis and the crystal x-axis (fig. 1) is 43.6” for the 5q model in both the corrugated and flat (a = 4.260 A) cases, and is 43.8” for the flat (a = 4.138 A) optimized case. For the CP model the value is 44.7” for the corrugated, flat (a = 4.260 A) and flat (a = 4.137 A) optimized models. These B values are very similar to those found earlier by Fuselier et al. [12], using two different potentials for N,-N, interactions. The value for the commensurate structure is very similar to that found by Bruch [17] using model Xl of ref. [21] for N,-N, and the Steele model [18] for N,-graphite. His optimized structure has a non-triangular centre-of-mass arrangement which leads to a larger a and smaller b lattice constant in the unit cell, but with almost no difference in the area per unit cell. This distortion changes the angle of tilt with respect to the glide directions by about 9” in the expected direction; for the Sq-Steele combination used at Bruch’s optimized lattice constants the angle differs from this by a degree, and the energy is lower than in the triangular system by almost 2%. The similarity between the results for flat and corrugated surface confirms the relative unimportance of the corrugation in determining the orientational structure [17]. The agreement in the lattice constants for the optimized, triangular 5q and CP monolayers on the flat surface is surprising, but also possibly coincidental. No correction [17] was applied during the optimization for the effects of three-body forces [32], either among the molecules, or between the molecules and the surface. The former is smaller for adlayers than for bulk solids, because the number of triangles is greatly reduced; the likely effect of including these terms is a small increase in the lattice constants. The latter, usually represented by the McLachlan [33] formula, leads to weakening of the attraction between mole-

G. Cardini,

S. F. O’Shea

/ Nitrogen

physisorbed

on graphite

237

cules at the same height above the surface, and hence to an expansion of the lattice. Studies of overlayers of xenon on graphite (221 have shown that the direct contribution of the McLachlan terms to the dynamical properties are almost negligible, and the main effect arises from the changes in force constant resulting from the relaxation of the lattice constant [29]. The lattice constants used here (a = 4.138 A and a = 4.260 A) almost certainly bracket the true optimum value, and the uncertainty associated with the interaction parameters [17] does not appear to justify a more elaborate determination, particular in the light of the neglect of corrugation effects. We have considered the dynamical consequences of a non-triangular centre-of-mass system, and give a brief account of the results below. 3. I. Lattice dynamics As mentioned earlier, the flat surface case has a slightly surprising dynamical symmetry in the harmonic approximation. At the r point the site group is and for every non-zero wave vector in the C,,, and the unit cell group D,,, reduced zone the group of the wave vector is C,. For every wavevector motions in-plane and out-of-plane belong to distinct symmetries, and at the r point librations and translations are further distinguished because of inversion symmetry. In the discussion we use the labels translational and librational to describe vibrations, but it should be clear this is correct only at the r point. Near r the mixing is weak, and the designation is still useful, but near the zone boundaries the mixing is strong, and the description inaccurate. In the case of the corrugated surface the site symmetry consists only of the identity, even at the r point. However, if the planes in the interior of the graphite substrate are approximated by a 10-4 potential [18], the site group is C, at the r point; the usual approximation of the surface potential by the leading terms in the Fourier expansion is equivalent to this treatment as far as symmetry is concerned. The corrugation, by lowering the symmetry so that branch crossings are eliminated, has a striking effect on the appearance of the dispersion curves but very little effect on the frequencies. Were it possible to study the vibrations by infrared and Raman spectroscopy there would also be important consequences in terms of selection rules. The phonon dispersion curves for the CP model are shown in fig. 2. Part (a) gives the curves for the flat surface (a = 4.260 A) with the out-of-plane branches drawn as solid lines and the in-plane branches as dashed lines. Branches of different types can cross, but those of a given type cannot, although they approach so close that it sometimes appears as though they do. The out-of-plane motions show little dispersion, being, for the most part, confined to the region from 50 to 60 cm -I. Such dispersion as there is occurs in the highest and lowest of the four branches, those which are predominantly librational in character and for which the frequency increases as the wave-

238

G. Card&

S.F. O’Shea / Nitrogen physisorbed

on graphite

Fig. 2. The dispersion curves for the CP-Steele combination (see text). Part (a), bottom, the flat, commensurate (a = 4.26 A) case, part (b) with the corrugated, commensurate part (c) with the flat, relaxed (a = 4.137 A) case. In part (a), the out-of-plane branches as solid lines, and the in-plane branches as dashed lines.

deals with case, and are drawn

length decreases. By way of contrast, the highest frequency in-plane branches, which are also librational, decrease relatively sharply in frequency as the wavelength is shortened for vibrations propagating in the r + X direction. This is easy to understand in terms of the molecular motions involved; at r the librating molecules tend to clash and at X they avoid one another more effectively. Along the r + Y direction the dispersion of the higher frequency in-plane modes, which again are strongly librational in character, can similarly be understood in terms of the changing contacts between neighbouring molecules; in this case the relative phases of the two molecules in the unit cell account for the fact that one branch softens as the wavelength shortens while the other stiffens. The two branches originating at I’ between 30 and 40 cm-’ are due mainly to translational motions, and show strong coupling between the two in-plane directions for the two molecules in the unit cell. The near-equivalence of these directions and the presence of avoided crossings suggest that

G. Card&

SF. O’Shea / Nitrogen physisorbed on graphite

239

these motions are likely to be anharmonic. At r the two in-plane acoustic modes have zero frequencies, and represent uniform translations of the rigid overlayer over the flat substrate. In fig. 2b the dispersion curves for the corrugated surface show that these translations are hindered by the corrugation, so that these modes have a limiting frequency of about 6 cm-‘. The low frequency of these’modes shows the weakness of the restraining potential, and the small splitting between the models is a reflection of the near-equivalence of small translations in the two in-plane directions. Apart from this feature, and the absence of branch crossings, due to the lowering of symmetry, the dispersion curves are virtually unchanged by the introduction of the corrugation. The avoided crossings between in-plane and out-of-plane motions inevitably lead to mixing between these motions, but an inspection of the eigenvectors shows that the mixing is confined to very narrow regions of q space, at least along the special directions. This weak coupling should not cause a great increase in the anharmonic character of the out-of-plane modes. It is worth noting that the corrugation has not shifted the frequency of the out-of-plane branches by more than a few percent anywhere, and always down in frequency. The Lennard-Jones model of molecule-surface interaction has little flexibility, and the relative strengths of the q = 0, flat surface, and q # 0, corrugation, terms can be varied only by changing LT.To study the consequences of increasing the corrugation we have varied IS and e to achieve the same overall binding energy but different levels of corrugation. $&eater corrugation can be achieved only by decreasing a, and with CI= 3.06 A and e/k: = 37.5 the layer sits at 3.115 A above the surface, with a classical binding energy of 1126 K; if the q = 0 term alone is used the height becomes 3.024 A and the binding energy 1098 K. Comparing the phonons for flat and corrugated cases shows the expected behaviour: the in-plane modes are virtually unchanged, and the out-of-plane frequencies are slightly higher (= 2 cm-‘) in the flat case. The corrugation terms have little effect on the curvature of the well for the motions in-plane, for which the force constant is dominated by mol~ule-mol~ule interactions, but they do change the height above the surface slightly. This is because the height dependence of the weak corrugation and strong flat surface terms is different, and the result is a modest change in the force constant for vertical motions. Even if the corrugation of the graphite proves to be substantially larger than in the Steele model 1181,as has been suggested to be case for the rare gases 134-363, the dynamics of the overlayer is not likely to be greatly modified unless the 4 = 0 part is also changed significantly. Optimization of the lattice constant for the flat surface results in a = 4.137 A. Out-of-plane translational modes are again slightly softened, but near the I‘ point out-of-plane librational modes are noticably (4-6 cm-‘) softened and the splitting between them increased. This is to be expected since the splitting is due to the coupling between the two molecules in the unit cell. For the in-plane modes, the compression of the lattice leads to a general expansion of

G. Cardint, SF.

240

O’Shea / Nitrogen physisorbed on graphite

the frequency scale. The splitting of the two high-frequency branches is not increased by the compression of the lattice, but the splitting of the middle two branches increases more than average. This overall picture is consistent with the notion that the out-of-plane motions are determined in large part by the molecule-surface potential while the motions in-plane are dominated by the molecule-molecule interactions. The dispersion curves for the 5q model are shown in fig. 3. A comparison of the curves for the flat (a = 4.26 A) model (fig. 3a) with those for the CP model (fig. 2a) shows that the general features are quite similar. The out-of-plane branches, especially that which begins as r as a pure translation, are almost unchanged from those in the CP model; the librational branch is slightly stiffened. The two high-frequency in-plane branches, pure librations at the r point, are systematically softer for 5q, except along the X -+ S direction where they are unchanged. The splittings along r - X and Y - F are only marginally increased. The two central in-plane branches are very similar to those of the CP model, with slightly increased splittings along Y + F. The two lowest branches show similar initial slopes along r + X and r --$Y, indicating the 80 60

40 20 80 60 40 20 80 60 40

_.d=L‘-\

/-.,

-=-=T,

20

1. p.

-

/I

0

I‘

,‘,f2, ’

____---.

x

Fig. 3. The dispersion of fig. 2.

;=-

--__--

s

K, Y

curves for the Sq-Steele

I combination

(see text), following

the conventions

G. Cardini, S.F. O’Shea / Nitrogen physisorbed on graphite

241

speed of sound is similar in both models, but with 5q these branches do not soften as much approaching the X point and remain consistently higher until the Y point. To investigate the effects of non-triangular centres-of-mass, we calculated the frequencies for a structure based on Bruch’s [17] cell constants and molecular orientation reoptimized for this combination of interaction potentials. The out-of-plane modes are changed only slightly, and so as to lower the dispersion. The in-plane modes suffer greater changes, mainly a softening of the highest frequency modes and a stiffening of the low frequency ones. The overall character of the results remains unchanged: a narrow dense band of out-of-plane modes is embedded in a broad band of in-plane modes which, at shorter wavelengths, are characterized by strong mixing of the coordinates. The effects of corrugation (fig. 3b) and of optimization of the lattice constant (a = 4.137 A) are essentially similar to those observed with the CP model.

0.20

0.15

1

> .E

E!

4 .lO

n 0.05

0.00

0.15 > .z

C

0.10

B 0.05

0.00 20

40

60

Wavenumbers Fig. 4. The one-phonon density of states for the Sq-Steele combination (see text) with a lattice spacing of a = 4.26 A. Part (a), bottom, is for the flat surface, and the dashed line details the density due solely to the out-of-plane modes. Part (b), top, is for the corrugated surface.

242

G. Cardini, S. F. O’Sheu / Nitrogen phymorhed

on gmphrte

3.2. Density of stutes The one-phonon density of states for the 5q model on a flat surface (fig. 4a) displays a considerable amount of detail. The sharp peaks below 40 cm ’ are due to the in-plane vibrations; these show some dispersion near the F point, where they are mainly translational in character, but tend to be flatter for shorter wavelengths and these plateaus give peaks in the density of states. The peaks between 40 and 50 cm-‘, the shoulder on the low frequency side of the peak near 60 cm-‘, and the plateau near the high-frequency cutoff are also due to motions in-plane. The in-plane librations have lower frequencies and overlap the frequencies of the out-of-plane branches, leading to a substantial peak in the density of states between 50 and 70 cm-’ (fig. 4) and the contributions due to the different motions are shown separately in fig. 4a. As shown by the dashed curve, the out-of-plane motions, both translations and librations, are localized between 52 and 62 cm- ‘. The data for the corrugated surface (figure 4b) is very similar, with the major difference being the presence of a low frequency cutoff at 6 cm ‘. This similarity emphasizes the relative unimportance of the corrugation in determining the dynamics of the overlayer for a given lattice constant. Integral thermodynamic properties derived from the density of states are likely to be similarly insensitive to the corrugation at all but the very lowest temperatures.

3. Molecular dynamics Gross properties, such as the order parameters, are in quantitative agreement with those Talbot et al. [16], and the energy per molecule is also consistent with theirs when allowance is made for the differences between the models. When analysing the mean squared displacements normal to the surface it became apparent that the results are anomalously large at 5 K. While investigating the origin of this results we discovered what we feel is a novel observation in simulations of “realistic” models: the systems is exhibiting recurrence behaviour of the type first discovered by Fermi, Pasta and Ulam (FPU) [37]. Equipartition of energy is not achieved [38], but instead the degree of excitation of different modes varies in time. This can be seen by examining characteristic coordinates of the motion as a function of time. In fig. 5a the mean height of the configuration above the surface at each timestep is plotted is evident, with the high for the entire 16400 step run. A definite “beating” frequency (- 60 cm-‘) oscillation due to the usual out-of-plane vibration being modulated in amplitude over thousands of cycles by anharmonic couplings with other modes. In the one-phonon approximation this coordinate represents the q = 0 out-of-plane translation model, and equivalent coordinates can be constructed for other modes, including those having non-zero q values.

G. Cardini, S. F. O’Shea / Nitrogen physisorbed

f

on graphite

243

3.335 3.325

e 1.05 d

-0.5 -1.5

C

0.0

! ‘ll”

‘Ill#J

steps Fig. 5. Collective coordinates as a function of time, given as multiples of the timestep (5 X IO-” s). The coordinates are (a) the mean molecular height above the surface, (b) the difference in mean heights of the two sublattices above the surface, (c) the mean out-of-plane librationat coordinate (see text), (d) the difference between the mean out-of-plane tibrationai coordinates for the two sublattices, and (e) the mean in-plane librational coordinate. These all refer to the run at 5 K, and (f) gives the mean height above the surface during the run at 17 K.

The antisymmetric out-of-plane translational model, where the two molecules in the unit cell move in opposite directions, yields a similar diagram (figure 5b) but the amplitude modulations are not in phase with those of the symmetric out-of-plane translation. Using the sine of the angle of tilt out-of-plane (or, equivalently, the cosine of the polar angle 8) as the librational coordinate for each molecule we also constructed the corresponding out-of-plane collective hbrational coordinates, and observe similar behaviour, although the structure is not as marked in these cases (figs. 5c and 5d). Further evidence of the FPU phenomenon can be obtained from the power spectrum of these collective coordinates by partitioning the run into smaller pieces (4096 steps here). The frequencies of the peaks in the spectra are essentially unchanged but the

amplitudes vary substantially from segment to segment. In view of the care that was taken to prepare a properly thermalized starting configuration this result is quite unexpected. While it is still possible that the procedure used was inadequate, this is unlikely since energy exchange occurs between modes, as evidenced by the amplitude fluctuations. One possible interpretation is that the phenomenon results from a type behaviour well known in the engineering literature, the excitation of a narrow band response in a resonant system by broad band noise [39]. The low dispersion in the out-of-plane modes discussed earlier means they can act as a narrow pass filter [40]. The recurrence behaviour is not exhibited by the in-plane librational modes. The molecular coordinate is the tangent of the angle between the projection of the bond in the surface plane and the glide axis. We have constructed collective coordinates using the average value obtained by summing over both sublattices (fig. 5e) or taking the difference between the sublattices (not shown). Each gives rise to noisy behaviour in time, shows no systematic amplitude modulation, and gives

0

20

40

cm-’

60

0

20

40

cm-’

60

G. Cardini, S. F. O’Shea / Nitrogen physisorbed on graphite

245

3

2-

l-kL!zz 0 20

1 40

cm-’

60

0

20

40

60

cm-’

Fig. 6. The dynamical structure factor. S(9, w), for the 5 K, left, and 17 K, right, runs. We report the data for the 9 9 values available along the r + X direction and the 7 kl values along r ---)Y. In the former case, the kl value is given by n(2~)/(9a), where n is the index plotted in the margin, and in the latter case by n(Zs)/(76). For each direction, the intensity scale is set by the larger peak among the 5 I: spectra, and the 17 K spectra are reported on the same scale for the r -+ Y direction, but are scaled up by a factor of 2 along r + X.

a Fourier transform exhibiting the strongest peaks in the positions expected from the phonon spectra, but with considerable lower intensity structure at lower frequency. This is consistent with greater anharmonic character and larger intermode coupling for motions in-plane, making thermalization much more likely here. The in-plane and out-of-plane motions are only very weakly coupled in the harmonic approximation and this behaviour appears to persist in the anharmonic system at low temperatures. At higher temperatures, as seen in the run at 17 K, the amplitude modulation is no longer as clearly evident. Instead, for example, the out-of-plane translation mode is modulated (fig. 5e) by a low-frequency motion which is probably the 6 cm-’ translation of the registered overlayer in potential well due to the corrugation. It could also be due in part to a difference band between this and other out-of-plane motions. Besides this 6 cm- ’ modulation there is some suggestion of amplitude modulation, but a more careful study will be needed to verify it. A detailed analysis of the FPU phenomenon is a major undertaking in itself, and would lead us away

G. Cardini, S. F. O’Shea / Nitrogen physisorbed on graphite

246

from the goal of the present work. In any case, such a study requires very accurate integration of the equations of motion because the trajectories can be degraded and the recurrence lost through the accumulation of numrical errors, whether due to rounding errors or the use of low-order integration schemes. We are working to implement a better integration technique, and to develop methods of analysis which will lead to a clearer understanding of the energetics of the system and its relationship to these dynamical features. In spite of the existence of the FPU recurrences, we can still achieve part of our goal of studying the growth of anharmonic effects in the vibrational structure with increasing temperature. The phonon Iinewidths and intensities are not reliable, varying as they do over the run at 5 K, but the line positions should be, within the resolution afforded by the total simulation length. A cautious estimate makes this about a wavenumber at 5 K and two wavenumbers at 17 K. The out-of-plane motions yield clean, well resolved phonons, while the in-plane motions are very strongly mixed and give noisy spectra. For the latter, S(q, o) (fig. 6) shows lots of structure and strong mixing among the

3

2

1 0

20

40

cm-’

sb

0

20

40

cm-’

60

G. Cardini,

0

20

40

cm-’

60

S. F. O’Shea / Nitrogen physisorbed on graphite

0

20

40

247

60

cm-’

Fig. 7. The one-phonon approximation to S(q, o) for the in-plane librations (see text). The conventions of fig. 6 are followed, including the resealing of intensity of the 17 K spectra for the r * X direction.

modes at larger q values, and the one-phonon spectra (fig. 7) have much of the structure seen in S(q, w) itself, although with different intensity distributions. Except near the zone boundary at X, the in-plane librations contribute very little intensity to S(q, w). In the one-phonon approximation the librational frequencies agree well with those obtained from lattice dynamics and also show the weak dispersion observed along r + X (fig. 7) but along r + Y (above) the peaks are barely distinguishable from the noise even in the one-phonon approximation, and cannot be analysed. The longitudinal and transverse phonons (not shown) account for the other major peaks in S(q, w), and both spectra show some structure below 10 cm-‘; it is most pronounced near the r point, where both have the expected peaks in the region of 6-8 cm-‘. Although both have strong peaks at frequencies which are consistent with the lattice dynamics, many of these are common to both spectra, as one would expect on the basis of the low symmetry and anharmonicity. Comparing S(q, w) for 17 K with that for 5 K, we see that the higher frequency features generally appear weaker, and that a strong feature at 6 cm-’ now occurs for

G. Cardini, SE

248

O’Shea / Nmogen physisorbed

on graphite

all q values. On an absolute scale the intensity at higher frequencies is somewhat reduced, and the peaks are broadened with a noticeable growth in the multiphonon background. The librational peaks (fig. 7) show a greater loss of intensity, but do not appear to soften; since the temperature is not more than two thirds of the transition temperature this is probably to be expected. The out-of-plane spectra (figs. 8 and 9) offer sharp contrasts to S(q, a), since they both consist of intense bands between 45 and 60 cm-‘, with nothing at other frequencies. Near the zone boundary the translational band (fig. 8) contains two narrowly separated peaks which can be resolved at 5 K although not at 17 K, although longer sampling should make it possible there also because well defined shoulders are evident. As mentioned earlier, the intensity of the two peaks changes over the course of the run at 5 K due to the energy exchange among the vibrational modes. The translational frequencies are essentially undispersed and agree in detail with the lattice dynamical results. The shoulders below 50 and at 60 cm-’ are due to coupling with the out-of-plane librational models. At 17 K the peak widths are no greater, and

9

8 I---

-

6---

0

20

40

cm-’

60

0

20

40

cm-’

60

G. Cardini,

S. F. O’Shea / Nitrogen ph.vsisorbed on graphite

249

6

0

20

40 cm-’

Fig. 8. The power spectra resealing of intensities.

60

0

40

20

60

cm-’ of the out-of-plane

translational

motions.

As in fig. 6, but without

may even be slightly smaller, although the data do not permit us to determine this reliably. A low frequency component is evident in the 17 K spectra, appearing as a doublet; however, the lower peak is at the threshold of detection for the length of run used (8192 steps). The harmonic librational bands for out-of-plane motion are more widely spaced than the translational bands and this is reflected in the MD data (fig. 8). Near the r point two peaks are evident, but three are sometimes resolved for intermediate q values due to mixing with the translation bands which lie between them. There is less evidence of coupling in the spectra from the 17 K run, so that the lines are noticeably narrower for most q values.

4. Conclusions The combination of lattice and molecular dynamics has given some useful insight into the dynamics of the registered monolayer of nitrogen on graphite.

G. Cardmi, S. F. O’Shea / Nitrogen physisorbed on graphtte

250

At low temperature, 5 K, the agreement between lattice and molecular dynamics is excellent, and continues to be good at 17 K. No substantial reduction in frequency with increasing temperature is observed, and this is due to the role of surface corrugation in the system. The corrugation serves to establish and maintain a lattice spacing which is expanded with respect to the optimum value determined by the intermolecular potential, and therefore increasing the temperature does not lead to an increase in lattice constant. The usual reduction of frequency with increasing temperature, which is due mainly to the reduction of force constants because of increasing intermolecular spacing, is not observed. The anharmonicity can be seen in the finite width of the phonon peaks, and the growth of low frequency intensity and of the multiphonon background with temperature. It is also responsible for the weak coupling which gives rise to the FPU phenomenon in the out-of-plane motions at 5 K. The coupling of in-plane and out-of-plane motions is weak, and seems to be detectable only in the presence of a 6 cm -’ peak in almost all spectra at 17 K. The motions in-plane are determined almost entirely by the intermolecular

8

0

20

40 cm-’

60

0

20

40 cm-’

60

G. Cardini, SF.

0

20

40

60

O’Shea / Nitrogen physisorbed on graphite

0

20

cm-’ Fig. 9. The power spectra of the out-of-plane of intensities.

40

251

60

cm-’ librational

motions.

As in fig. 6, but without

resealing

potential, except that the corrugation plays a vital role in determining the lattice constant and the band gap below 6 cm-‘. On the basis of the lattice dynamic calculations for the overlayer, little distinction can be made between the CP and 5q models. There is little prospect of experimental determination of in-plane dispersion curves. The scattering of atoms off the surface layers involves reversing the component of momentum normal to the surface, a process that is much more strongly coupled to the vibrational modes polarized normal to the surface than to those polarized in-plane. In any event, since major questions still exist concerning the two-body part of the molecule-molecule potential [41] we do not anticipate being able to contribute usefully to the discussion of surface mediation and other three body forces even if frequencies for in-plane motions were available. The out-of-plane motions are dominated by the gas-surface potential, and it is here that atomic beam spectroscopy offers immediate progress. One obvious difficulty which can be anticipated on the basis of our calculations is that the out-of-plane motions are essentially undispersed, and the translations and librations have similar frequencies and

252

G. Cardini, S. F. O’Shea / Nitrogen physisorbed on graphite

intensities in the one-phonon approximation. The resolution of the beam experiments may not be sufficient to separate the peaks, as indeed our resolution is not for some q values. However, the lack of dispersion is due in part, at least, to the use of the rigid substrate appro~mation, if one is to judge on the basis of the results obtained for rare gases on multilayers of graphite v41.

The FPU recurrence observed here may be an artefact of the assumptions (classical mechanics, rigid substrate, etc.) used here, but it may also be observable in the experimental system. We are continuing our investigation of it, and hope to consider which observable properties might best distinguish between FPU and non-FPU behaviour. If this behaviour is due to the presence of narrowly dispersed and near degenerate vibrations, then many adsorbed systems may exhibit similar characteristics,

Our interest in this system has been stimulated by continuing contacts with Sam Fain, Mike Klein, Jim Morrison and Bob Thomas, and Steven Sibener’s beam experiments for xenon on silver set us going on this particular project; Giacinto Stoles convinced us the beam experiments for nitrogen on graphite are feasible. M. Keramat Ali provided very helpful comments on FPU behaviour. This work was made possible by a grant from NSERC Canada, and by computer resources generously allocated by the University of Lethbridge. The staff of the computer centre went out of their way to be helpful, as always.

References [I] [2] [3] [4] [5] [6] f7] [8] [9] [lo] [ll] [12] [13] [14] [15]

S. Ross and J.P. Ohvier, On Physical Adsorption (Interscience, New York, 1964). W.A. Steele, The Interaction of Gases with Solid Surfaces (Pergamon, Oxford, 1974). O.E. Viiches, Ann. Rev. Phys. Chem. 31 (1980) 463. C. Boekel, J.P. Coulomb and N. DuPont-Pavlovsky, Surface Sci. 116 (1981) 369. B.F. Mason and B.R. Williams, Phys. Rev. Letters 46 (1981) 1138. K. Gibson and S.J. Sibener, to be published. J. Eckert, W.D. Ellenston, J.B. Hastings and L. Passell, Phys. Rev. Letters 43 (1979) 1329. T.T. Chung and J.G. Dash, Surface Sci 66 (1977) 559. A.D. Migone, H.K. Kim, M.H.W. Chan. J. Talbot, D.J. Tildesley and W.A. Steele, Phys. Rev. Letters 51 (1983) 192. M.H.W. Chan, A.D. Migone, K.D. Miner and Z.R. Li, Phys. Rev. B30 (1984) 2681. R.D. Diehl and S.C. Fain, Jr. Surface Sci. 125 (1983) 116. C.R. Fuselier, N.S. Gillis and J.C. Raich, Solid State Commun. 24 (1978) 747. A.B. Harris and A.J. Berlinsky, Can. J. Phys. 57 (1979) 1852. S.F. O’Shea and M.L. Klein, Chem. Phys. Letters 66 (1979) 381. O.G. Mouritsen and A.J. Berlinsky, Phys. Rev. Letters 48 (1982) 181.

G. Cardini, S. F. O’Shea / Nitrogen physisorbed on graphite [16] [17] (181 [19] [ZO] [21] [22] [23] [24] [25] [26] 1271 [28] [29] [30] [31] [32] [33] [34] [35] (361 [37] [38] (391 [40] (411

253

J. Talbot, D.J. Tildesley and W.A. Steele, Mol. Phys. 51 (1984) 1331. L.W. Bruch, J. Chem. Phys. 79 (1983) 3148. W.A. Steele, J. Phys. (Paris) 38 (1977) C4-61. P.S.Y. Chung and J.G. Powles, Mol. Phys. 32 (1976) 1383. C.S. Murthy, S.F. O’Shea and I.R. McDonald, Mol. Phys. 50 (1983) 531. C.S. Murthy, K. Singer, M.L. Klein and J.R. McDonald, Mol. Phys. 41 (1980) 1387. G. Cardini and S.F. O’Shea, Phys. Rev. B15, in press. W.A. Steele, Surface Sci. 36 (1973) 317. F.W. de Wette, B. Firey, E. de Rouffignac and G.P. Alldredge, Phys. Rev. B28 (1983) 4744. S. Califano, V. Schettino and N. Neto, Lattice Dynamics of Molecular Crystals (Springer, Berlin, 1981). M.L. Klein, in: Computer Modelling of Matter, ACS Symposium Series No. 86, Ed. P. Lykos (American Chemical Society, New York, 1978) p. 94. L. Van Hove, Phys. Rev. 95 (1954) 249. J.J. Weis and M.L. Klein, J. Chem. Phys. 63 (1975) 2869. G. Jacucci and M.L. Klein, Solid State Commun. 32 (1979) 437. G. Cardini, SF. O’Shea, M. Marchese and M.L. Klein, Phys. Rev. Letters, submitted. G. Ciccotti, M. Ferrario and J.P. Ryckaert, Mol. Phys. 47 (1982) 1253. H. Margenau and N.R. Kestner, Theory of Intermolecular Forces (Pergamon, Oxford, 1969). A.D. McLachlan, Mol. Phys. 7 (1964) 381. G. Vidali and M.W. Cole, Phys. Rev. B29 (1984) 6736. R.J. Gooding, B. Joos and B. Bergerson, Phys. Rev. B27 (1983) 7669. B. Joos, B. Bergerson and M.L. Klein, Phys. Rev. B28 (1983) 7219. E. Fermi, J.R. Pasta and S.M. Ulam, Los Alamos, Report LA-1940, May 1955; reprinted in E. Fermi, Collected Works, Vol. II (University of Chicago Press, Chicago, 1965). N.J. Zabusky, in: Mathematical Models in Physical Sciences, Eds. S. Drobot and P.A. Viebrock (Prentice-Hall, 1962) p. 99. D.E. Newland, Random Vibrations and Spectra1 Analysis (Longman, London, 1975). L. Brillouin, Propagation of Waves in Periodic Structures (Dover, New York, 1953). M.S.H. Ling and M. Rigby, Mol. Phys. 51 (1984) 855.