JOURNAL OF CHEMICAL PHYSICS

VOLUME 112, NUMBER 2

8 JANUARY 2000

Diffusion anomaly from analytical formula Subir Sarkara) School of Physical Sciences, Jawaharlal Nehru University, New Delhi, 110 067 India

A. V. Anil Kumarb) Solid State and Structural Chemistry Unit, Indian Institute of Science, Bangalore, 560 012 India

Subramanian Yashonathc) Solid State and Structural Chemistry Unit, Indian Institute of Science, Bangalore, 560 012 India and Supercomputer Education and Research Centre, Indian Institute of Science, Bangalore, 560 012 India

共Received 28 June 1999; accepted 21 September 1999兲 An analytic expression is derived for the diffusion coefficient of a sorbate in a crystalline porous solid with bottlenecks. The diffusion coefficients obtained from the analytic expression is found to agree well with the molecular dynamics results. It is also found to reproduce the temperature dependence of the levitation effect for zeolites Y and A. The present calculations provide a strong theoretical support for the levitation effect obtained so far purely from molecular dynamics calculations. © 2000 American Institute of Physics. 关S0021-9606共99兲50547-3兴

as, e.g., aluminophophates 共VPI-5兲 also exhibit such anomaly.15 Further, it has been shown that the peak in D occurs when the levitation parameter ␥ defined as14

Understanding diffusion processes in porous media is a subject of fundamental as well as technological importance.1–3 Examples of such processes in nature are ion diffusion in biomembranes, reverse osmosis, molecular sieving—all of which depend on diffusion properties of molecules in porous media. There exist a large number of experimental and computational studies of diffusion in porous media in the literature which are listed in the review by Bates et al.2 More recently, there have been some studies on diffusion in zeolites.4–10 Some of the computational approaches have yielded significant insights into issues such as interaction between framework vibrations and diffusing species,6,7 single file diffusion,8,9 use of transition state model among others.10 Barrer11 discusses several theoretical approaches to understand diffusion in zeolites. But most of these relate to the macroscopic or Fickian diffusivity. There are few analytic treatments of the problem of diffusion, in particular self-diffusion in porous media to date. Many of these studies eventually resort to use of Ising or lattice gas models for the calculation of the diffusion properties.12,13 Recent studies of dependence of diffusion coefficient D, on sorbate diameter ␴ , have shown the existence of two regimes: 共i兲 a linear regime in which D⬀(1/␴ 2 ) and 共ii兲 an anomalous regime in which D exhibits a peak.14 It has been shown that for spherical Lennard-Jones systems interacting via purely dispersive and repulsive interactions, such a peak is found irrespective of the topological and geometrical features of the void network thereby supporting the view that the effect is universal.15 In view of this it is likely that this effect will have a role to play in all situations where diffusion takes place in a confined media. Indeed, it has been found that spherical sorbate in hosts other than zeolites such

␥⫽

共1兲

approaches unity. Here ␴ sz is the sorbate-porous solid Lennard-Jones interaction parameter. In the case of zeolites, this is the sorbate-oxygen parameter. ␴ w is the diameter of the bottleneck viz., the window defined as the center of mass–center of mass distance between diagonally placed oxygens. For zeolites such as Y and A this is the dimension of the 12-ring and the 8-ring, respectively, which have approximately dimensions of 8 Å and 4.5 Å. When ␥ approaches unity, the sorbate diameter is comparable to the window diameter. Thus, Eq. 共1兲 says that when this happens, a peak in the diffusivity is to be expected. Derouane and co-workers16 have suggested that under certain conditions diffusivity would increase based on their theoretical-cumcomputational approach. However, there are as yet no experimental evidence or rigorous theoretical treatments in support of the existence of such a peak. In this work, we propose an analytic treatment which provides an expression for the diffusion coefficient for sorbates in porous systems consisting of three dimensional networks with bottlenecks. Interestingly the expression reproduces the above anomaly in the diffusivity. Results obtained here are compared with those of long molecular dynamics 共MD兲 simulations in three zeolites viz., A, Y, and ␳ . Zeolites are aluminosilicates with extended three dimensional network of corner shared TO4 共T⫽Si,Al兲 tetrahedra which give rise to voids of different sizes and shapes. For the purpose of derivation we shall assume that there are large voids with narrow windows interconnecting these voids. Examples of this are several: A, Y, and ␳ zeolites. Figure 1 shows two cages along with the interconnecting window for these three zeolites. The diameter of the cages is in the range

a兲

Electronic mail: [email protected] Electronic mail: [email protected] c兲 Also at Condensed Matter Theory Unit, Jawaharlal Nehru Center for Advanced Scientific Research, Bangalore 560012, India. Electronic mail: [email protected] b兲

0021-9606/2000/112(2)/965/5/$17.00

2 7/6. ␴ sz ␴w

965

© 2000 American Institute of Physics

Downloaded 04 Jan 2005 to 130.102.2.60. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

966

J. Chem. Phys., Vol. 112, No. 2, 8 January 2000

Sarkar, Kumar, and Yashonath

FIG. 2. Schematic representation of a neck ⍀ in a zeolite.

representation of one such neck ⍀. From elementary kinetic theory, the number of atoms hitting the neck per unit area around the point sជ on the neck is ␾ (sជ ) which is given by18

␾ 共 sជ 兲 ⫽

FIG. 1. Two neighboring cages with the interconnecting 共a兲 8-ring window of zeolite A, 共b兲 12-ring window of zeolite Y, and 共c兲 8-ring window of zeolite ␳ .

11–13 Å while the window diameter varies between 4.5 to 7.5 Å. We shall present the derivation for a situation with cubic symmetry which will lead to isotropic diffusion. If P(nជ c ,t) is the probability for the atom to be in the cage labeled by the vector nជ c ⫽(n cx ,n cy ,n cz ) at time t, then the dynamics of P(nជ c ,t) is governed by

⳵ P 共 nជ c ,t 兲 ⫽␪ ⳵t

冋兺 ␣ជ

册

P 共 ␣ជ ,t 兲 ⫺q P 共 nជ c ,t 兲 ,

共2兲

where the sum is over all the immediate neighbors ␣ជ of nជ c and ␪ is the rate of loss of probability through any bottleneck normalized with respect to the probability of being in that cell. Here q is the number of exits for a given cage. q⫽4 for zeolite Y and 6 for zeolites A and ␳ . ␪ can only depend on the local concentration of sorbates and we will demonstrate that it approaches a finite limit as the concentration of sorbates goes to zero. Please note that 关 兺 ␣ជ P( ␣ ជ ,t)⫺q P(nជ c ,t)兴 is nothing but the discrete representation of a 2 “ 2 P where a is the lattice constant and “ 2 is the Laplacian operator. This immediately leads to the equation

⳵P ⫽␪a2 “2P ⳵t

共3兲

from which we read off the diffusion coefficient D to be ␪ a 2 . What remains now is to evaluate ␪ . Earlier molecular dynamics simulations have shown that the rate determining step for diffusion is the passage through this window.17 Therefore we do it by assuming a situation of quasiequilibrium and by applying some elementary results of kinetic theory of gases. The pore within a cell is a connected volume of space denoted by R with generally more than one bottleneck which act as exit points. Let us consider one such neck ⍀ and find out what is the loss of probability through it. Figure 2 is a

冑

k BT n 共 sជ 兲 , 2␲m

共4兲

where T is the temperature, k B is the Boltzmann constant, m is the mass of the atom, and n(sជ ) is the number density of atoms around the point sជ . Please note that the equation of local flux given in Eq. 共4兲 does not require the atoms to be noninteracting. In other words, Eq. 共4兲 should be valid for low as well as high concentration of sorbates. Thus the total flux through the neck is given by

冕␾ជ ⍀

共 s 兲 dA⫽

冑 冕 k BT 2␲m

⍀

n 共 sជ 兲 dA.

共5兲

Remember that ␪ is the ratio of the rate of loss of probability per window or bottleneck of the unit cell divided by the probability to be within that cell. In terms of the number of atoms this translates into the flux of atoms per window or bottleneck of the cell divided by the number of atoms within the unit cell. Thus

冑k B T/2␲ m ␪⫽

冕

R

冕

⍀

n 共 sជ 兲 dA 共6兲

,

n 共 rជ 兲 dV

where n(rជ ) denotes the number density at a point rជ within the ␣ cage. The integral in the denominator is actually a volume integral over the accessible part of the cage, R, i.e., the pore space. Now we use the assumption of quasiequilibrium and the low concentration limit of the sorbates to replace n(sជ ) and n(rជ ) by ␦ e ⫺V(sជ )/k B T and ␦ e ⫺V(rជ )/k B T , respectively, where ␦ is a constant of proportionality which will eventually cancel itself out between the numerator and denominator. Here we have used the low concentration limit. V(rជ ) 关or V(sជ )兴 is the potential experienced by the diffusing atom due to the interaction with the zeolite at the point rជ 共or sជ ). Therefore,

冑k B T/2␲ m ␪⫽

冕

R

冕

⍀

e ⫺V(sជ )/k B T dA

e ⫺V(rជ )/k B T dV

.

共7兲

The diffusion coefficient is given by this expression for ␪ multiplied by a 2 where a is the smallest repeating unit which

Downloaded 04 Jan 2005 to 130.102.2.60. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

J. Chem. Phys., Vol. 112, No. 2, 8 January 2000

is usually the lattice constant. In the case of the zeolites considered here the value of a is dimension of the supercage or the ␣ cage. Several comments are in order regarding the derivation. Please note that it is possible to obtain an expression similar to Eq. 共7兲 from transition state theory 共TST兲.10 However, the present derivation gives several insights on the present problem which are not obtained if we follow the TST. These are 共a兲 we have demonstrated that, in the limit of concentration of sorbates going to zero, ␪ has a well-defined value which is determined by the interaction with the backbone of the zeolite. Since Eq. 共3兲 constitutes the description of diffusive motion with 具 兩 xជ (t) 2 兩 典 proportional to t, we have provided a theoretical basis for expecting that this should be seen at arbitrarily large time scales, even beyond what one can usally go up to in a MD integration—which is often rather limited. 共b兲 Since the motion governed by Eq. 共2兲 is diffusive only on a time scale on which the atom has gone through several unit cells, this sets a lower limit on the time scale in MD below which the motion cannot be expected to be diffusive. This observation is important for a numerical as well as experimental determination of the diffusion constant from 具 兩 xជ (t) 兩 2 典 vs t data. In PFG-NMR measurements and also scattering experiments, one should ensure that the sorbate has diffused over a distance several times a. Finally 共c兲 what is the criterion for quasiequilibrium to exist? This will hold provided the atom undergoes many collisions with the wall of the pore before leaking out to a neighboring cell. Typical time scale for collision with the wall ( ␶ c ) is a/ v (T) where v (T) is the characteristic speed at temperature T. Time scale ( ␶ d ) for leaking to a neighboring cell is 1/␪ . Thus

␶c Ⰶ1 ␶d implying that D a 冑3k B T/m

Ⰶ1.

Even though Eq. 共7兲 represents a simple analytical formula for the diffusivity, the rather complex structure of the backbone of the zeolite leads to a very complicated potential function. As a result, although it would be ideal to compute the diffusivity from Eq. 共7兲 analytically, we can do it only numerically at present. Dependence of self-diffusion coefficient D on sorbate diameter ␴ have been obtained from Eq. 共7兲 for zeolites A, Y, and ␳ which have a unit cell composition of Na32Ca32Al96Si96O384 , Na32Al48Si144O384, and Cs1.15共SiAl兲48O96D, respectively. In order to verify if Eq. 共7兲 predicts the values of self-diffusivity D accurately, we will then compare them with those obtained from MD simulations. MD results for zeolite Y are already available from the literature 共see Refs. 14 and 22兲. For zeolites A and ␳ no such results exist in the literature. Therefore, MD simulations under different conditions have also been carried out here on zeolites A and ␳ . The crystal structures of zeolites required for evaluating V(sជ ) and V(rជ ) in Eq. 共7兲 as well for MD simulations were taken from x-ray diffraction studies.19–21 Details of sorbate–

Diffusion anomaly from analytical formula

967

TABLE I. Details of MD simulations on zeolite A and ␳ . N s and N z are the number of sorbate and zeolite atoms, respectively, in the simulation; L is the duration of the MD simulation. System size is the number of crystallographic unit cells taken in each direction. Zeolite system and size

a 共Å兲

Ns

Nz

L(ns)

A (2⫻2⫻2) ␳ 共4⫻4⫻4兲

24.555 14.601

64 64

5120 9280

3 3

sorbate and sorbate–zeolite intermolecular interactions for zeolite A and Y are as given in Ref. 22. The interaction parameters for zeolite ␳ are the same as those for A and Y. The details of MD runs on zeolite A and ␳ are given in Table I. Integrals in Eq. 共7兲 were evaluated at a grid size of 0.1 Å. Calculations at several other grid sizes, viz., 0.2 Å and 0.05 Å indicated that a 0.1 Å grid size was adequate and the integrals in Eq. 共7兲 had indeed converged. The limits of the integrals were chosen as follows: 共i兲 zeolite A: all points within a sphere of 6.0 Å radius centered at the cage center for the volume integral and all points within a circle of radius of 3.0 Å centered at the window center for the surface integral 共ii兲 zeolite Y: the cutoff here are 6.0 Å and 4.5 Å, respectively, for the volume and surface integrals. 共iii兲 Zeolite ␳ : the cutoff here are 7.3 Å and 3.0 Å. Note that the value of the integral does not depend crucially on the cutoff so long as a reasonable estimate of cage or window radius has been chosen. This is because V(rជ ) or V(sជ ) increases rapidly when grid points are close to the wall or overlapping the atoms of the zeolite leading to zero or negligible contribution to the integral. The values of D so obtained are plotted in Figs. 3 and 4. Along with these, points from MD simulations are also shown for comparison. Figure 3 shows the diffusion coefficient as a function of 1/␴ 2 where ␴ is the Lennard-Jones diameter of the diffusing species for zeolites A, Y, and ␳ . It is seen that values obtained from Eq. 共7兲 indicated by 〫 in the figure agrees well with those obtained from the MD calculations indicated by ⫻ in the figure for all the three zeolites. The agreement of the diffusion coefficient is not only good for small ␴ but also for large ␴ where D peaks. Thus the levitation effect or diffusion anomaly is well reproduced by Eq. 共7兲. As stated earlier 共also see Ref. 14兲, the peak in the diffusion coefficient appears as ␥ →1. This is because, at or near ␥ ⬇1, the dimension of the bottleneck is close to the minimum in the sorbate–zeolite interaction in the LennardJones potential. As a result, the energy barrier for passage through the bottleneck is smallest when ␥ ⬇1. Thus, even though the sorbate is now relatively larger in size, it can get past the bottleneck more easily as compared to smaller sized sorbate which encounters a higher energetic barrier. Another way of interpreting the findings is as follows: when the sorbate is small in comparison with the void dimension, the sorbate is strongly attracted from the surface of the host nearest to it. As a result the sorbate is strongly bound within the porous solid at preferred positions which are generally referred to as the physisorption sites. When the levitation parameter approaches unity, the sorbate dimension is comparable to the void dimension. A consequence of this

Downloaded 04 Jan 2005 to 130.102.2.60. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

968

J. Chem. Phys., Vol. 112, No. 2, 8 January 2000

FIG. 3. Diffusion coefficient D vs 1/␴ 2 for Lennard-Jones sorbates in 共a兲 zeolite A at T⫽140 K, 共b兲 zeolite Y at T⫽150 K, and 共c兲 zeolite ␳ at T⫽100 K. Values obtained from expression 共7兲 (〫) are compared with the MD results (⫻). MD results were taken from Ref. 21 for zeolite Y. MD results on zeolites A and ␳ are from the present study 共see Table I兲.

is equal attraction of the sorbate from opposite surfaces of the porous host which essentially cancel each other. The net force on the sorbate is therefore negligible. The sorbate now behaves more like a free particle and less like a bound particle even though it is still confined within the porous media. As this effect might be of considerable importance both from fundamental as well as industrial points of view, it is gratifying that expression 共7兲 is able to reproduce the levitation effect. In order to probe a little further, we wanted to see if the recently reported MD results on temperature dependence of the levitation effect can also be reproduced by Eq. 共7兲. Previous MD simulations on zeolite Y suggest that the peak disappears at higher temperatures.22 Is this disappearance of the peak predicted by Eq. 共7兲? Figure 4共a兲 shows a plot of D against 1/␴ 2 for four different temperatures. It is evident that the temperature dependence is also reproduced by the analytical expression. In particular, Eq. 共7兲 predicts that the peak observed at low temperatures disappears at higher temperatures in agreement with the MD results. Figure 4共b兲 shows the temperature dependence of the levitation effect in zeolite A from Eq. 共7兲. MD results for zeolite A 共details of run are in Table I兲 are also plotted in the same figure. It is evident that the agreement between the values of self-diffusion coefficient from Eq. 共7兲 and the MD results are good. One of the interesting applications of the levitation effect is for separation of mixtures. Normal separation in industries

Sarkar, Kumar, and Yashonath

FIG. 4. D vs 1/␴ 2 for Lennard-Jones sorbates in 共a兲 zeolite Y at four different temperatures (T⫽100K, 150 K, 200 K, 300 K兲 and 共b兲 zeolite A at three different temperatures (T⫽100 K, 140 K, 200 K兲. Values obtained from expression 共7兲 (〫) along with MD results (⫻) for comparison are shown. MD results on zeolite Y are from Ref. 21 while that on zeolite A is from the calculations carried out by us in the present study.

as well as in laboratories seem to employ only the linear regime where the activation energy for the sorbates are large. This necessitates use of higher temperatures. If one were to make use of the anomalous regime reported here, good separation at lower temperatures can be achieved since the sorbate whose dimension lies in the anomalous regime will have lower activation energy. Diffusion anomaly can also achieve selective separation of just one component from the mixture by using a porous solid whose dimension is comparable to the dimension of the component to be separated. It is therefore clear that the analytical expression derived here actually can yield a good estimate of the diffusion coefficient which is otherwise obtained from rather long and tedious MD runs. The effort in the evaluation of expression 共7兲 is at least an order of magnitude less as compared to MD simulations. However, the assumptions that have been made in deriving 共7兲 may not always be satisfied. In particular, the assumption that the bottleneck is narrower than the cage diameter is not valid in all zeolites. Hence, in zeolitic structures where this assumption is not satisfied, e.g., zeolite ZSM-5 or the aluminophosphate VPI-5, expression 共7兲 may not give the correct estimate of the diffusion coefficient. The authors are thankful to the Council of Scientific and Industrial Research, New Delhi for partial support of this research. The authors thank Department of Science and

Downloaded 04 Jan 2005 to 130.102.2.60. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp

J. Chem. Phys., Vol. 112, No. 2, 8 January 2000

Technology for project assistantship to one of us 共A.V.A.K.兲. S.S. thanks the Indian Institute of Science and the Jawaharlal Nehru Centre for Advanced Scientific Research, Bangalore for supporting the visit during which a part of this work was done. Contribution No. 1414 from the Solid State and Structural Chemistry Unit. J. Karger and D. M. Ruthven, Diffusion in Solids 共Wiley, New York, 1992兲. 2 S. P. Bates and R. A. vanSanten, Adv. Catal. 42, 1 共1998兲. 3 P. Demontis and G. B. Suffritti, Chem. Rev. 97, 2845 共1997兲. 4 R. V. Jasra, N. Choudhary, and S. G. T. Bhat, Ind. Eng. Chem. Res. 35, 4221 共1997兲. 5 F. D. Magalhaes, R. L. Laurence, and Wm. Curtis Conner, J. Phys. Chem. B 102, 2317 共1998兲. 6 G. Schrimpf, M. Schlenkrich, J. Brickmann, and P. Boff, J. Phys. Chem. 96, 7404 共1992兲. 7 P. Demontis, G. B. Suffritti, and A. Tilocca, J. Chem. Phys. 105, 5586 共1996兲. 1

Diffusion anomaly from analytical formula

969

K. Hahn, J. Ka¨rger, and V. Kukla, Phys. Rev. Lett. 76, 2762 共1996兲. K. Hahn and J. Ka¨rger, J. Phys. Chem. 100, 316 共1996兲. 10 R. L. June, A. T. Bell, and D. N. Theodorou, J. Phys. Chem. 95, 8866 共1991兲. 11 R. M. Barrer, Zeolites and Clay Minerals as Sorbents and Molecular Sieves 共Academic Press, London, 1978兲. 12 C. Saravanan, F. Jousse, and S. M. Auerbach, Phys. Rev. Lett. 80, 5754 共1998兲. 13 S. Y. Bhide and S. Yashonath, J. Chem. Phys. 111, 1658 共1999兲. 14 S. Yashonath and P. Santikary, J. Phys. Chem. 98, 6368 共1994兲. 15 Sanjoy Bandyopadhyay and S. Yashonath, J. Phys. Chem. 99, 4286 共1995兲. 16 E. G. Derouane, J-M. Andre, and A. A. Lucas, J. Catal. 110, 58 共1988兲. 17 S. Yashonath and P. Santikary, J. Phys. Chem. 97, 3849 共1993兲. 18 F. Reif, Fundmentals of Statistical and Thermal Physics 共McGraw–Hill, New York, 1985兲 p. 273. 19 J. J. Pluth and J. V. Smith, J. Am. Chem. Soc. 102, 4704 共1980兲. 20 A. N. Fitch, H. Jobic, and A. Renouprez, J. Phys. Chem. 90, 1311 共1986兲. 21 J. B. Parise et al., J. Phys. Chem. 88, 2303 共1984兲. 22 S. Yashonath and C. Rajappa, Faraday Discuss. 106, 105 共1997兲. 8 9

Downloaded 04 Jan 2005 to 130.102.2.60. Redistribution subject to AIP license or copyright, see http://jcp.aip.org/jcp/copyright.jsp