13778
J . Phys. Chem. 1993,97, 13778-13787
Sorbate Properties and Cage-to-Cage Diffusion of Argon in NaCaA: A Molecular Dynamics Study+ Subramanian Yashonatb'$t*gand Prakriteswar Santikaryt Solid State and Structural Chemistry Unit and Supercomputer Education and Research Centre, Indian Institute of Science, Bangalore-56001 2, India Received: July 15, 1993; In Final Form: September 27, I993@
Detailed molecular dynamics simulations of argon in zeolite NaCaA are reported. Thermodynamic, structural, and dynamical properties of the sorbate as a function of temperature have been obtained. The properties calculated include various site-site radial distribution functions, different energy distribution functions, selfdiffusion coefficients, the power spectra, and properties relating to cage-to-cage diffusion. The results suggest that sorbate is delocalized above 300 K. Both modes of cage-to-cage diffusion-the surface-mediated and centralized diffusion-are associated with negative barrier heights. Surprisingly, rate of cage-to-cage diffusion is associated with negative and positive activation energies below and above 500 K. The observed differences in the behavior of the rate of cage-to-cage diffusion between Xe-NaY and Ar-NaCaA systems and the nature of the potential energy surface are discussed. Presence of sorbatezeolite interactions results in significant enhancement in the rate of cage-to-cage diffusion and rate of cage visits. It is shown that properties dependent on the long-time behavior such as the diffusion coefficient and the rate of cages visited exhibit the expected Arrhenius dependence on temperature.
1. Introduction Considerable interest has been generated during the past few years in the study of molecules adsorbed in zeolites. Computer simulation studies have proven to be highly useful in the study of these systems. Faujasites and silicalite (and its analogue ZSM5) have been, by far, the most widely studied structural types. Faujasites have relatively large cages interconnected through windows. Silicalite, in contrast, consists of intersecting straight and sinusoidalchannels. Properties of a number of small to large sized adsorbates sorbed in silicalite have been investigated by computer simulationduring the past few years.'" For example, June et a1.I have investigated the sorption of alkanes at low occupancy in silicalite. Demontis et a1.* have investigated methane in silicalite. Several simulations have been carried out on zeolite Y and its The grand canonical ensemble Monte Carlo simulation of Rowlinson and c o - w ~ r k e r son ~ , ~xenon and methane sorbed in zeolites X and Y yielded significant insight into the structure as well as thermodynamic aspects. Dynamical studies of adsorption of small molecules such as xenon and benzene have also been reported.l&l2 In comparison to investigationsof zeolite Y and silicalite, computer simulation investigations of sorption in zeolite A have been rather limited. There are, however, a few recent studies in the literature. Demontis et al.,13for example, have investigated structural and vibrational properties of anhydrous phase of Linde zeolite 4A. Kono and Takasaka14 have investigated the sorption properties of argon and nitrogen in dehydrated 4A by the Monte Carlo method. van-Tassel et al.15have reported Monte Carlo calculationsof Xe in NaA. Earlier '29XeNMR spectroscopic investigations of Xein NaA have shown that Xe cannot diffuse from one a-cage to another through the eight-ring window.16 In this work, we report molecular dynamics (MD) results of argon in NaCaA zeolite. We comparethermodynamic,structural, and dynamical properties obtained for the present system with the experimental results. Results have also been compared with the results of MD simulations on Xe in NaY zeolite. These Contribution No. 962 from the Solid State and Structural Chemistry Unit. *Solid State and Structural Chemistry Unit. s Supercomputer Education and Research Centre. Abstract published in Aduance ACS Absrrarrs, December 1, 1993.
0022-365419312097-13778$04.00/0
Figure 1. Two a-cages of A zeolite and the interconnecting eight-ring window of approximate diameter 4.5 A are shown. The diameter of the a-cage is approximateIy 11.4 A.
comparisonsare worthwhilesince zeolites A and Y have a-cages of similar size and differ only in the topology and the size of window interconnecting the a-cages.
2. Structure of Zeolite NaCaA The structure of Linde 4A reported in the literature was employed in the present work.17 The space group is Fmjc with a = 24.555 A and unit cell composition Na32Cas2A196Si960384. The sodiumand calcium atoms occupy positions close to the center of the six-ring windows. A unit cell of NaCaA consists of eight supercages. Thesesupercages are octahedrallyplaced with respect to each other and have a diameter of 11.4 A (see Figure 1). They are interconnected through an eight-ring window of 4.5-A diameter. All calculations were carried out with eight unit cells of zeolite NaCaA. Therefore, a total of 64 supercages were considered in the simulations. Below, the words supercage and cage are used interchangeably. 3. Intermolecular Potential Functions Guest-zeolite interactions were modeled in terms of shortrange (6-1 2) Lennard-Jones interactions, 6
12
(1) @az(raz)= -Aazf raz + Bazf raz where a = Ar and z = Si, 0, Na, Ca. The potential parameters reported by Kiselev and Pham Quang DuI8 have been employed 0 1993 American Chemical Society
MD Simulations of Ar in Zeolite
The Journal of Physical Chemistry. Vol. 97, No. 51, 1993 13779
TABLE I: Intermolecular Potential Parameters for Guest-Guest and Guest-Zeolite Interactions atom A (lo3W/mol, A6) B (lo6kJ/mol, AI2) guest-guest Ar-Ar guest-zeolite Ar-O Ar-Na Ar-Ca
6.264
9.849
4.701 1.641 18.390
4.173 2.261 27.421
in the present study. The interaction parameters between guest and Si atoms are taken to be zero. This is justified since the Si atom is surrounded by bulky oxygens, thereby preventing the close approach of the guest species, and the fact that Si has a small polarizability. Earlier MD simulationsusing these potential parameters have yielded reasonable values for the diffusion coefficient.ll Many-body induction interactions are expected to contribute less than 10% to the total interaction energy. These interactions are not included due to the enormous computer time required for their computation. However, earlier simulations on xenon in NaY suggest that even without their inclusion results are obtained which are in agreement with experiment. Since the contribution from induction interactions should be smaller for argon than for xenon, we expect the error to be even smaller for argon in NaCaA. In any case, the Kiselev model provides a well-tested effective potential. The potential parameters are listed in Table I. The guest-guest interactions were modeled in terms of (6-12) Lennard-Jones form
+
u(r) = -A/r6 B/r12 where the parameter A = 4sa6and B = 4td2.The values C A ~ - A ~ = 119.8 K and U A ~ - A = ~ 3.41 A were taken from the literat~re.1~ 4. MD Calculations All calculationswere carried out in the microcanonical ensemble at fixed (N,V,E). Cubic periodic boundary conditions were employed. The zeolite atoms were assumed to be rigid and were not included in the MD integration. The integration was carried out using the Verlet scheme. Calculations were carried out at seven different temperatures, 139,189,261,341,433,525,and 585 K, at a loading of c = 1 Ar/cage. A time step of 20 fs was found to yield good conservation of energy and momentum. Properties were calculated from configurationsstored at intervals of 0.10ps. Equilibration was performed over a duration of 100 ps during which velocities were scaled to obtain the desired temperature. The production runs were carried out for 600 ps at all the temperatures.
5. Results and Discussions 5.1. Structure and Dynamics of Sorbates. 5.1.1. Thermodynamic Properties. Guest-guest interaction energy,
(3) and guest-zeolite interaction energy, N Ugh
=
&
dar(raz)
(4)
a=l z= ,Na.Ca
and the total energy are listed in Table 11. We did not find any calorimetric measurements for Ar-NaCaA system. However, the heat of adsorption is 16.7 kJ/mol for methane at room temperature.20
TABLE II: Equilibrium Thermodynamic Properties for Argon in NaCaA at Different Temperatures T (K) CUB) (kJ/mol) (Ugh) (kT/mol) (U)(kJ/mol) a t (kJ/mol) ~~
139 189 261 341 433 525 585
-0.41 -0.36 -0.31 -0.16 -0.15 -0.13 -0.12
-15.55 -14.52 -13.32 -1 1.96 -1 1.87 -11.06 -10.85
-15.96 -14.88 -13.63 -12.12 -12.02 -11.19 -10.97
17.1 1 16.45 15.79 14.95 15.61 15.42 15.83
This may be comparedwith 15.8kJ/mol for the heat of adsorption at 261 K obtained by us. 5.1.2. Structureand Energetics. Radial distribution functions (rdf) at different temperatures between the sorbate Ar and the zeolite atoms 0, Si, Na, and Ca are displayed in Figure 2. At the lowest temperature (139K) it can be seen that the peaks are well-defined, indicating that the sorbate is largely localized in the supercages of NaCaA. Figure 3 shows a plot of Ar-Ar rdf. An intense peak at 4.0 A is observed at 139 K. Around 261 K, even though the features of the rdfs are still well defined, the peak heights have decreased significantly. The well-defined peaks of various guest-host as well as guest-guest rdfs disappear by 341 K. On the basis of the rdfs, it appears that the localized to delocalized transition takes place somewhere between 260 and
340 K. The distribution functions for guest-host interaction energy, f(Ugh), are shown in Figure 4 at different temperatures. The distribution has significant intensity on the lower energy side at 139K. At higher temperatures, a peak with maxima near -15 kJ/mol and a shoulder near -1 1 kJ/mol areobserved. Increasing the temperature from 261 to 341 K changes the most intense peak inf(U,h) from -15 to -1 1 kJ/mol. This suggests that the sorbate delocalization is associated with the shift of the most intense peak in f(Ugh) to -1 1 kJ/mol from -15 kJ/mol. This implies that sorbate (whose interaction energy with the host lattice NaCaA is stronger than -15 kJ/mol, that is, lower than or equal to -15 kJ/mol) is localized. Methane in NaY also exhibited similarbehavior, showing a bimodal guest-host energy distribution curve, and a similar shift in the intensity was observed from low energy to high energy. Further, this shift was found to be associated with the localized to delocalized transition albeit at a somewhat lower temperature of -200 K.21 Figure 5 displays the single sorbate density distribution, p(r), inside the supercage as a function of the distance from the cage center at various temperatures. At 139 and 341 K, the density of Ar is mainly concentrated near the inner surface of the supercage. Above 341 K, the population of Ar begins to increase towards the cage center. It was shown in our earlier work on Xe in NaY that significant population near the cage center alters the mechanism of diffusion from one cage to another.22 We show in section 5.2 that the change in p(r) has a similar effect on cage-to-cage diffusion in zeolite NaCaA. 5.1.3. Dynamical Properties. Figure 6 shows the variation of mean square displacement, ( u 2 ( t ) ) ,as a function of time for various temperatures. Diffusion coefficients at various temperatures obtained from the Einstein relationf9 D = (u2(t))/2dit (6) valid for long times are listed in Table 111. Here di is the dimensionalityof the system,taken as 3 in the present case. These values may be compared with that of methane in NaCaA at 293 K and 6 molecules/cage obtained from PFG-NMR, for which results are available. Karger and Pfeiferz3 estimated the intracrystalline diffusion coefficient to be around 1.7 X le9m2/s and the long-range diffusion coefficient to be 2.2 X lo-' m2/s. The intracrystalline diffusion coefficient is obtained by measuring the mean square displacement over distances which are significantly smaller than the size of microcrystallites. The present
Yashonath et al.
13780 The Journal of Physical Chemistry, Vol. 97, No. 51, 1993 2.5
1 c
-139K
2 .o
0.5
t
189K
*
261K
0
341K
X
585K
- 139K
. I .
0.00
8
L
3.75
2.0
- 139K + 189K
(d)
(b)
*
3.00 1.5
m
0 X
261K 341K 585K
2.25 L
1.0
v
01
1.50
* 261K
0.5
0
X
0
I
I
4
8
r,
8
341K 585K
0.75
12
0
0
4
r , 8,
a
12
Figure 2. Radial distributionfunctionsbetween the guest argon and the zeolite atoms (a) 0, (b) Si, (c) Na, and (d) Ca shown at different temperatures. 4
0.4
- 139K
- 139K + 189K
3
n
*
261K
o
341K
X
585K
0.3
t
261K
*
3L1K
o
585K
Q a )
:2
W
0
1c
0.2
0.1
1 0.0 I
0
4
r, A
8
1:
Figure3. Radial distribution function for gueat-guest(Ar-Ar) interaction at different temperatures for an adsorbate concentration of c = 1 Ar/ cage.
result (0.65 X lo-* m2/s) is for a single crystal and therefore should be compared with the intracrystalline diffusioncoefficient of Karger and Pfeifer. Since CHI and Ar are of roughly the same size, and noting that the present results are for 1 Ar/cage, the agreement with the value reported by Karger and Pfeifer is reasonable. Results obtained here may be compared to the 0.7 X le8 m2/s for D at 393 K in our earlier calculation using the
-15.0
-10.0
-5.0
I
Ugh, kJ/mol Figure 4. Distribution of energy of interaction between the guest and the host, flu@),at different temperatures.
same set of potential parameter^.^^ The small discrepancy between the present results and our earlier results arise from the fact that calculations have been carried out on 64 argon atoms for 600 ps in the present instance as compared to 8 argon atoms for 200 ps in the earlier calculation. The present results therefore are considerably more accurate. Figure 7 displays an Arrhenius plot of the diffusion coefficient. The activation energy, E,, obtained by fitting all the seven points lying in the range 139-585 K is 1.8 kJ/mol. This may be compared with the activation energy for
The Journal of Physical Chemistry, Vol. 97, No. 51, 1993 13781
MD Simulations of Ar in Zeolite
- 18.0
1 5
-139K + 341K
*
1 2
o X
433K 525K 585K
- 18.5
f-
9
-19.0
I
v
a 6
- 19.5
3
-20.0 1
3
5
7
9
( 1 /Tx 1 03, K-’)
Figure 7. Arrhenius plot of diffusion coefficient with temperature.
Figure 5. Single sorbate density distribution, p(r), as a function of the distance r from the cage center inside the supercage at different temperatures.
6
I
-31 -3
’
I
I
0
3
J
6
2400 1800
observed earlier for Xe in Figure 8 shows a plot of In ( d ( t ) ) against In T and In ( uz(t) ) against 1 / T for short and long times, respectively. At short times (t I1.82 ps) the mean square displacement shows a linear dependence on temperature, and at long times an Arrhenius dependence is observed. Recently Kolb and AmraniZ6have observed similar type of behavior for argon diffusing in silicalite by MD simulation. Presentlywe are working on model systems to understand the factors which influence the crossover. Figure 9a displays the power spectra obtained by the Fourier transformation of the velocity autocorrelation function. At low temperature (1 39 K), two distinct peaks are observed: one around 45 cm-I and another at 14 cm-I. These two peaks are separated by a minima near 26 cm-l. Earlier investigations of xenon in NaY zeolite suggest significantcontributionto the high-frequency peak from the velocity autocorrelation function for the velocity component perpendicular to the inner surface of the cage. It was also found that motion parallel to the inner surface of the cage gives rise to significant contribution to the low-frequency peak. The parallel and perpendicular componentscan be obtained from the expressions
a,(?)
= b(t).P(t)
(7)
1200
where b(t) is the velocity and P(t) is the unit vector connecting thesorbate to thecenter of the cage. Two velocity autocorrelation functions, CII and cl, may be defined corresponding to these components:
600
-
n 0
100
200
300
400
t , p s --* Figure 6. Plot of the mean square displacement against time for different temperatures. Inset shows a plot In (u2(r)) vs In r.
TABLE 111: Diffusion Coefficients at Different Temperatures Obtained from Emtein’s Relationship (Eq 6) Using the Long-Time Slope of the Mean Square Displacement after the Crossover T (K) 139 189 261 341 433 525 585 D X 108 (m2/s) 0.30 0.50 0.65 0.83 0.89 0.91 1.08 xenon in NaY of 4.2 kJ/mol. The reason for the difference in the activation energy will be discussed in a later section. Inset of Figure 6 shows a log-log plot of mean square displacement against time. The curves clearly show a crossover from ballistic motion, characterized by G ( t ) p: t z , to diffusive motion, characterized by u2(t) 0: t . Such a crossover has been
Figure 9b displays the power spectra for the parallel and perpendicular components at two different temperatures, 139 and 585 K. The sorbate is highly localized at 139 K, whereas at 585 K it is highly delocalized. The first peak in the power spectrum at 139 K extends up to 26 cm-I, whereas the second peak extends from 26 to 80 cm-1. Table IV lists the intensity of the first (Z(1)) and second peaks (Z(2)) for the parallel and perpendicularvelocity components, at two different temperatures. At 139 K, the parallel and perpendicular componentscontribute in the same proportion to the first and second peaks. This is in contrast to the Xe in NaY system, where it was found that c l made a larger contribution to the high-frequency peak at low temperat~res.2~ At 585 K, the ratio I(2;cl)/Z(2;cll) 0.77 as
Yashonath et al.
13782 The Journal of Physical Chemistry, Vol. 97, No. 51, 1993 1.0
6.0
18.0
I
(a)
-1 3 9 K 4.5
x
t
261K
*
505K
vi
c
0
3.0
4J
-c 1.5
0.0 0
80
40
u,
120
1( 3
cn7-I
2.25
- Parallel
(b)
n (b)
1.50
0.0
lo
(139K)
t
Perpendicular ( 1 3 9 K )
0
Parallel ( 5 8 5 K )
*
Perpendicular (585)o
I
I
I
40
80
u,
120
1 60
cm-'
Figwe9. (a) Power spectraat differenttemperaturesfor argon in NaCaA
obtained by Fourier transformation of the velocity autocorrelation function. (b) The power spectra for the velocity component parallel and perpendicular to the cage surface at different temperatures.
Figure8. (a) Logarithmof mean squaredisplacementbefore theobserved crossover, (see inset of Figure 6). In ( u * ( t ) ) , is plotted against the logarithm of temperature, In T. The filled circles show the short-time mean square displacement before the crossover, and the open circles show the long-time mean square displacement after the crossover. The best least-squares line for the short-time mean square displacements is alsoshown. (b) Logarithmof meansquaredisplacement is plottedagainst reciprocal temperature. The best least-squares line for the long-time mean square displacement is also displayed. compared toZ(1;cL)/Z(lq) = 0.72, which means that the second peak (the high-frequency peak) has a relatively larger contribution from the perpendicular velocity component. Figure 9a clearly shows that with increase in temperature the low-frequency peak gains in intensity while simultaneously the high-frequency peak loses intensity. 5.2. Short-Time Behavior: Cage-to-Cage Diffusion. 5.2.1. Arrheniw Behavior. Earlier investigations have shown that cageto-cage diffusion in zeolites is an important subprocess which influences strongly the overall diffusion rate. Here we present various properties related to cage-to-cage diffusion of argon in NaCaA. Figure 10 shows the variation of the total energy (U) = (Ugh)+ ( Ugg) with d, the perpendicular distance between the sorbate and the plane of the eight-ring window. We know that the sorbate has to pass through the eight-ring window to migrate
TABLE Iv: Intensity of the First (I(1))and Second Peaks (I(2))in the Power Spectra Obtained from Fourier Transformation of the Parallel and Perpendicular Velocity Components at Low and High Temperatures 139 K 585 K CII
CI
I(1) 20.7 16.4
I(2) 32.0 23.5
~~
I(1) 34.9 25.2
I(2) 28.3 21.7
to another cage. Here, d is less than zero when sorbate is in the parent cage (that is, the cage from which the sorbate is diffusing to the neighboring cage), d = 0 when the sorbate is in the plane of the eight-ring window, and d is greater than zero when the sorbate is in the daughter cage (the cage to which the sorbate has diffused). Hence, d represents a typical reaction coordinate. It is clear from Figure 10that the barrier height for a sorbate passing through theeight-ring window is negative. Thismay becontrasted with the positive barrier height observed for diffusion across the 12-ringwindow in zeolite Y .22 The barrier height shown in Figure 10 is based on short time or local variation of potential energy. For a detailed discussion of some of these points, see ref 22. Inset of Figure 10 shows (Ugh)as a function of time. It is clear from the inset that the actual crossover from one cage to another takes only about 2 ps. Finally, we note that the observed behavior of (U)vs d arises from a change in the traversed trajectory with the temperature of the system.22 The behavior of the rate of cage-to-cage diffusion with
MD Simulations of Ar in Zeolite
The Journal of Physical Chemistry, Vol. 97,No. 51, 1993 13783
I
24.68
-10
high T
t -14 -16
-3 -2
.z4.98F0
2
1/T x lo3
--I
1 low T
s 24.88
1/T
xld +
I
I
/
I
/ 3
-16
, 1
-1 8 -6
24.68
1.0
139K
-3
0
3
K-'
I
I
1
5.0
6.0
7.0
4
high temperatures.
TABLE V Rate of Cageto-Cage Diffusions (k) and Fractional Contribution from Centralized (rp" < 3 A) and Surface-Mediated (P$'> 3 A) Modes at Different Temperatures before Correcting for Recrossings' fractional contribution from kc X c.d. mode s-m. mode T (K) (/sorbate/s) ($c 3 A) c p 3 4 0.04 0.12 0.29 0.40 0.56 0.60 0.68
b.0
1/T xl$,
d,a
7.08 (2720) 5.90 (2269) 4.60 (1760) 4.48 (1720) 4.51 (1644) 4.51 (1734) 4.58 (1760)
3.0
Figure 11. Arrhenius plot of kcwith temperature. Expanded plots below and above 500 K are shown separately in the inset. Two additional runs at 718 and 801 K were carried out above 500 K to confirm the trend at
F i p e 10. Temperature dependence of the variation of the total energy interaction, U,as a function of the distanceof thesorbatefrom the window plane, d. Inset shows the plot of Uph vs timejust before, during, and after the crossover in the interval -4 to +4ps. The time, t = 0, corresponds to the time when the sorbate is on the window plane.
139 189 26 1 341 433 525 585
2.0
0.96 0.88 0.7 1 0.60 0.44 0.40 0.32
a Total number of cage-to-cage diffusions during the 600 ps run is listed in parentheses. temperature would be interesting since for the present system a negative barrier height has been obtained. In Table V we list k,, the rate of cage-to-cage diffusion for several temperatures. It is seen that the rate decreases with temperature up to 500 K and then increases. Figure 11 shows an Arrhenius plot (In k, vs 1/ r ) for the whole range of temperatures. Inset shows the same plot below and above 500 K, respectively. The activation energy from the slope is -0.95 kJ/mol (SO0 K), respectively. This transition from a negative to a positive temperature dependence of k, is surprising. The factors responsible for the peculiar observed dependence of k, on temperature will be discussed in the next section. The necessity of correcting for recrossings by means of dynamical correction factor has been pointed out by Chandler and c o - w o r k e r ~ . ~According ~ * ~ ~ to this formalism, the dynamical correction factor@t) corrects for the recrossings whenever rcorr,the time for the recrossing, is significantly less than the cage residence time, 7,. This formalism has been proposed for canonical ensemble simulation^.^^^^^ Since we have carried out calculations in the microcanonicalensemble, we have to use an alternative criterion to eliminate recrossings. As can be seen from the inset of Figure 10, diffusion across the eight-ring window occurs within 2 ps. Hence, the time taken for
TABLE VI: Rate of Cage-to-Cage Crossovers (ke)and Fractional Contribution from Centralized (rp"< 3 A) aod Surface-Mediated (r"' > 3 A) Modes at Different Temperatures after elimination of AU Recrossing Taking Place within 1 ps of Each Other fractional contribution from c.d. mode s-m. mode k, X T (K) (/sorbate/s) (rf:c3A) ($>>A) 139 189 26 1 34 1 433 525 585
6.17 5.04 4.15 4.09 4.00 4.21 4.34
0.04 0.14 0.25 0.43 0.59 0.63 0.66
0.96 0.86 0.75 0.57 0.41 0.37 0.34
a sorbate on the dividing surface to enter into the daughter cage is about 1 ps. We have therefore calculated k, after eliminating all crossings which occur within 1 ps of each other. Table VI lists the crossings at various temperatures after incorporating this correction. The activation energy obtained from the Arrhenius plot of the data is -0.85 and 1.44 kJ/mol below and above 500 K, respectively. This is not very different from the values obtained for uncorrected k,. 5.2.2. Characteristics of Cage-to-Cage Diffusion. In Figure 12a we show the temperature dependence of the distribution of interaction energy of the sorbate at the window, ew = Ugh at d = 0. It isseenthatthemaximaofthedistributionshiftincreasingly towards higher energy. In Figure 12b we show the distribution of the distance of the sorbate from the window center at d = 0. It is seen that the maxima innraw)shift to higher distances, Le., away from the window center. This may be compared withflew) andflr,,) distributions in Xe-NaY system, where no perceptible shift was observed. The maximum in the case of zeolite A lies at around 0.5 A, which is close to the window center, whereas in the case of zeolite Y the maximum was at 1.3 A. This may be partly attributed to the relatively narrow (4.5-A diameter) size of the 8-ring window, as compared to Nay, where 12-ring window has a diameter of 8 A. The width of the distribution is also larger in Xe-NaY system, which again reflects the larger diameter of the window in Y zeolite. The observed shift of the maxima in flew)andflr,,) may be attributed to the shallow minimum in the potential energy curve when plotted as function of raw. Figure 12c shows the temperature dependence of the variation of the distance of the adsorbate from the parent cage center, r, and the
Yashonath et al.
The Journal of Physical Chemistry, Vol. 97, No. 51, 1993
13784
- ‘39K 0.30 n
+
261K
*
341K
o
585K
(a 1
I
(a)
-
rm P
38
-10
0
E \-12
a;
7 Y
W ‘c
J-
0.15
-14
0.0 -
.18.0
-14.0
-10.0
- 1 61 -6
-2.0
-6.0
ew,kJ/mol
I 3
6
d.a
- 101
5
I
1
-139K
(b) 4
t
261K
*
341K
Q
585K
3
A 3
1 0
I -3
-0
- 1 2 t
7
I
u
LO W
+
2
-
1
1
-3
I 0 d,
-
-1 3 9 K
(C) 1
2-
+
261K
*
341K 5 a 5 ~
9-
6-
~
3-
-
0
6
A
Figure 13. Temperaturedependence of the variation of Uwith d for (a) the centralized diffusion mode, $‘ < 3 A, and (b) the surface-tomediated mode, rp” > 3 A.
0
1
I 3
3
6
9
p ‘
,s
12
15
Figure 12. (a) Distribution function, flew), is shown as function of temperature. Here e, is the guest-host potential energy; e, = U, at d ~Oduringcage-to-cagediffusion.(b) Distributionofthedistancebetween argon migrating from one cage to another and the window center, raw, at d = 0 when the sorbate is in the plane of the eight-ring window as a function of temperature. (c) Temperature dependence of the variation ofrpversusrdduringcage-to-cagediffusion. Herer,andrdare respectively the distance between the sorbate and center of the parent and daughter
cages. daughter cage center, rd during cage-to-cage diffusion. It is seen that at high temperatures, sorbates approach the cage centers more closely. A similar behavior was observed for Xe in Nay. 5.2.3. Centralized and Surface-Medimted Diffusion Modes. Results on Xe in NaY have shown that cage to cage diffusion
proceeds, broadly speaking, via two different modes. They are surface-mediated (s-m.) and centralized diffusion (c.d.) modes. The two modes are distinguished by r: > 3 A and r: < 3 A, respectively. Here r: is the minimum distance between the sorbate and the center of the parent cage during the interval starting 4 ps before the sorbate has diffused across to the neighboringcage and until the sorbate has actually diffused into the neighboringcage, (-4,Ops). In the s-m. diffusion ($> 3 A), the sorbate keeps close to the inner surface of the a-cage, prior to diffusion to the neighboring cage. In the c.d. mechanism (r: < 3 A), the sorbate is close to the center of the cage prior to crossing over to the neighboring cage. We have analyzed the trajectories during (-4,+4ps) of cage-to-cagediffusion events to obtain details about the mechanism of diffusionof Ar in NaCaA. Figure 13 shows the variation of the total energy (U)with d for $‘ < 3 A and r: > 3 A. The barrier height for centralized diffusion is found to be negative. Negative barrier height was also found for this mode of diffusion for Xe in N a y . Surfacemediated mode is associated with negative barrier height, as is seeninFigure 13. Thisresult isatvariancewithwhatwasobserved for Xe in NaY where positive barrier heights were obtained for s-m. mode. As we shall see, this unexpected result is responsible for most of the observed differences in the behavior of kc and other properties between Xe-NaY and Ar-NaCaA systems.Table VI1 lists the barrier heights for c.d. and s-m. modes. These have been obtained from Figure 13. These results are consistent with the lower activation energy of 1.8 kJ/mol obtained from the Arrheniusplot of diffusion coefficient for Ar-NaCaA as compared to 4.2 kJ/mol for X e N a Y system. Single sorbate density profile, p ( r ) , in the a-cage, where r is the distance from the cage center, shows a peak around 4.0 A (see
MD Simulations of Ar in Zeolite
The Journal of Physical Chemistry. Vol. 97, No. 51, 1993 13785
TABLE Vn: Barrier Heights for Cage-to-Cage and for Centdzed 8nd Surface-Mediated Diffusions Obtained from a Plot of Potential Energy Against d
T (K)
U. (kJ/mol)
139 189 26 1 341 43 3 525 585
-0.9 1 -1.25 -1.82 -2.05 -2.29 -2.45 -2.57
u.($ < 3 A)
u. (P > 3 A)
(kJ/mol)
d/mol)
-1.6 -2.05 -2.11 -2.28 -2.45 -2.62 -2.91
-0.92 -1.08 -1.42 -1.54 -1.65 -1.82 -1.94
Figure 5). From Figure 5,we note that the population is negligible near the center at low temperatures. Hence, cage-to-cage migration takes place mainly via s-m. mode at low temperatures. Cohen de Lara and Kahn have obtained from neutron diffraction a similar profile for ~ ( r )They . ~ also ~ find at higher temperatures p ( r ) gains intensity near the cage center. With increase in temperature, the population near the cage center becomes appreciable, leading to a significant contribution from the c.d. mode to the overall rate of cage-to-cagediffusion. These changes in the mode of diffusion with temperature are similar to those observed for Xe in Nay. The difference, however, is the fact that the barrier height is negativefor s-m. mode for Ar in NaCaA. This is responsible for the negative activation energy at low temperatures obtained from the plot of In k, vs 1/T (Figure 11). Since the barrier for the c.d. mode is also negative, one expects the activation energy in the plot of In k, vs 1/ T to remain negative. The change from negative to positive activation energy around 500 K in the Arrhenius plot shown in Figure 11 is therefore surprising. Figure 14a displays schematically the variation of the potential energy for the s-m. and the c.d. modes of diffusion. From this figure it is clear that the sorbate can diffuse from one cage to another without any activation energy via the s-m. mode. This is, however, not the case for the c.d. mode since the sorbate has to be activated to higher energy from the region near the inner surface of the cage to that of the cage center which is higher in potential energy. Figure 14b shows a schematic plot of the variation of potential energy for c.d. mode. The activation energy required for a sorbate via c.d. mode is indicated by a vertical arrow in the figure. If this interpretation is correct, the activation energy should be negative when the predominant mechanism for cage-to-cage diffusion is s-m. mode and positive when the predominant mechanism for migration from one cage to another is the c.d. mode. Table V, which lists the values for the rate of cage-to-cage diffusion, shows that the fraction of cage-to-cage diffusions via the s-m. mode is higher below 500 K. By 525 K, the c.d. mode has become the predominant mode of cage-to-cage diffusion, contributing about 60%of all cage-to-cage crossovers. Thus, the results of the analysispresented in Table V are consistent with the interpretation we have proposed. It is worthwhile to compare the potential energy surface for Xe-NaY and Ar-NaCaA systems. Variation of the potential energy for Xe-NaY system is shown schematically in Figure 14c. In the case of Xe-Nay, the potential energy in the vicinity of the inner surface of the a-cage, e,, is lower than the energy at the window, e,, which is lower than the energy at the cage center, e,: e, < e, < e,. In the case of Ar-NaCaA, e, < e, < e,. The principal differences in the behavior of k, with temperature between Xe-NaY and Ar-NaCaA may be largely attributed to this difference in the underlying potential energy surface. Recent work suggests that apart from geometrical factors the sorbate-zeolite interactions play a significant role in cage-tocage diffusion.30 In order to investigate the importance of the sorbate-zeolite interactions in cage-to-cage diffusion of argon in NaCaA, we have carried out a series of MD runs at several temperatures in the absence of sorbatezeolite dispersion interactions. More specifically, the parameter A, was taken as zero in eq 1 for a = Ar and z = 0,Na, and Ca. The results for the
window
cage
window A i - NaCa A
I
I
I
I
window
cage
winQw
window
caqe
w indvu
1
I
I
c.d.
I
I
cage
window
window
caw
window
I
I
I
1
1
window
window
+
cage
I
window Figure 14. (a) Schematicillustration of the variation of potential energy for centralized (P < 3 A) and surface-mediated (rp”> 3 A) diffusion modes for Ar in NkaA. Theverticalarrow indicates the excitationthat
is essential before the sorbate can diffuse across the eight-ring window via centralizeddiffusion mode. (b) The variation of the potentialenergy for a sorbate near the inner surface of the a-cage and diffusing via the centralized diffusion mode is shown. (c) Variation of potential energy for the centralized and surface-mediateddiffusion modes for xenon in NaY is shown.
TABLE VIIk Rate of Cage-to-Cage Crossovers (4) 8nd Fractional Contribution from Centralized (1,” < 3 A) a d Surface-Mediated (T > 3 A) Mode at Different Temperatures in tbe Absence of Dispersion Interactions (A.z = 0 in Eq 1) between the Guest and Zeolite Atoms’ fractional contribution from k, X 1 P 0 c.d. mode s-m. mode T (K) (/sorbate/s) ((33‘4) ( p 3 4 3 A) Modes at Different Temperatures fractional contribution from ~
kv X T (K) 139 189 26 1 34 1 433 525 585
2
~~
(/sorbate/s)
c.d. mode (r: < 3 A)
1.75 2.11 2.78 2.94 3.1 1 3.33 3.58
0.18 0.33 0.47 0.62 0.76 0.82 0.87
s-m. mode (C’3A) 0.82 0.67 0.53 0.38 0.24 0.18 0.13
I
I
I
I
23.98
-. ../ Y L
-
23.681
23.381
This leads to an enhancement in the overall diffusion coefficient due to increase in the value of di in eq 6. The value of di lies between 2 and 3 in the presence of sorbate-zeolite dispersion interactions. The enhancement thus obtained can at most be 1.5 times thevalue when the dispersion interactions areabsent. Thus, this alone cannot explain the observed increase. A large contribution to the increase in the rate of cage-to-cage diffusion arises from the fact that the velocity component perpendicular to the inner surface of the cage cannot lead to cage-to-cage crossover. It is the component of velocity that is parallel to the inner surface of the cage that makes significant contribution to the cage-to-cage diffusion, and this is the major component of the velocity when the dispersion interaction is switched on. This is evident from Figure 9 and data presented in Table IV. Yet another important factor responsible for increase in k, is the negative barrier for diffusion across the eight-ring window in the presence of the dispersion interaction. This leads to reduction in E, in the Arrhenius expression, leading to an increase in k,. We have seen that the diffusion coefficient exhibits an Arrhenius dependenceon temperature over the whole range of temperatures investigated (see Figure 7). This may be contrasted with the Arrhenius plot of k, in Figure 11. The differencein their behavior as a function of temperature arises from the fact that D is dependent on the long-time behavior whereas k, is dependent on the short-time behavior. Here, the short time refers to times shorter than the crossover time of approximately 1.8 ps obtained from the plot of In ( u * ( t ) )against In t , depicted in Figure 6. Long time obviously refers to times longer than the crossover time. The system follows the Arrhenius law at long times. This is clear from Figure 8b. Rate of cage-to-cage diffusion, on the other hand, is a property that is dependent on the details in the vicinity of the eight-ring window such as the potential energy surface, etc. In order to confirm this, we define another property, k,, which we refer to as rate of cage visits. This is calculated by counting only those cage-to-cage diffusions, i to j, which are not preceded or followed by diffusion from j to i. Table IX shows therateofcagevisits (k,) at different temperatures. This property (k,) should be dependent on long-time behavior since the sorbate has to diffuse through the interior of the cage before it can diffuse to a neighboring cage. Figure 15 shows a plot of In k, against 1/T.All the calculated points fall on a straight line. This indicates that non-Arrhenius behavior is indeed restricted to properties such as k,, which are dependent on short-time behavior. 6. Conclusions
Results on argon in zeolite NaCaA indicate that the argon atoms are localized below 300 K. At higher temperatures, argon atoms are delocalized and diffuse from one a-cage to another. Cage-to-cagediffusions may be said to proceed via two alternative mechanisms-the surface-mediated and the centralized diffusion modes. These findings are similar to those reported for xenon in N a y . Faujasites or Y zeolites have a-cages of similar size, but the cages are interconnected tetrahedrally and the windows
1
’.,
I
I
I
3
5
7
9
1 / T x l 03,K-’
Figure IS. Arrhenius plot of the rate of cage visits k, with temperature.
are larger in diameter. It is shown that the smaller eight-ring window in zeolite A leads to several unexpected results. First, it is found that the rate of cage-to-cage diffusions, k,, is significantly higher for argon in zeolite NaCaA even though the window in A is significantly narrower. Further, it is found that the barrier for both surface-mediated and centralized diffusion modes is negative. The activation energy obtained from the Arrhenius plot of k, is found to be negative below 500 K. Above this temperature the activation energy is positive. This change from negative to positive activation energy around 500 K is attributed to the change in the principal mode of cage-to-cage diffusion from surface-mediated to the centralizeddiffusion mode. This also justifies the analysisof cage-to-cagediffusions in terms of these two principalmodes of diffusion from one cage to another. The non-Arrhenius behavior of the k, is attributed to its dependence on the short-time (52 ps) behavior. An Arrhenius plot of the diffusion coefficient yields a value of 1.8 kJ/mol for the activation energy which may be compared with 4.2 kl/mol for Xe in Nay. It is shown that the expected Arrhenius behavior is recovered when properties dependent on long-time behavior, such as the diffusion coefficient and the rate of cage visits, are investigated. It is found that the cage-to-cage diffusions are either zero or negligible in the absence of sorbate-zeolitedispersion interactions. The sorbate-zeolite dispersion interactions enhance the rate of cage-to-cage diffusions by a few orders of magnitude, This enhancement is attributed to three factors. First, decrease in the dimensionality for the sorbate motion from 3 in the absence of sorbate-zeolite dispersion interactions to some value between 2 and 3 leads to an increase in D as well as in k,. Second, the predominance of surface-mediated mode in the presence of sorbate-zeolite dispersion interactions and, last, the negative barrier height for diffusion across the eight-ring window lead to an increase in k, as well as in k,. These results and those on xenon in zeolite Y indicate that the cage residence times for argon in NaCaA are smaller than those for xenon in Nay. It should be possible to verify these results experimentally.
Acknowledgment. One of the authors (P.S.)gratefully acknowledges the Council of Scientific and Industrial Research (C.S.I.R.), New Delhi, India, for the award of a fellowship. References and Notes (1) June, R. L.; Bell, T. A.;Thwdorou, D. N. J. Phys. Chem. 1990,91, 1508.
(2) Demontis. P.;Suffriti, G. B.: Fois. E. S.:Quarticri.S. J . Phys. Chcm. 1992, 96, 1482. (3) Pickett, S. D.; Nowak, A. K.; Thomas, J. M.; Peterson, B.K.; Swift, J. F. P.: Chcetham. A. K.: den Ouden. C. J. J.: Smit. B.: Post, M. F. M. J. Phys. &em. 1990, 94, 1233. (4) Zikanova, A.; Bulow, M.; Schlodder, H. Zeolites 1987. 7, 11s.
MD Simulations of Ar in Zeolite ( 5 ) Titiloye, J. 0.; Parker, S.C.; Stone, F.S.;Catlow, C. R. A. J. Phys. Chem. 1991.95, 4038. (6) Nowak, A. K.; den Ouden, C. J. J.; Pickett, S.D.; Smit, B. J. Phys. Chem. 1991,95,848. (7) Rowlinson, J. S.;Woods,G. B. Physfcu A 1990, 164, 117. (8) Woods,G. B.; Rowlinson, J. S.J. Chem. Soe., Faraday Trans. 2 1989, 85, 765. (9) Schrimpf, G.; Schlenkrich, M.; brickmann, J.; Boff, P. J. Phys. Chem. 1992,96, 7404. (IO) Demontis, P.; Yashonath, S.;Klein, M. L. J . Phys. Chem. 1989.93, 5016. (11) Yashonath, S.Chem. Phys. Lett. 1991, 177, 54. (12) Yashonath. S . J. Phys. Chem. 1991, 95, 5877. (13) Demonth, P.; Suffriti, G. B.; Quartieri. S.;Fois, E. S.;Gamba, A. J . Phys. Chem. 1988,92,867. (14) Kono, H.; Takasaka, A. J . Phys. Chem. 1987.91, 4044. (15) van Tassel. P. R.; Davis, T. H.; McCormick, A. V. Mol. Phys. 1991, 73, 1107. (16) Chmelka, B.F.;Raftary, D.; McCormick, A. V.;de Menorval, L. C.; Levine, R. D.; Pines. A. Phys. Reo. Lett. 1991, 66, 580.
The Journal of Physical Chemistry, Vol. 97, No. 51, 1993 13787 (17) Pluth, J. J.; Smith, J. V. J. Am. Chem. Soc. 1980, 102,4704. (1 8) Kiselev, A. V.; Du,P. Q. J . Chem.Soc., Furuduy Trans. 2 1981,77,1, (19) Allen, M. P.; Tildesley, D. J. Computer Simulation of Liquids; Clarendon Press: Oxford, U.K., 1987. (20) Cohen De Lara, E.; Kahn, R.; Goulay, A. M. J. Chem. Phys. 1989, 90, 7482. (21) Yashonath.S.;Thomas,J. M.;Nowak,A.K.;Cheetham,A.K.Nurure 1988, 331, 601. (22) Yashonath, S.;Santikary, P. J. Phys. Chem. 1993, 97, 3849. (23) Karger, J.; Pfeifer, H. Zeolites 1987, 7, 90. (24) Santikary, P.; Yashonath, S . J. Chem. Soc., Faruduy Trans. 1992, 88, 1063. (25) Santikary, P.; Yashonath, S.;Ananthakrishna, G. J . Phys. Chem. 1992,96, 10469. (26) Amrani, E. S.;Kolb, M. J. Chem. Phys. 1993, 98, 1509. (27) Chandler, D. J . Chem. Phys. 1978,68, 2959. (28) Montgomery, J. A.; Chandler, D.; Berne, B. J. J. Chem. Phys. 1979, 70, 4056. (29) Cohen De Lara, E.; Kahn, R. J . Phys. 1981,42, 1029. (30) Yashonath, S.; Santikary, P. Mol. Phys. 1993, 78, 1.
