- Email: [email protected]
Modelling diffusion in zeolites with cellular automata
Zeolites and Related Materials: Trends, Targets and Challenges Proceedings of 4th International FEZA Conference A. Gédéon, P. Massiani and F. Babonneau (Editors) © 2008 Elsevier B.V. All rights reserved.
701
Modelling diffusion in zeolites with cellular automata Pierfranco Demontis, Federico G. Pazzona, Giuseppe B. Suffritti Università degli Studi di Sassari, via Vienna 2, I-07100 Sassari, Italy
Abstract Exploiting the analogy between the partitioned spatial structure of zeolites and the space-discrete nature of Cellular Automata (CA), we developed a probabilistic Cellular Automaton (CA) able to capture both static and transport properties of diffusing species in zeolites at a coarse-grained level. Our model uses cells to represent pores. It focuses on the sensitivity of each cell upon its instantaneous occupancy to mimic the particleframework and particle-particle interactions of adsorbates into real zeolite pores. It makes use of local partition functions and kinetic barriers to build up a simple and fast evolution rule allowing our model to reproduce data such as adsorption isotherms, global and local distributions of occupancies inside of the zeolite pores, diffusivities, correlations in space and time, etc. from experimental and/or atomistic simulation. The local and parallel nature of our CA, together with its drastically reduced number of degrees of freedom makes it a powerful tool to enlarge the space-time scales of numerical simulations of diffusion in zeolites.
Keywords: Cellular Automata, Zeolites, Diffusion, Adsorption Isotherms
1. Introduction When reduced to its essential constituents, the network of connected channels and cages which build up the three-dimensional framework of a zeolite can be conveniently represented as a set of structured lattice points (cells) exchanging adsorbed particles according to well defined local rules. These are exactly the basic ingredients of Cellular Automata, parallel algorithms providing a space-time discrete environment in which physical systems can be modelled in a reductionistic approach, in order to cover large scales of space and time. [1] We constructed a CA to model intercage diffusion and static properties of simple molecules adsorbed in zeolites. [2-4]
2. The Model The framework of the zeolite of interest is represented in our CA via a network of M cells with the same topology. Along with it, each cell is able to exchange particles with its own first-neighboring cells. The total number of particles, N, and the temperature, T, are held constant in each simulation. The coarse-grained properties of the diffusing species are modeled through the characteristic cell parameters: each cell is made up of a finite number of exit sites Kex (allowing intercell particle transfers) and of inner sites Kin (not allowing such transfers). Their sum, K = Kex + Kin, defines the capacity of each cell due to mutual exclusion among particles. We are assuming the adsorbate dimension to
702
P. Demontis et al.
be such that only one molecule at time can cross the window connecting two adjacent cells. In Fig. 1 a single cell is pictured which represents the automaton representation of a cage of an LTA-type zeolite. As can be seen, the 6 exit sites, which corresponds to the locations near of the windows of the real LTA zeolite, are characterized by a well defined structure (which imprints the topology to the entire lattice as shown in Fig. 2) while the inner sites have not a well-defined arrangement and, independently of their number Kin, their ensemble can be pictured as a sphere. Moreover, Fig. 2 suggests that adjacent exit sites connect adjacent cells. That means that a molecule will be able to migrate from cell to cell (that is, from pore to pore) only through jumps between adjacent exit sites. A molecule adsorbed in a site contributes with the energy Eq.1 :
f α (n) = ε α + φ(n) Figure 1. a CA unit
to the total energy of the cell. In Eq. 1, α = ex, in depending on cell (cubic symmetry). the type of site, εex and εin are the framework parameters, that is they represent the adsorption potential of an exit and an inner site, respectively, and φ (n) is the mean-field interaction that each particle “feels” because of the presence of other n–1 adsorbed molecules in the same cell. This is an essential ingredient of our model. Such a strict locality of interactions (interactions are restricted to each cell) is the heart of the model. As we shall briefly illustrate, it allows to construct a global evolution rule through the composition of local rules applied synchronously on all the cells or on small groups of connected cells in which the lattice can be partitioned. The system evolves in descrete time steps. At each time step, the occupied sites of each cell define its Figure 2. some connected cells: a small 3D portion of the instantaneous configuration, and the partial occupancies (i.e. the number of exit and inner sites occupied, nex and nin) define a level with effective energy Eq.2 : F (n ex , nin ) = n ex f ex (n) + nin f in (n)
The local partition function is defined as: Eq.3 : Q(n ,T ) =
¦
nex + nin = n
§ K ex ·§ K in · − F ( nex ,nin ) / kBT ¸¸¨¨ ¨¨ ¸¸e © n ex ¹© nin ¹
where kB is the Boltzmann’s constant and the sum runs over all values of the partial occupancies satisfying a total cell occupancy equal to n. Such a function is computable since the maximum occupancy K is usually not large. The weights of the levels contained in the local partition function play a central role in the evolution rule.
703
Modelling diffusion in zeolites with cellular automata
Normalized Self-Diffusivity
The model evolves in time through a randomization-propagation scheme[5]. At each time step the following sequence of operation is applied: 12
1.50 T = 300 K Mean Field Simulation
10 8
T = 1300 K Mean Field Simulation
1.25 1.00
0.8
6
0.75
0.6
4
0.50
0.4
2
0.25
0.2
0 0
2
4
6
8
0.00 10 12 14 16 0
Loading
(a)
2
4
6
8
10 12 14 16
Loading
(b)
T = 2500 K Mean Field Simulation
1.0
0.0 0
2
4
6
8
10 12 14 16
Loading
(c)
Figure 3. diffusion profiles at different temperatures obtained by setting φ(n)= 0 and εex–εin= 10 kJ mol–1, to show an example of the variety of possible kinetic behaviors of the model.
1) A randomization operator changes the internal configuration of all cells without varying their respective occupancies. This operator works on all the cell simultaneously. 2) A propagation operator allows molecules in the exit sites to try to jump to the adjacent exit sites, therefore migrating to the adjacent cell. This operator is blocksynchronous:[6] that is, the lattice is dynamically partitioned into blocks (each block is a pair of adjacent cells), and different partitions are alternated in such a way as to allow each cell to communicate with all its 6 neighboring cells without introducing any extra correlation or memory effect; on each lattice partition, the propagation operator works simultaneously on all blocks. Details about randomization and propagation can be found in Refs. [2-4] Such evolution rule neglect completely the time-scale of internal migration processes, while it focuses only on the inter-cell migration time-scale. Nevertheless, a track of the internal state is kept by means of a differentiation Δf(n) = fex(n) fin(n) between the energy parameters of exit and inner sites (see Eq. 1), which determines the thermodynamic accessibility of the exit sites and therefore rules the tendency of a cell to retain or to release molecules outside. Together with Kex and Kin, the parameters fex(n) and fin(n) determine the shape of the adsorption isotherm. Kinetic effects are modelled in first instance through the modeling of Δf(n), and may be refined with the use of kinetic barriers at the interface between neighboring cells. For simplicity, in the following applications we will assume only thermodynamic barriers to rule the jumps on the lattice.
3. Some Applications Without any interaction assumption, the model exhibits a temperature-dependence both in thermodynamic and kinetic properties. With varying the temperature, the model covers four of the five diffusion profiles found by Kärger et al. through PFG-NMR experiments [7] (Fig. 3a: types V at low loadings and VI at intermediate-high loadings; Fig. 3b: type II; Fig. 3c: type I). With the introduction of the interactions through the adjustable parameter φ(n) it is possible to fit experimental isotherms to obtain some informations about the local behaviour of each cell. In Fig. 4 we reported the results of the fitting of a set of experimental data from Jameson et al. [8] about the adsorption of Xe in CaA zeolite (which is an LTA type, the same of the cells pictured in Figs. 1 and 2, with maximum occupancy K = 8). The resulting adsorption isotherm is shown in Fig.4a,
704
P. Demontis et al.
1.2
-40 -45 -50 -55
Mean Field Simulation
1.0 0.8 Interaction Parameter
Our model (Mean Field) Our model (Simulation) Jameson et al.
Normalized Self-Diffusivity
Chemical Potential / kJ mol–1
-35
0.6 0.4 0.2
0.2 0.0 -0.2 -0.4 1 2 3 4 5 6 7 8
-60
0
1
2
3
4
5
6 –1
Loading / molec cage
7
8
(a)
0.0
Occupancy
0
1
2
3
4
5
6 –1
Loading / molec cage
7
8
(b)
Figure 4. (a) Adsorption isotherm of Xe in CaA zeolite, and (b) resulting diffusion profile and interaction potential (inset) together with a Lennard-Jones 10-12 fit.
in direct comparison with experimental reference data, while in Fig. 4b the resulting diffusion profile is shown. The framework parameters have been set as εex= –44 kJmol–1 and εin= –50 kJ mol–1, and the resulting interaction parameter φ(n) is reported in the inset of Fig. 4b, together with a Lennard-Jones 12-10 fitting curve of the same kind of the one proposed by Ayappa [9]: ª§ σ ·12 § σ ·10 º n ¸ − ¨¨ ¸¸ », with r (n) = Eq.5 : φ LJ (n) = ε LJ «¨¨ (r2 − r1 ) + r1 r n ( ) K «© r (n) ¸¹ » © ¹ LJ ¬ ¼ As can be seen, the fit works well for low occupancies. For higher occupancies, repulsion effects become important. Cages close to saturation tend to modify their effective volume to accommodate the adsorbed molecules and partially compensate the repulsion effects. In this region the fit is not good. However, more efforts should be made to connect in a more rigorous way the observed behaviour and the trend of the interaction parameter.
References [1] B. Chopard, M. Droz, Cellular Automata Modelling of Physical Systems, 1998, Cambridge University Press, Cambridge, England. [2] P. Demontis, F.G. Pazzona, G.B. Suffritti, J. Phys. Chem. B, 110 (2006) 13554 [3] P. Demontis, F.G. Pazzona, G.B. Suffritti, J. Chem. Phys, . 126 (2007) 194709 [4] P. Demontis, F.G. Pazzona, G.B. Suffritti, J. Chem. Phys., 126 (2007) 194710 [5] J.-P. Boon, D. Dab, R. Kapral, A.T. Lawniczak, Phys. Rep, . 273 (1996) 55 [6] T. Toffoli, N. Margolus, Cellular Automata Machines: A New Environment for Modeling, 1998, Cambridge University Press, Cambridge, England. [7] J. Kärger, D.M. Ruthven, Diffusion in Zeolites and Other Micorporous Materials, 1992, Wiley, New York. [8] C. Jameson, A.K. Jameson, R. Gerald II, A.C. de Dios, J. Chem. Phys,. 96 (1990) 1690 [9] K.G. Ayappa, J. Chem. Phys,. 111 (1999) 4736
701
Modelling diffusion in zeolites with cellular automata Pierfranco Demontis, Federico G. Pazzona, Giuseppe B. Suffritti Università degli Studi di Sassari, via Vienna 2, I-07100 Sassari, Italy
Abstract Exploiting the analogy between the partitioned spatial structure of zeolites and the space-discrete nature of Cellular Automata (CA), we developed a probabilistic Cellular Automaton (CA) able to capture both static and transport properties of diffusing species in zeolites at a coarse-grained level. Our model uses cells to represent pores. It focuses on the sensitivity of each cell upon its instantaneous occupancy to mimic the particleframework and particle-particle interactions of adsorbates into real zeolite pores. It makes use of local partition functions and kinetic barriers to build up a simple and fast evolution rule allowing our model to reproduce data such as adsorption isotherms, global and local distributions of occupancies inside of the zeolite pores, diffusivities, correlations in space and time, etc. from experimental and/or atomistic simulation. The local and parallel nature of our CA, together with its drastically reduced number of degrees of freedom makes it a powerful tool to enlarge the space-time scales of numerical simulations of diffusion in zeolites.
Keywords: Cellular Automata, Zeolites, Diffusion, Adsorption Isotherms
1. Introduction When reduced to its essential constituents, the network of connected channels and cages which build up the three-dimensional framework of a zeolite can be conveniently represented as a set of structured lattice points (cells) exchanging adsorbed particles according to well defined local rules. These are exactly the basic ingredients of Cellular Automata, parallel algorithms providing a space-time discrete environment in which physical systems can be modelled in a reductionistic approach, in order to cover large scales of space and time. [1] We constructed a CA to model intercage diffusion and static properties of simple molecules adsorbed in zeolites. [2-4]
2. The Model The framework of the zeolite of interest is represented in our CA via a network of M cells with the same topology. Along with it, each cell is able to exchange particles with its own first-neighboring cells. The total number of particles, N, and the temperature, T, are held constant in each simulation. The coarse-grained properties of the diffusing species are modeled through the characteristic cell parameters: each cell is made up of a finite number of exit sites Kex (allowing intercell particle transfers) and of inner sites Kin (not allowing such transfers). Their sum, K = Kex + Kin, defines the capacity of each cell due to mutual exclusion among particles. We are assuming the adsorbate dimension to
702
P. Demontis et al.
be such that only one molecule at time can cross the window connecting two adjacent cells. In Fig. 1 a single cell is pictured which represents the automaton representation of a cage of an LTA-type zeolite. As can be seen, the 6 exit sites, which corresponds to the locations near of the windows of the real LTA zeolite, are characterized by a well defined structure (which imprints the topology to the entire lattice as shown in Fig. 2) while the inner sites have not a well-defined arrangement and, independently of their number Kin, their ensemble can be pictured as a sphere. Moreover, Fig. 2 suggests that adjacent exit sites connect adjacent cells. That means that a molecule will be able to migrate from cell to cell (that is, from pore to pore) only through jumps between adjacent exit sites. A molecule adsorbed in a site contributes with the energy Eq.1 :
f α (n) = ε α + φ(n) Figure 1. a CA unit
to the total energy of the cell. In Eq. 1, α = ex, in depending on cell (cubic symmetry). the type of site, εex and εin are the framework parameters, that is they represent the adsorption potential of an exit and an inner site, respectively, and φ (n) is the mean-field interaction that each particle “feels” because of the presence of other n–1 adsorbed molecules in the same cell. This is an essential ingredient of our model. Such a strict locality of interactions (interactions are restricted to each cell) is the heart of the model. As we shall briefly illustrate, it allows to construct a global evolution rule through the composition of local rules applied synchronously on all the cells or on small groups of connected cells in which the lattice can be partitioned. The system evolves in descrete time steps. At each time step, the occupied sites of each cell define its Figure 2. some connected cells: a small 3D portion of the instantaneous configuration, and the partial occupancies (i.e. the number of exit and inner sites occupied, nex and nin) define a level with effective energy Eq.2 : F (n ex , nin ) = n ex f ex (n) + nin f in (n)
The local partition function is defined as: Eq.3 : Q(n ,T ) =
¦
nex + nin = n
§ K ex ·§ K in · − F ( nex ,nin ) / kBT ¸¸¨¨ ¨¨ ¸¸e © n ex ¹© nin ¹
where kB is the Boltzmann’s constant and the sum runs over all values of the partial occupancies satisfying a total cell occupancy equal to n. Such a function is computable since the maximum occupancy K is usually not large. The weights of the levels contained in the local partition function play a central role in the evolution rule.
703
Modelling diffusion in zeolites with cellular automata
Normalized Self-Diffusivity
The model evolves in time through a randomization-propagation scheme[5]. At each time step the following sequence of operation is applied: 12
1.50 T = 300 K Mean Field Simulation
10 8
T = 1300 K Mean Field Simulation
1.25 1.00
0.8
6
0.75
0.6
4
0.50
0.4
2
0.25
0.2
0 0
2
4
6
8
0.00 10 12 14 16 0
Loading
(a)
2
4
6
8
10 12 14 16
Loading
(b)
T = 2500 K Mean Field Simulation
1.0
0.0 0
2
4
6
8
10 12 14 16
Loading
(c)
Figure 3. diffusion profiles at different temperatures obtained by setting φ(n)= 0 and εex–εin= 10 kJ mol–1, to show an example of the variety of possible kinetic behaviors of the model.
1) A randomization operator changes the internal configuration of all cells without varying their respective occupancies. This operator works on all the cell simultaneously. 2) A propagation operator allows molecules in the exit sites to try to jump to the adjacent exit sites, therefore migrating to the adjacent cell. This operator is blocksynchronous:[6] that is, the lattice is dynamically partitioned into blocks (each block is a pair of adjacent cells), and different partitions are alternated in such a way as to allow each cell to communicate with all its 6 neighboring cells without introducing any extra correlation or memory effect; on each lattice partition, the propagation operator works simultaneously on all blocks. Details about randomization and propagation can be found in Refs. [2-4] Such evolution rule neglect completely the time-scale of internal migration processes, while it focuses only on the inter-cell migration time-scale. Nevertheless, a track of the internal state is kept by means of a differentiation Δf(n) = fex(n) fin(n) between the energy parameters of exit and inner sites (see Eq. 1), which determines the thermodynamic accessibility of the exit sites and therefore rules the tendency of a cell to retain or to release molecules outside. Together with Kex and Kin, the parameters fex(n) and fin(n) determine the shape of the adsorption isotherm. Kinetic effects are modelled in first instance through the modeling of Δf(n), and may be refined with the use of kinetic barriers at the interface between neighboring cells. For simplicity, in the following applications we will assume only thermodynamic barriers to rule the jumps on the lattice.
3. Some Applications Without any interaction assumption, the model exhibits a temperature-dependence both in thermodynamic and kinetic properties. With varying the temperature, the model covers four of the five diffusion profiles found by Kärger et al. through PFG-NMR experiments [7] (Fig. 3a: types V at low loadings and VI at intermediate-high loadings; Fig. 3b: type II; Fig. 3c: type I). With the introduction of the interactions through the adjustable parameter φ(n) it is possible to fit experimental isotherms to obtain some informations about the local behaviour of each cell. In Fig. 4 we reported the results of the fitting of a set of experimental data from Jameson et al. [8] about the adsorption of Xe in CaA zeolite (which is an LTA type, the same of the cells pictured in Figs. 1 and 2, with maximum occupancy K = 8). The resulting adsorption isotherm is shown in Fig.4a,
704
P. Demontis et al.
1.2
-40 -45 -50 -55
Mean Field Simulation
1.0 0.8 Interaction Parameter
Our model (Mean Field) Our model (Simulation) Jameson et al.
Normalized Self-Diffusivity
Chemical Potential / kJ mol–1
-35
0.6 0.4 0.2
0.2 0.0 -0.2 -0.4 1 2 3 4 5 6 7 8
-60
0
1
2
3
4
5
6 –1
Loading / molec cage
7
8
(a)
0.0
Occupancy
0
1
2
3
4
5
6 –1
Loading / molec cage
7
8
(b)
Figure 4. (a) Adsorption isotherm of Xe in CaA zeolite, and (b) resulting diffusion profile and interaction potential (inset) together with a Lennard-Jones 10-12 fit.
in direct comparison with experimental reference data, while in Fig. 4b the resulting diffusion profile is shown. The framework parameters have been set as εex= –44 kJmol–1 and εin= –50 kJ mol–1, and the resulting interaction parameter φ(n) is reported in the inset of Fig. 4b, together with a Lennard-Jones 12-10 fitting curve of the same kind of the one proposed by Ayappa [9]: ª§ σ ·12 § σ ·10 º n ¸ − ¨¨ ¸¸ », with r (n) = Eq.5 : φ LJ (n) = ε LJ «¨¨ (r2 − r1 ) + r1 r n ( ) K «© r (n) ¸¹ » © ¹ LJ ¬ ¼ As can be seen, the fit works well for low occupancies. For higher occupancies, repulsion effects become important. Cages close to saturation tend to modify their effective volume to accommodate the adsorbed molecules and partially compensate the repulsion effects. In this region the fit is not good. However, more efforts should be made to connect in a more rigorous way the observed behaviour and the trend of the interaction parameter.
References [1] B. Chopard, M. Droz, Cellular Automata Modelling of Physical Systems, 1998, Cambridge University Press, Cambridge, England. [2] P. Demontis, F.G. Pazzona, G.B. Suffritti, J. Phys. Chem. B, 110 (2006) 13554 [3] P. Demontis, F.G. Pazzona, G.B. Suffritti, J. Chem. Phys, . 126 (2007) 194709 [4] P. Demontis, F.G. Pazzona, G.B. Suffritti, J. Chem. Phys., 126 (2007) 194710 [5] J.-P. Boon, D. Dab, R. Kapral, A.T. Lawniczak, Phys. Rep, . 273 (1996) 55 [6] T. Toffoli, N. Margolus, Cellular Automata Machines: A New Environment for Modeling, 1998, Cambridge University Press, Cambridge, England. [7] J. Kärger, D.M. Ruthven, Diffusion in Zeolites and Other Micorporous Materials, 1992, Wiley, New York. [8] C. Jameson, A.K. Jameson, R. Gerald II, A.C. de Dios, J. Chem. Phys,. 96 (1990) 1690 [9] K.G. Ayappa, J. Chem. Phys,. 111 (1999) 4736