J. Phys. Chem. 1995, 99, 10932-10941
10932
Computer Simulation of the Structure, Energetics, and Diffusion Properties of p-Xylene in Zeolite Na-Y G. Schrimpf, B. Tavitian, and D. Espinat* Department of Applied Physico-Chemistry and Analysis, Institut Franqais du Pitrole, 1 et 4 Avenue de Bois Priau, BP311, 92506 Rueil-Malmaison Cedex, France Received: January 27, 1995; In Final Form: May 2, 1995@
The properties of p-xylene adsorbed in zeolite Na-Y are studied by molecular simulations. Constrained reaction coordinate minimizations are used to calculate the minimum-energy path for the diffusion. The predicted adsorption site inside the supercage is found to be in good agreement with neutron diffraction results. Extensive molecular dynamics simulations are performed in order to study the behavior of p-xylene as a function of loading at temperatures above 500 K. As important information for the interpretation of experimental data, the distribution of the molecules over the supercages is obtained. Due to the very long simulation runs the diffusion coefficients of p-xylene in the model system are determined with high accuracy. The diffusion coefficients obtained for low sorbate concentration are close to the largest experimental values, which scatter over more than 2 orders of magnitude. It is found that the diffusivity decreases drastically with increased sorbate loading. Furthermore, the microscopic details of the diffusion process are analyzed by the investigation of site and cage residence times and by orientation correlation functions. The simulations reveal that site-to-site migration inside a supercage is often associated with inversion of the orientation of the aromatic plane.
1. Introduction
The adsorption of aromatics in zeolites is of great current interest in the context of separation and catalysis processes such as alkylation, transalkylation, and isomerization of m-xylene into p-xylene. For instance, the separation of p-xylene from C8 aromatics (a mixture of the xylene isomers and ethylbenzene) is performed on the industrial scale using selective adsorption on zeolitic molecular sieves.' The separation characteristics of a zeolitic sorbent depends on its cationic form, the degree of the cation exchange, and the coverage with aromatics. Synthetic faujasite-type zeolites in their sodium form are commonly used and are converted into particularly selective sorbents by cation exchange, thereby improving their pardmeta selectivity.* A large number of microscopic and macroscopic measurements, such as crystal log rap hi^^-^ and i ~ ~ f r a r estudies d ~ ? ~ and thermodynamiclo,' and diffusion coefficient I 2 - l 5 measurements, have been carried out on aromatics adsorbed in faujasitetype zeolites in order to understand their sorption and diffusion properties. The need to gain a detailed insight into the behavior of zeolite-sorbate systems on the molecular scale has also inspired many computational chemistry studies. A recent overview of the works in this field and the different techniques used has been given by Catlow.I6 Among them, the molecular dynamics (MD) technique gives access to time-dependent properties and has been already often used to study the diffusional behavior of sorbates. Only little work, however, has been undertaken until now to investigate the static and dynamic properties of aromatic compounds in faujasite-type zeolites by means of MD simulations."~'~ This is mostly due to the fact that the low sorbate mobility of these systems requires very long simulation times in order to achieve proper sampling of the configurational space. In this work, we present the results of molecular simulation studies of p-xylene adsorbed in zeolite Na-Y. In order to
* To whom correspondence @
should be addressed. Abstract published in Advance ACS Abstracts, June 15, 1995.
0022-3654/95/2099-10932$09.00/0
systematically localize the adsorption sites and to obtain the minimum-energy path for the diffusion of the guest molecules, we use constrained reaction coordinate minimizations. MD simulations are performed to study the sorbate behavior at different concentrations and temperatures. The simulations are performed at temperatures above 500 K to enable the accurate investigation of the diffusion process on the time scale of the simulations. Static properties determined include the spatial distribution of the sorbate molecules inside the supercages as well as their distribution over the supercages. The dynamic properties studied include the self-diffusivity of the sorbate molecules, residence times in characteristic regions of the zeolite framework, and the rotational decorrelation of the sorbate molecules. 2. The Model
The framework of zeolite Na-Y, whose structure is depicted in Figure 1, consists of cuboctahedral sodalite cages linked together in a tetrahedral arrangement by six-membered rings of oxygen atoms to form large cavities, called supercages. The supercages are interconnected by windows that are formed by rings consisting of 12 oxygen atoms and 12 tetrahedral S i l d atoms. The letter T will be used to indicate tetrahedrally coordinated atoms. One unit cell of zeolite Na-Y is used to construct the simulation box. As determined by Fitch and cow o r k e r ~ the , ~ unit cell is cubic with space group Fdgm and a lattice parameter a = 24.85 A. The unit cell contains eight supercages. In the model system, all extraframework cation sites SI and SII are fully occupied by sodium ions. This corresponds to a SUA1 ratio of 3/1 and the composition T1920384Nus. In order to reduce the computing time, we make the approximation that the zeolite atoms are fixed at their crystallographic positions. Comparative MD studies with a rigid and a flexible framework modelI9 have recently been performed for xenon20 and p-xylene2' in zeolite Na-Y. They have shown that the diffusional properties of the sorbates, in which we are 0 1995 American Chemical Society
p-Xylene in Zeolite Na-Y
J. Phys. Chem., Vol. 99,No. 27, 1995 10933 TABLE 1: Charges (in Elementary Charge Units) Assigned to the Force Centers forcecenter charge
T
0
Na
C,
H
CH3
1'2 -0'7
o'8
-0.146 (C,-H) -0.078 (C,-CH3)
0.153
0.064
TABLE 2: Lennard-Jones Parameters for the Xylene-Zeolite and Xylene-Xylene Interactions interaction pair
c,-0 C,-Na H-0 H-Na CH3-O CH3-Na
Figure 1. Sketch of a supercage of zeolite Na-Y with examples for the sodium positions. The SI sites are located at the centers of the hexagonal prism, and the positions of the SI]sites are in front of the 6-T-rings inside the supercages. The T atoms are located at the vertices, and the oxygen atoms occupy positions close to the centers of the lines.
most interested here, are not significantly influenced by the framework vibrations, In our model, p-xylene is assumed to be a rigid, planar molecule, consisting of carbon atoms, hydrogen atoms, and methyl groups. Carbon-carbon, carbon-hydrogen, and carbonmethyl bond lengths are fixed at 1.40, 1.10, and 1.48 A. The carbon-carbon-carbon, carbon-carbon-hydrogen, and carboncarbon-methyl angles are set to 120". The zeolite and the sorbates are assumed to interact via a pairwise-additive potential between atoms of the adsorbed aromatic molecules and atoms of the zeolite. The site-site interactions are modeled with a Lennard-Jones plus Coulomb potential,
where i and j indicate atoms of the sorbate molecule and the zeolite, respectively, and rij is the distance between these atoms. 6, and a , are the Lennard-Jones parameters, and qi and q, are the partial charges of the sites. The interaction parameters are taken from the literature. For the interactions of the methyl group we adopt the parameters used in a simulation study of methane and butane in silicalite?2 whereas the partial charges of the zeolite atoms and the parameters for the Lennard-Jones interactions of carbon and hydrogen are taken from a work of Demontis et al.,I7who studied the behavior of benzene in zeolite Na-Y. The parameters for the methyl-sodium interaction are derived by applying the Lorentz-Berthelot combination rules.23 The Lennard-Jones interactions between the guests and the T atoms are neglected. This approximation is justified by the fact that the SUA1 atoms are well-shielded from the sorbates by the tetrahedrally coordinated oxygens. For the sorbate-sorbate interactions we again use a LennardJones plus Coulomb potential. Interactions of the methyl group are taken from Goodbody et al.22 The carbon and hydrogen interactions are modeled according to Yashonath et al.24 The force field used in ref 24 gives reasonable predictions for the properties of liquid and solid benzene. The charge distribution for p-xylene is obtained by AM1 calculations. We find that the partial charges resulting from AM1 calculations for benzene are by a factor of 1.25 larger than the ones used in the benzenebenzene interaction model of Yashonath et al.24 In order to be consistent with this model for the aromatic interactions, the charges obtained from the AM1 calculations of p-xylene are scaled by the factor of 1.25.
mol)
a@)
interaction pair
1.069 0.175 0.691 0.132 1.108 0.308
2.998 3.328 2.589 2.765 3.214 3.460
C,-C, C,-H C,-CH3 H-H H-CH3 CH3-CH3
6
(u/mol)
a(&
0.3980 0.1458 0.6997 0.0534 0.2563 1.2300
3.466 3.188 3.598 2.910 3.320 3.730
The Lennard-Jones interactions are truncated at a distance rcut,and the potentials are modified to shifedforce potentials Vsf such that the forces and the functions themselves become zero at the cutoff distance:
for rI/2 rcut It has been recently shown by Dufner et al.25 that this potential truncation function can also be used for the calculation of the Coulombic interactions in the present system. They found that this technique reproduces the electrostatic field of faujasites with a good approximation as compared to the computationally time consuming Ewald summation technique.26 A cutoff radius of 10 8, is used for the Lennard-Jones and the Coulomb interactions. The partial charges are given in Table 1, the LennardJones parameters are summarized in Table 2.
3. Molecular Dynamics Simulations The behavior of p-xylene was investigated for 8, 16, and 24 molecules in the simulated system, corresponding to mean loadings of 1, 2, and 3 molecules per supercage ( molecules/ sc), respectively. The simulations were performed in a hightemperature region between 500 and IO00 K, where the mobility of p-xylene is sufficiently large in order to be investigated on the time scale accessible to molecular dynamics simulations. Even in this temperature range, runs of up to 100 ns real time turned out to be necessary to examine the long-range diffusion. For a loading of 3 molecules/sc we investigated the system additionally at 1400 K in order to check whether the diffusion coefficients follow a temperature dependence according to the Anhenius law. The equations of motion were integrated by means of the Verlet alg~rithm?~ and the periodicity of the system was taken into account by applying the usual minimum image convention.28 A test run at the temperature of 1000 K in the microcanonic ( W E ) ensemble showed that a time step of less than 1 fs is necessary in order to obtain good energy conservation during the long simulation runs required to examine the diffusional properties of the sorbate (up to 100 ns). In order to overcome the problem of energy drift and to enable the use of larger integration steps, we decided to perform the simulations in the canonic (NVT) ensemble. This was done by coupling the system to a stochastic thermostat, recently developed by Kast and c o - w o r k e r ~ . ~The ~ thermostat acts as a viscous
Schrimpf et al.
10934 J. Phys. Chem., Vol. 99, No. 27, 1995 TABLE 3: Simulation Parameters of the Runs with Calculated Energetic Data and Diffusion Coefficients"
I
$lane of 12-T-jing
j
-50
1
500
1
650
1 1
800
2
1000 650
2
800
2
1000 650
3 3 3 3
800
lo00 1400
507 100 -4.2 656 40 -3.8 816 25 -3.4 1023 20 -3.0 671 35 -7.2 826 15 -7.2 1023 10 -6.9 646 80 -10.9 806 25 -10.8 1007 20 -10.6 1411 5 -10.1
-73.7 -67.6 -60.4 -52.6 -67.9 -62.5 -55.6 -70.2 -65.9 -60.9 -52.4
-77.8 -71.3 -63.8 -55.6 -75.2 -69.7 -62.5 -81.1 -76.7 -71.6 -62.4
-82.1 -76.8 -70.6 -64.1 -80.7 -76.5 -71.0 -86.4 -83.4 -79.9 -74.2
0.90
4.50 12.00 30.00 1.80 7.20 15.80 0.08 0.40 1.40 6.70
Loadings in average number of molecules per supercage, temperatures in K, energies in kJ/mol, length of the run (tsim) in nanoseconds, m2 s-]. Ttemis the temperature self-diffusion coefficients (D)in of the stochastic thermostat, and Tsysis the actual temperature of the system. a
,E -60 a x
-70
c
-80 -90 -10
-8
-6
-4 -2 0 distance [A]
4. Results
4.1. Adsorption Sites and Energetics. Energy minimizations under constraints were used for the systematic search for local minima. Figure 2 shows the minimum guest-host interaction energy as a function of the distance between the center-of-mass of p-xylene and the plane of the window connecting two supercages. For each distance, the energy minimization was started from 10 randomly chosen center-ofmass positions and orientations for p-xylene. Because of the symmetry of the zeolite, the energy function is mirror symmetric with respect to the plane of the window. Only two crystallographically different types of energetic minima are observed. The minima at -8.8, -4.2, and 4.2 correspond to the principal adsorption sites (C-sites) located inside the supercages in front of sodium ions at positions SII(see Figure 3, bottom). At this site the center of the aromatic lane is located on the C3-axis of the crystal at a distance 2.7 from the sodium ion. The aromatic plane is oriented parallel to the plane of the 6-Tring, and the methyl-methyl vector of p-xylene is located in
K
4
6
Figure 2. Minimum energy of adsorption Uads for p-xylene in zeolite Na-Y as a function of the distance between the center-of-mass of the molecule and the plane of the window.
W-site
medium which reduces the mobility of the sorbate. Therefore, runs using different coupling constants (relative mass of the particles of the thermostat) were performed in order to investigate the influence of the coupling constant on the diffusion coefficient. We found that a coupling constant of 0.1 reduces the diffusivity by a factor 2-3 compared to the W E result. A relative mass of 0.001 yielded the same value for the diffusion coefficient as the NVE simulation and was therefore adopted for all production runs. Test simulations were also performed using different integration times. We found that the integration time can be increased up to 5 fs without influencing the investigated sorbate properties. As a consequence of the weak coupling of the sorbate molecules to the thermostat and the relatively large integration step (5 fs), however, the temperature of the system is in general slightly higher than the temperature of the thermostat. Starting from a uniform distribution of the sorbates over the supercages, the systems were equilibrated for about 1 ns. During the production runs, the positions of the sorbate molecules were stored in intervals of 1 ps for later analysis. Table 3 contains the simulation parameters of the single runs. Also shown are the temperatures of the thermostat, which are used in the text to label the single runs. The MD simulations were performed on a Fujitsu VP 2400/10. The CPU time requirement for 1 ns of simulation time ranges from 40 to 110 min, depending on the sorbate loading.
2
N
Figure 3. Orientation of p-xylene at the C-site inside the supercage (bottom) and at the W-site (top) resulting from the energy minimizations. The sodium ions at the SIIsites are represented by large gray spheres.
the mirror plane of the zeolite. These results are in full agreement with the neutron diffraction data obtained by Czjzek and co-workers for a sorbate concentration of 1.1 molecules/ S C . ~ The above diffraction experiment showed that the major part (88%) of the molecules have the orientation computed by us (orientation A). For 12% of the molecules, however, it was found experimentally that the methyl-methyl vector forms a 30" angle with the mirror plane of the crystal (orientation B). The latter orientation could not be reproduced by energy minimizations with a single molecule. In order to check whether orientation B can be a consequence of steric hindrance between the sorbates, we performed energy minimizations at higher loadings. Indeed, several of the energy minimizations performed with three and four molecules placed in one supercage yielded sorbate configurations with orientation B occupied. We conclude, therefore, that orientation B might not be a stable orientation at the adsorption site, but rather a consequence of cluster formation between p-xylene molecules. This hypothesis is supported by the distribution functions for the supercage occupancy of the sorbate, obtained by the molecular dynamics simulations discussed below. A second type of adsorption site is found in the center of the window (W-site), as depicted at the top of Figure 3. The guesthost interaction energy at the W-site (-58 kT/mol) is much more unfavorable than at the C-site (-89 kJ/mol). An occupation
p-Xylene in Zeolite Na-Y
J. Phys. Chem., Vol. 99, No. 27, 1995 10935
of the window is therefore to be expected only for high sorbate loadings, as observed experimentally for benzene in zeolite NaY.3 The aromatic ring of benzene is oriented coplanar with the plane of the window at the W-site. Such an orientation is not possible for p-xylene. The repulsive interactions between the atoms forming the 12-T-ring and the methyl groups result in an orientation where the planes of the aromatic ring and the window form an angle of about 30". As can be seen from Figure 2, the C-sites are separated from each other by high energetic barriers. The minimum energy barrier for the migration between two adsorption sites inside a supercage (sites at -8.2 and -4.2 A) is 35 kJ/mol, whereas the one for the intercage diffusion is 37 kJ/mol. The energy of adsorption U d Sresults from the time average of the potential energies (V). In the present system Uads is composed of the guest-guest interaction energy Vgg and the guest-host interaction energy Vgh:
35 30
1
-
I
a) 1 moVsc
25
0,
20 15 10
5
n 1
1.5
2
2.5
r [AI
3
3.5
25
b) Ta8OOK
I
4
4.5
4
4.5
.. .
: '.
20 15
(3)
t 0, 10
Table 3 shows the calculated energies of adsorption u a d s with their single contributions (V,,) and (Vgh) and the heats of adsorption Hads, resulting from H a d s = U a d s - RT. We find that for a given temperature the energy of adsorption decreases with the loading. This is a well known property of zeolitesorbate systems, usually caused by the increase of the sorbatesorbate interactions when the loading becomes higher. The interesting phenomenon characteristic for the present system is that not only (Vgg) but also (Vgh) decreases upon increased sorbate concentration. An explanation of this phenomenon can be provided by the analysis of the radial distribution functions g(r) of the center-of-mass of p-xylene with respect to the centers of the supercages, defined as (4) where Nscis the total number of supercages in the system, Qcms is the macroscopic density of the sorbate molecules, and nscems is the number of center-of-supercage / center-of-mass of p-xylene pairs in the system. Figure 4a shows the functions g ( r ) for a loading of 1 molecules/sc at different temperatures. At 500 K the sorbate molecules are confined to regions close to the adsorption sites, which are 3.4 A away from the center of the supercage. At higher temperatures the distribution functions are more extended toward shorter distances. However, the probability density for short distances remains very small, even at 1000 K. This shows that the region located near to the cage center is energetically very unfavorable for p-xylene. The fact that (Vgh) decreases with increasing sorbate concentration can be explained by the evolution of the shape of the distribution function g ( r ) for different loadings. From Figure 4b it can be seen that the amplitude of g(r) for short distances decreases with increasing loading, which means that the occupation of the energetically unfavorable region near the center of the cage is prevented by the mutual steric hindrance of the sorbate molecules. This result also explains why the decrease of (Vgh) with loading becomes more pronounced for higher temperatures, where the molecules have an increased tendency to sample the energetically unfavorable region close to the center of the supercage. In order to compare the energetic properties of the model system with experimental data, we performed simulations covering a real time duration of 500 ps in the zero coverage limit (24 molecules in the system and no sorbate-sorbate interaction) at temperatures 300 and 420 K. The heat of
5 0 1
1.5
2
2.5
r [AI
3
3.5
Figure 4. Radial distribution functions g(r) of the center-of-mass of p-xylene with respect to the center of the supercage: (a) for 1 molecules/ sc at different temperatures; (b) for 1, 2, and 3 moleculeslsc at 800 K.
adsorption calculated by us at 420 K is -80 kJ/mol. This compares well with the experimental result of Ruthven and Goddard,'O who obtained H a d s = -77 kJ/mol by applying the van't Hoff equation to the measured values of the Henry constants for the adsorption isotherms at 403 and 443 K. Recent calorimetric measurements of the enthalpies of adsorption at room temperature by Simonot-Grange and c o - w ~ r k e ryielded s~~ H a d s = -85 k J h " at the zero coverage limit. This is again in good agreement with our simulation result (Hads = -83 kJ/mol at 300 K). However, both much higher and much lower values for Hads, ranging from 50 kJ/mol at 523 K (Bellat et al.3') to 113 kJ/mol at 423 K (Santacesaria et aL3*), have also been reported in the literature. 4.2. Supercage Occupation. The interpretation of experimental data for zeolite-sorbate systems is often complicated by the fact that the distribution of the sorbate in the inner void space of the host framework, e.g., the neighborhood relation for the sorbate molecules, is unknown. This problem occurs, for instance, when the X-ray or neutron diffraction experiments are being interpreted with a view to localize adsorption sites for molecules in zeolites. Various sites or various molecular orientations are frequently observed, and the interpretation of these results depends on the underlying assumptions regarding the sorbate distribution with respect to the zeolite cavities. Detailed local information pertinent to the sorbate distribution is not easily obtained experimentally, although the multiple quantum NMR turned out to be a promissing tool to gain insight into the heterogeneity of the sorbate distribution. Chmelka et al.33have applied this technique to the study of hexamethylbenzene in zeolite Na-Y at 573 K. Their proton multiple quantum NMR results indicate that for macroscopic loadings of 1 and 2 molecules/sc the sorbate molecules are uniformly distributed over the supercages. Figure 5 shows the distribution of supercage occupation for p-xylene obtained from our simulations. The histograms are calculated by assigning a molecule to the supercage with the smallest center-of-cagekenter-of-mass distance. For mean
Schrimpf et al.
10936 J. Phys. Chem., Vol. 99, No. 27, 1995 Drobabilitv
7
(%I
6000
1
50
l
I
I
I
I
I
I
5000
2!
40
I
1 mollsc, T=500 K
average loading: 1 molecule / cage
c.
30 20 10
E
4000
-0n % 2 4U
3000
E
0 1
0
2
3
4
E
2000 1000
0 0
10
20
30
40 50 60 time I ns
70 80
90 100
Figure 6. Ensemble average of the mean-square displacement of p-xylene in zeolite Na-Y for a loading of 1 molecules/sc at 500 K. TABLE 4: Fraction of Occupied Window Regions (See Text for Details) 0
80
1
2
3
4
loading (molecules/sc)
1 1
1000 K
2 2 2
1400 K
~
average loading: 3 molecules / cage
20
3 3 3 3
0 0
1
2 3 molecules / cage
4
Figure 5. Distribution of occupation of the supercages at different temperatures for mean loadings of 1, 2, and 3 molecules/sc.
loadings of 1 and 2 molecules/sc the supercage occupancies are inhomogeneous, while for the highest loading studied the opposite behavior is observed: nearly every supercage is occupied by three molecules. The center-of-mass criterion for assinging the molecules to supercages, however, is not sufficient to provide a detailed insight into the mutual arrangement of the guests. For example, a occupancy of 4 can result from four guest molecules, each one located completely inside the supercage, as well as from four molecules, all located close to one of the four windows of the supercage. In order to be able to differentiate between those situations, we also adopted an additional criterion for the cage occupation: a molecule is assigned to a given supercage only if all of its aromatic carbons are located inside the supercage. Otherwise, the molecule is assigned to the window region. Using the latter criterion, the resulting supercage occupancies are found to be almost identical with the ones depicted in Figure 5. This is a consequence of the very small occupancies of the window regions (Table 4). For the average loadings 1 and 2 molecules/sc, the fraction of supercages which are occupied by three molecules increases with decreasing temperature. It is likely that this trend would continue also for temperatures lower than the ones studied here. This hypothesis is supported by energy minimizations of clusters containing between one and four molecules placed in the same supercage. The results show that the energy of adsorption u& decreases with the number of molecules up to a cluster size of three (Table 5). This is caused by the increase of the interactions between the sorbates. The guest-host interactions remain unaffected up to a cluster size of 3. Insertion of a fourth molecule leads to significant molecular rearrangements, resulting in less favorable guest-host interactions. However, this
fraction of occupied window regions (%)
500 650 800 000 650 800 000 650 800 000 400
1
800 K
40
(K)
1
0650 K
60
Ttem
0.1 0.5 1.3 2.2
0.5 1.4 2.6 0.1 0.3 1.o 3.5
TABLE 5: Results of the Energy Minimizations Performed for Different Numbers of Molecules inside a Supercage no. of molecules 1 2 3 4
vgh
(kJ/mol)
-89.0 -88.9 -88.8 -84.2
V,, (kJ/mol)
0.0 -4.6 -8.2 -12.4
Uads
(kJ/mol)
-89.0 -93.5 -97.0 -96.6
increase of interaction energy is compensated by the equivalent decrease of the guest-guest interaction energy. As a consequence, the energies of adsorption for three and four molecules are nearly the same. At low temperatures the energetic contribution dominates the adsorption behavior and the entropic effects become less important. From our energy minimization results one can therefore conclude that at lower temperatures the tendency for the aggregation of three or even four molecules inside a supercage increases. Diffraction studies on zeolite sorbate systems are often performed at very low temperatures around 5 K. We predict that under these conditions, even for low average loadings of 1 molecules/sc, a large fraction of supercages would be occupied by more than two p-xylene molecules. 4.3. Coefficients of Self-Diffusion. The dynamic property of most interest in the study of zeolite-sorbate systems is the diffusion coefficient of guest molecules. The coefficients of self-diffusion were computed from the mean-square center-ofmass displacement as a function of time. The long simulation times used enabled us to precisely determine the long-range diffusivity. Figure 6 shows the ensemble averaged mean-square displacement resulting from the simulation with 1 molecules/ sc at 500 K. The value 5000 A2,reached at 100 ns, corresponds to the root mean-square displacement of about three unit cells of the zeolite. The diffusion coefficient D can be calculated from the slope of the mean-square displacement according to
p-Xylene in Zeolite Na-Y loo[
,
I
J. Phys. Chem., Vol. 99, No. 27, 1995 10937 ,
I
,
,
1
I
I
1
i [
100
0
7
1
8 .
I
I t
t
E
uptake rate
1 moVsc 2moUsc 3 moVsc
0.001 0.6 0.8
1
1
+ 0
1.2 1.4 1.6 1.8 lOOOK/T
2
2.2 2.4 2.6
Figure 7. Arrhenius plot for the calculated self-diffusion coefficients of p-xylene in zeolite Na-Y at different sorbate loadings. Also shown are the results of PFG-NMR measurementsi4 at a loading of 0.9 moleculeshc and the results of uptake rate measurement^.^^ Uptake rate diffusivities were found to be nearly independent of the sorbate concentration.
TABLE 6: Activation Energies Ea and Preexponential Factors DOfor the Self-Diffusion of p-Xylene in Zeolite Na-Y Determined from the Mean-Square Displacements of the Centers-of-Massof the Molecules loading (molecules/sc) DO mz s-I) E, (kJ/mol) 1 9.4 29 2 9.9 34 3 3.3 45 ~
the Einstein relation:34
D = (Ir(0) - r(r)12)/6r
(5)
The average indicated by the symbols ( ) was carried out over all p-xylene molecules and multiple time origins. The computed self-diffusivities are given in Table 3. Figure 7 shows the diffusion coefficients as a function of the inverse temperature. For all loadings studied, the temperature dependence of D is found to follow the Arrhenius law:
D = Doexp(-E,lRT)
(6)
The activation energies Ea and the preexponential factors DO (Table 6 ) were determined from least squares fits. For the average loading 1 molecules/sc the computed activation energy is 29 kJ/mol, which is inside the range of experimental diffisivity data obtained using various macroscopic and microscopic techniques. For example, Ea = 32 kJ/mol results from the PFG-NMR measurements on low silica Na-Y (SUAl = 1.8) by Germanus et al.,I4 while uptake rate measurements of Goddard and R u t h ~ e non~ natural ~ faujasite (SUA1 = 2.3) are consistent with Ea = 26 Wmol. The dependency of the activation energy on the Si/Al ratio of the zeolite remains an open question, since controversial results have been reported. 1 4 , 3 5 3 It can be seen from Figure 7 that the absolute values of the diffusion coefficients obtained for 1 molecules/sc are in close agreement with the PFG-NMR results of Germanus et al. for 0.9 molecules/sc; at 500 K our value for D is about 3 times lower than the experimental one, On the other hand, the uptake rate diffusivities are over 2 orders of magnitude smaller than the PFG-NMR ones. The reasons for this discrepancy, generally found for C8 aromatics in faujasite-type zeolites, is still a matter of d i s c u s ~ i o n . ~Furthermore, ~-~~ the macroscopic and micro-
1'
0.6
0.8
1
1.2
1.4
1000 K / T
1.6
1.8
I
2
Figure 8. Arrhenius plot for the decorrelation times for the out-ofplane rotation (tout) and in-plane rotation (r,J at different sorbate loadings.
scopic techniques lead to different results concerning the concentration dependence of the diffusion coefficients. The PFG-NMR data suggest a rather rapid decline of diffusivity as saturation is approached, whereas the uptake rate diffusivities show no evidence of any such trend. Again, the molecular dynamics results are in agreement with the PFG-NMR data as we find a strong loading dependence of the diffusivity. However, the simulations predict that the activation energies should increase significantly with loading, while the corresponding PFG-NMR results exhibit a slightly decreasing trend. 4.4. Rotational Dynamics. The rotational dynamics of p-xylene was investigated by studying orientation autocorrelation functions c(t), defined as
4)= (U(O).U(O)
(7)
where u is the unit vector describing the molecular orientation in space. We investigated two kinds of motions, namely, the rotation in the aromatic plane around the C6-axis of the aromatic system (in-plane rotation) and the rotation of the aromatic plane itself (out-of-plane rotation). The in-plane rotation was studied by following the time evolution of the unit vector pointing from one methyl group to another, whereas the out-of-plane rotation was studied using the unit vector normal to the aromatic plane. Both types of time correlation functions show an approximate simple exponential decay, c(t) = exp(-tlz), and can, therefore, be characterized by a single time constant, the decorrelation time z. Figure 8 presents the calculated time constants for the in-plane and out-of-plane rotations (z,,, and rout,respectively) as a function of the inverse temperature. As can be seen in the figure, both types of decorrelation time follow an Arrheniustype temperature dependence of the form z = zo exp(E,lRT)
(8)
with activation energies Ea from 11 to 12 !d/mol for the inplane rotation and from 29 to 32 H/mol for the out-of-plane rotation. For all loadings the in-plane decorrelations are much faster than the out-of-plane ones (about 1 order of magnitude at 650 K). Furthermore, their low activation energies indicate that these rotations are not strongly coupled to the translational
Schrimpf et al.
10938 J. Phys. Chem., Vol. 99, No. 27, 1995 diffusion of p-xylene. It will be shown in the next chapter that the time constants Ti, are always smaller than the residence times at the adsorption sites; e.g., C6-axis rotation is fast compared to migration between adsorption sites. The decorrelation of the in-plane rotation becomes slower as the sorbate loading is increased. This shows that molecules located at adjacent sites inside a supercage interact significantly and thus reduce each other's rotational mobility in the aromatic plane. The blocking of this rotation has been previously observed by Czjzek and co-worker~,~~ who studied the behavior of p-xylene in zeolite (Na,Yb)-Y using quasi-elastic neutron scattering on the time scale of lo-" s. In further agreement with the scattering experiment is the fact that up to a temperature of 460 K the out-of-plane rotation of p-xylene could not be observed. The molecular dynamics results extrapolated to 460 K indicate that the out-of-plane rotation can be observed on a time scale of s, which is 2 orders of magnitude larger than the one used in the experiment. As mentioned above, the decorrelation is a much slower process for the out-of-plane than for the in-plane rotation. The decorrelation times routare also considerably longer than the average residence times at the adsorption sites (see next chapter). This indicates that the out-of-plane rotation is mainly caused by the diffusional migration of the sorbate. The fact that the diffusion of sorbate molecules between sites is accompanied by orientational decorrelation can be used to derive the diffusion coefficients of the guest molecules from the line shape analysis of 2H-NMR data. This technique has been applied, for example, by Burmeister et al.36and Bull et a1.'* to investigate the mobility of aromatics in faujasite-type zeolites. They have interpreted the NMR decorrelation time as the mean site residence time z and applied the EinsteinSmoluchowski relation for the 3-dimensional diffusion in a cubic lattice,
D = 1216r
(9)
to calculate the diffusion coefficients D . Burmeister et al. assumed a constant jump length 1 of 11 A, while Bull et al. used a jump length of 5 A. The former value corresponds to the distance between the centers of the supercages, and the latter length corresponds roughly to the distance between neigboring adsorption sites (C-sites) for aromatics in faujasites. The question arises as to whether for the present system a constant jump length 1 can be assumed in eq 9 in order to calculate the diffusion coefficients from the decorrelation times of the out-of-plane rotations. Indeed, for a loading of 1 molecules/sc we find a jump length of 4.6 A, which is independent of the temperature. This length is close to the one used by Bull et al. However, assuming this jump length for higher loadings would overestimate the diffusivities and, moreover, would not reproduce the correct activation energies for the diffusion. For example, a jump length of 4.6 A predicts that the diffusion coefficient does not change when the loading is increased from 1 to 2 moleculeslsc. This is in contrast to the diffusivities obtained from the mean-square displacements of the guests. Actually, the jump length is found to become smaller at higher loadings and also temperature dependent. For 2 molecules/sc 1 increases from 2.8 A at 650 K to 3.5 8, at 1000 K, and for 3 molecules/sc 1 increases from 1.1 at 650 K to 1.9 A at 1400 K. The reasons for the variation of the jump length can be generally understood by the following considerations. The longrange diffusion results from the migration between adsorption sites located in different supercages. The decorrelation of the out-of-plane rotation, however, results from both: from the
Figure 9. Sketch of the T-atom framework of two neighboring supercages illustrating the translational motion of the center-of-mass of p-xylene. Shown is a single molecule trajectory of duration 2.2 ns, taken from the simulation at 500 K with a loading of 1 molecules/sc. The sodium ions at the SIIsites are represented by large gray spheres.
migration between adsorption sites inside a supercage and from the migration between supercages. The jump length I will therefore depend on the ratio of intercage to intracage migration. For low loadings, where the interactions between the sorbate molecules are not important, the energetic barriers for the intercage and intracage migration are roughly the same. This can be deduced from the energy profile in Figure 2 and is also reflected in the fact that the activation energies for the out-ofplane decorrelation and for the diffusion have the same value (29 kJ/mol). As a result, the ratio between intercage and intracage diffusion does not depend on temperature. This means that the fraction of jumps which contribute to the long-range diffusion remains constant, and as a result, the jump length is temperature-independent. Increasing the loading up to 2 molecules/sc causes a reduction of the jump length due to the decrease of the ratio of intercage to intracage jumps. Moreover, this ratio, and therefore the jump length 1, is now temperaturedependent, because the energetic barriers for the intercage and intracage migration differ. While the activation energy for the out-of-plane decorrelation remains 29 kJ/mol, the activation energy for the diffusion is now higher by about 5 kJ/mol. Thus, an increase of the temperature leads to an increase of the intercage to intracage jump ratio and therefore to an increase of the jump length 1. The explanation for the variation of the jump length for the loading of 3 molecules/sc can be done along the same lines. 4.5. Residence Times. The diffusion coefficients of sorbates quantify the long-range mobility but provide no information about the mechanism of the diffusion process. An example of the translational motion of the center-of-mass for p-xylene is given in Figure 9. The figure shows a single molecule trajectory of duration 2.2 ns, taken from the simulation at 500 K with a loading of 1 molecules/sc. It can be seen that the molecule is strongly localized in front of the sodium ions at sites SII. Apparently, the residence times at the adsorption sites are much longer than the time spans needed for the change of adsorption site: e.g., the diffusion has features characteristic for a hopping process. In order to gain a more quantitative insight into the diffusion mechanism, we calculated the distribution function e(t)for the residence time at the adsorption site (C-site) as well as inside the supercage. The functions Q(t)are normalized as
p-Xylene in Zeolite Na-Y
J. Phys. Chem., VoE. 99, No. 27, 1995 10939
p-Xylene was assigned to the appropriate regions using the nearest distance criterion for the center-of-mass of the molecule with respect to the sodium ions at sites SIIand with respect to the centers of the supercages. The calculated mean residence times at the adsorption sites (tads) and the mean supercage residence times (tsc)are listed in Table 7 and shown in Figure 10 as a function of the inverse temperature. For the concentrations of 1 and 2 molecules/sc we find that the mean residence times approximately follow the Arrhenius-type temperature dependence. The activation energies for both (tads) and (tSC)are about 29 kJ/mol at these sorbate loadings. At the loading of 3 molecules/sc the temperature dependence is significantly less pronounced for (tads) than for (tsc).In addition, the activation energy seems to decrease with decreasing temperature for the loading of 3 molecules/sc. The effect of the mutual interactions of the guest molecules on their translational mobility can be easily observed in the diagram. Due to steric hindrance both (fads) and (tsc)increase with increased loading. The mobility inside a supercage relative to the mobility between supercages is quantified by the ratio (tsc>/ (rads). For example, at 800 K and for loadings of 1,2, and 3 molecules/sc the ratio is 9, 11, and 44,respectively. These values demonstrate that with increasing sorbate concentration the intracage migration process becomes faster compared to intercage migration. This trend was observed for all conditions studied here, although the exact ratio (&)/(tads) for an individual loading depends on temperature. It is furthermore interesting to consider the integrals Z(t) of the distribution functions e(t),
I( t ) =
hr@ (t') df
TABLE 7: Calculated Values of the In-Plane and Out-of-Plane Decorrelation Time, q,, and rmt,Respectively, the Mean Residence Times inside the Supercages and at the Adsorption Sites, (tW)and (tab),Respectively, and the Mean Lifetime of the Orientation of the Aromatic Plane with Respect to the Center of the Supercage, (t,,l.,>n Tsys (K)
loading (molecules/sc) 1
507 656 816 1023 671 826 1023 646 806 1007 1411
2 3
tin
tout
(PS) 7.8 4.2 2.6 1.7 5.4 3.6 2.6 8.7 6.0
(PSI
4.5
3.0
(rSJ (PSI 750.0 171.3 63.1 29.6 324.8 113.4 57.1 3501.0 1186.4 321.0 76.8
390.3 85.4 30.1 12.0 73.4 28.1 12.1 240.6 54.3 27.0 9.1
(rads)
(tplane)
(PS)
(PSI 50.1 16.5 6.9 3.4 18.6 8.8 5.1 28.6 17.8 11.5 6.0
113.3 22.5 7.3 3.1 27.4 10.0
4.8 56.0 27.2 13.6 5.4
See text for definition of the single properties.
,.
c
K . g
,
i
I
100
'P
E
3 10
which give the fraction of residence times smaller than or equal to a given time t. Figure 1l a shows the integrals fads(t) of the distribution of the adsorption site residence times for a temperature of 800 K. The distribution of the residence times is broad, and the times reach values up to about 100 ps. It is worth noting that the fraction of events with residence times smaller than 1 ps is very large (30-40%). Information about the influence of the mutual interactions of the guests on the diffusional behavior inside a supercage can be gained by analyzing the temporal sequence of the adsorption sites which are occupied by the single molecules. As a first approximation, the intracage diffusion of p-xylene can be considered as a Markow process consisting of jumps between adsorption sites: e.g., it is assumed that the probability for a jump from one site to another depends only on the actual position of the molecule and is therefore independent of its past. We now consider the hopping sequence A B C, where A, B, and C denote adsorption sites inside the same supercage. Because a supercage contains four adsorption sites in a tetrahedral arrangement, the probability that a particle returns to the initial site (C = A) is 1/3. However, this is valid only for a single molecule inside a supercage, as long as neighboring adsorption sites are not occupied by other particles. For higher loadings, the retum probability is expected to depend on the rate of molecular rearrangements inside the supercage after a jump event A * B and will therefore also be a function of the residence time at site B. Figure I l b shows the probability for the return event eret(t) versus the residence time at site B. For the lowest concentration studied, we indeed find the theoretically estimated value of 1/3 for e&), which is constant for all residence times. At a loading of 2 molecules/sc the return
-
1 ; ,2.*, b':- ; 1 ,
,
1 mollsc
Q
3 mollsc
1 0.6
0.8
1
1.2 1.4 1000 K / T
1.6
1.8
2
Figure 10. Arrhenius plot for the mean residence times inside the supercages, ( f S c ) (dotted lines), and at the adsorption sites, lines), at different loadings.
(fads)
(solid
probability is significantly increased for residence times shorter than 20 ps. For 3 molecules/sc the residence time dependence is very pronounced: the return probabilities are larger than 70% for residence times shorter than 10 ps. The time span necessary for complete molecular randomization after a jump event is found to be of the order of 100 ps. The integrals Is&) of the distribution of the supercage residence times are depicted in Figure 1IC for the temperature of 800 K. The residence times range from 1 ps to 1 ns for lower loadings and can be as large as 10 ns for 3 molecules/sc. In spite of the long mean residence times, there is a significant fraction of events with residence times smaller than 1 ps. In analogy to our previous analysis of the migration between adsorption sites, we have considered the cage-to-cage hopping sequence A B C and calculated the probability e & that a particle returns to the initial supercage (C = A). The results, presented in Figure lld, can be interpreted as follows. The centers of the supercages of faujasite-type zeolites form a diamond lattice; e.g., the supercages are interconnected in a tetrahedral arrangement. From the assumption that the longrange diffusion is a Markov process of jumps between neigh.-)
10940 J. Phys. Chem., Vol. 99, No. 27, 1995
Schrimpf et al.
........ ........ . I
1
I
'......I
.
2/3
113
-..........
1 molkc 2 moVsc
----
-
3 mol/sc .....
2/3
1 13
0
1
10
1000
100 t/Ps
Figure 11. Integrals [(t) of the distributions of the residence times and the correspondingretum probabilities e&) for the migration between the regions (see text) at 800 K for different loadings. Left: for the adsorption sites; right: for the supercages. boring supercages, it follows, therefore, that the return probability is 1/4. This value is reached for all concentrations in the limit of long residence time. However, for residence times shorter than 100 ps the retum probability is much higher, even for the lowest loading studied. This result is caused by' three reasons. First, for higher loadings, the molecules entering supercage B might be immediately scattered back into supercage A by molecules occupying supercage B. Second, the molecules entering supercage B are preferentially adsorbed at one of the three sites that are located next to the window connecting the supercages A and B, whereas the fourth site can be occupied only after the time of the order of the mean adsorption site residence time. Therefore, if the residence time of the molecule in supercage B is not sufficiently long for the intracage positional randomization, its retum probability to cage A is momentarily increased. Third, a small percentage of the molecules leaving supercage A may become trapped inside the shallow local minimum in the window region (Figure 2). Their vibrational motion inside the energy well is registered by the geometrical criterion (also commonly used in other simulation studies) as a series of frequent A B A events. Although the number of molecules actually trapped may be very small, there may be a large number of those events per every single molecule, which results in the marked increase of the nominal return probability. Therefore, the mean supercage residence times (tsc)as calculated here cannot be used in a straightforward manner to derive the diffusion coefficients by a simple model that is assuming a Markow process of jumps between supercages. 4.6. Intrasupercage Migration and Molecular Reorientation. The decorrelation times of the out-of-plane rotation tout are longer than the mean residence times at the adsorption sites (tads) and shorter than the mean residence times inside the supercages (tSc>.This is observed for all the conditions studied here (see Table 7). Therefore, out-of-plane decorrelation must arise to a great extent from site-to-site migration inside the supercages. In order to gain more information about the molecular reorientation process inside the supercage, we analyzed the function
-
+
30) = ~,(~).~,,,(~)
(12)
for the single molecule trajectories. In eq 12 uno,,,, denotes the
I
I
1
g 0 -1
I
0
200
400
600
800
time I ps
I
1000 1200 1400
Figure 12. Bottom: the function s ( t ) calculated from a trajectory of a single molecule at 500 K and the loading of 1 molecules/sc. The molecule does not leave its supercage during the time span shown. Top: the instantaneous adsorption site of the molecule, indicated by the height of the solid line. unit vector normal to the aromatic plane of the sorbate and u ~denotes ~ the ~ unit- vector ~ pointing ~ from the instantaneous center-of-mass of p-xylene to the center of the supercage. The function s(t) reaches values of 1 or - 1, when the aromatic plane is oriented coplanar to the intemal surface of the supercage. The minimum energy configuration of p-xylene at the C-site corresponds to such an orientation. The function s(t) enables one to distinguish between the two possible mechanisms for the intrasupercage migration. Either the aromatic ring remains in contact with the wall of the cage during the site-to-site migration, or the migration is accompanied by detachment of the aromatic plane from the wall followed by the inversion of the aromatic plane with respect to the center of the supercage. The latter process occurs when s(t) changes sign. The mean lifetimes of the aromatic plane orientation (tplme),calculated from the number of sign changes of the function s(t), are listed in Table 7. In Figure 12 we compare the molecule orientation as given by the function s(t) calculated from a trajectory of a single molecule at the loading of 1 moleculesk at 500 K, with the temporal sequence of the occupied adsorption sites (indicated by the numbers 1, 2, 3, and 4) inside a single supercage. As can be seen from the figure, most site-to-site migrations are accompanied by inversion of the orientation of the aromatic
p-Xylene in Zeolite Na-Y plane with respect to the center of the supercage. One exception is observed at the time of 540 ps. Here, s(t) changes sign, but the orientations at the initial and final adsorption site are identical. Strong fluctuations of the molecule orientation, leading to short-time inversion of the sign of s(t), occur also while the molecule is attached to a single adsorption site (for times 900,1080,1180, and 1330 ps in Figure 12). The behavior of the function s(t) gets more complex at increased loadings and at elevated temperatures. However, from the visual inspection of single trajectories the following general statements can be made: (i) complete inversion of the aromatic plane within the region of an adsorption site is never observed, (ii) most of the site-to-site migrations are accompanied by inversion of the aromatic plane, and (iii) the fraction of site-to-site migrations without inversion gets larger with increasing sorbate concentration. The latter phenomenon is thought to be due to the sterical requirement for the inversion process.
5. Conclusions The properties of p-xylene adsorbed in zeolite Na-Y were studied by molecular simulations. Constrained reaction coordinate minimizations were used to calculate the minimumenergy path for the diffusion. The predicted adsorption site inside the supercage was found to be in good agreement with neutron diffraction results. Extensive molecular dynamics simulations were performed in order to study the behavior of p-xylene for mean loadings of 1, 2, and 3 molecules per supercage at temperatures above 500 K. The distribution of the molecules over the supercages was obtained, giving important information for the interpretation of experimental data. For mean loadings of 1 and 2 molecules/sc the supercage occupancies were found to be far from homogeneous. The guest-host interaction energy was found to decrease with increased loading. This phenonemon could be explained in terms of the spatial distribution of p-xylene inside the supercages. The calculated heats of adsorption are in close agreement with experimental data. Due to very long simulation runs the diffusion coefficients of p-xy!ene could be determined with accuracy. The diffusion coefficients obtained for low occupancy are in close agreement with the PFG-NMR results. It was found that the diffusivity decreases drastically with increased sorbate loading. The dynamics of the in-plane and out-of-plane rotation was investigated using orientation autocorrelation functions. The in-plane rotation is fast compared to the mean residence time at the adsorption site, while the out-of-plane rotation is strongly coupled to the migration between the adsorption sites. The microscopic details of the diffusion process were studied by means of site and cage residence times. Furthermore, the correlation between the intrasupercage migration and molecular reorientation was investigated. The simulations reveal that siteto-site migration inside a supercage is frequently accompanied by the reorientation of the aromatic plane with respect to the center of the supercage.
Acknowledgment. We thank Dr. A. Radlinsky for carefully reading the manuscript.
J. Phys. Chem., Vol. 99,No. 27, 1995 10941 References and Notes (1) Neuzil, R. W. U.S. Patent 3,558,730, 1971. Neuzil, R. W. U.S. Patent 3,558,732, 1971. (2) Neuzil, R. W. US.Patent 3,626,020, 1971. (3) Fitch, A. N.; Jobic, H.; Renouprez, A. J. J. Phys. Chem. 1986,90, 1311. (4) Czjzek, M.; Vogt, T.; Fuess, H. Zeolites 1991, 11, 832. (5) Czjzek, M.; Fuess, H.; Vogt, T. J. Phys. Chem. 1991, 95, 5255. (6) Mellot, C.; Espinat, D.; Rebours, B.; Baerlocher, C.; Fischer, P. Catal. Lett. 1994, 27, 159. (7) Klein, H.; Fuess, H. In Zeolites and Related Microporous Materials: State of the Art 1994; Studies in Surface Science and Catalysis; Weitkamp, J., Karge, H. G., Pfeifer, H., Halderich, W., Eds.; Elsevier: Amsterdam, 1994. (8) Lian Su, B.; Manoli, J. M.; Potvin, C.; Barthomeuf, D. J. Chem. Soc., Faraday Trans. 1993, 89, 857. (9) Barthomeuf, D.; de Mallmann, A. Ind. Eng. Chem. Res. 1990, 7, 1437. (10) Ruthven, D. M.: Goddard, M. Zeolites 1986. 6. 275. (1 1) Bellat, J. P.; Simonot-Grange, H. M.; Jullian, S. C. R. Acad. Sci. 1993, 316, 1363. (12) Kiirger, J.; Pfeifer, H. Zeolites 1989, 9, 267. Kiirger, J.; Ruthven, D. M. Zeolites 1987, 7, 90. (13) Eic, M.; Goddard, M. V.; Ruthven, D. M. Zeolites 1988, 8, 258. (14) Germanus, A.; KSirger, J.; Pfeifer, H.; Samulevic, N. N.; Zdanov, S. P. Zeolites 1985, 5, 91. (15) Jobic, H.; Bee. M.; Kiirger, J.; Pfeifer, H.; Caro, J. J. Chem. Soc., Chem. Commun. 1990, 341. (16) Catlow, C. R. A. Modelling of Structure and Reactiviry in Zeolites; Academic Press: London, 1992. (17) Demontis, P.; Yashonath, S.; Klein, M. L. J. Phys. Chem. 1989, 93, 5016. (18) Bull, L. M.; Henson, N. J.; Cheetham, A. K.; Newsam, J. M.; Heyes, S. J. J. Phys. Chem. 1993, 97, 11776. (19) Schrimpf, G.; Schlenkrich, M.; Brickmann, J.; Bopp, P. J. Phys. Chem. 1992, 96, 7404. (20) Schrimpf, G.; Brickmann, J. Accepted for publication in J. Comput.Aided Mater. Des. (21) Mosell, T.;Schrimpf, G.; Brickmann, J. Manuscript in preparation. (22) Goodbody, S. J.; Wantanabe, K.; MacGowan, D.; Walton, J. P. R. B.; Quirke, N. J. Chem. Soc., Faraday Trans. 1991, 87, 1951. (23) Fischer, W.; Brickmann, J. Ber. Bunsen-Ges. Phys. Chem. 1982, 86, 650. (24) Yashonath, S.; Price, S. L.; McDonald, I. R. Mol. Phys. 1983, 64, 361. (25) Dufner, H.; Schlenkrich, M.; Brickmann, J. Submitted for publication in J. Compur. Chem. (26) Ewald, P. P. Ann. Phys. 64,1921, 253. (27) Verlet, L. Phys. Rev. 1967, 159, 98. (28) Allen, M. P.; Tildesley, D. J. Computer Simulations of Liquids, Claredon Press: Oxford, 1987. (29) Kast, S. M.; Nicklas, K.; Biir, H. J.; Brickmann, J. J. Chem. Phys. 1994, 100, 566. (30) Simonot-Grange, M. H.; Bertrand, 0.;Silverdier, E.; Bellat, J. P.; Saulin, C. Submitted to J. Ther. Anal. (31) Bellat, J. P.; Simonot-Grange, M. H.; Jullian, S. C. R. Acad. Sci. 1993, 316, 1363. (32) Santacesaria, E.; Gelosa, D.; Picenoni, D.; Danise, P. J. Colloid Inreface Sci. 1984, 98, 467. (33) Chmelka, B. F.; Pearson, J. G.; Liu, S. B.; Ryoo, R.; de Menorval, L. C.; Pines, A. J. Phys. Chem. 1991, 95, 303. (34) McQuarrie, D. A. Statistical Mechanics; Harper and Row: New York, 1976. (35) Goddard, M.; Ruthven, D. M. Zeolites 1986, 6, 283. (36) Burmeister, R.; Schwarz, H.; Boddenberg, B. Ber. Bunsen-Ges. Phys. Chem. 1989, 93, 1309. (37) Kiirger, J.; Ruthven, D. M. Zeolites 1989, 9, 267. (38) Kiirger, J.; Ruthven, D. M. Diffusion in Zeolites and Other Microporous Solids; John Wiley & Sons: New York, 1992. (39) (39) Czjzek, M.; Jobic, H.; Bee, M. J. Chem. Soc., Faraday Trans. 1991.87, 3459. JP9502868
