Langmuir 2007, 23, 6431-6437
6431
Examining the Accuracy of Ideal Adsorbed Solution Theory without Curve-Fitting Using Transition Matrix Monte Carlo Simulations Haibin Chen and David S. Sholl* Department of Chemical Engineering, Carnegie Mellon UniVersity, Pittsburgh, PennsylVania 15213 ReceiVed February 7, 2007. In Final Form: March 19, 2007 Ideal adsorbed solution theory (IAST) is a well-known approach to predicting multicomponent adsorption isotherms in microporous materials from experimental or simulation data for single-component adsorption. A limitation in practical applications of IAST is that useful calculations often require extrapolation of fitted single-component isotherms beyond the range for which data are available. We introduce a molecular simulation approach in which the intrinsic accuracy of IAST can be examined in a context that avoids any need to perform curve fitting with single-component data. Our approach is based on using transition matrix Monte Carlo to define single-component adsorption isotherms for arbitrary bulk-phase pressures from a single simulation. We apply our approach to several light gas mixtures in silica zeolites and a carbon nanotube to examine the intrinsic accuracy of IAST for these model systems.
1. Introduction The adsorption of gases inside microporous materials is the basis of many industrial processes for gas purification and separation.1 Although multicomponent adsorption can be assessed experimentally, performing experiments that systematically explore the full range of potential operating conditions for a material of interest is time consuming. Performing accurate experiments to characterize single-component adsorption, by contrast, is relatively routine. As a result, a significant amount of work has been done on methods to predict multicomponent adsorption from single-component data. Among these methods, thermodynamic approaches that do not rely on a particular physical model for the single-component isotherms are especially useful. The best known of these models is ideal adsorbed solution theory (IAST),2 which has been successful for multicomponent adsorption equilibria in a wide variety of systems. Because it treats each adsorbed mixture as being ideal, IAST is not accurate in all situations. It has been pointed out that IAST frequently becomes less accurate at high densities of the adsorbed species, even for mixtures that are relatively ideal.3-7 One approach to dealing with this nonideality is to approximate the activity coefficients of the adsorbed mixture, via real adsorbed solution theory (RAST).8,9 Other approaches to this challenging problem include vacancy solution theory,10 the use of virial excess mixing coefficients,11-13 and extended isotherms suitable for multicomponent systems.14,15 * To whom correspondence should be addressed. E-mail: sholl@ andrew.cmu.edu. (1) Yang, R. T. Gas separation by adsorption processes; Imperial College Press: London, 1997. (2) Myers, A. L.; Prausnitz, J. M. AIChE J. 1965, 11, 121. (3) Dunne, J.; Myers, A. L. Chem. Eng. Sci. 1994, 49, 2941. (4) Murthi, M.; Snurr, R. Q. Langmuir 2004, 20, 2489. (5) Goj, A.; Sholl, D. S.; Akten, E. D.; Kohen, D. J. Phys. Chem. B 2002, 106, 8367. (6) Akten, E. D.; Siriwardane, R.; Sholl, D. S. Energy Fuels 2003, 17, 977. (7) Hu, X.; Do, D. D. AIChE J. 1995, 46, 1585. (8) Costa, E.; Sotelo, J. L.; Calleja, G.; Marro´n, C. AIChE J. 1981, 27, 5. (9) Talu, O.; Zwiebe, I. AIChE J. 1986, 32, 1263. (10) Ding, L. P.; Bhatia, S. K. AIChE J. 2002, 48, 1938. (11) Appel, W. S.; LeVan, M. D.; Finn, J. E. Ind. Eng. Chem. Res. 1998, 37, 4774. (12) Qi, N.; LeVan, M. D. Ind. Eng. Chem. Res. 2005, 44, 3726. (13) Qi, N.; LeVan, M. D. Ind. Eng. Chem. Res. 2005, 44, 3733. (14) Bai, R.; Yang, R. T. J. Colloid Interface Sci. 2001, 239, 296. (15) Nguyen, C.; Do, D. D. Langmuir 2001, 17, 1552.
Although IAST does not rely on any specific model for the adsorption of individual species, practical application of IAST does of course require that the functional form of the singlecomponent isotherms be specified. In practice, these isotherms are fitted to experimental (or simulation) data obtained over a finite range of pressures. An intrinsic feature of IAST is that to make predictions about multicomponent adsorption it is frequently necessary to use information from the single-component isotherms from pressures outside the range for which data is available. The main aim of this paper is to introduce a calculation method for molecular simulations of adsorption where this seemingly inevitable uncertainty due to extrapolation of the singlecomponent isotherms is avoided. Molecular simulations of adsorption can play a useful role in understanding multicomponent adsorption, in large part because obtaining multicomponent adsorption is considerably more straightforward in simulations than in experiment.16-19 The usual approach to simulating adsorption in a well-defined microporous material is grand canonical Monte Carlo (GCMC) simulation.20 GCMC is conceptually attractive because it closely mimics experimental measurements of adsorption, but it is not the only simulation method that can reliably describe adsorption. Transition matrix Monte Carlo (TMMC) is an alternative sampling method that can increase the numerical efficiency of standard GCMC approaches.21,22 In this paper, we use an algorithm that combines IAST with TMMC and histogram reweighting (HR) methods. This algorithm applies IAST without requiring any curve fitting. The results of using this approach provide insight into the intrinsic accuracy of IAST that have not previously been possible. The details of our algorithms are presented in Section 2. In Section 3, the predictions of our TMMC-based IAST methods are compared with conventional approaches to using IAST based on curvefitting and to standard binary GCMC simulations for several adsorbed mixtures. We conclude in Section 4 by summarizing our findings. (16) Fuchs, A. H.; Cheetham, A. K. J. Phys. Chem. B 2001, 105, 7375. (17) Heuchel, M.; Snurr, R. Q.; Buss, E. Langmuir 1997, 13, 6795. (18) Krishna, R.; van Baten, J. M. J. Phys. Chem. B 2005, 109, 6386. (19) Vlugt, T. J. H.; Krishna, R.; Smit, B. J. Phys. Chem. B 1999, 103, 1102. (20) Frenkel, D.; Smit, B. Understanding Molecular Simulations: From Algorithms to Applications; Academic Press: London, 1996. (21) Errington, J. R. J. Chem. Phys. 2003, 118, 9915. (22) Chen, H.; Sholl, D. S. Langmuir 2006, 22, 709.
10.1021/la700351c CCC: $37.00 © 2007 American Chemical Society Published on Web 04/19/2007
6432 Langmuir, Vol. 23, No. 11, 2007
Chen and Sholl
2. Methods and Simulation Details 2.1. Ideal Adsorbed Solution Theory. IAST is a well-known technique for predicting the adsorption equilibria for components in a multicomponent mixture using only data from the purecomponent equilibria at the same temperature on the same adsorbent. Myers and Prausnitz2 originally derived the theory for a twodimensional homogeneous adsorbed phase with a temperatureinvariant area that is equally accessible to all components. The theory was later extended to treat a three-dimensional adsorbed phase.23,24 For an N-component system, the Gibbs isotherm3 is N
dψ )
∑n d ln f i
(1)
i
i)1
where ψ is the grand potential, fi is the fugacity of component i, and ni is the adsorbed amount for pure component i in equilibrium with a bulk phase with the specified pure-component fugacity. Equation 1 may be integrated as follows
∫
f 0i
0
n0i (t) d ln t ) ψ0i ) ψ
(2)
where the superscript “0” indicates a property of the pure component evaluated in the mixture at temperature T and grand potential ψ. The final equality in this expression indicates that the grand potential must be the same for each adsorbed component. Multicomponent vapor-adsorption phase equilibria at total pressure P can then be written as fi ) Pφiyi ) f 0i γixi
(3)
where yi and xi are the vapor-phase and adsorbed-phase mole fractions of species i, respectively, φi is the fugacity coefficient of component i in the multicomponent vapor phase, and γi is the activity coefficient of component i in the adsorbed phase. Equations 2 and 3 can be coupled to yield the total loading via 1/nt )
∑[x /n (f )] + RT∑(∂ ln γ /∂ζ) i
0 i
0 i
i
i
T,x
(4)
i
Here, ζ is the grand potential density, defined as ζ ) -RTψ/V. The individual component loadings are then simply ni ) xint
(5)
IAST assumes that each activity coefficient is equal to unity for all conditions, so eq 3 simplifies to Pyi ) P0i xi
(6)
and the second term on the right side of eq 4 vanishes. This expression is written in a form suitable when the vapor phase can be treated as ideal. To use IAST, it is of course necessary to specify the form of the single-component isotherms. In practice, this is done by fitting a continuous function to a discrete set of adsorption data.25-28 This approach creates two possible sources of uncertainty in the accuracy of IAST. First, if the functional form selected for the single component isotherms does not precisely fit the available data, then IAST cannot precisely describe the adsorbed mixture. This type of uncertainty is relatively easy to control by testing multiple functional forms when fitting the available data. A more subtle but equally important source of uncertainty is that eq 2 frequently requires integration of the pure-component isotherms to fugacities outside the range of available (23) Karavias, F.; Myers, A. L. Mol. Simul. 1991, 8, 51. (24) Vuong, T.; Monson, P. A. Adsorption 1999, 5, 295. (25) Zhu, W.; Hrabanek, P.; Gora, L.; Kapteijn, F.; Moulijn, J. A. Ind. Eng. Chem. Res. 2006, 45, 767. (26) Li, S.; Falconer, J. L.; Noble, R. D.; Krishna, R. Ind. Eng. Chem. Res., in press. (27) Skoulidas, A. I.; Sholl, D. S.; Krishna, R. Langmuir 2003, 19, 7977. (28) Skoulidas, A. I.; Sholl, D. S. J. Phys. Chem. B 2001, 105, 3151.
data, even when considering multicomponent mixtures with total pressures similar to those for which pure-component data is available. This difficulty is particularly severe in examples with large adsorption selectivity during integration of the weakly adsorbed components in eq 2. This observation suggests that curve fitting inevitably introduces errors into IAST predictions that go beyond those associated with the physical approximations in the theory and that no obvious means is currently available to assess these errors. To directly evaluate the physical assumptions embodied in IAST, it is necessary to have a means to measure the adsorbed phase activity coefficients. This can be accomplished provided that a means of establishing the true multicomponent adsorption isotherm is available. Following Dunne and Myers,3 we can back-calculate the activity coefficients from binary mixture adsorption isotherms. First, ψ is calculated from ψ(P,y1) ) ψ(P,0) +
∫
y1
y1)0
(n1/y1 - n2/y2) dy1
(7)
This value of ψ then can be inserted into eq 2 for each component to determine the corresponding f0i . By assuming that the vapor-phase fugacity coefficients, φi, are unity, the activity coefficient γi can then be determined from eq 3. In this paper we use this approach only to estimate the adsorbed activity coefficients. To this end, we used GCMC to calculate binary adsorption isotherms for species of interest at a range of compositions for total gas-phase pressures from 0 to 100 bar. The dual-site Langmuir isotherm used to fit each pure-component adsorption isotherm and an empirical function that contains 16 parameters29,30 was used to fit the full range of binary adsorption data. The principal advantage of this method is that it allows the activity coefficients to be determined analytically from method outlined above. 2.2. Transition Matrix Monte Carlo and Histogram Reweighting. If an atomically detailed model of adsorbates in a microporous adsorbent is available, then several molecular simulation methods are available to compute pure-component or mixture adsorption isotherms.20,31 The best known method is GCMC, a method that is entirely analogous to experimental measurements in which adsorption is characterized separately at a discrete set of state points. Highly efficient sampling schemes have been developed that allow GCMC to be used even for adsorbates with many degrees of freedom at high densities.17,20 An attractive feature of GCMC is that it is relatively straightforward to use it to compute binary adsorption data at state points of interest.5,17,32-36 A relatively recently developed alternative to GCMC is to sample the adsorption isotherm using TMMC.21,37-43 Smith et al. introduced the transition probability method to estimate the weights of a multicanonical data set.37,38 Recent applications of TMMC have included calculations of interfacial tension in Ising models,39,40 calculation of single-component liquid-vapor coexistence,21 and calculation of binary fluid properties.43-45 We have used the TMMC methods developed in this earlier work to compute pure-component and binary (29) Skoulidas, A. I.; Bowen, T. C.; Doelling, C. M.; Falconer, J. L.; Noble, R. D.; Sholl, D. S. J. Membrane Sci. 2003, 227. (30) Chen, H.; Sholl, D. S. J. Membrane Sci. 2006, 269, 152. (31) Allen, M. P.; Tildesley, D. J. Computer simulation of liquids; Oxford University Press: Oxford, 1987. (32) Clark, L. A.; Gupta, A.; Snurr, R. Q. J. Phys. Chem. B 1998, 102, 6720. (33) Yang, J.-H.; Clark, L. A.; Ray, G. J.; Kim, Y. J.; Du, H.; Snurr, R. Q. J. Phys. Chem. B 2001, 105, 4698. (34) Schuring, D.; Koriabkina, A. O.; de Jong, A. M.; Smit, B.; van Santen, R. A. J. Phys. Chem. B 2001, 105, 7690. (35) Schenk, M.; Vidal, S. L.; Vlugt, T. J. H.; Smit, B.; Krishna, R. Langmuir 2001, 17, 1558. (36) Challa, S. R.; Sholl, D. S.; Johnson, J. K. J. Chem. Phys. 2002, 116, 814. (37) Smith, G. R.; Bruce, A. D. J. Phys. A 1995, 28, 6623. (38) Smith, G. R.; Bruce, A. D. Phys. ReV. E 1996, 53, 6530. (39) Fitzgerald, M.; Picard, R. R.; Silver, R. N. Europhys. Lett. 1999, 46, 282. (40) Fitzgerald, M.; Picard, R. R.; Silver, R. N. J. Stat. Phys. 2000, 98, 321. (41) Wang, J.-S.; Tay, T. K.; Swendsen, R. H. Phys. ReV. Lett. 1999, 82, 476. (42) Wang, J.-S.; Swendsen, R. H. J. Stat. Phys. 2002, 106, 245. (43) Shen, V. K.; Errington, J. R. J. Chem. Phys. 2005, 122, 064508. (44) Singh, J. K.; Errington, J. R. J. Phys. Chem. B 2006, 110, 1369. (45) Shen, V. K.; Errington, J. R. J. Chem. Phys. 2006, 124, 024721.
Examining the Accuracy of IAST without CurVe-Fitting
Langmuir, Vol. 23, No. 11, 2007 6433
adsorption isotherms for microporous adsorbents.22 Here, we outline the principal ideas needed to use TMMC for single-component adsorption calculations. The central idea in using TMMC is to define the system of interest in terms of macrostates, S, and microstates, s. In the case of adsorption, the number of particles of the adsorbed species provides a natural set of macrostates, while the detailed configuration of the adsorbed particles defines a microstate. During a TMMC simulation, a collection matrix that contains information about all attempted transitions (not just the accepted transitions as in a standard GCMC simulation) between macrostates is updated after each Monte Carlo event. The net macrostate probability (particle number probability distribution) Π(N) can then be determined by detailed balance from the accumulated collection matrix.21,30,43 In practice, the TMMC simulation involves defining a chemical potential reference state, µ0. The macrostate probability associated with other values of the chemical potential is given in a simple way by HR:41 ln Π(N;µ) ) ln Π(N;µ0) + βN(µ - µ0)
(8)
Once the macrostate probability is known, the adsorption isotherm is defined by Nmax
∑ Π(N,µ)N
〈N〉 )
N)1
Nmax
(9)
∑ Π(N,µ)
N)1
Here Nmax denotes the particle number above which the macroscopic probability distribution for adsorbed phase is negligibly small. For the purposes of this paper, the crucial feature of eq 9 is that it defines the pure-component adsorption isotherm for arbitrary chemical potentials without curVe fitting once Π(N) is known accurately. We can therefore use pure-component isotherms derived in this way to apply IAST for the first time in a way that avoids the uncertainties inherent in all previous applications of this theory that relied on curve fitting of isotherms over a finite range of data. In the next section we pursue this idea by comparing IAST predictions based on TMMC-derived pure-component adsorption isotherms with binary adsorption data generated using GCMC. Because binary adsorption data can be generated using GCMC with a high level of accuracy, the results from these calculations should be considered as defining the “true” mixture adsorption in the model adsorbents we consider.
3. Results We present here five examples of using IAST to predict binary adsorption in ordered nanoporous materials: CH4/CF4 and CH4/ CO2 in a defect-free (10,10) single-walled carbon nanotube (SWNT), CH4/CF4 and CH4/CO2 in silicalite, and CH4/CO2 in the all silica form of zeolite DDR. Each of these materials is relatively simple from a chemical point of view, but they differ in the characteristics of their pore geometry. Carbon nanotubes define very smooth cylindrical pores. The smoothness of these pores causes molecular diffusion to be extremely rapid in carbon nanotubes, an observation with interesting implications for fabricating membranes from these materials.30,46-51 Silicalite, a zeolite that has been studied with a wide variety of experimental (46) Skoulidas, A. I.; Ackerman, D. M.; Johnson, J. K.; Sholl, D. S. Phys. ReV. Lett. 2002, 89, 18501. (47) Chen, H.; Sholl, D. S. J. Am. Chem. Soc. 2004, 126, 7778. (48) Holt, J. K.; Park, H. G.; Wang, Y.; Stadermann, M.; Artyukhin, A. B.; Grigoropoulos, C. P.; Noy, A.; Bakajin, O. Science 2006, 312, 1034. (49) Hinds, B. J.; Chopra, N.; Rantell, T.; Andrews, R.; Gavalas, V.; Bachas, L. G. Science 2003, 303, 62. (50) Majumder, M.; Chopra, N.; Andrews, R.; Hinds, B. J. Nature 2005, 438, 44. (51) Newsome, D. A.; Sholl, D. S. Nano Lett. 2006, 6, 2150.
Table 1. Maximum Particle Number Per Simulation Box in TMMC Simulations, Nmax, and the Corresponding Values in Different Adsorbent Structures carbon nanotube
silicalite
DDR
Nmax molecule/nm Nmax molecule/uc Nmax molecule/uc CH4 CF4 CO2
65 34 150
13.2 6.9 30.5
171 131 181
21.4 16.4 22.6
127 N/A 192
21.2 N/A 32.0
and theoretical methods,52,53 has a three-dimensional channel structure. The pores in DDR are made up of cages that are separated by relatively narrow windows that make up a twodimensional pore network. Because of the narrow windows in DDR, it is of interest for potential applications in light gas separations.54-56 Most of the adsorbent structures and interatomic potentials were taken from our previous work on single- and binarycomponent adsorption in these materials.22,46,47,57,58 The interatomic potentials for molecules adsorbed in silicalite and DDR were assumed to be identical. Interactions of CO2 with DDR were described using the potential derived by Makrodimitris et al. for CO2 adsorption in silicalite.59 The structure of DDR was defined using the atomic coordinates measured experimentally by Gies (neglecting the template molecules also present in this experimental crystal structure).54 In our calculations with CO2, partial charges on Si and O atoms of +2e and -1e, respectively, were used. Each adsorbent was assumed to be rigid, so no relaxation of the adsorbent atoms due to adsorption was considered. All calculations were performed at 298 K. Simulations of the SWNT used an isolated nanotube 4.43 nm long with periodic boundary conditions along the pore axis. For silicalite and DDR, the simulation volumes were 2 × 2 × 2 and 3 × 2 × 1 unit cells, respectively. The quadrupole-quadrupole interactions between pairs of CO2 molecules decay quickly as the molecules’ separation is increased. As a result, the Coulomb interactions between CO2 molecules were handled using a simple radial cutoff.58 Coulomb interactions between CO2 molecules and atoms in the silica zeolite’s frameworks were precomputed on a fine spatial grid using a large brute force calculation. The potential energy from these interactions at arbitrary positions was then defined via interpolation from the precomputed grid.5 Unless otherwise stated, the single-component TMMC adsorption isotherms discussed below were calculated using a reference chemical potential corresponding to an ideal gas pressure of 100 bar. More simulation details can be found in our previous work.22 All pressures reported below are given in terms of ideal gas pressures, so vapor-phase fugacity corrections would be needed to relate these pressures to real pressures for some of the high-pressure conditions we consider. In our GCMC simulations of binary adsorption, 5 × 106 MC events were used to equilibrate the system and data was then collected during another 5 × 106 MC events for each bulk gas composition. Nmax was determined by continuing calculations as N0 was increased stepwise until N0 ln(Π(N0)) < -100 when Π is normalized such that ∑N)0 Π(N) ) 1. More details regarding Nmax are listed in Table 1. (52) Ka¨rger, J.; Ruthven, D. Diffusion in Zeolites and Other Microporous Materials; John Wiley and Sons: New York, 1992. (53) Keil, F. J.; Krishna, R.; Coppens, M. O. ReV. Chem. Eng. 2000, 16, 71. (54) Gies, H. Z. Kristallographie 1986, 175, 93. (55) Tomita, T.; Nakayama, K.; Sakai, H. Micro. Meso. Mater. 2004, 68, 71. (56) Krishna, R.; van Baten, J. M.; Garcia-Perez, E.; Calero, S. Ind. Eng. Chem. Res., in press. (57) Skoulidas, A. I.; Sholl, D. S. J. Phys. Chem. B 2002, 106, 5058. (58) Skoulidas, A. I.; Sholl, D. S.; Johnson, J. K. J. Chem. Phys. 2006, 124, 054708. (59) Makrodimitris, K.; Papadopoulos, G. K.; Theodorou, D. N. J. Phys. Chem. B 2001, 105, 777.
6434 Langmuir, Vol. 23, No. 11, 2007
Chen and Sholl
Figure 1. Contour plots of |(ci,GCMC - ci,IAST)|/ci,GCMC calculated as described in the text for binary adsorption of CH4/CF4 mixtures in silicalite at 298 K. The left plot shows the result for CH4, while the right plot shows the result for CF4. In each case, the results are shown as a function of the total bulk-phase fugacity and the bulk-phase mole fraction of CH4.
Figure 1 shows |(ci,GCMC - ci,IAST)|/ci,GCMC (expressed as a percentage) for CH4/CF4 mixtures adsorbed in silicalite at 298 K. Here, ci,IAST is the adsorbed amount of species i determined from applying IAST using the pure-component isotherms calculated using TMMC and ci,GCMC is the adsorbed amount of species i determined from our GCMC mixture simulations. These data have been presented previously.22 To generate these data, GCMC simulations were performed with the state points spanning the full range of possible adsorbate compositions from 0 to 100 bar. Specifically, GCMC was used for total pressures of 0.1, 0.2, ..., 0.9, 1.0, 2.0, ..., 9.0, 10, 20, ..., 100 bar with CH4 mole fractions of 0.01, 0.05, 0.1, 0.15, 0.2, ..., 0.85, 0.90, 0.95, and 0.99 at each total pressure. It can be seen from Figure 1 that IAST is quite accurate in this system for many bulk-phase conditions, giving relative errors in the adsorbed amounts of both species of less than 5%. One exception to this observation occurs when the bulk phase is rich in one component and total pressure is higher than roughly 20-30 bar. In this case, the IAST predictions are in error by as much as 16-18%, particularly for the more weakly adsorbed species, CH4. Similar results have been observed previously from experimental and simulation studies of this adsorbed mixture by Heuchel et al.17 Figure 2 shows CH4/CF4 selectivity for adsorption from equimolar bulk mixtures in the SWNT and silicalite. The mixture adsorption data needed to define this selectivity was calculated in three independent ways. First, adsorption data at discrete state points were calculated using binary GCMC. These results are shown using unfilled symbols in Figure 2. In both materials the more strongly adsorbing CF4 is selectively adsorbed, but the adsorption selectivity decreases as the bulk-phase pressure increases. At very high pressures, CH4 is selectively adsorbed in both materials rather than CF4. Similar entropy-driven changes in adsorption selectivity are well known in other nanoporous materials.60 The second approach used in Figure 2 was a conventional application of IAST based on fitting curves to the pure-component isotherms. Specifically, a dual-site Langmuir isotherm27,61 was fitted to pure-component data with fugacities between 0 and 100 bar. This range of conditions is intended to mimic the limited range of pressures that can be routinely accessed experimentally. In applying IAST as shown in Figure 2, it is necessary to use results from the fitted single-component isotherms that lie outside the range of the data used to fit the isotherms. For CH4/CF4 in silicalite at a total pressure of 100 and 500 bar, for example, it (60) Krishna, R.; Smit, B.; Calero, S. Chem. Soc. ReV. 2002, 31, 185. (61) Krishna, R.; Paschek, D. Ind. Eng. Chem. Res. 2000, 39, 2618.
Figure 2. Selectivity for adsorption equilibria with 50:50 gas-phase CH4/CF4 mixtures adsorbed in silicalite and (10,10) SWNT at 298 K as computed using IAST-TMMC-HR (solid curves), a traditional application of IAST (dashed curves), and binary GCMC (symbols). The dilute loading correction (dash dot curve) was calculated from IAST-TMMC-HR using pure-component adsorption isotherms calculated with TMMC using a low-pressure reference state, as described in the text.
is necessary to integrate eq 2 up to fugacities of 117 (87) and 492 (508) bar for CH4 (CF4), respectively. In the (10,10) carbon nanotube at the same conditions, it is necessary to integrate eq 2 up to fugacities of 86 (119) and 329 (1036) bar for CH4 (CF4), respectively. It can be seen from Figure 2 that the IAST predictions based on this approach are not accurate in the limit of low pressures. This is a limitation of the fitted isotherms, not a limitation of IAST, since IAST is exact in this limit and the selectivity approaches the ratio of the Henry’s constants of the two adsorbing species.2,5,17 Because we fitted our isotherms in order to achieve the best fit for the complete set of pure-component adsorption data the resulting curves are slightly inaccurate in predicting the Henry’s constants for each species. In principle, this problem can be avoided by accurately calculating (or measuring) the Henry’s constant and constraining the fitted isotherms to include this value. When doing simulations of the type we have performed here, MC methods exist to rapidly calculate Henry’s constants.62 In a third approach to this problem, we used IAST combined with the single-component isotherms calculated using TMMC (62) Maginn, E. J.; Bell, A. T.; Theodorou, D. N. J. Phys. Chem. 1995, 99, 2057.
Examining the Accuracy of IAST without CurVe-Fitting
and HR. That is, the macrostate probability Π(N) was determined for each component by TMMC simulations and pure-component adsorption isotherms, n0i (fi), were then calculated at any state of interest by applying eq 9. Equation 2 can then be integrated numerically by the knowledge of single component isotherm n0i (fi). We denote these results by IAST-TMMC-HR. The agreement between IAST-TMMC-HR and the binary GCMC results in the Henry’s law region in Figure 2 is not perfect, particularly for the (10,10) SWNT. Recall that in TMMC simulations, the isotherm is generated from simulations performed at a reference chemical potential, which in these examples was chosen to be equivalent to a bulk-phase pressure of 100 bar. Because the macrostate probabilities at dilute loadings can differ by orders of magnitude from the macrostate probability at the reference state, the macrostates most relevant to the dilute limit are not always sampled with high accuracy even after use of biased sampling.21,22 This means that the Henry’s constants for single-component adsorption are computed less precisely than the isotherms at nondilute loadings. For the purposes of understanding the accuracy of IAST this is of little importance, since as we have already discussed, IAST is highly accurate in this limit. It is possible to explore the convergence of the IASTTMMC-HR results by recomputing the single-component isotherms with a different reference state. Figure 2 shows results computed via this method for the (10,10) carbon nanotube using a reference state of 3 bar for CH4 and 1 bar for CF4. This calculation yields good agreement with the IAST result in the limit of low pressures, as expected. Two conclusions can be drawn from Figure 2. First, IAST is remarkably accurate for both adsorbed mixtures in both materials over a broad range of pressures. That is, the assumption of ideality in the adsorbed mixtures is accurate in these examples. This conclusion is also clear for silicalite from Figure 1, although this figure also highlights the some regimes where the accuracy of IAST is less than with an equimolar gas phase. Independently, Figure 2 indicates that the dual-site Langmuir isotherm used to fit the single-component isotherms accurately captures the characteristics of these isotherms even in the range of dense loadings where no data were used to fit the curves. Figure 3 shows CO2/CH4 adsorption and selectivity from equimolar bulk phases in a (10,10) SWNT, silicalite, and DDR. In each case, results from binary GCMC are compared with results from IAST based on fitting isotherms to single-component data and to IAST based on single-component isotherms derived from TMMC. For the curve fitting approach, single-component GCMC data spanning pressures from 0 to 100 bar were used to determine the fitting parameters for a dual-site Langmuir isotherm. The TMMC-based results for CO2/CH4 adsorption in the (10,10) carbon nanotube show that IAST is highly accurate for this mixture over the full range of pressures shown in Figure 3a. At total pressures larger than ∼10 bar, the TMMC-based IAST result slightly overestimates the amount of CH4 that is adsorbed, resulting in a slight underestimation of the adsorption selectivity. It is interesting to note that when IAST is used with the standard curve-fitting approach, the result matches the actual selectivity very accurately in the same pressure range. A surprising observation that can be taken from comparing these two applications of IAST is that the agreement for the more traditional approach is due to a fortuitous cancellation of errors between the small inaccuracy inherent in IAST and the small inaccuracy inherent in using the fitted isotherms to describe the true singlecomponent isotherms in this pressure range. As with the example in Figure 2, there is a systematic discrepancy between the IAST results based on curve fitting and the true results in the limit of
Langmuir, Vol. 23, No. 11, 2007 6435
Figure 3. Binary adsorption isotherms and adsorption selectivities for 50:50 gas-phase CO2/CH4 mixtures in equilibrium with (a) a (10,10) SWNT, (b) silicalite, and (c) DDR at 298 K as computed using IAST-TMMC-HR (solid curves), IAST-curve-fitting (dashed curves), and binary GCMC simulations (symbols). The dilute loading correction (dash-dot curve) in (b) was calculated using IAST-TMMCHR with a reference state of 1 bar for CO2 and 3 bar for CH4.
low pressures. As stated above, this arises from imprecision in describing the Henry’s constants of the two adsorbing species with an isotherm that is fitted to a full range of pressure values, not an inaccuracy of IAST. We reiterate that this inaccuracy can be avoided if the fitted isotherms are constrained to have Henry’s
6436 Langmuir, Vol. 23, No. 11, 2007
constants in agreement with the true results. This phenomenon also appears in the other examples presented below, and we will not discuss it separately for these examples. The adsorption of CO2/CH4 mixtures in silicalite, as shown in Figure 3b, is the first example we have shown where IAST is not accurate for equimolar bulk phases. In particular, IAST overpredicts the adsorbed amount of CH4 for total pressures larger than ∼10 bar, leading to a systematic underestimate of the adsorption selectivity. As with the carbon nanotube, the TMMCbased IAST result, which we view as the most precise test of the intrinsic accuracy of IAST, shows somewhat stronger deviations from the binary GCMC results than the outcome involving curve fitting of the single-component isotherms. In this example, unlike the results shown in Figure 3a, the TMMC-based results show a small inaccuracy in the limit of low pressures. This inaccuracy occurred for the same reasons as in the case of CH4/CF4 mixture discussed above, namely, limited sampling of configurations with very low adsorbate densities when using a relatively highpressure reference state for the TMMC simulation. This inaccuracy can be reduced by using a lower pressure reference state, as shown in Figure 3b by the results from IAST calculations with a reference state of 1 bar for CO2 and 3 bar for CH4 The results for adsorption of CO2/CH4 mixtures in DDR are shown in Figure 3c. In this case it is clear from the TMMC-based IAST calculations that IAST substantially overestimates the adsorption selectivity of this mixture for total pressures larger than ∼1 bar. This occurs because IAST simultaneously overestimates the adsorbed amount of CO2, the more strongly adsorbing species, and underestimates the adsorbed amount of CH4. When the predicted selectivities are relatively large, as they are in this example, applying IAST requires use of the single-component adsorption isotherms at pressures considerably exceeding the total pressures of interest in the mixture. Calculating the IAST results at a total pressure in the mixture of 100 bar, for example, required single-component pressures of 54 and 689 bar in the curve fitting method and 54 and 850 bar in the IASTTMMC-HR method for CO2 and CH4, respectively. In light of this, it is somewhat surprising that the IAST results based on curve fitting of the single-component isotherms is in close agreement with the more precise TMMC-based results for the elevated pressures shown in Figure 3c. This observation suggests that the dual-site Langmuir isotherm accurately captures the details of CO2 and CH4 adsorption in this material even for dense adsorbate loadings. In order to confirm that the limitations observed on the accuracy of IAST in the examples above stemmed from the nonideality of the adsorbed mixtures, we estimated the activity coefficients of the adsorbed species for two examples. Using binary GCMC data, the activity coefficients can be calculated from eqs 7, 2, and 3. Equation 7 can be integrated either numerically or analytically. For simplicity, we fitted our binary GCMC adsorption data to empirical functions (eqs 1 and 2 in ref 30) that allowed eq 7 to be solved analytically. We used an analytical solution of eq 2 based on fitting the single component isotherm with dual-site Langmuir function. The activity coefficients calculated in this way are shown in Figure 4. The activity coefficients calculated from the CH4/CF4 mixture adsorbed in silicalite are close to unity at all pressures considered in Figure 4. That is, this adsorbed mixture is close to being ideal. This is consistent with the good performance of IAST for this mixture (cf. Figure 1). The activity coefficients for the CO2/CH4 mixture adsorbed in DDR shown in Figure 4 deviate considerably from unity for total pressures larger than ∼1 bar. This is consistent with our previous observation that IAST shows systematic deviations from
Chen and Sholl
Figure 4. Activity coefficients for 50:50 gas-phase CH4/CF4 mixture (dashed curves) adsorbed in silicalite and for 50:50 gas-phase CH4/ CO2 mixture (solid curves) adsorbed in DDR. The dilute loading correction (dash dot curve) for the adsorbed mixture in DDR was predicted by numerically integrating eq 2 using a reference state of 1 bar for CO2 and 3 bar for CH4.
the properties of the actual adsorbed mixture in this pressure range. The DDR results in Figure 4 highlight the somewhat approximate nature of our calculated activity coefficients because they do not approach unity at low pressures. This physically incorrect behavior is a result of imperfect curve fitting in this regime. We recalculated these activity coefficients using the TMMC-generated single-component isotherms generated using reference states of 1 bar for CO2 and 3 bar for CH4; the results from these calculations are shown as dashed-dotted curves in Figure 4, and it can be seen that these activity coefficients have the correct behavior in the low-pressure limit. We reiterate that the main aim of calculating these activity coefficients was to demonstrate that it was the nonideality of the adsorbed mixtures in DDR that lead to the imprecise predictions of IAST for these mixtures. If calculations of this kind were to be used within the framework of real adsorbed solution theory, then the accuracy of the calculated activity coefficients would be more important than it is in making the qualitative comparison we have made above.
4. Conclusion Our calculations provide a perspective on the intrinsic accuracy of IAST that has not previously been available. The key property of our approach is that TMMC is used to generate the singlecomponent isotherms in a manner that defines the full adsorption isotherm with no need for curve fitting to finite data sets. In the examples we studied, the traditional approach to IAST involving curve fitting for the single-component isotherms gave results that were in general in good agreement with the more precise TMMC-based approach, although there were instances where inaccuracies intrinsic to the curve-fitting procedures masked inaccuracies in the predictions of IAST because of fortuitous cancellation of errors. We hope it is clear that the methods we have introduced cannot be used when using IAST based on experimental data; these methods can only be used in conjunction with molecular simulations of adsorption. Nevertheless, we feel that these methods can give a useful alternative to standard GCMC methods for understanding the adsorption behavior of atomically detailed models of microporous adsorbents. Our TMMC-based IAST methods can be used with approximately the same level of computational effort as using IAST in conjunction with traditional
Examining the Accuracy of IAST without CurVe-Fitting
single-component GCMC simulations but with lower intrinsic numerical uncertainty. The TMMC approach therefore seems attractive in any instance where simulations would be used to generate data for use in IAST calculations. Moreover, we have previously described extensions of the TMMC methods used here that are applicable to directly computing the full multicomponent adsorption isotherm, albeit with considerably greater computational effort than is needed to compute the singlecomponent results. These methods can therefore form the basis for a hierarchical simulation approach to studying new materials63,64 in which IAST is used to initially screen the materials for (63) Sarkisov, L.; Duren, T.; Snurr, R. Q. Mol. Phys. 2004, 102, 211. (64) Sholl, D. S. Acc. Chem. Res. 2006, 39, 403.
Langmuir, Vol. 23, No. 11, 2007 6437
multicomponent adsorption and more precise calculations can then be performed where the accuracy of these IAST-based predictions must be established. If IAST is shown to be inaccurate for a particular system of interest, accurate mixture isotherms can be computed using standard GCMC methods or variants on the TMMC approach.22 The latter method has the advantage that it yields information on the mixture isotherm at arbitrary state points, not just a discrete set of state points as is available with standard GCMC. Acknowledgment. This work was partially supported by the National Energy Technology Laboratory and by the National Science Foundation (CTS 0406855). LA700351C
