Molecular dynamics simulations of rare gas matrix deposition

C@-nical

Physics Cbmbl Physics 186 (1994)

t85-204

Molecular dynamics simulations of rare gas matix deposition Dspormrcn ofPItyhal

ZtuchamaFraenkel, Yehuda Haas * ChemWy and the FcvLar Centerfor Light Induced Procewes, The Hebmw University of Jeru.wlem, Jerusalem, Ismel

Received 28 F&muy 1994; in final form 29 April 1994

Mdecular dynamics simulations of the tmpping of molecules in a rare gas .matfix based on a fee lattice is presented. In distktion with most previous simulations, a site’s sWcture is not obtained by a guess and optimize method, but is found by allowing the lattice to grow by adding atoms from the vapor onto the pqntial field of a “seed” crystal. Simple pairwise potentials are used to simulate the growth of ram gas matrices containing some small guests, including the diatomic molecules N2, NBt and Br,. changea in the structun and energy of the solid arc cslculated, ss well as the vibrational frequencies of the trapped mokcoles, yielding the m&ix shift_ Novel site geometries, in which the host lattice atoms are displaced from their

rquil&iunl crystal positions, sre ohserved for guest moleculesthat can stronglydistort the lattice becauseof asymmetryor size mismatches.In these cases, rahd cooling of the incomingatom is essentialto the formationof several distinct trappingsites.

1.IntrQductlon Matrix isolation was originally developed mainly to study the properties of reactive intermediates [ 11,but in recent years it is used increasingly also to study the spectroscopy of stable molecules and their reactions. Rare gas hosts are of particular interest, since their chemical inertness makes them the simplest “solvents” and thus an ideal model for comparmg gas phase energetics and dynamics to condensed phase behavior. In many cases, the interaction with the neighboring atoms or molecules is weak enough, so that energy levels and thus spectral features can be compared with calculations performed on the isolated molecule. However, different micro-environments do lead to differences in some properties of the trapped species, such as spectral shifts, radiative or non-radiive decay rates and possibly also chemical reactivity, grouped together under-the term “site effects” [Z-4]. * Coming

author.

0301-0104/94/$07.00 0 1994 F&evier Science B.V. All rights reuerved SSD10301-0104(94)00160-C

It has become increasingly evident in the last few yerusthatsuchweakinter&4mscanbequite@ortant when the dynamics of the trapped species are investigated. Studies of photodissociation [ 51, energy transfer [ 61, diffusion [7] (for a general discussion of diffusion in lattices see Ref. [ 8 I), bimolecular reactions [ 91, charge transfer [ lo] and other dynamic processes appear to depend to a large extent on the detailed environment in which the process is taking place, namely, on the properties of the sites. These, in turn, are largely determined at the preparation stage of the matrix which is known to be of crucial importance in determining its properties. ‘Ihus, the dynamics of matrix trapped species can be controlled, or at least influenced, by experimental parameters such as the rate of deposition, the guest/host ratio, the temperature during deposition and annealing procedures. The increasing awareness of the role of site structure and energetics has also led to a growing number of attempts to simulate them theoretically starting, as much as possible, from basic princi-

186

R. Fraenkel, Y. Haas / Chemical Physics 186 (I!%%) 185-204

ples. In this respect, rare gas matrices are characterized by the simplest, and probably best known, interaction potentials, making them the natural choice for bridging the gap between gas and condensed phase behavior. Much of our current understanding of the dynamics of matrix isolated species arises from theoretical simulations using Monte Carlo (MC) [ 111, molecular dynamics (MD) [ 121 or other methods [ 131. In these studies, the standard procedure is to assume that the guest molecule occupies a well defined site in a crystal lattice. The structure of these sites is usually obtained by a numerical search method, in which empirical or calculated potentials are used to determine the potential energy of the molecule located in the site. The computation is started by a physically reasonable guess (for instance, a substitutional site), and the final configuration is arrived at by an energy minimization procedure, involving small changes in the coordinates of the lattice atoms near to the guest. The difficulty in predicting the existence, let alone the nature of several distinct sites in a given matrix has been often alluded to [2,15-171. The presence of several distinct sites is quite common in matrix isolation spectroscopy; both atoms [ 141 and molecules [ l-4,15-1 81 are believed to be trapped in well defined locations. The present work was motivated partly by recent experimental studies [ 19,201 that revealed the presence of several distinct matrix sites in a highly diluted matrix of anthracene in argon. It was suggested that these sites are formed during the deposition, and reflect energetically favorable loci of similar stabilization capacity for accommodating anthracene molecules in solid argon environment. In view of the large size of anthracene compared to an argon atom, one might expect severe disruption of the host matrix structure, and a large number of different sites [ 11. However, only three major sites were actually observed by absorption spectroscopy, along with some others with smaller intensity [ 191. In a more recent study [20], Wolf and Hohlneicher reported 11 different sites for anthracene in an argon matrix. In an attempt to formulate a reasonable model that will account for this finding, which is typical of many matrix isolated systems, we discuss a method that simulates the time history of trapping of molecules in a rare gas matrix with as little a priori assumptions as possible.

The model should be able to account for the following experimental observations: - Any molecule can be trapped in a rare gas matrix, but doped crystals can be grown only if the certain compatibility conditions are met. In particular, large guests cannot be incorporated into the crystal lattice of a small rare gas atom. Thus, a matrix can be used to prepare non-equilibrium solid solutions [ 211, while a crystal cannot. - The properties of the matrix are determined to a large extent by the deposition conditions. Typically, high temperatures (of the order of T,/3, where T,,, is the melting point of the RG solid) lead to formation of optically clear matrices; a strong tendency to form aggregates necessitates the use of highly dilute samples (host/guest ratios of 10000 or more) to ensure trapping of isolated single molecules. At lower temperatures, larger guest/host ratios can be used, while still maintaining a low aggregate concentration. Lowering the deposition temperature results in a less ordered matrix - more defects, such as vacancies, are observed and the optical quality is inferior. - For any molecular species (monomer, aggregates) several distinct trapping sites are often observed in a matrix. Annealing can sometimes reduce their number, but usually does not lead to a single site. The number of sites tends to increase at lower deposition temperatures. In contrast to other methods (such as Monte Carlo) the molecular dynamics (MD) approach enables one to follow the real time evolution of a matrix by deposition of atoms from the gas phase onto a cold surface, including motion of atoms in and near the trapping site. For instance, the vibrational frequency of a trapped diatomic molecule can be simply obtained from the plot of the internuclear distance as a function of time. The method thus provides a straightforward way of calculating the changes in the characteristic frequencies of the trapped molecules with respect to the gas phase values, and comparison with the experimentally observed matrix shift provides a means to test the applicability of the potentials used in the simulation. We are not aware of previous such MD work, not even aimed at the more modest goal of simulating the formation of a pure rare gas solid crystal by gas deposition onto a cold surface. Studies of related processes, such as the simulation of the growth of Lennard-Jones fee crystals from the melt [ 221, or of molecular beam

R Fmunkl,Y.Haas/ChanicalPkysics186(1994)185-204

growth of a Lenmud-Jones crystal @ase [23], have been repotad. ation is b the calcubrdiomof classiealtr@ee&riesofgas atomsimpingingona cold seed Qrystalof a well defined stn%hue. we apply itfirsatotheclpasaaf’apura~,witbouta~stmolecule, in order to test its capability to reprodwell knowntrends in mat&deposition technology. Simple pair wise Lennar&Joaes petentials allow the use of a desktopworkstationtosimulatethegrowthafamatrix consisting of several hundred atoms. Some interesting demils of the deposition. process are revealed, such as thetr@eetoryof nearthe!#uface,andthe Pnal p&ion. The effect time need&d‘for of the mat& rempenature on the s-m ofthe newly dqxhted~yers is consistent with experimental results, and provides insight into the rules of thumb widely used. Next, we consider the incorporation of guest molecules in the matrix and show how several sites con be formed. Within the limits of the p&-wise potential approximation. the effect of the guest’s size, its interaction go&k&d as well as the matrix' temperature and the cooling rates are discussed. A preliminary ace&t of this we& was given in Reti. [24,2$], A further motivation for this work arises from many experiments on molecules embedded in clusters generated in supersonic jets [26,27]. Such clusters are often considered as intern&ate species, bridging the gapbetween isolated molecules and salvatedones. One of the pmperties~that has been meaatues in the jet, and correlatedwith bulk (i.e., matrix) behavior, is the spectral shift observed for v&rational or alee@mic transitions. The bulk vahre used in these comparisons may sometimes depend on the preparation conditions, and one may wish to know~more ‘about the intimate environment of the guest in the matrix. In thii paper we restrict attuntion to vibrational shifts of diatomic molecules only.

The stable form of rare gas crystals at low temperatnresisanf~stn&ue (28,29], whibmregaschrsters, even big ones, are icosahed&or doak&e&d [30]. MD or MC elustmsimulations, using pair-wise Lennard-Jones @otentials, have been able to reproduce the

187

experimentals&u~[31].The~parameters neededfor these simulations are derived f&m experic m~~on~isotvapoursf2~,29].It.hrs~n shownthattheuseofthesepotentialsalonedoesnot leadtothegrowthofaf&crysUlfromthegasphase [32] - the inclusion of three-body in@aetions is required. Nevertheless, pair-wise potentials do in fact faithfullyreprodua&tbemajorpestdshegotemialquite wellasfaresshortrangeintemctionsareoo~.In an attempt to bane& from the simplicity of pair-wise potential modeling, and yet retain tha bc structure as far as possible, the simulation starts with Ia fee crystal teroplataasa“~‘.~nniaaiir;ig~~Ivyallorwing astrtlrmdr~gasatonsis,~~,~ma300Entsr ervoir, to impinge on thjs seed which, is held at a low temperature. The final crystal structure is therefore dewrmined only by the ~roteatialsand#he snuting, conditions (a similar method was used in the MI&Esimulations mentioned above [ 23 I ) . This ‘ganxudm avoids the intemating and impormnt issue of the formation m~smofthefi~tstsecE,aindEoEUsagonthd)quedan of crystal growth and the cond&us required to form a stable structure. It is noted that in principle, the crystal may grow in, werent shapes and forms and need not repreducetheseedstruotum.Tbisapproa&isakinto some cluster growth calculations [ 33 3, but has not yet been reported for matrices, as far as we know. Most of our actual computations were carried out fer the specific case of argon as the host, and we shall use this example for concreteness. Theseedwasfarm(ldbyatrsmginganumberofargon atoms in an fee unit cell (a = 5.32 A at 25 R), and duplicating it few times, depending on the problem at hand, as detailed belbw. The initial unit bell size for any desii temperauuu was set according to litemture values ‘[1] and each atom was given a random initial velocity chosen fi+oma Maxwell~l%&zmtum distribution at that temperattue. The positions and velocities of the atoms were computed by numerically solving Newton’s equations of motioa using a velocity Verlet integrator [34,35]. After the new coordinates wore obtained, they were used to recalcnlate the potentials and hence the new forces, accelerations and velocities, and the procedure was repeated as long as needed Comerva&n of the total energy and, time reversal reprodu&&tywereusedasoriteriatoeheckthevalidity of the procedm. Energy was found to bateonserved to better than 10” over a few thousand ‘time steps if

188

R, Fraenkel, Y, Haus / Chemical Physics 186 (1994) 185-204

they were kept sufficiently small ( 1 or 2 fs). Under these conditions, the crystal structure was found to be stable and self sustaining at all temperatures up to the melting point. The temperature of the system was defined as T=2&/3Nk, where Er =O.SCmiU? is the total kinetic energy, N is the total number of atoms and k is Boltzmann’s constant. Periodic boundary conditions [34] were imposed to simulate an infinite layer in two dimensions (X and y, along the 100 and 010 planes) so that the crystal was allowed to grow (or evaporate) only in the f z direction (001 plane). The 001 plane was chosen since it has an {ab} repeat in an fee crystal, and thus its reproduction leads to an fee structure. Experimental evidence [ 1] indicates that the first layer of argon on any solid support is usually the dense 111 plane. Crowing an fee crystal on this plane requires an { abc} repeat, while an (ab) repeat leads to a hcp structure. A pairwise potential fails to stabilize the fee structure over the hcp one - in fact, the latter is more stable by about 0.01% (ch. 9 in Ref. [ 281, p. 562). Since we are interested at this point in fee sites only, the study of 111 plane growth is deferred to a later stage. The potential function for the interaction between non-bonded atoms is given by the additive pairwise Lennard-Jones expression V,=C~E[(~/R~)‘~-(U/R~)~],

(1)

u

where R, is the distance between atoms i and i. The potential was truncated at R = 2.7~ (9.18 8, for Ar-Ar interaction) and the seed crystal lateral dimensions were at least twice this value. The intramolecularpotential between the atoms comprising a diatomic molecule AB, was represented by a Morse potential (Eq. (2)) where D is the dissociation energy, R the interatomic distance, and R, the equilibrium distance in the gas phase. V*_,=D(l-exp[

-P(R-R,)]}‘.

(2)

The total potential energy is thus given by VToy. = V&B + CVpur-A + CV&_B + cv,+Ar

,

(3)

where VA_B is the Morse potential for AB, and the other terms are the pairwise Lennard-Jones potentials. Forces were computed from the potential gradient, and Newton’s equations of motion were used to calculate the

resulting accelerations (and velocities) throughout this work, which is thus restricted to the classical mechanics approximation. A simulation was started by constructing the seed crystal as described above. After the initial assignment of velocities, the atoms were allowed to move under the potential field for about 1000 time steps, a period during which the system approached equilibrium, with each atom oscillating around its equilibrium position. Deposition was begun after the crystal seed was equilibrated at the chosen temperature. Initially, the kinetic energy of the atom to be deposited was selected at random from a Maxwell-Boltzmann distribution maintained at 300 K. Occasionally, the initial direction and velocity of the in-coming atom were such, that its trajectory did not lead it to the vicinity of the matrix. In these cases the calculation was terminated after a pre-determined time interval, and a new atom was introduced. Experience showed that if the ad-atom was allowed to start with the mean kinetic energy characteristic of 300 K in the z direction only, (that is, vertical to the surface) computer time was saved and the result appeared to be indistinguishable from that obtained by the random energy and direction sampling method. The gas phase atom trajectory was begun at a zdistance of 2.6~ (8.84 A for argon) from the seed surface while the x and y coordinates were selected at random. The equations of motion of all atoms (seed crystal atoms plus the new one) were then solved using the procedure described above. The temperature of the system was calculated after each step, and as expected, a gradual heating of the solid phase was observed. The total kinetic energy of the matrix was kept constant by periodically cooling it during the deposition, normally every 1000 time steps. Cooling (or warming) was simulated by resealing the velocity of each atom by (Til TA If29 where Ti is the initial lattice temperature and Tf the final one. This procedure mimicked the cooling of the matrix by the cryostat, as expected when rapid energy exchange with the cold sample holder takes place. The cooling of the ad-atom was simulated for two extreme limits: slow cooling, in which the atom was allowed to equilibrate with the matrix solely under the action of the potential field, a process that lasted typically several thousand time steps. In this “slow cooling” case, a new atom’s trajectory was begun at a fixed time (10 ps) after reaching the surface. Rapid or fast

cooling was simulated by equilibrating the atom’s temperature with that of the matrix, as soon as it reached a predetermined vertical distance from the surface (z,,& . When this scheme ukgaused, the atom’s velocity was adjusted to conform to the matrix’ temperature when the vertical distance was smaller or equal to O&J (a is the lattice constant, ha = 3.192 A for argon at 25 K). Since throughout the deposition, the uppermost surface of the lattice could not be descr&d as a horizontal plane, cooling was performed v&en the vertical distance of the ad-atom from at least two surface atoms was smaller than zmi,,.The ,choice of & as O&J was made by examining the effect of the value of k,, on the efficiency of cooling the ad-atom. At d&distance, cooling was found to be most efficient. Temperature equilibration was performed by assigning to the adatoms a velocity selected at random from the distribution characterizing the matrix’ ~$mpi$rature. At this point another 300 K atom was released from the reservoir and the procedure repeated. The structure of the trapRing sites, as calculated from the final exact location of each atom in the resulting solid phase was determined by averaging the coordinates of the surrounding atoms over a period of 10 ps, at 1 K. Cooling the system lto 1 K was performed in order to eliminate unce&n@ies dne to ‘&cc random motion of the atoms at more elevated temperatures. Annealing, when appli+, was simulated by gradually warming the matrix to the desired temperature, and then cycling the system for ld time steps at that temperature. The annealing cycle was completed by gradual retooling to the Ilnal temperature. When a guest molecule vvasdded, the same procedure was used, with the appropriate LJ’potential parameters. The parameters used in this work are listed in Table 1. AI1 calculations were performed on a Silicon Graphics Indigo R4OtXlwork stat&.

3. ReamIts 3.1. Pure argon matrices The trajectory of a singIe atom from its initial position to its final equil&ium position can ‘be vlsualized in several ways. In the present ease the interesting parameters am the vertical &placement (here, the z coordinate) from the matrix’ surface and the velocity. Fig.

Table1 PmwWers used in the simuhtions Molecularpair

AFAr Ar-Xe Ar-sF~ Ar-N Ar-Br (NW Ar-Br (Br2) Molecule

N2 Br2 NBr

Leaaard-JolIesp

a(A)

s/k(K)

3.40 3.65 4.5 3.09 3.46 3.51

120 178 153 56 164 164

Morsepammekrs D (cm-‘)

B (A-‘)

R. (A)

77752 24557 25335

2.7316 1.586 1.586

1.097 2.2809 1.769

1 shows the velocity of the ad-atom for two matrix tempemmres, and under different cooling protocols. The initial acceleration towards the surface is clearly seen in all cases, followed by a collision with the surface, and eventually cooling down to the lattice’s temperature. When the deposition is carried out at a low temperature and rapid cooling is applied, the ad-atom’s velocity is seen to change abruptly on impact (approximately 4 ps after the start). In all other cases, the velocity change is much more gradual. The bottom panels show the velocity profile of a bulk argon atom, for comparison. Fig. 2 shows the results of typical depositions carried out at 1 K and at 25 K (slow cooling). In the low temperatum case, the matrix contains many vacancies and its structure is far less ordered than in the higher one. At 25 K the crystal grows in an orderly fashion, reproducing fairlyaccurately the se&s structure. Since the structures shown are snapshots taken at random, they depict the ins&nto~ ow position of each atom. Therefore, deviations from the “ideal” structure, caused by thermal motiou are more pronounced at 25 K. The 1 K matrix, in contrast, shows very nearly the exact equilibrium location of each atom, since thermal motion is pm&ally arrested at this temperature. At25K,thenewlydepo&edatomisobservedto almost invariably end up on top of a “hole” formed by four neighboring atoms on the seed crystal’s surface, thus starting a new layer displaced by 0% in the pos-

R. Fraenkel Y.Haas /Chemical Physics 184 (1994) I S-204

25K

cooled atom

8OOv J

uncooled atom

600-

I

uncooled atom

. .

1’

lattice atom

‘“:I

“““““““” 7

0

4

8

12

16

.,.

0

8.3.

4

8

I.

12

16

0

time, ps Fig.1. Velocity profilesobservedat some typical depositionconditionsof argonas a functionof time. The left panel refersto a 5 K matrix,and the right hand one, to a 25 K one. The top partshows the velocity of the ad-atomfor rapidcooling conditions on impact (see text for an exact definition), the middle one the velocity of the ad-atomfor the case of slow cooling and the bottom one, the velocity of a bulk atom.

itive z direction. The crystal is found to grow layer by layer in an (ab} repeat sequence, maintaining, in this case, the fee structure. Occasionally, during the run, atoms were observed to accumulate in a restricted section and start the formation of steps. The “steps”, however, were never found to grow significantly, and were seen to collapse after a few more atoms hit the surface. The gradual growth of the solid phase may be visualized by taking snapshots at certain time intervals. Fig. 3 shows some snapshots taken for a deposition temperature of 5 K and using the rapid cooling mode. The matrix is seen to grow in an irregular pattern and to

contain many vacancies. This appears to be due to interception of newly arrived atoms by the previously deposited atoms, on its way to the surface - once trapped, the low temperature prevents any further motion. Annealing at 35 K followed by cooling back to 1 K led to the creation of a relatively well ordered crystal, but some vacancies still remain (Fig. 4). The simulated time duration needed to add 200 atoms to the crystal varied due to different starting conditions, but was typically about 2 ns. This rate corresponds to a deposition rate of about 3 X 106 mmol/s cm2 - many orders of magnitude larger than real deposition rates. The running time on a Silicon Graphics

A~~~~~f~~s~~~~~a~by using Letmar+Jones :wters characterixmg SF,, treated85 a sphere. Inthja case the sintuktion led to a bravely different sides: argon atoms tend to cluster around the guest molecules, and an ahnost.dendritic growth is observed for some time. Thus, four layers of argon atoma form near the gnest at- .svhiie ~Of~~~Of~~~~~C~S~~~ uncovered. This situation b expected to lead to a highly

Deposition at SK (fmt cooling)

Ind&o ‘I&IWOworh-station for reaching this size was

36 hours. 3.2. Depositionof guest molecules The s~ulation procedure was the same as for pure argon, except that the first species to be deposited on the template crystal was the guest molecule, to allow coverage by a maximum number of rare gas atoms. The disti~~on between a template atom and a deposited argon atom disappeared during the deposition process - as will be seen for all individuaksystems studied, the final structure around the guest bore no “memory” of the borderline between template and added atoms. Thus, the early deposition of the guest atom did not bias the final site’s structure in any observable way. The ~~~mction betwekn argon and guest atoms was assumed to be represented by a LJ potential, with the parameters listed in Table 1. 3.2.1.5’pkerkaL gtwts Xenon was chosen to s~~~~e case of a “small” perturbation of the potential. The resulting crystal lattice is very similar to that obtained, under simihir conditions, with pure argon. The guest atom was found to

R. Fraenkel, Y. Haas /Chemical Physics 186 (1994) 185-204

192

distorted lattice, riddled with faults. Simulating an even larger guest (CT=5.7 A, same E as for SF,) resulted in complete separation; the guest molecule was continuously pushed away from the growing argon matrix, and did not incorporate in it (Fig. 5). In a macroscopic sample this situation would probably lead to complete separation of the guest from the rare gas crystal. Experimentally, incorporation of large molecules into a rare gas crystal is much more difficult than into matrices l-361. The case of SF, was studied in some more detail [ 241. It was found that deposition at 25 K led to the formation of a rather symmetric site, in which one guest molecule replaces four argon atoms, two each from adjacent planes parallel to the 001 plane. The guest is located at the center of a tetrahedron formed by these four atoms, i.e. its center is midway between two argon atom layers. (The deposition was carried out over a template of four layers of argon atoms, and the center of mass of SF, is located between layers 5 and 6). The two dimensional projection of this four substitutional site is reproduced in Fig. 6a: it is composed of twelve equidistant nearest neighbors, each of which is slightly displaced away from the central molecule, with respect to the perfect lattice positions, forming a truncated tetrahedron (see also Fig. 2 in Ref. [ 241). When deposition is carried out at a low temperature and with rapid cooling of the ad-atom, a different situation was

Fig. 5. The result of co-deposition of a large molecule (a= 5.7 A) with argon. The argon atoms are pushed away from the guest and two “phases” are formed.

observed. Figs. 6b and 6c show some typical results; the number of atoms removed from the site in which the guest finds itself varied between 4 and 7, and they were located around it in a non-symmetric fashion. In some cases, the regular pattern of the subsequently deposited layers was disturbed, and a dislocation could be seen to develop. Annealing of such matrices at 35

1K deposition before annealing O

25K deposition after annealing

A

Fig. 4. The effect of annealing on a distortedmatrix.L.&tpanel: a 1 K deposited matrix,before and a&x annealingat 35 K. Right panel: a 25 K deposited matrix,with the same numberof atoms. (The snapshotswere taken at& cooling the matricesto 1 K.)

R. Fmenkel, Y.Haas/Chemical Phytics I&5(1994) MS-204

SQ

a)Z.%

cl=

inargon

a0

0.0

193

distortion of the lattice. The two argon atoms lying near the molecule along the moleculai’ruiis are essentially atthesamelocation~asinppureargoncrystal;The position of only four argon atoms are seen to change beyond,‘the thermal fluctuations: these era the atoms nearest to the N-N bond, symmehicaIly situated in the plane perpendicular to the molecular axis (one of them is m&ted by 1 in Fig. 7), as shown in Table 2 and Fig. 7. The vibrational frequency of the trapped Nz molecule can be calculated from the Fourier transform of the oscillatory motion of the N+N distance about the equilibrium position. It is found that the classical oscillation frequency (w,) is 2359.3 cm-’ in the gas phase and 2360.1 f 0.5 cm- ’ in the SS site obtained in the simulation. This result was obtained at several different deposition temperatures and cooling conditions. S~hurath’s group [3,37-391 reported an extensive study of the NBr molecuIe, discussing the data in terms of several distinct sites. Single or double substitution sites were considered, the former calculated to he more stable than the latter, by a large margin. In thepresent simulations [ 251, the potemials used were the same’as in Refs. [ 37391 (omWng the induetive cWribution, cf. Table 1) , and &position was simtil& as described above. It was found that the most frequekly formed site (30 out of 45 runs) at any temperature and uuder

F’ig.6. ‘IXe fitmchmof an argon matrix containhlg SF, deposited

WKte# diftkrentcoaditioas. only the layers closest to the guest molecule am shown: The left fmmz shows atoms in the nesmst layers (5 and 6), andthe rightoee &hows atomsin the nextnearestlayers (4 and 7). (the axfa ale in A.ngsqm units).Largecircles: atoms belonging to the lower laya in each tiame. Large sqllslcs: atoms belonging to the hi&u layer in each kame. Small symbols: the position of the atom in a perfectargoncrystal.(a) Deposition at 25 K, slowcooling.(b. c) Deposition at 5 K. fast cooling. The results of two independent runs am shown, to demonstmte the different

pa&us obtainWeundertheseconditions’in separateins. K led to a much more ordered lattice, the dislocation disappeared, but some vacancies were seen to persist. For example, a six-atom substitutional site, after annealing, becomes more symmetric, but does not revert to the apparently more stable four-atom substitutional site. 3.2.2. Non-spherical guests All simulations of deposition of a nitrogen (Nz) molecule in,argon perf& thus far (21 Gmrlation runs, under differentcunditions), led invariably to the formation of a single atom substitutional site, with little

Fig.7.Them

ofthctmppingsikofN2insrgon.Themolede (black symbols) occupies a single snbstitntionaisite; Aqon atoms whose position changed significantly by the preseace of tlic guest adownkgfey.TtIecshowthc~~~~fthe change. For nnmedcal dues see Table 2.

R. Fraenbl,

194

Y. Haas/ Chemical Physics 186 (1994) 185-204

Table 2 Displacementof argonatoms caused by trappedmolecules (in lob3 A) *

Molecule

Atom

N

Ax

Ay

AZ

N2b

AR

1

4

92

ncc

-72

-97*8

NBr d

1 2 3 4 5 6

1 2 2 4 2 2

nc -186 nc -40 -60 nc

-124 60 100 28 -63 -79

nc nc 253 34 13 -80

-124 178&2 235fl5 35fll 80f7 105f5

Br*C

1 2 3 4 5 6

2 8 4 4 8 4

75 53 nc -89 37 -130

-72 -73 25 nc 33 -105

nc -26 nc 85 nc 263

-108f5 -96i9 55flO 129f 17 23f5 320*9

’ Ax, Ay and AZare the displacement of the atoms from their equilibriumposition in a perfect fee lattice. The signs signify the directions of the change for the numbered atoms only, along the axes. The magnitudeof the change for equivalent atoms (defined in the text) is the same, directions depend on the location of each atom with respect to the molecule. AR is t _ change in the radial distance between the atom and the origin or the coordinatesystem, fixed to coincide with the center of mass of the molecule; Its magnitudeand sign are the same for all equivalent atoms (a minus sign means movementtowards the molecule). The errorbarsrepresentthe standarddeviation.N is the numberof equivalent atoms in each group. bFor the atom numberingcode, see Fig. 7. ’ nc - negligible change (less than 20 r& . ’ For the atom numberingcode, see Fig. 9. ’ For the atom numberingcode, see Fig. 12.

all cooling conditions was a singly substituted site (Fig. 8a). The bromine atom occupied an argon atom site, and the nitrogen atom was located in an interstitial site between two 111 planes, along one of the sides of a unit cell - essentially reproducing the structure proposed in Ref. [ 381. The nearest argon atoms are pushed away from the guest molecule, again in excellent agreement with Ref. [ 381. Since we have ignored the inductive term to the potential, its contribution to the determination of the trapping site structure appears to be minor. Fig. 9 shows the changes in the argon atom positions nearest to the trapped molecule in the site, termed the SSl site, and Table 2 lists the numerical values. As seen from the figure, the argon atoms whose positions are significantly changed by the guest molecule,

are located in a cylindrical symmetry around the NBr molecule. The numbered atoms in the figure represent motions of a group of atoms with the same symmetry. The largest shift is undergone by atom 3 and its equivalent, and the next largest by the two type 2 atoms. The stress induced by these relocations in the matrix, along with the force exerted by type 5 and 6 atoms, slightly compresses the N-Br bond. The four type 4 atoms, radially located with respect to the molecular axis, are only slightly shifted from their natural positions. Fig. 10 compares the changes in the N-Br distance as a function of time for NBr trapped in the SSl site (a) at 25 K with gaseous NBr (c) and displays also the corresponding Fourier transforms. The N-Br distance in the gas phase reproduces the experimental value of R,= 1.769, whereas in the SSl site, the bond

b. SS2 site

Fig. 8. The structureof different NBr trapping sites in an argon matrix,a plot of layers 5 and 6. The NBr molecule is located in layer 5. The symbols repmsentingthe argon atoms are explained in the caption to Fig. 6. The black circle representsthe Br atom and the nitrogenatom is denoted by N. (a) SSl site, (b) SS2 site, (c) MS site.

be 696.5*0,7 an-‘. Aamealing at 40 K lid to the convelsionofthaYiesrifEgtos$“1’ones‘~ngthe lattertobemorestabh!. Dqn&#ien iE, lower Uqreratures an& ‘Uder rapid Cookingcomhti~namwrpDted irtn$nn cases, in th0qpsarance’of m&i@ su (M@s)idw. Th4,$rnOtwe of one of these sitiisshown in @g. f@c,und Rig. 18 shows the metton ah&&e MLTJ,dManceaJ a Mction~oftime, and theko&t&dT mti. The vibr&ional fMc@Miiubw~cu4culti’fortheMl:sit&e varied in the rauge691~ cm- ‘, namely a blue&ifi of 0 to 3 cm-l with mqect to the ‘gauphake. Expel-imentally, two ntajor trappkig sit& erthibi~g a srmtll retf s?r@ with respect to the gas plMu weae’~obsk%ved in addition tothe SS site mention&, +I&$ w&e&own to be metastable, since anneaQng at 37 IK~&Ito a dramatic rcduct&r in t&it popslatloa; They were ~tentativtly assigned in R&r. [37-891 to d0uble substitutional CDS) sites. ‘The present ~simulation was not able to reproduce the ex@riarlenb9lI&& oftbi+ rod shifk The rklabon of thb tx&dWi ‘%%?I&es tr) tlie experimental ones is’discuss.kc!in Section 44.

&ctedby tlk

in Table

2.Thenllmkingcodeiscxplailiedintheiext. is seen to be “squeezed” a bit to an equihbrium .distance of 1.765 A.The vibrati~pmtio~~ of tiqe molecule is seen to be strongly modulated by the thermal motion of the lattice atoms, but the frequency, as &rived from the Fourier transform, is only slightly changed. We find a frequency of 691 cmrl for free NBr, compared to the ,eqnuimental value of 691.75 cm-’ and for the trapped molecule, 700.5 f 1 cm-‘, implying a matrix shift for this site of 8.6 cm- ‘. The experimental shift found for the most stable single sub stituted (SS) site (#5 in Ref. [39]), is 8.35 cm-‘. In addition to the previously described site, another SS site was sometimes observed (termed SS2), in which the Br atom also occupied the argon lattice site, and the nitrogen atom was roughly in an interstitial position, padEel to a 111 plane (Fig. 8b). The argon atom along the molecule axis, (near the nitrogen), is pushed away from its previous lattice position, while some argon atoms in the second shell of the NBr are missing. This site was obtained in six runs (out of 45)) when rapid coohng of the ad-atoms was used. ‘Ihe vibrational frequency of NBr in this site was found to

~,.‘oi-_.--,u 693

b.MS. 2% I+

,.‘oi I 00.511.522J30 lime, ps

l.‘&.i

zoo4006m~

frequency,cm-1

Fig. 10. Tim evolution of the N-Br distance at 25 K. (a) NBS in a SSl site. (b) NBrinaMSsi~. (c)Free (gasphm) molecule,Note theSligtdy8hostcriIWWomicdistaaceofti~~~in tllcSS1sim.aztdrhet3llhmd~oflowfnaqueacy*. cwadbytlm8lIo&y~mahfx.

1%

R. Fraenkel, Y. Haas / Chemical Physics 186 (1994) 185-204

The relative potential energy of the different sites could be readily calculated once their structure was determined. The total potential energy was calculated at 1 K, a temperature at which the contribution of the kinetic energy is minimized. The calculation was performed for matrices containing the same number of atoms, by averaging the total potential energy over a period of 10 ps. The results are summart‘xed in Table 3 and shown graphically in Fig. 11. The effect of the annealing procedure was checlced by calculating the energies of some chosen matrices that were deposited under rapid cooling conditions before and after annealing. It was found that annealing stabilized the matrices by lO-20%, making their energies about equal to the matrices prepared by the slow cooling method, or, at times, even somewhat lower. Twenty two deposition runs were performed with Brz as a guest molecule, using the potential parameters of Schurath [ 41 and of Herman and Beme [ 461. The experimental vibrational shift is in this case to the red, possibly indicating a weakening of the Br-Br bond energy and increasing of the bond length. It was found that a red shift is indeed obtained, albeit with a smaller magnitude than the experimentally observed one: the frequency of the trapped Brz was 324.3f0.2 cm-‘, while that of the isolated molecule was 327.5 f 0.1 cm- ’ (the experimental red shift is 7 cm- ’ ) . The bond length is increased slightly, from 2.8209 to 2.832 A. As seen from Fig. 12, the diatomic molecule fits into a double substitutional site, in which it replaces two nearest neighbors argon atoms. The argon atoms placed next to the bromine ones along the bond’s axis (atom 1 and its equivalent on the other side of the molecule in Fig. 12 and Table 2) are slightly shifted inwards, leading to a pulling force that causes the slight elon-

gation of the Br-Br bond. The eight type 2 atoms also contribute to that effect. The largest geometry change, however, is found for the four argon atoms nearest to the molecular axis and displaced radially from it (the four type 6 atoms in the figure), that are pushed OUCwards from the molecule. Numerical values are listed in Table 2 for the matrix atoms nearest to the guest. In contrast with the cases of N2 and NBr, atoms in the second and even third layer away from the guest are appreciably displaced ( = 50 n&) from their equilibrium position in the pure argon crystal. It was noted by Pratt and Chandler [47], that since the Ar-Br potential depends on the polarizability of the molecule, which in turn depends on the Br-Br separation, a LJ potential with fixed parameters cannot be realistic. Herman and Beme [46] proposed a modified potential, in which E for the Br-Ar interaction depends on the Br-Br distance and showed that it leads to a red shift upon solvation in liquid argon at 300 K. However, using this potential in our simulation did not change the geometry of the site., but the calculated vibrational frequency of the trapped molecule shifted slightly to the blue. Beme’s proposed correction was found to lead to a red shift for bromine in liquid argon at 300 K, a result reproduced in our calculations. Herman and Beme found that increasing the pressure led to a blue shift - a physically plausible result. Thus, using the actual geometry of an argon crystal, which is equivalent to the high pressure conditions in Beme’s model, we calculate essentially the same blue shift as Ref. [46]. However, experiment shows that the shift is to the red. It thus appears that in the more rigid solid structure, the empirical proposed correction is insufficient.

Table 3 Total potentialenergy of matricescontainingNBr, results of differentdepositionconditionsl Deposition temperature(K)

cwing~mode

nb

EnergyC(X103cm-‘)

UC ( x 103cm-‘)

1 1 5 5 13 13 25 25

slow fast slow fast slow fast slow fast

3 4 8 13 3 4 5 5

- 268.6 - 257.5 -271.4 - 262.0 -271.8 - 266.4 -271.7 - 267.6

2.3 0.9 1.5 4.0 0.4 4.1 1.7 1.6

’ 420 argon atoms in each matrix. ’ n is the numberof depositions. c The averageenergy for n depositions. CT is the standarddeviation.

R Fnsenkel, Y. Haas / Chemical Physics I86

(1994)185-204

197

that the pmpertiw of the trapped molecules, such as

-2.7: 10’ -2.6: 10’ -2.6; Id energy, 12

-2.5: 10’ -2.;

cm-l

,

I

10 8

n

6 4

theirspectra,arenotsignificantlycsasegad.Ttremethod has also been used to study reaotion naeehauisms, mostly by identifying reactive intermediatas. In this case the trappad molecules are N~hangedby the reaction, and one has to consider the effect of the matrix on the dynamius ofthe reaction. Since translation& and often also rotutional, motions are severely restricted in a matrix, the properties of the site in which the reacting molecubs are located may affect the ‘rate and the outcome of the reaction. It has therefore become increasingly evident, that proper interpretation of trapped molecules’ dynamics and reactivity require detailed knowledge of the site’s properties [ 4Q,413. The simulations were primarily aimed at reproducing the prqerties of the immediate neighborhood of the trapped species. The main, purpose was to predict the characteristic fea@res of different trapping sites, and to help in accounting for the spectral and dynamic properties of molecules trapped in’them. Therefore, the model cannot be applicable to properties related to larger scale dimensions, such as light saWzing. Nevertheless, insight can be g&xi on some wch features as in the case of diffusion of atoms in tha matrix. Several mechanisms have been propqn&, the simplest one being vacancy migration, whieh is essentially the only

2 0 -2.8’10’

-2.74 10’ -2.68’ 10’ -2.62’ Id energy,

-2.5; 10’

-2.5’10’

cm-’

Fig. 11. (a) A bistogmm-wingthe number of depositionsresulting in NBr containingmatria3swith a given total potentialenergy, pm~bydi~rentmethdds.~h~isu)ocm~‘wl&.The~l .codeontbeupperrigbtd4mtotestbedepositiontempcatam(K)and thecoolingm~-f~(fast)cootiagis~byfc,slow coolingbysc.(b)The~~stogram,aaangedtoshowthewmber of depositions leading to a given site stmctum.

4. Dlseusslon The main application of matrix isolation is the trapping of a given species in a solid, inert environment, so that its ptqerties can be studied at leisure. It is assumed that the interaction with the sunounding atoms is small compared to intramolecular forces, and

Fig. 12. llle sauctme ofaBr,tmppiagsite.Titecolorcodefortbe [email protected] tbelugsnMlberofalpllatWllt3WhocspaaidonC~.SiglliIicfUltly.eScOmpRdtomatike+conEaioingN,aIldNBr.

198

R. Fraenkel, Y. Haas / Chemical Physics 186 (I 994) 185-204

effective one in a solid crystal. In a real matrix, the contribution of motion along grain boundaries, for instance, cannot be neglected. The model provides some insight as to the formation of vacancies, and predicts that a matrix prepared at low temperature and under rapid cooling conditions, contains many vacancies. It is thus expected to allow more facile diffusion than one prepared under self annealing conditions, in agreement with experience. Ideally, the only required input consists of the interaction potentials, and experimentally controllable parameters such as the support’s temperature and the gas pressure and flow rate. For simplicity, we haved used simple Lennard-Jones potentials for the Ar-Ar interaction, but more realistic ones are used in the literature [42]. Comparison of the results with experiment should help in assessing the validity of the potentials used, and also to introduce modifications, if necessary. 4.1. Matrix growth kinetics If equilibrium conditions can be maintained during deposition, the thermodynamically most stable configuration should prevail. In practice, one finds that the site distribution can be altered by changing the deposition temperature and rate, and that metastable sites can be generated. It is thus accepted that matrix isolation is often kineticaEZycontrolled. We found that deposition in an argon matrix at 25 K leads primarily to a single site, that appears to be the most stable one, as checked by annealing and by energy calculations (see below). Other sites were found to be formed only when more extreme conditions: cooling to 13,5 or, better, 1 K was required, along with the rapid cooling procedure. By following the course of the deposition, it could be established that the new sites could be generated, provided the motion of the incoming atom was restricted as soon as it hits the surface. Otherwise, local heating due to the energy transfer from the incoming atom led to local annealing and thus to the energetically most favorable site. As Fig. 1 shows, if the ad-atom is allowed to move solely under the action of the potential when the matrix was held at 5 K, its velocity profile is similar to that of an atom hitting a 25 K surface. Indeed, under these slow cooling conditions, the matrices generated were very similar. This point is also illustrated by the energy calculations shown in Table 3 - the total energy of matrices prepared under slow cooling con-

ditions was very nearly the same at all temperatures, the difference between 1 and 25 K depositions being 1.1%. In contrast, under rapid cooling conditions, the difference is 4%. It thus turns out that in order to reproduce the experimental observation that deposition at low temperatures is known to result in a much more porous matrix [ 411, rapid cooling must be used. It was found that in this case the atom was largely restricted to stick to the location reached randomly; This trend led to the creation of point defects, in which semi-dendritic, rather than epitaxial, growth was observed. In the presence of a foreign guest, several different trapping sites could be generated when the guest did not fit neatly into a substitutional site. The physical reason for the more realistic result obtained by the rapid cooling mode may be traced to the fact that at low temperatures, all matrix atoms are essentially at the bottom of the potential well. The incoming atom loses much of its kinetic energy upon collision with the stationary atom, and takes its place, being thus effectively frozen into place. At higher temperatures, the matrix’ atoms are most often found, on the average, at the classical turning points and a collision with the incoming atom sets both in motion around the collision location. The rapid cooling mode is thus “natural”, in a classical system, only for systems at or very near to the absolute zero. Experiment shows, however, that disordered matrices can be prepared also at higher temperatures ( 10-15 K for argon). This can be accounted for by recalling that at these temperatures most of the atoms occupy the lowest vibrational level of the lattice, and as any quantum oscillator are found mostly in the equilibrium position and not near the classical turning points. Collision with the incoming atom can excite this atom to a higher energy level, trapping the incoming one in the zero point energy level. The rapid cooling mode can be viewed as a manifestation of the quantum nature of the system, which affects the dynamics to a larger extent at low temperatures. The restricted number of simulations performed so far limits the statistical validity of our conclusions, but it appears that the number of distinct sites was always small, even under extreme cooling conditions. This result, that agrees with experimental observations, is explained by the fact that in the immediate vicinity of the randomly selected impact point, there is only a small number of preferred sites, that can be reached

R Fmenkel,I’. Ii-/

by very small 4Y4mph&

ChemicalPhysics1% (1991) M-204

motions. EvEonitoFing the

motionoftheincoluingatomwoBnd,thatonatrivalat thesulfaoe,itsmotionisrestrietedtotbsvicin&yofthe 1 siteiteneouW&fhdanditsneamstneighhors - long distance motions were not observed Ieven at2&ILTheleasamigkmamoftlmatomieposiCionsto reach the final oonfarmation is&em a huge number drai#domsmall motiotrs, and in no case was , sewi to move 424~1tily anarmorany

acrwadiitsw32wlaFgeasalkWhh~ter. ARnealhngk&y$im~byti~~‘metlhod.It may tveused to i&W@ the more s&Mesites, but we foundthatitclblanot~al~~ci~,~~,wbern a “porous” matrix was@qared at a low @mpuratum, annealing was seen to remove most, but not all of the vacanoiee. Once a vaMutcy is ttqpedby an orderly array of ‘host atoms, each in its eqmMtrlum ~pasition, small amplitude motions are unable to remove it. ‘Ihus, a matrix deposited under conditions favoring the creation of vacancies, will contain many ofthem even after annealing. Acomparisonwiththerealgrowthroteofamatrix is of some intorest. A typical deposition rate is 2-3 mmolJh cm3 - about nine orders of magnitude (!) slower&an the simnl~rate, W fastest experimental rates are obtained in the single crystal gWwth experiments, and are of the order of 100 mmol/h cm* - still much slower than the simulated rate. Tlkis result has been noted in other MD 5imuM0ns [X3], ml *flets the~thatinonlctlr~~ve~coanputar;tlmo,cyrwr,rsleusiets the a&atoms one after tlte other in ra@id 5uwasion, very close to the surface. In teal systems, deposition rates are determinudby d&fusion and by mass and heat tmnqort rates. The low thermal WMIWMQ of solid rare gases has to be ,talren intoac&otmt when a macro scopic matrix is formed. B we allowed for essentially intXantaneonscooling,mabingrapiclsuccession of ad-atoms possible. 4.2. Site sts-mure The sbwtwe of a trapping siteshas often been estimated by allowing the guest molecule to occupy aspace

n&d slmulatlons, that the meat stable shes are those, t&&h intxudrm&W~tiest d&&on into the latthk A weakness ofthi@ method is that it may fail to

199

detect possible energy minima that involve a simultaneous change in the positions of several atoms. The approach presented here, offers a way to circumvent this limitation and thus could help to reveal less intuitive structures. It should also be useful in recovering local minima, that may be kinetically stabilized. Therefore, less personal bins is involved, and the method should be limited only by the accuracy and validity of the potentials. lhe difficulty in growing a rare gas crystal containing an impurity is for by the model. An illuminating case is mokule is apparin an armn lattice, ently too big to be and indeed, a crystal (in d&in&ion with a matrix) could not be grown [ 16,171. In xenon, the cavity formed is just large enough, and a dihrte crystal can be grown. ‘Ihis interpretation is subwantiatedby our calculations [ 241. Analysis of the high resolution IR spectrum of SF, in argon led to the suggestion of a tetrahedral site formed by the removal d four atoms “though it is difficult to display” (Ref. [ 161, p. 3206). As we have shown, the tet&edWsitekereakednaturallybylasinga(ppr;opljat&Upot&laials.~~~found inthesimula&mconWinstheC~axesthatwereproposed to explain splittings. SitloG a real lSFBmolecule is not spherical, it may occupy the site in some orientations that lead to slightly dif%t%nt~speetralshifts. For instance, a stable tetrahedral site aau be visualixed, in which the C3 axis of SF, coinbides with the axis connecting the centers of the triangles with the centers of the hexagons. The C, axes of the trnncated tetrahedron coincide with the thme Cs axes of SF6 ‘(Fig. 13) and the three F atoms, when viewed down the C3 axis, are seentoberotatedby6@wi&mqnWtotbethreeargon atoms. A high energy conf~&mtion is erqected if the SF, molecule were rotated atom& this axis so itbat the threeFatomswouldliecl~to~~~~s. A subsequent suggestion [ 21, that o the atom, substitutional site is present inthe matrix, is not Wbatantiated by our results, In fae& no such site could be found in any of the runs up to now. Another experimental observation concerning SF,

was the cuppeazancs d &&hwly bwul p&s ‘&atdisappear upon rvmnreraling [ 161. The authors t,qpsst6d that them are due to motto&n&e SFd in a loose, relatively unstable site, aad prrqokl ‘a fiwcatam substitutioud

situ @hich collapses to tlm more stable four atom substitutional site on’annealimg. ‘Such sites are

200

R. Fraenkel, Y. Haas / Chemical Physics 186 (1994) 185-204

commonly observed site in the simulations by a large margin. Mote extensive calculations, that will include other site geometries such as hcp, may be required in order to simulate the formation of the red shifted metastable sites. For both symmetric diatomic molecules, N2 and Br,, only one predominant site was found: the former appears to fit in snugly into a single substitutional site, while the latter, found in a double substitutional one, leads to significant changes in the location of only the four argon atoms lying near the molecular axis and perpendicular to it; These atoms are not expected to affect the molecular frequency. The formation of only one structure, for these axially symmetric molecules even under extreme cooling conditions, probably relates to the fact, that less stable sites are formed if removal of an argon atom releases stress in a molecular bond. Fig. 13. A possible cotiguration for SF, in a tetmhedmlsite, formed by the 12 argon atom truncatedtetrahedron.

4.3. Sitestability

indeed formed in our simulation, when deposition is carried out onto a cold matrix and using rapid cooling. It is noted that in the cited experiments, deposition was carried out at the relatively low temperature ( 10 K) . LeRoy and coworkers [ 431 studied the structure and the infrared spectra of SF&r,, clusters (II = l-100). and were able to calculate the spectral shift using the anisotropic potential of Pat et al. [ 441. The same potentials may be used in a calculation involving matrix isolation. However, Perera and Amar [45] used LJ potentials to generate clusters, and their results for 35 K are quite similar to those of LeRoy et al. at 30 K. Thus, the gross features of the interaction are faithfully reproduced with the simpler spheric potentials, at least at relatively high temperatures. The case of NBr also illustrates the main theme of this work - a site geometry “naturally” evolved from the model. We found that under the conditions of the simulations, DS (double substitutional) sites rarely form, as expected for a site of much higher energy compared to a SS one. Two SS sites were obtained, for which the stretching frequencies are 696.5 and 700.5 cm- ‘. Seranski et al. [ 391 assign the most stable site for which the NBr stretch frequency is shifted by 8.35 cm-’ to the blue of the gas phase value, to a SS site. This result is substantiated by our model: we find a SS site shifted to the blue by 8.6 cm- ‘, which is the most

The Lennard-Jones potentials used in this work proved to be satisfactory for the purpose of simulating the trapping of guests in a fee lattice. As mentioned in Section 2, these potentials will not lead to the formation of a fee crystal in a free growing system, or, when any other lattice plane except the 001 one is used as the deposition plane. As shown by our simulations, by using the 001 plane as a template, both the structure and the energy of an argon crystal were reasonably well reproduced. The calculated energy of the argon crystal was 82.4 mevlatom, close to the experimental value (Ref. [ 281, ch. 4, p. 242). Introduction of an impurity into the solid can lead to considerable distortions in the perfect structure, and one of the main motivations was the attempt to clarify the formation mechanism of different sites, and to determine their relative stability. It is thus of interest to compare the energetics of different sites, or, more precisely, the relative energies of the crystallites containing different trapping sites. We illustrate this point by reference to NBr containing matrices, in which the deposition simulation revealed several different sites. We have calculated the total potential energy of the matrices obtained, using Eq. (3)) the summation being made over all the argon atoms. In each case the number of argon atoms was 420 but the structure of the matrices differed considerably. The results for 45 depositions

R Fraenkel, Y. Htws/Ck&d

aresummaWd4nTable3andinF~. W the toWIemgy for the s-s S$l s$te.was~y.~%

ll.&wasfound conWing the

4ttc to the inkractio and not necessarily between the argon atoms and the guest. Thus if one calculates only the &EmsN&zing directly to NBr (i.e. omitting the last term fkom Eq. (311, ~o~~~~

isthemain~conof the m#rix, At h&her motion of the atoms is large enieqb tQ@at&a@ mo&of tJbW?B*g.

~,~D10Elt,~lesit&:,(SSl~ihis’oase),isthelone that confkrms be& with that S-W, a&bough it involves local stmss in the vicinity of the guest and in the. htr~~&c&~ ,Wd: of the d&pm& m&e&e. w&Xl the ~~~ti~~ ~o~~~ ,@eFmitf&Wenergy exchange this is the mpst probable ~~~~ as shoarn by ti fact that qt $J $4cx m&r slow cooling, most,of the dep~itions ended in ita formation.

Physics I86 f1994~IS-204

201

202

R. Fraenkel, Y. Haas/ Chemical Physics 186 (1994) 185-204

vibrational parameters w, and W,X,by fitting a Morse potential to the effective potential (cf. Eq. ( 1) ) in the direction of the molecular axis. It could be shown, that with the parameters used in the calculation, only the SS site was stable. The DS site was used as an example of a metastable site that may be present in the matrix, although it cannot be formed under equilibrium conditions. The MD method, in contrast with the methods discussed above, allows the attainment of the final matrix structure, without any explicit initial guess. As in the other methods, its ability to reproduce experimental data depends on the accuracy of the potentials. However, for a given potential it has the advantage of a simpler physical interpretation, and provides a straightforward means to study dynamic processes, such as a chemical reaction. In the present calculation, the vibrational spectral shift for a diatomic molecule AB is obtained simply from the Fourier transform of the time dependent separation between A and B. In the classical limit, all the information needed to calculate this motion can be incorporated into the potential and thus other parameters, such as the dipole moment of the molecule and its polar&ability are not required. Thus, the shift provides a simple means to check the validity of the potentials used, by comparison with experiment. The results show that in an argon matrix, the shifts can be related to forces along the vibrational coordinate. Considering only the most stable sites, one finds Nz is situated in a site that exerts no strain along the bond, and the frequency shift is negligible. The site in which Br, is found exerts a force that pulls the two atoms outwards, loosening the bond, and a red shift is obtained. NBr is in a tight one, the bond is squeezed somewhat, and the vibrational frequency increases. It has been noted by many authors [ 541, that simple two body interaction potentials do not lead to the correct prediction of frequency shifts, particularly if they are to the red, and that three body forces must be explicitly considered. Mirsky and Manz [ 551 showed that an optimized geometry based on pairwise potentials similar to the ones used here ‘(they used a Buckingham potential instead of a Lennard-Jones one) could not reproduce the experimental red shift for CO in an argon matrix. An empirical relation, in which the intermolecular potential between the argon atom and the atoms in the molecule was allowed to depend linearly on x = R - &, was found to lead to the desired agreement

with reasonable parameters. At the same time, the small three body correction did not affect the trapping site’s geometry or stabilization. Such empirical relations can thus be a means to obtain intermolecular distances and force constants for trapped molecules [ 5O,H].

5. Conclusions The simulation has been shown to reproduce some salient features of rare gas matrices, such as the fee crystal structure, formation of several trapping sites and the more porous nature of matrices deposited at low temperatures. Quantitative calculations of the site’s stability can be made, and the spectral shifts can be evaluated. The model appears to successfully reproduce blue spectral shifts in the vibrational spectra, but not the red shifts. This deficiency is probably due to an under-estimation of the attractive intramolecularpotential with respect to the attractive part of the intermolecular one. Several improvements of the potentials may be used to remedy this situation: one is by adding three body interactions; Another, by using more realistic pairwise potentials, particularly in the attractive part. Such potentials can be based on empirical expressions, or on ab initio calculations. Nevertheless, the simple Lennard-Jones pairwise potential was found to satisfactorily reproduce the argon crystal structure (to within 0.01 A) with the starting template. It is therefore believed that the gross features of the site structures obtained will not be substantially changed by the use of more accurate and sophisticated potentials. As mentioned in the Introduction, hcp trapping sites are frequently observed in rare gas matrices, even though the stable form of the pure crystal is fee [28,29,56]. They are believed to be part of threedimensional hcp domains, that may form during the vapor deposition process. The method can be extended to include simulations of the formation of these sites by using a 111 crystal face as the growing side of the template crystal. The rapid cooling mode, which was introduced in order to obtain better agreement with experiment, turns out to be a better approximation than the slow cooling one. It was found to reproduce the “open” structure that is known from experiment to be formed upon rapid deposition at very low temperatures. When guest mol-

R. Fraedd, Y. HaadChcmical Physics186 (1994) 185-2Od

eculeswtxepreslesrt,&iamcdewasfkmdtobeessential in promoting the formation of less stable trapping sites. In practi&, lattice illl$&&tiotts lklay al&o‘be introduced by the sim+neou& grixvth of ,ii uktrix from mllnr nurFzaatisn*m, reauM+&!, #J$dPP~*Ps andlr#%wmi~.-tke~~~~~‘Gallbe used to simulatMhatIJrooessalso, fiseds,ti b&6%M&ed

by furthercomp6ter b@etimen&.

Achroosllsdgemeat We thank Prof-r RB. Gerber for many enlightening discussions, and Iko&%sorR. Naaman for a delightful converq+qon. The Farlqs Ce&er for Light InducedProcesses is supporkd by the Miia Gesellschaft fik die ForaohlmgmbH, Muni& Gemlany.

References [ I] HE. Hallam, ed., Vibrational svpy

of trapped species (Wiley, London, 1973); B. Meyer, Low temperature spectroecopy (Elsevier, Amsterdam, 1983); ML. Klein, ed., Springer series in clnunical physics, Vol. 4. Inert gases (Springer, Berlin, 1984) ; L. Andmws and M. Moskovits, eds., Chemistry and physics of matrix-isolated species (NorthHolland, Amsterdam, 1987). [2] B.I. Swanson and L.H. Jones, in: Vibrational spectra and structme. Vol. 12, ed. JR. Durig (Blsevier, Ams&dam, 1983) PP. 2-67.. [3] A.C.Becker,K.-P.LodemarmandU.Schurath,J.Chem.Phys. 87 (1987) 6266. [4] J. Langen, K.-P. Lodemann and U. Schurath, Chem. Phys. 112 (1987) 393. [5] R. Schriever, M. Chergui and N. Schwentner, J. Chem. Phys. 93 (1989) 9206; E.I. Tarasova, A.M. Ratner, V.M. Stepanko, LYa. Fugol, M. Chergui, R. S&lever and N. Schwenmer, J. Chem. Phys. 98 (1993) 7786. [6] R. Alimi, R.B. Gerber and V.A. Apkarian, J. Chem. Phys. 89 (1988) 174. [7] D. LaBrake and B. Weitx, Chem. Phys. Letters 211 (1993) 430. [8] A.D. LeClaim, in: Physical chemistry, an advanced tmatise, Vol. 10, ed. W. Jest (Academic Press, New York, 1970) pp. 261-331. [ 91 H. Fmi, in: Vibrational spectra and structure. Vol. 20. ed. J.R. Durig (Elsevier, Amstetdam, 1992) pp. 2-65; H. Frei andG.C. Pimentel. in: Springer series in chemical physics, Vol. 4. Inert gases, ed. M.L. Klein (Springer, Berlin, 1984) p. 139; Ann. Rev. Phys. Chem. 36 ( 1985) 491; R.N. Peru& Chem. Rev. 85 (1985) 77.97.

203

[ 101 M. Fs$ardo and V.A. A&a+. J. Chem, phys. 85 (1986) 5660; 89 (1988) 4124,54102. 111I ME. Fjardo. J. Ckem. Pkys. 98 ( 1993) 119. [12]D.TbinuPsldHJ.Bnulda~BJ.BanAJ.~.83 (1985) 230; L.M. Raff. J. cboso, Plrya. 93 (1990) 3160; R. Allmi,A.Br&snnnnndR.B,Gailmr,J.ahmn.I1Rys.91(1989) 1611;R.Alimi.A.V.ApkarianandR.B.Gu&~u,J.C!hen~Phys. 92 (1990) 35511 R. Al&@ R.B. Gesb nnd A.V. Apk&an, Phys. Rev. I&em 66 (1991) 1295. [ 131 S. f%&nnkuv R. Viarut&1. Primelo and Id. B&et, J. Mol. Struct. 161(1987) 265; L.C. BallingandJ.JwHgl#, JhChem. P&s. 81 (1984) 675; J.S. Win& J. C&&t. Phys. 94 (1991) 5275. 1141 J. Homes and J. S&&r, C&em. Phys. 74 (1983) 433; JJ. Ww andL.C. Balling, J. Clam, phys. 73 ( ldaB0)944; J.C. ~,RL.Moway,b~.~&W.Kim.P.N.San;rtzand L. AndPaws, J. c&ml Pltys. 74 (1981) 6349. [ 151 R. Gunde, H.J. Keller, T.-K. Ha and H.H. Gu&luud. J. Phys. Cbm. 95 (1991) 28(12;V.J.Lang and J.S. Wkm, 1. Chem. Phys. 94 (1991) 92To: P. Bresse& amm. phys. 172 (1993) 315. [la] B.I.SwansonandL.H.Jones, J.CkemPhys.74 (1981) 3205. [17] L.&JonesbdB.1.Swansen,J.C!bm.Phys.79(1983) 1516. [ 181 R. Gmuk. P. Felder and ‘H&H. G&h&d. Chem. Phys. 64 (1982) 313. [ 191 R. Fraenkel, U. Samuni, Y. Haas and B. Dick, Chem. Phys. Letters 203 (1993) 523. [20] J. Wolf and G. Hohlneicher, Chem. Phys. 181 (1994) 185. [21] M. Ptnger, B. Asmussen and C.J. Carlyle, J. Chem. Phys. 89 (1988) 1030. [ 221 E. Burke, J.Q. Bmughton and G.H. Gilmer. J. Chem. Phys. 100 (1994) 247. [23] B. Dodson. Crit. Rev. Solid State Mater. Sci. 16 (1990) 115; M. Schneider, Phys. Rev. Letters 55 ( 1985) 604; S.M. Paik and S. Das !%uma,Phys. Rev. B 39 (1989) 1224. [24] R. Fraenkel and Y. Haas, J. Chem. Phys. 100 (1994) 4323. [25] R.FraenkelandY.Haas.Chem.Phys.Lettem22O(1994) 185. [26] J. Jortner, Z. Physik D 24 ( 1992) 247. [27] G. S&es, ed, The chemical physics of atomic and molecular clusters (North-Holland, Amsterdam, 1990). [28] M.L. Klein and T.A. Venables, Ram gas solids (Academic Press, New York, 1976). [29] GL. Pollack, Rev. Mod. Phys. 36 (1964) 748. [30] J. Farges, M.F. de Femudy, B. Raoult and G. To&et, Advan. Chem. Phys. 70 (1988) 45; G. Torchet. J. Fargea, M.F. de Feraudy and B. Raoult, in: The chemical physics of atomic and molecular clusters, ed. G. Stoles ( North-Hollaud.Amsterdam, 1990) p. 513. 1311 B. Raoult, J. Farges. M.F. de Femudy and G. Torchet. Phil. Mag. B 60 (1989) 881; B.J. van der Waal, J. Chem. Phys. 98 (1993) 4909. [32] T.KiharaandS.Koba,J.Phys.Scc.Japan7(1952)348;B.W. van der Waal. Phys. Rev. Letters 67 (1991) 3263. [33] Q. Chang and V. Buch, J. Clnun. Phys. 92 (1990) 5004. 1341 M.P. Allen and D.J. Tlldesley. Computer simulation of liquids (Oxford Univ. Press, New York, 1987).

204

R Fraenkel, Y. Haas /Chemical Physics 186 (1994) 185-204

[35] L. Verlet, Phys. Rev. 159 (1967) 98; W.C. Swope, H.C. Andersen, P.H. Berens and K.R. Wilson, J. Chem. Phys. 76 (1982) 637. [36] B.A. Swanson, L.H. Jones, J.A. Ekberg and H.A. Fry, Chem. Phys. Letters 126 (1986) 455. [37] AC. Becker and LJ.Schurath, Ber. Bunsenges. Physik. Chem. 91 (1987) 1238. [38] M. Winter, K. Seranski and U. Schurath, Chem. Phys. 159 (1992) 235. [39] K. Seranski, M. Winter and U. Schurath, Chem. Phys. 159 (1992) 247. [40] W. Rudnick, R. Hansel. H. Nahme and N. Schwentner, Phys. Stat. Sol. (a) 87 (1985) 319; W.G. Lawrence and A.V. Apkatian, J. Chem. Phys. 97 (1992) 2224. [41] W. Langel, W. Schuller, E. K&ringer, W.H. Fleger and H.J. Lauter, J. Che.m. Phys. 89 (1988) 1741; E. Kn&inger, W. Schuller and W. Langel, Faraday Discussions Chem. Sot. 86 (1988) 285. [42] R.A. AZ&, in: Inert gases, ed. M.L. Klein (Springer, Berlin, 1984) p. 5; R.A. Aziz and M.J. Slaman, Mol. Phys. 58 (1986) 679. [43] D. Bischenauer and R.J. LeRoy, J. Chem. Phys. 88 (1988) 2898;M.A.KmeticandR.J.LeRoy,J.Chem.Phys.895(1991) 6271; D.J. Chartrand, J.C. Shelley and R.J. LeRoy, J. Phys. Chem. 95 (1991) 8310.

[44] R. Pac. J.J. Valenthdand J.B. Cross, J. Chem. Phys. 77 (1982) 5486. [45] L. Pemra and F.G. Amar, J. Chem. Phys. 95 (1990) 4884. [46] M.F. Herman and B.J. Berne. J. Che.m.Phys. 78 (1983) 4103. [47] L.R. Prattand D. Chandler, J. Chem. Phys. 72 ( 1980) 4054. [48] A.D. Buckhrgham, Proc. Roy. Sot. A 248 (1958) 169; A.J. Barnes, in: Vibrational spectroscopy of trapped species, ed. HE. Hallam (Wiley, New York, 1973) ch. 4. [49] G.C. Pimentel and S.W. Charles, Pure Appl. Cbem. 7 ( 1963) 111. [50] H. Friedman and S. Kimel, J. Chem. Phys. 43 ( 1965) 230, 3925. [51] A. War&l and S. Lifson, J. Chem. Phys. 53 (1970) 582. [52] S. Lifson and A. Wamhel, J. Cbem Phys. 49 (1968) 5116. [53] A. Warshel and M. Levitt, QCPE 11.247 (1974) - QCPE/PI 247: a program for consistent force field evaluation of equilibrium nuclear configurations and vibrational fmquencies of molecules. [54] A.I. Krylova, A.V. Nemukhin and N.F. Stepanov, J. Mol. Struct. 262 (1992) 55. [55] K. Mirsky. Chem. Phys. 46 (1980) 445; K. Mirsky and J. Manz, Chem. Phys. 46 (1980) 457. [56] A.C. Becker and U. Schurath,Chem. Phys. Letters 142 (1987) 313; C. Blindauer, M. Winter, 0. Sild, G. Jansen, B. Hess and U. Schurath,J. Phys. Chem. 97 (1993) 10002.