pubs.acs.org/Langmuir © 2011 American Chemical Society
Atomistic Simulations of CO2 and N2 within Cage-Type Silica Zeolites Lindsey Madison, Henry Heitzer, Colin Russell, and Daniela Kohen* Chemistry Department, Carleton College, Northfield, Minnesota 55057, United States Received October 21, 2010. Revised Manuscript Received December 9, 2010 The behavior of CO2 and N2, both as single components and as binary mixtures, in two cage-type silica zeolites was studied using atomistic simulations. The zeolites considered, ITQ-3 and paradigm cage-type zeolite ZK4 (the all-silica analog of LTA), were chosen so that the principles illustrated can be generalized to other adsorbent/adsorbate systems with similar topology and types of interactions. N2 was chosen both because of the potential uses of N2/CO2 separations and because it differs from CO2 most significantly in the magnitude of its Coulombic interactions with zeolites. Despite similarities between N2 and CO2 diffusion in other materials, we show here that the diffusion of CO2 within cage-type zeolites is dominated by an energy barrier to diffusion located at the entrance to the narrow channels connecting larger cages. This barrier originates in Coulombic interactions between zeolites and CO2’s quadrupole and results in welldefined orientations for the diffusing molecules. Furthermore, CO2’s favorable electrostatic interactions with the zeolite framework result in preferential binding in the windows between cages. N2’s behavior, in contrast, is more consistent with that of molecules previously studied. Our analysis suggests that CO2’s behavior might be common for adsorbates with quadrupoles that interact strongly with a material that has narrow windows between cages.
1. Introduction Microporous adsorbents that can be used for CO2 capture have received a lot of attention given their potential use in reducing the emission of greenhouse gases.1 In particular, zeolites and metalorganic framework materials seem to be promising candidates that have potential applications in separation processes, catalysis, and gas storage.2-4 Characterization of behavior on the molecular level is essential to a rational search for materials that are well suited for a given application and may also facilitate process design. Furthermore, recent efforts at screening large numbers of materials3-6 are rooted in our detailed knowledge of diffusion and adsorption within nanometer-scale pores. Given that diffusion and adsorption control performance of materials in many applications,7-9 understanding how carbon dioxide and other gases’ behavior is determined by the topology of the pore may prove critical to selecting appropriate materials for specific applications. Molecular simulations have played a useful role in approaching this issue because it is possible to use these methods to probe the properties of diffusing molecules on the atomic length scales that define the differences between materials.10-14 A particularly important use of molecular simula(1) D’alessandro, D. M.; SmitB.; Long, J. R. Angew. Chem., Int. Ed. 2010, 49, 6058. (2) Choi, S.; Drese, J. H.; Jones, C. W. ChemSusChem 2009, 2, 796. (3) Keskin, S.; Sholl, D. S. Langmuir 2009, 25, 11786. (4) Yazaydin, A. O.; Snurr, R. Q.; Park, T. H.; Koh, K.; Liu, J.; LeVan, M. D.; Benin, A. I.; Jakubczak, P.; Lanuza, M.; Galloway, D. B.; Low, J. J.; Willis, R. R. J. Am. Chem. Soc. 2009, 131, 18198. (5) Krishna, R.; van Baten, J. M. Chem. Eng. J. 2007, 133, 121. (6) Liu, B.; Smit, B. Langmuir 2009, 25, 5918. (7) Yang, R. T. Gas Separation by Adsorption Processes; Imperial College Press: London, 1997. (8) K€arger, J.; Ruthven, D. M. Diffusion in Zeolites and Other Microporous Solids; Wiley: New York, 1992. (9) Keil, F. J.; Krishna, R.; Coppens, M. O. Rev. Chem. Eng. 2000, 16, 71. (10) Heuchel, M.; Snurr, R. Q.; Buss, E. Langmuir 1997, 13, 6795. (11) Yang, J. H.; Clark, L. A.; Ray, G. J.; Kim, Y. J.; Du, H.; Snurr, R. Q. J. Phys. Chem. B 2001, 105, 4698. (12) Krishna, R.; Vlugt, T. J. H.; Smit, B. Chem. Eng. Sci. 1999, 54, 1751. (13) Dubbeldam, D.; Calero, S.; Vlugt, T. J. H.; Krishna, R.; Maesen, T. L. M.; Beerdsen, E.; Smit, B. Phys. Rev. Lett. 2004, 93. (14) Smit, B.; Maesen, T. L. M. Chem. Rev. 2008, 108, 4125.
1954 DOI: 10.1021/la104245c
tions is to examine the properties of adsorbed mixtures.15-17 Performing and interpreting experiments to characterize mixture diffusion in microporous materials in a detailed way is challenging. In molecular simulation, in contrast, treating the mixture diffusion is only slightly more involved than treating species adsorbed as single components, at least for mixtures with two or three components. In particular, the diffusion of CO2 and other gases in cage-type materials possessing pores that are commensurate in size with the adsorbate are being increasingly studied because they present peculiar behaviors that might increase their promise as separation materials.5,18,19 In this article, we use molecular simulations to examine the behavior of CO2 and N2 carefully in cage-type pure-silica zeolites. We have previously described the adsorption and diffusion properties of CO2 and N2 in three all-silica zeolites by assessing the influence of connectivity and pore shape on their behavior. Silicalite, ITQ-7, and ITQ-320,21 have identical chemical compositions but different pore structures: silicalite and ITQ-7 have a 3D network of intersecting channels, and ITQ-3 has a cagelike structure. In these three materials, CO2 diffuses more slowly than N2 but otherwise the behavior of these gases within ITQ-7 and silicalite is quite different than within ITQ-3. In ITQ-7 and silicalite, the loading dependence of diffusion is very similar for CO2 and N2, the apparent activation energies for the diffusion of each adsorbate are similar in the two materials, and the diffusion properties of adsorbed mixtures can be easily understood. In contrast, none of these are true within ITQ-3. Our analysis suggested21 that this behavior might be common for adsorbates that interact strongly with the walls of the narrow passages that connect large cavities in cage-type zeolites, as is the case for CO2 (15) Sholl, D. S. Acc. Chem. Res. 2006, 39, 403. (16) Sanborn, M. J.; Snurr, R. Q. AIChE J. 2001, 47, 2032. (17) Sanborn, M. J.; Snurr, R. Q. Sep. Purif. Technol. 2000, 20, 1. (18) Jee, S. E.; Sholl, D. S. J. Am. Chem. Soc. 2009, 131, 7896. (19) Krishna, R.; van Baten, J. M. Sep. Purif. Technol. 2008, 61, 414. (20) Goj, A.; Sholl, D. S.; Akten, E. D.; Kohen, D. J. Phys. Chem. B 2002, 106, 8367. (21) Selassie, D.; Davis, D.; Dahlin, J.; Feise, E.; Haman, G.; Sholl, D. S.; Kohen, D. J. Phys. Chem. C 2008, 112, 16521.
Published on Web 01/13/2011
Langmuir 2011, 27(5), 1954–1963
Madison et al.
Figure 1. Structure (top) and schematic illustration of the topology (bottom) of ZK4 and ITQ-3. For ITQ-3, both pictures display two sets of three interconnected cages. For ZK4, in the top image an R cage is in the center surrounded by β cages and in the bottom image four R cages are shown.
but not N2. To investigate this finding, we have chosen to continue our study of the roots of CO2’s and N2’s behavior within ITQ-3 and also to examine their behavior within the paradigm cage-type zeolite, ZK4 (the all-silica analog of LTA). ZK4’s simple structure and chemical composition offer hope that the principles illustrated here can be generalized to other adsorbent/ adsorbate systems with similar topology and type of interactions. Thus, the work presented here complements earlier studies of the adsorption and diffusion of CO2 and N2 in silica zeolites20,21 but focuses on the behavior of adsorbates in two siliceous zeolites with narrow windows connecting large cages, ZK4 and ITQ-3. Figure 1 shows the structures and schematic illustrations of the pore topology of these two zeolites. ZK4 is the all-silica analogue of the Linde type A (LTA) zeolite, and it has been extensively studied with molecular simulations (e.g., refs 22-25). The structure of LTA is well known; it has a 3D channel system composed of large cavities (R cages) that are about 10 A˚ in diameter. Each R cage is connected to six other R cages by eight-ring pores of about 4.1 A˚ in diameter. The ZK4 unit cell is cubic, its length is 24.555 A˚,24 and it contains eight smaller β cages that surround an R cage. ITQ-3 (structure type ITE), however, has a 1D pore system with small channels of about 4 A˚ in diameter made up of eightmembered rings that open to larger cavities. A second straight channel runs through the material but is too narrow to accommodate guest molecules. Figure 2 shows pores’ and adsorbates’ dimensions; the fact that these are commensurate in size underscores the importance of the molecular characterization of behavior within these systems. The Figure also shows that the biggest difference between the two adsorbates studied here is the magnitude of their quadrupole moment. (22) Fritzsche, S.; Haberlandt, R.; Hofmann, G.; Karger, J.; Heinzinger, K.; Wolfsberg, M. Chem. Phys. Lett. 1997, 265, 253. (23) Schuring, A.; Auerbach, S. M.; Fritzsche, S.; Haberlandt, R. J. Chem. Phys. 2002, 116, 10890. (24) Demontis, P.; Fenu, L. A.; Suffritti, G. B. J. Phys. Chem. B 2005, 109, 18081. (25) Beerdsen, E.; Dubbeldam, D.; Smit, B. J. Phys. Chem. B 2006, 110, 22754.
Langmuir 2011, 27(5), 1954–1963
Article
Figure 2. (a) Window dimensions (in angstroms).22(b) Cartoon showing approximate molecular dimensions (in angstroms) estimated using the Lennard-Jones parameters and bond lengths. Experimental quadrupole moments are also displayed.
The remainder of the article is structured as follows. In section 2, we describe our computational methods. Section 3 presents our Results and Discussion. Section 3.1 details our single-component and mixture diffusion results, and our predictions for binary selfdiffusion coefficients from single-component data are presented in section 3.2. These two sections set the stage for section 3.3, where we examine the molecular causes of the unusual behavior within cage-type materials. We conclude in section 4.
2. Methods and Models The atomistic models and simulation methods used in this article are those described and validated in our previous work.20,21 They are summarized below. Interactions between adsorbed molecules and the zeolite framework atoms are modeled following the method of Makrodimitris et al.26 Each O atom in the zeolite is assumed to interact with each site on the adsorbed molecules through both a Lennard-Jones potential and electrostatic interactions, but only electrostatic interactions are considered between adsorbates and Si atoms in the zeolite. We have used the parameters reported by Makrodimitris et al.26 as model LJCB.JBTLC because this model is successful in reproducing experimental isotherms for CO2 and N2 on silicalite. An important assumption in our work is that potential parameters can be transferred directly to model the same molecules adsorbed in ITQ-3 and ZK4. This approach seems plausible because each material has precisely the same composition, but given the absence of experimental data in these two materials, we are unable to assess the validity of this assumption directly. Zeolites are relatively stiff; therefore, as is the case in most atomistic studies of adsorption in zeolites, we assume that the adsorbent can be treated as a rigid framework with atomic positions fixed by the observed crystallographic data.27 In the case of ZK4, we also assume the presence of an extra atom in the center of the smaller β cage28 to block these experimentally inaccessible regions of the zeolite. (26) Makrodimitris, K.; Papadopoulos, G. K.; Theodorou, D. N. J. Phys. Chem. B 2001, 105, 777. (27) Fuchs, A. H.; Cheetham, A. K. J. Phys. Chem. B 2001, 105, 7375. (28) Bates, S. P.; van Well, W. J. M.; van Santen, R. A.; Smit, B. J. Am. Chem. Soc. 1996, 118, 6753.
DOI: 10.1021/la104245c
1955
Article
Madison et al.
Figure 4. Single-component adsorption isotherms for CO2 (filled symbols) and N2 (unfilled symbols) at 298 K in ZK4 (diamonds) and ITQ-3 (triangles) as computed from GCMC. Table 1. Diffusion Coefficients at Infinite Dilution
ZK4 ITQ-3
Figure 3. Diffusion coefficients computed from our single-component EMD calculations for CO2 (filled symbols) and N2 (unfilled symbols) in ZK4 (top) and ITQ-3 (bottom) at 298 K. Squares correspond to self-diffusion coefficients, and circles correspond to corrected diffusion coefficients. Carbon dioxide was modeled as a linear triatomic, and nitrogen was modeled as a diatomic molecule with fixed bond lengths and bond angles. Both species were modeled by charged LennardJones (LJ) centers using the TraPPE force field developed by Potoff and Siepmann.29 These potentials quantitatively reproduce the vapor-liquid equilibria of the neat systems and their mixtures. LJ parameters for the unlike-pair interactions are calculated with the Lorentz-Berthelot combining rules.30 Host/guest and guest/ guest interactions were modeled by the sum of short-range interaction terms using Lennard-Jones potentials and Coulombic terms accounting for the long-range interaction between the quadrupole moment of the guests and the zeolite’s electrostatic field. Equilibrium molecular dynamics (EMD) simulations were performed using methods that are very similar to those described by Skoulidas and Sholl.31,32 Well-developed methods were used to examine both self-diffusion and the diffusion coefficients that are relevant to net mass transfer. Self-diffusion coefficients (Dis) and corrected diffusion coefficients (Dio) for each adsorbed component i were calculated from trajectories obtained using EMD and the corresponding Einstein relationships.15 Each EMD simulation started with a grand canonical Monte Carlo (GCMC) pre-equilibration (2.5 106 Monte Carlo moves) (29) Potoff, J. J.; Siepmann, J. I. AIChE J. 2001, 47, 1676. (30) Allen, M. P.; Tildesley, D. J. Computer Simulations of Liquids; Oxford Science Publications: Oxford, UK, 1994. (31) Skoulidas, A. I.; Sholl, D. S. J. Phys. Chem. B 2001, 105, 3151. (32) Skoulidas, A. I.; Sholl, D. S. J. Phys. Chem. B 2002, 106, 5058.
1956 DOI: 10.1021/la104245c
DsN2(0) (10-8 m2/s)
DsCO2(0) (10-8 m2/s)
1.04 0.188
0.52 0.05
followed by EMD equilibration (2.5 105 MD steps). After equilibration, production runs of 2 106 MD steps were performed using a Nose-Hoover thermostat to keep the temperature constant. The time step in all simulations was 1 fs. Computing corrected diffusion coefficients require averaging over multiple calculations; therefore, at least 10 independent runs were performed for each of the results shown. In addition, to ensure equilibration and convergence in our EMD simulations, the number of runs was increased in those calculations with a small number of molecules per unit cell. In all cases reported here, it was confirmed that, after equilibration, the observed mean square displacements were linear in time and thus the systems were experiencing normal diffusion. Single-component and mixture adsorption simulations were carried out using GCMC.27,33 In our GCMC simulations, four types of trial moves were used: attempts to insert a new molecule in the zeolite, attempts to delete an existing molecule from the zeolite, and attempts to translate or rotate an adsorbed molecule. Typically, simulations were run for about 5 106 Monte Carlo moves, and the first half of these moves were used for equilibration and were not included in the averaging. The chemical potential of the bulk gas phase in our single-component calculations was related to the gas pressure using the ideal gas equation of state.
3. Results and Discussion 3.1. Single-Component and Mixture Diffusion. The diffusion coefficients computed from our single-component EMD calculations for CO2 and N2 in ZK4 at 298 K are shown in Figure 3. Each set of diffusion coefficients is shown as a function of the fractional loading, θ, of the adsorbed species. θ is calculated using saturation loading values derived from Langmuir fits to the single-component adsorption isotherms shown in Figure 4. Table 1 lists the diffusion coefficients at infinite dilution. Figure 3 shows that DsCO2(θ) decays much faster than DsN2(θ) within both ZK4 and ITQ-3. This is in contrast to behavior seen in materials with intersecting channels, such as silicalite and ITQ-7, where DsCO2 (θ) and DsN2(θ) decay at very similar rates.21 As with (33) Frenkel, D.; Smit, B. Understanding Molecular Simulation: From Algorithms to Applications; Academic Press: London, 1996.
Langmuir 2011, 27(5), 1954–1963
Madison et al.
Article
Table 2. Apparent Activation Energy at Infinite Dilution
ZK4 ITQ-3
ZK4 ITQ-3
EaN2 (kJ/mol)
EaCO2 (kJ/mol)
0.07 1.1
8.2 9.8
EaN2 (K)
EaCO2 (K)
9 138
990 1174
Figure 5. Correlation between apparent activation energy and diffusivity. Results are presented for CO2 (filled symbols) and N2 (open symbols) on ITQ-3 (blue triangles), ZK4 (orange diamonds), silicalite (black circles), and ITQ-7 (green squares). The black line corresponds to a linear fit to the data for CO2 in all materials but for N2 only in silicalite and ITQ-7.
the self-diffusion coefficient, the behavior of the corrected diffusion coefficients in ZK4 is similar to that of ITQ-3 and quite different from that of the other two materials. The loading dependence of the corrected diffusivity of both CO2 and N2 is very close to the self-diffusion coefficients in ZK4 and ITQ-3, but not in silicalite and ITQ-7. This is especially true for CO2, for which these are essentially identical. It is often assumed that the apparent activation energies of selfdiffusion correlate well with diffusion rates. In our previous work, it was noticed that this does not seem to be the case within ITQ-3. To investigate this further, apparent activation energies were calculated within ZK4 by computing Ds(0) at 200, 298, 350, and 400 K and then fitting these results to the usual Arrhenius equation. The resulting activation energies are listed in Table 2. Figure 5 shows e-Ea/kT versus Ds(0) for ZK4, ITQ-3, ITQ-7, and silicalite.21 Note how Ds(0) correlates loosely with the apparent activation energies for CO2 in all materials, but for N2, this correlation is seen only in silicalite and ITQ-7. The linear fit to all data except for N2 in ZK4 and ITQ-3 underscores this fact. (The magnitude of Ea is not able to account for the low diffusivity of N2 in cage-type zeolites.) This is another indication in our results that the behavior within materials with cages is quite interesting. We examine the roots of this behavior in section 3.3. Mixture simulations allow us to study how DsN2 depends on CO2 loading and vice versa. Figure 6 shows the dependence of the computed self-diffusion coefficients of CO2 and N2 for both pure components and equimolar mixtures as a function of the total fractional loading. In both ZK4 and ITQ-3, N2 diffuses much more slowly in the mixture than as a pure species, whereas CO2 diffuses almost as fast in the mixture as a pure species (making comparisons at the same total loading). That is, in these materials, Langmuir 2011, 27(5), 1954–1963
Figure 6. Self-diffusion coefficients of CO2 (red) and N2 (blue) computed from our EMD calculations for both the pure components and equimolar mixtures in ZK4 (top) and ITQ-3 (bottom) at 298 K versus the total fractional loading. Filled symbols indicate pure-component calculations, and unfilled symbols denote equimolar mixtures.
CO2 molecules hinder diffusion more than N2 molecules do; later in this article, it will be shown that this is a result of segregated adsorption. This is not true in the case of silicalite and ITQ-7, but it has been reported19,34 for CO2 in DDR and other zeolites with cages and narrow pores. 3.2. Predicting Binary Self-Diffusion Coefficients from Single-Component Data. One useful way to analyze the results for mixture diffusion is to ask whether the results for adsorbed mixtures can be predicted using only single-component information. Krishna and Paschek35 developed a useful correlation that predicts mixture self-diffusivities from single-component information. This correlation has been shown to give accurate results for mixture diffusion in silicalite and ITQ-7,21 carbon nanotubes, and metal-organic framework materials.36-38 Figure 7 presents a scatter plot of the predictions from Krishna and Paschek’s correlation for all loadings/compositions for N2 and CO2 for ITQ-3, ITQ-7, silicalite,21 and ZK4 versus the results obtained in the MD simulations. It is apparent in Figure 7 that Krishna and Paschek’s formulation predicts the behavior of the mixture relatively accurately in ITQ-7 and silicalite but not in (34) (35) (36) (37) (38)
Krishna, R.; van Baten, J. M. Chem. Phys. Lett. 2007, 446, 344. Krishna, R.; Paschek, D. Phys. Chem. Chem. Phys. 2002, 4, 1891. Krishna, R.; van Baten, J. M. Ind. Eng. Chem. Res. 2006, 45, 2084. Skoulidas, A. I.; Sholl, D. S.; Krishna, R. Langmuir 2003, 19, 7977. Keskin, S.; Liu, J. C.; Johnson, J. K.; Sholl, D. S. Langmuir 2008, 24, 8254.
DOI: 10.1021/la104245c
1957
Article
Madison et al.
Figure 7. Parity plot between normalized self-diffusion results from the binary MD calculations and those predicted using Krishna and Paschek’s formulation. Results are presented for CO2 (filled symbols) and N2 (unfilled symbols) on silicalite (black circles), ITQ-3 (blue triangles), ITQ-7 (green squares), and ZK4 (orange diamonds).
ITQ-3 and ZK4. Within the context of this formulation, the failure can be traced back21 to the fact that single-component MD simulations showed that DsCO2(θ) and DoCO2(θ) were remarkably close to each other. This is the case within both cage-type materials (Figure 7) and indicates that molecular jumps are minimally correlated.35,37 Similar DsCO2(θ) and DoCO2(θ) values result in very large CO2/N2 exchange coefficients, which implies vanishing correlations between CO2 and N2 molecules. However, this is not the behavior exhibited by these adsorbates within ZK4 and ITQ-3 where CO2 blocks N2, so it is not surprising that the formulation is not as accurate in these materials as it is in ITQ-7, silicalite, and many others.36-38 An underlying assumption of Krishna and Paschek’s formulation is that the sites occupied by the adsorbates in the mixture are equivalent. Moreover, the formulation is not designed to take into account situations in which molecules of one species hinder the diffusion of molecules of the other species in an asymmetric way. That CO2 and N2 within both cage-type materials studied do not occupy the same sites and that there is segregated adsorption occurring within these materials will become clear in the next section. 3.3. Molecular Understanding of CO2 and N2 Diffusion in Materials with Cage-Type Zeolites. In this section, preferred sites of adsorption for each species and the free energy and potential energy profiles for diffusion that these molecules experience within each material will be investigated. This will clarify the origin of the following behaviors within cage-type materials ITQ-3 and ZK4: (a) the pure CO2 and N2 self-diffusion loading dependence that differ significantly from one another; (b) the apparent activation energy of single-component N2 self-diffusion that does not correlate well with its self-diffusion coefficient; (c) the presence of CO2 in a mixture slows the diffusion of both CO2 and N2 more than the presence of N2; and (d) the poor performance of a correlation for predicting mixture self-diffusion that works well for many other systems. Note that none of these are true in silicalite or ITQ-7, both of which have a 3D network of intersecting channels. Furthermore, the role of adsorbate rotations will also be investigated because it supplements our understanding of linear-molecule behavior within cage-type zeolites. Behaviors a-d ultimately originate in the fact that within the cage-type materials investigated here CO2 and N2 adsorb in sites that differ qualitatively. Figure 8 shows the probability for CO2 1958 DOI: 10.1021/la104245c
Figure 8. CO2 center of mass (in red) and N2 center of mass (in blue) probabilities along relevant coordinates at 298 K at infinite dilution. Shaded regions correspond to narrow channels between cages.
and N2 at 298 K and infinite dilution along the relevant 1D coordinate. These probabilities were obtained by integrating 3D normalized histograms of the locations of molecular centers of mass collected during EMD simulations, and equivalent results are obtained if the data is collected in a GCMC simulation. The integrations were performed over a series of planes perpendicular to a reaction coordinate representative of diffusion in each material. In ITQ-3, this trivially corresponds to the y direction, and for ZK4, we chose to use a radial coordinate that originates at the center of the R cage. Note that in the later case the plotted probabilities were obtained by dividing the number that resulted after integration by the volume of the concentric shell at the corresponding radius. Figure 8 shows that the preferred sites are clearly not the same for CO2 and N2 within each material. Moreover, for CO2 in both ITQ-3 and ZK4 there is a preferred site of adsorption within the narrow channel connecting the cages. In contrast, there is a lack of N2 probability at those locations. Segregated adsorption such as that described here has been reported not only in our previous work within ITQ-3 but also in other zeolites with narrow windows connecting larger cages19,34 or in zeolites with independent pore networks.39 This suggests that there are a number of zeolite structures for which the properties we have seen here may be relevant. In what follows, we explore the reasons and consequences of this behavior. (39) Sant, M.; Leyssale, J. M.; Papadopoulos, G. K.; Theodorou, D. N. J. Phys. Chem. B 2009, 113, 13761.
Langmuir 2011, 27(5), 1954–1963
Madison et al.
Article
Figure 10. Free-energy profiles within ZK4 for the CO2 center of mass (in red) and the N2 center of mass (in blue) at 298 K.
Figure 9. CO2 center of mass (in red) and N2 center of mass (in blue) probability along relevant coordinates at 298 K at infinite dilution when only dispersive forces are included. Shaded regions correspond to narrow channels between cages.
The plots in Figure 8 reveal the reason underlying the poor performance of the correlation approach examined above (Figure 7) for predicting mixture self-diffusion. This correlation is based on the assumption that both adsorbing species compete for similar sites in the adsorbent and that mixing on these sites occurs in a simple way. For N2/CO2 mixtures in these cage-type materials, this assumption is clearly not correct because the two species prefer quite different adsorption sites. Examining the distribution of molecules also reveals why CO2 molecules hinder the motion of N2 molecules more than the reverse: by sitting in the narrow passage between cages, CO2 blocks the motion of other molecules. The presence of CO2 hinders the ability of N2 to cross the narrow necks and thus to diffuse from cage to cage. It is illuminating to examine the relative importance of dispersive versus electrostatic forces in determining the behavior of these systems. To probe this issue, we obtained probability plots like the ones shown but with partial charges on each diffusing molecule set artificially to zero. That is, we performed EMD calculations using the same potential parameters as for Figure 8 but including only dispersion interactions. The results are shown in Figure 9. The change in the behavior of nitrogen in either material is minimal, consistent with the small quadrupole moment of this molecule. However, the change in the behavior of carbon dioxide is substantial. In the absence of Coulombic interactions, N2 and CO2 behave almost identically. Given that the biggest difference between these molecules (Figure 2) is the magnitude of their quadrupole moment, this result is not surprising. Langmuir 2011, 27(5), 1954–1963
However, it reinforces the notion that when the interactions of the adsorbate with the zeolite are strong, this neck between cages is no longer a transition state for diffusion as it is often assumed but rather the location of a preferred site. Note that when studying the diffusion of methane in LTA, Fritzsche and co-workers22 found that CH4’s preferred locations were either in the channel or in front of the channel depending on the strength of the parameters used to describe the interactions. In the CH4 case, however, the change in potential parameters did not involve Coulombic interactions. The fact that CO2’s large quadrupole results in its distinctive behavior compared to that of other similarly sized molecules within cage-type zeolites underscores the limits of using notions such as the degree of confinement. In this context, the degree of confinement is typically defined as the ratio of a characteristic size of a molecule with respect to a characteristic pore, or window, of a zeolite.40 The degree of confinement has proven to be a useful metric to gain insight into the microscopic mechanisms of diffusion within porous materials,37,40 but here we show that CO2’s behavior is dictated by the size of its quadrupole rather than by its dimensions. Figures 8 and 9 also underscore the fact that the topology of ITQ-3 is more complex than that of ZK4. ZK4 has symmetrical cages, but ITQ-3 does not. This was one of the reasons behind choosing ZK4 as a model material for our studies. In the rest of this article, plots will be shown only for ZK4 with the understanding that similar conclusions can be drawn from the results for ITQ-3. The behavior within cage-type materials can be further understood in terms of barriers to diffusion. Figure 10 shows freeenergy profiles for pure CO2 and N2 at infinite dilution along the radial diffusion coordinate for ZK4. Free-energy profiles are calculated as βF(R) = -ln ÆP(R)æ,41 where P(R) represents the probabilities shown in Figure 8. Figure 10 shows that for N2 the largest free-energy barrier lies inside the connecting channel whereas for CO2 there is a barrier at the entrance to this channel. Note how for both molecules there is a barrier at the center of the cage (R=0 A˚). The relatively large distance between the zeolite wall and adsorbates around R ≈ 0 A˚ results in weak interactions and the ensuing barrier to diffusion at the center of the cage. The different loading dependences of CO2’s and N2’s diffusion coefficients within each material (Figure 3) can be understood in (40) Krishna, R.; van Baten, J. M. Microporous Mesoporous Mater. 2008, 109, 91. (41) Dubbeldam, D.; Beerdsen, E.; Vlugt, T. J. H.; Smit, B. J. Chem. Phys. 2005, 122, 224712.
DOI: 10.1021/la104245c
1959
Article
Madison et al. Table 3. Free-Energy Barriers to Diffusiona ΔF (infinite dilution) (K)
ITQ-7 CO2 ITQ-7 N2 silicalite N2 (average) silicalite CO2 (average) ZK4 N2 ITQ-3 N2 ITQ-3 CO2 ZK4 CO2 a Values for silicalite and ITQ-721 have been context for our discussion. (See the text.)
535 559 644 718 1004 1075 1458 1928 given to provide a
Figure 12. Three-dimensional free-energy map for CO2 within ZK4 at 298 K and infinite dilution. The ZK4 structure is shown in dark gold; the regions where it can be seen correspond to the β cages that have essentially zero probability of containing a diffusing molecule. The arrow points to the barrier located at the entrance of the narrow channel.
Figure 11. Free-energy profiles (F(R) þ factor) and U(R) for the CO2 center of mass (in red, above) and the N2 center of mass (in blue, below) along the radial coordinate at 298 K. The factor has been chosen in each case so that the entropy is always arbitrarily positive.
terms of the location of the barrier. Smit and co-workers42,43 have used free-energy profiles to understand diffusion in nanoporous materials for the case when the barrier is located in the narrow connector channel, as is the case for N2. For those systems, it is not unusual to find that, as in the cases discussed in this work, the diffusion coefficients do not change much with loading because two effects oppose each other. On one hand, as the loading increases, the crowding of the probable locations within cages increases. This results in a higher free energy compared to that experienced by an isolated molecule at these locations with the consequent decrease in the height of the barrier to diffusion. (42) Beerdsen, E.; Dubbeldam, D.; Smit, B. Phys. Rev. Lett. 2006, 96, 044501. (43) Beerdsen, E.; Dubbeldam, D.; Smit, B. Phys. Rev. Lett. 2005, 95, 164505.
1960 DOI: 10.1021/la104245c
On the other hand, as the loading increases, the collisions became more frequent and thus crossing the barrier becomes more difficult. These two effects result in diffusion coefficients that do not change much with loading. This contrasts with the previously21 examined case of CO2 within ITQ-3, a case when the barrier is not located in the narrow connector channel. In this situation, barriers to CO2 diffusion increase substantially with loading corresponding to the pronounced loading dependence of CO2 diffusion coefficients seen in Figure 3. It is also useful to note how the size of the barriers differs among materials (Table 3), where ΔF values have been calculated as Fmax(q) - Fmin(q) and q is the reaction coordinate representative of diffusion in each material. The relatively large free-energy barrier for diffusion within ITQ-3 and ZK4 could be the reason for the similarity between Ds(θ) and Do(θ) within cage-type materials. Large barriers to diffusion probably result in weakly correlated molecular jumps between adsorption sites and thus a limited amount of collective motion. The nature of the free-energy barriers to diffusion can be further understood by examining the relative contribution of energetic and entropic effects. The energy contribution can be isolated by calculating Z Uðx, y, qÞe - Uðx, y, qÞ=kT hUðqÞi ¼ x, y where U is the Boltzmann averaged potential energy of 250 molecules with centers of mass at x, y, and q and random orientations. The entropic contribution can be inferred by using F = U - TS. Figure 11 shows F(R) (at infinite dilution) and U(R) at 298 K for ZK4. Notice that UN2(R) is less negative than UCO2(R), which can be explained by the relatively strong Coulombic interactions between CO2’s quadrupole and the zeolite. These plots show that within ZK4 the free-energy barrier to N2 diffusion is purely entropic in nature. N2 has a potential energy minimum at the same location as its free energy maximum (the window between cages). The existence of entropic barriers to diffusion has been described previously for the diffusion of methane in ZK4-zeolite.23,25 Langmuir 2011, 27(5), 1954–1963
Madison et al.
Article
Figure 13. CO2’s (top) and N2’s (bottom) free energy within ZK4 as a function of both R and j.
The fact that the free-energy barrier is entropic for N2 explains the lack of correlation between the apparent activation energy and the diffusion rates discussed earlier (Figure 5). The correlation that was proposed captures the energy contribution to the barrier but not the Arrhenius pre-exponential factor, which takes into account entropic effects. It is not surprising then that when a barrier is entropic in nature, as for N2 in cage-type materials, the correlation is less meaningful. Figure 11 also shows that the free-energy barrier to CO2 diffusion is energetic in nature. This corresponds nicely to the fact that for both cage-type materials there is a good correlation between e-Ea/kT and Ds(0). It is intuitive why CO2 has a freeenergy minimum inside the narrow channel: the Coulombic Langmuir 2011, 27(5), 1954–1963
interactions between the close walls and adsorbate are attractive enough to overcome the unfavorable loss of entropy that is the result of the decrease in available volume. What is not clear is why there are free-energy barriers located at the entrances of the channels. The location of these barriers for ZK4 (around R =4.7 A˚) is seen very clearly in Figure 12. Figure 12 shows a 3D free-energy map for CO2 within ZK4 at 298 K and infinite dilution. The area shown is exactly the same as that shown in the ZK4 structure in Figure 1 (top), and the black regions correspond to the β cages that have essentially zero probability of containing a diffusing molecule. This 3D map highlights the fact that CO2 molecules move along the walls, that the narrow channels are the most stable locations, and that the center of the R cage is the least-stable DOI: 10.1021/la104245c
1961
Article
Madison et al.
the coordinates of the center of mass at the final time and B ro represents the coordinates of the center of mass at the initial time. The angular velocity, ωB, was calculated from the angular momentum, LB, and the moment of inertia, I, by the following equation ! L ! ¼ ω ¼ I
3 P i¼1
mi ! r i ð! r i -! v COM Þ 3 P i¼1
Figure 14. Translational (solid lines) and angular (dashed lines) velocity autocorrelation functions for CO2 (red) and N2 (blue) in ZK4.
location (barring the β cages that are not accessible to diffusing molecules and thus are excluded from this analysis). The defined landscape shown in Figure 12 also suggests that CO2 rotates or reorients itself as it diffuses. To understand further the origin of the barrier at the entrance of the channel, we have examined the dependence of the free energy on the orientation of molecules. The orientation of a molecule, j, has been defined as the angle between the diffusion coordinate and a vector defined by the location of the terminal atoms of the molecule. Thus, j = 0 corresponds to a molecule that is parallel to the diffusion coordinate, and j = 90 corresponds to one that is perpendicular. Figure 13 shows CO2’s (top) and N2’s (bottom) free energy within ZK4 as a function of both R and j. First, notice how the variations are much more pronounced for CO2 than for N2, which indicates that orientation effects are not as important for N2 as for CO2. For both, however, it is unfavorable for the molecules to be perpendicular to the diffusion coordinate at large radii (in the channel), but at lower radii (center of the cage), the orientations are more varied. Note how Figure 13 shows that the lowest-free-energy path for a CO2 molecule to cross the threshold of the channel involves its reorientation. In this case but not for N2, it appears that there is a well-defined path. Specifically, as the CO2 molecule crosses the free-energy barrier when it enters the channel, it must reorient from a perpendicular to a parallel orientation; to leave the channel, it must do the reverse. This is the likely cause of the free-energy barrier at the entrance of the channels that is seen in the free-energy profile for CO2 in ZK4. The fact that when Coulombic interactions are ignored (Figure 9) this barrier disappears suggests that the need for reorientation originates in the unfavorable oxygen-oxygen Coulombic interactions that occur between CO2 and the ring of oxygens at the entrance of the channel. The examination of translational and angular velocity auto correlation functions allows us to gain an additional appreciation of the difference in the behavior of CO2 and N2. The normalized autocorrelation value for time t was calculated for both angular and translational velocity as44 M X N ! Æ a i ðtj Þ 3 ! a i ðtj þ tÞæ 1 X C! ! ðtÞ ¼ a a MN j ¼ 1 i ¼ 1 Æ! a i ðtj Þ 3 ! a i ðtj Þæ
where the vector aB is either B v COM for translational velocity or ωB r o where B r f represents for angular velocity. Note that B v COM=rBf - B (44) Fritzsche, S.; Osotchan, T.; Schuring, A.; Hannongbua, S.; Karger, J. Chem. Phys. Lett. 2005, 411, 423.
1962 DOI: 10.1021/la104245c
m i ð! r i 3! r iÞ
where mi is the mass of atom i, B r i is the coordinate of atom i, and v COM is the relative velocity of atom i compared to the center B ri-B of mass. The summation over N includes every molecule in the simulation, and the summation over M includes autocorrelations over multiple relaxation periods. The results of these calculations are shown in Figure 14. As was noted earlier in this article, the qualitative features discussed here are also present in the case of ITQ-3. The CO2 traces in this Figure show that both the translational velocity autocorrelation function, Cvv, and the angular velocity autocorrelation function, Cww, oscillate about zero and eventually reach a relaxation state where the correlation value is zero. Cvv and Cww do not oscillate with the same frequency; therefore, it is not clear if the two motions are coupled to each another. These oscillations may be indicative of the wiggling motion of CO2. A visual inspection of CO2 in ZK4’s narrow channel revealed that CO2 spends a lot of time wiggling; this wiggling might be better described as frustrated rotations that molecules experience as they move until they reach an orientation that is appropriate to cross the free-energy barrier. A comparison of CO2 and N2’s traces reveals that the oscillations are much more prevalent in CO2. Furthermore, when the autocorrelation functions are calculated in the absence of Coulombic forces, the CO2 traces are no longer oscillatory and behave much more like N2 (results not shown). These results suggest that frustrated rotations are important in the behavior of linear molecules with strong quadrupoles when they are adsorbed within cage-type zeolites.
4. Concluding Remarks We have used molecular simulations to examine the behavior of N2 and CO2 in cage-type siliceous zeolites with the goal of further documenting and understanding previously reported interesting behavior. Despite similarities between N2 and CO2 diffusion in other materials, we show here that the diffusion of CO2 within two cage-type zeolites is dominated by an energy barrier to diffusion located at the entrance to the narrow channels connecting larger cages. This barrier originates in Coulombic interactions between the zeolites and CO2’s quadrupole and results in well-defined orientations for the diffusing molecules. Furthermore, CO2’s favorable electrostatic interactions with the zeolite framework result in preferential binding in the windows between cages. This contrasts with N2’s behavior, which is much more in line with that of molecules previously studied. N2 diffusion is dominated by the existence of an entropic barrier to diffusion located in the windows between cages, and reorientation does not play an important role in diffusion. To our knowledge, this is the first study in which the distinctive behavior of CO2 within cage-type zeolites is shown to be rooted in its large quadrupole and is characterized in terms of rotational behavior. The differences between N2 and CO2 in all silica cage-type zeolites have interesting implications for understanding and predicting the properties of CO2/nonpolar adsorbate mixtures. Langmuir 2011, 27(5), 1954–1963
Madison et al.
Within cage-type materials, these species will likely not directly compete for the same adsorption sites, CO2 might block the diffusion of the other species, “simple” mixture theories would have limited applicability, and estimated activation energies might not be a useful predictor of diffusion rates. Furthermore, the fact that CO2’s peculiar behavior is a result of its quadrupole suggests that studying materials with a more complex chemical composition (and thus a more intense electric field) but similar
Langmuir 2011, 27(5), 1954–1963
Article
pore structure to the ones studied here could yield microporous materials uniquely suited to practical CO2 separations. Acknowledgment. We gratefully acknowledge the National Science Foundation (CHE-0520704), the Petroleum Research Fund (PRF no. 45141-B5), Carleton College, and the Howard Hughes Medical Institute for computing resources and stipend support to carry out this research.
DOI: 10.1021/la104245c
1963
