Occupancy Dependency of Maxwell-Stefan Diffusivities in Ordered Crystalline Microporous Materials.

Molecular dynamics simulation data for a variety of binary guest mixtures (H2/CO2, Ne/CO2, CH4/CO2, CO2/N2, H2/CH4, H2/Ar, CH4/Ar, Ar/Kr, Ne/Ar, CH4/C2H6, CH4/C3H8, C2H6C3H8, CH4/nC4H10, and CH4/nC5H11) in zeolites (MFI, BEA, ISV, FAU (all-silica), NaY, NaX, LTA, CHA, DDR) and metal-organic frameworks (MOFs) (IRMOF-1, CuBTC, MgMOF-74) show that the Maxwell-Stefan (M-S) diffusivities, Đ 1, Đ 2, Đ 12, are strongly dependent on the molar loadings. The main aim of this article is to develop a fundamental basis for describing the loading dependence of M-S diffusivities. Using the ideal adsorbed solution theory, a thermodynamically rigorous definition of the occupancy, θ, is derived; this serves as a convenient proxy for the spreading pressure, π, and provides the correct metric to describe the loading dependence of diffusivities. Configurational-bias Monte Carlo simulations of the unary adsorption isotherms are used for the calculation of the spreading pressure, π, and occupancy, θ. The M-S diffusivity, Đ i , of either constituent in binary mixtures has the same value as that for unary diffusion, provided the comparison is made at the same θ. Furthermore, compared at the same value of θ, the M-S diffusivity Đ i of any component in a mixture does not depend on it partner species. The Đ i versus θ dependence is amenable to simple interpretation using lattice-models. The degree of correlations, defined by the ratio Đ 1/Đ 12, that characterizes mixture diffusion shows a linear increase with occupancy θ, implying that correlations become increasingly important as pore saturation conditions are approached.

Introduction

Ordered crystalline microporous materials such as zeolites (alumino-silicates), metal–organic frameworks (MOFs), and zeolitic imidazolate frameworks have wide applications as catalysts, adsorbents, and as perm-selective layers in membrane separations.[1−8] The design and development of catalytic and separation processes requires reliable and accurate models to describe intracrystalline diffusion of mixtures of guest molecules.[3,7,9] Intracrystalline diffusion of reactants and products invariably exert a strong influence on the conversion and selectivity of catalyzed reactions.[10−12] For mixture separations in a fixed-bed adsorber, intraparticle diffusion limitations cause distended breakthrough characteristics and usually lead to diminished separation effectiveness.[13,14] Diffusional effects may become strong enough to over-ride the influence of mixture adsorption equilibrium and become the prime driver in fixed-bed separations.[14−18] The selectivities in membrane separations are governed by a combination of mixture adsorption equilibrium and mixture diffusion characteristics.[4,19−21]

It is widely recognized that the most convenient and practical approach to modeling n-component mixture diffusion is to adopt the Maxwell–Stefan (M–S) formulation that relates the intracrystalline molar fluxes N to the chemical potential gradients[6,7,9]where R is the gas constant (=8.314 J mol–1 K–1), ρ represents the framework density of the microporous crystalline material, and the component loadings q are defined in terms of moles per kg of framework. The x in eq are the component mole fractions of the adsorbed phase within the micropores

The Đ characterize species i–wall interactions in the broadest sense. The Đ are exchange coefficients representing interaction between components i with component j. At the molecular level, the Đ reflect how the facility for transport of species i correlates with that of species j. Conformity with the Onsager reciprocal relations demands the symmetry constraint

Specifically, for a binary mixture, that is n = 2, the M–S eq can be re-written to evaluate the fluxes N explicitly by defining a matrix [Λ]

Combining eq with 4 we derive the following explicit expression for calculation of the elements of the 2 × 2 dimensional square matrix [Λ]

The ratios Đ1/Đ12, and Đ2/Đ12 quantify the degrees of correlation. The magnitude of Đ1, relative to that of Đ12, determines the extent to which the flux of species 1 is influenced by the chemical potential gradient of species 2. The larger the degree of correlation, Đ1/Đ12, the stronger is the influence of diffusional “coupling”. Generally speaking, the more-strongly-adsorbed-tardier partner species will have the effect of slowing down the less-strongly-adsorbed-more-mobile partner in the mixture.

The elements of [Λ] cannot be determined directly from experimental measurements. However, Λ are directly accessible from molecular dynamics (MD) simulations[22] by monitoring the individual molecular displacements

In this expression, n and n represent the number of molecules of species i and j, respectively, and r(t) is the position of molecule l of species i at any time t. The three M–S diffusivities Đ1, Đ2, and Đ12 can be backed out from the MD-simulated values of Λ for the binary mixture; detailed procedures are provided in the Supporting Information accompanying this publication. As illustration, Figure presents the M–S diffusivities for four different mixture/host combinations: CO2/N2 in MFI, CH4/Ar in BEA, Ne/Ar in CHA, and CH4/C2H6 in IRMOF-1. It is noteworthy that the M–S diffusivities are not constant but vary with the total mixture loading qt = q1 + q2.

Figure 1

M–S diffusivities Đ1, Đ2, and Đ12, backed out from MD-simulated values of Λ, for equimolar (q1 = q2) binary mixtures: (a) CO2(1)/N2(2) in MFI, (b) CH4(1)/Ar(2) in BEA, (c) Ne(1)/Ar(2) in CHA, and (d) CH4(1)/C2H6(2) in IRMOF-1. The x-axis is the total mixture loading qt = q1 + q2. Also plotted are the MD simulation data for the corresponding unary diffusivities.

Formally speaking, the M–S eqs and 4 serve only to define the M–S diffusivities Đ1, Đ2, and Đ12; for practical applications, we need reliable procedures for estimation of these diffusivities. An important, persuasive advantage of the M–S formulation is that the M–S diffusivities Đ1 and Đ2 for mixture diffusion may be identified with the corresponding M–S diffusivities for unary diffusion that are more easily accessible from either experiments or MD simulations.[23] To test this hypothesis, the MD-simulated values of the unary diffusivities are also plotted in Figure . For the four sets, there is reasonably good agreement between the unary diffusivities and the corresponding values in the mixture.

For binary mixtures of guest constituents that have significantly large differences in saturation capacities, the agreement between the M–S diffusivities in the mixture is not as good as portrayed in Figure , especially as pore saturation conditions are approached. This is evidenced by the data presented in Figure for CH4/H2 in MFI, CH4/C3H8 in BEA, CH4/C2H6 in NaY zeolite, CH4/nC4H10 in IRMOF-1, CO2/H2 in IRMOF-1, and CO2/H2 in CuBTC. The departures between two sets of data on the M–S diffusivities Đ1 and Đ2 plotted in Figure stem from the fact that the comparisons on the basis of total molar loadings are not based on a sound fundamental footing, as we shall demonstrate in this article.

Figure 2

M–S diffusivities Đ1 and Đ2, backed out from MD-simulated values of Λ, for equimolar (q1 = q2) binary mixtures: (a) CH4(1)/H2(2) in MFI, (b) CH4(1)/C3H8(2) in BEA, (c) CH4(1)/C2H6(2) in NaY zeolite, (d) CH4(1)/nC4H10(2) in IRMOF-1, (e) CO2(1)/H2(2) in IRMOF-1, and (f) CO2(1)/H2(2) in CuBTC. The x-axis is the total mixture loading qt = q1 + q2. Also plotted are the MD simulation data for the corresponding unary diffusivities.

This article has fourfold objectives. First, using the ideal adsorbed solution theory (IAST) of Myers and Prausnitz,[24] we develop arguments to demonstrate that comparisons of the diffusivities in the mixture with the constituent unary diffusivities need to be based on equality of spreading pressures, and not the total molar loadings. Second, we derive an expression for the occupancy, θ, as a function of the spreading pressure; the derived θ serves as a convenient and practical proxy for the spreading pressures, and is the appropriate parameter to be used as x-axes in Figures and 2. The third objective is to show that unary M–S diffusivities, Đ1 and Đ2, when compared at the same occupancy θ are in good agreement with those determined from the MD simulations for binary mixtures, not just for the data in Figures and 2 but for a wide variety of guest mixtures (H2/CO2, Ne/CO2, CH4/CO2, CO2/N2, H2/CH4, H2/Ar, CH4/Ar, Ar/Kr, Ne/Ar, CH4/C2H6, CH4/C3H8, C2H6C3H8, CH4/nC4H10, and CH4/nC5H11) in zeolites (MFI, BEA, ISV, FAU (all-silica), NaY, NaX, LTA, CHA, DDR) and MOFs (IRMOF-1, CuBTC, MgMOF-74). The fourth objective is to show that degrees of correlations, Đ1/Đ12, and Đ2/Đ12, are linearly dependent on the occupancy θ.

The Supporting Information accompanying this publication provides (a) structural details for zeolites and MOFs considered and analyzed in this article, (b) configurational-bias Monte Carlo (CBMC) simulation methodology,[23,25,26] (c) MD simulation methodology,[23] (d) CBMC simulation data of the unary adsorption isotherms, along with dual-Langmuir–Freundlich data fits, (e) detailed derivation of the IAST calculation procedures for the spreading pressure, and its proxy θ, using the unary adsorption isotherms determined from CBMC simulations, (f) MD simulation data sets for unary and binary mixture diffusion for each mixture/host combination (a total of 70 data sets), and (g) procedures for estimation of the degrees of correlation for mixture diffusion.

Thermodynamics of Mixture Adsorption

The thermodynamics of mixture adsorption has an important bearing on the diffusion characteristics within microporous crystalline host materials because the guest constituent molecules exist entirely in the adsorbed phase. The Gibbs adsorption equation[3] in differential form is[27,28]

In eq , A represents the surface area per kg of framework, q is the molar loading, μ is the molar chemical potential, and π is the spreading pressure.

At phase equilibrium, equating the component chemical potentials, μ, in the adsorbed phase and in the bulk fluid phase mixture, we write

Briefly, the basic equation of IAST of Myers and Prausnitz[24] is the analogue of Raoult’s law for vapor–liquid equilibrium, that iswhere x is the mole fraction in the adsorbed phase defined by eq , and P0 is the pressure for sorption of every component i, which yields the same spreading pressure, π, for each of the pure components, as that for the mixturewhere q0(f) is the pure component adsorption isotherm. The units of , also called the adsorption potential,[29] are mol kg–1. Eq suggests that the fundamentally correct procedure for comparing the unary M–S diffusivities and those representing the mixture diffusion characteristics must be done on the basis of equal adsorption potentials, that is a proxy for the spreading pressure.

For the simplest scenario in which the binary mixture is made up of components, whose unary isotherms are described by the 1-site Langmuir isotherm, with equal saturation capacitieswe derive the following expression for the mixture occupancy (detailed derivations are provided in the Supporting Information)

For each of the guest/host combinations investigated in this study, CBMC simulations[23,25,26,30] of the unary adsorption isotherms were performed in order to determine the unary isotherms. In every case, the unary isotherm characteristics required use of the more general dual-Langmuir–Freundlich model to describe the unary isotherms

The dual-site Langmuir–Freundlich model fit parameters for every guest/host combination is tabulated in the Supporting Information accompanying this publication. Analytic integration of eq , in conjunction with eq , yields

As illustration, Figure a presents IAST calculations of the adsorption potential plotted as a function of the molar loading for equimolar (q1 = q2) binary CO2(1)/H2(2) mixtures in IRMOF-1 at 300 K. For molar loadings lower than 25 mol kg–1, the value of is the same for each component as for the mixture. However, for molar loadings >25 mol kg–1, the equality in the spreading pressures as demanded by eq can only be achieved at different molar loadings of the unary components and the mixture. Indeed, if the MD data for the M–S diffusivities are plotted as a function of , the two sets of M–S diffusivities are in good agreement with each other; see Figure b. Comparison of Figure e and 3b underscores the need for a proper thermodynamic comparison yardstick for diffusivities.

Figure 3

(a) Adsorption potential plotted as a function of the molar loading for equimolar (q1 = q2) binary CO2(1)/H2(2) mixtures in IRMOF-1 at 300 K. (b) Comparison of the unary M–S diffusivities Đ1 and Đ2 with those backed out from mixture MD simulations, plotted as a function of the adsorption potential.

From Figure a, it is to be noted that the adsorption potential increases exponentially with the molar loading as the pores become increasingly saturated. Consequently, it is much more convenient in practice to compare the diffusivities on the basis of occupancy, θ, defined by the following generalization of eq

Eq degenerates to eq for the 1-site Langmuir isotherm; the occupancy θ is to be viewed as a convenient, and practical, proxy for the spreading pressure, π.

M–S Diffusivities as a Function of Occupancy

The same set of MD simulation data in Figure are plotted in Figure with the occupancy θ as x-axes. In each of the six guest/host combinations, there is much closer agreement between the unary diffusivities and those characterizing mixture diffusion. Comparisons analogous to those presented in Figure are presented in the Supporting Information for 70 different mixture/host combinations. The same conclusions drawn from Figure hold in most, but not all, of these cases. There are two scenarios in which the M–S diffusivity in the mixture deviates to a significant extent from the corresponding unary M–S diffusivity. For water/methanol and water/ethanol diffusion in MFI and FAU zeolites, the M–S diffusivities of either guest molecule in the mixture are significantly lower than the corresponding unary diffusivity because of the molecular clustering caused by hydrogen bonding.[31] Similar departures between unary M–S diffusivities and those characterizing mixture diffusion may also be expected for highly polar guest molecules such CHN, CH2N2, CH2O, and C2H4O. For diffusion of CO2-bearing mixtures in cage-type zeolites such as LTA, DDR, and ERI, CO2 gets preferentially, and strongly, adsorbed at the narrow windows of these zeolites, hindering the diffusion of partner molecules. As a consequence, the M–S diffusivity of the partner molecule falls significantly below the corresponding value of the unary M–S diffusivity. Detailed analysis and explanation of the hindering effects caused by segregated adsorption effects are provided in earlier works.[32−34]

Figure 4

M–S diffusivities Đ1 and Đ2 backed out from MD-simulated values of Λ, for equimolar (q1 = q2) binary mixtures: (a) CH4(1)/H2(2) in MFI, (b) CH4(1)/C3H8(2) in BEA, (c) CH4(1)/C2H6(2) in NaY zeolite, (d) CH4(1)/nC4H10(2) in IRMOF-1, (e) CO2(1)/H2(2) in IRMOF-1, and (f) CO2(1)/H2(2) in CuBTC. The x-axis is the occupancy θ defined by eq .

A further, not fully appreciated, advantage of the M–S formulation is that the M–S diffusivity of any species in a mixture is also not influenced by the choice of the partner molecules. To illustrate this, Figure provides data on the M–S diffusivity of CO2, Đ, determined from MD simulation data for diffusion of a variety of binary mixtures of CO2 and different partner species in six different host materials (MFI, FAU, LTA, CHA, IRMOF-1, CuBTC). For any host material, we note that the diffusivity of CO2 in a binary mixture is practically independent of the partner species. Furthermore, when compared at the same occupancy, θ, the values of Đ are nearly the same in the mixture as those determined for unary diffusion, indicated by the open symbols in Figure . Similar conclusions hold for the M–S diffusivity of CH4 in binary mixtures containing different partner species, in six different host materials (FAU, NaY, NaX, BEA, IRMOF-1, CuBTC); see Figure .

Figure 5

M–S diffusivity, Đ, of CO2-determined MD simulation data for diffusion of a variety of equimolar (q1 = q2) binary mixtures of CO2 and different partner species in (a) MFI, (b) FAU (all silica), (c) LTA, (d) CHA, (e) IRMOF-1, and (f) CuBTC. The x-axes represent the fractional θ defined by eq . Also shown in open symbols are the MD simulations of Đ, for unary CO2 diffusion.

Figure 6

M–S diffusivity, Đ, of CH4-determined MD simulation data for diffusion of a variety of equimolar (q1 = q2) binary mixtures of CH4 and different partner species in (a) FAU (all silica), (b) NaY (48 Al), (c) NaX (86 Al), (d) BEA, (e) IRMOF-1, and (f) CuBTC. The x-axes represent the fractional θ defined by eq . Also shown in open symbols are the MD simulations of Đ, for unary CH4 diffusion.

Figure presents the data on the M–S diffusivity, Đ, of Ar determined from MD simulations for diffusion of a variety of binary mixtures of Ar and different partner species in MFI, FAU, and IRMOF-1. The M–S diffusivity of Ar is the same whether it diffuses on its own or in the presence of any other partner molecule.

Figure 7

M–S diffusivity, Đ, of Ar-determined MD simulation data for diffusion of a variety of equimolar (q1 = q2) binary mixtures of Ar and different partner species in (a) MFI, (b) FAU-Si, and (c) IRMOF-1. The x-axes represent the fractional θ defined by eq . Also shown in open symbols are the MD simulations of Đ, for unary Ar diffusion.

Use of the generalized definition of occupancy θ (determined using eq ) as a comparison metric also allows a simpler description of the occupancy dependence of the M–S diffusivities; for example, the M–S diffusivity of CO2 in MFI, FAU, and CHA (see Figure a,b,d), CH4 in FAU, NaY, NaX, and BEA zeolites (see Figure a–d), and Ar in FAU (see Figure b) conform reasonably well with a simple lattice model in which the hopping frequency of molecular jumps is proportional to the number of unoccupied sites.

The success of the simple model in these cases is directly ascribable to the fact that the occupancy defined by eq takes proper account of all of the isotherm characteristics, such as inflections, that influence diffusivities.

For the other guest/host combinations, the Đ versus θ dependences are more complicated and require models such as that of Reed and Ehrlich[35] that account for molecule–molecule interactions.[36−38]

Degree of Correlations

Figure shows MD simulation data for the degree of correlations, Đ1/Đ12, for diffusion of equimolar binary mixtures H2/CO2, N2/CO2, CH4/Ar, Ne/Ar, CH4/C2H6, and CH4/C3H8 in a variety of host materials. For any guest/host combination, Đ1/Đ12 is seen to increase linearly as the pore occupancy increases; correlation effects are enhanced as the micropores become increasingly occupied. The degree of correlations is weakest in cage-type structures such as CHA, DDR, ERI, and LTA that have narrow eight-ring windows in the 3.6–4.2 Å size range. In such structures, the windows allow the intercage hopping of only one molecule at any given instant of time; consequently, the jumps are practically uncorrelated.[39] On the other hand, correlations are strongest in one-dimensional channel structures (e.g., BTP-COF, MgMOF-74, NiMOF-74), intersecting channels (e.g., MFI, BEA, ISV), and “open” structures (e.g., IRMOF-1, CuBTC, FAU, NaY, NaX) consisting of large cages separated by wide windows.[39] Procedures for estimation of the degree of correlations are provided in the Supporting Information.

Figure 8

MD simulation data for the degree of correlations, Đ1/Đ12, for diffusion of equimolar (q1 = q2) binary mixtures (a) H2/CO2, (b) N2/CO2, (c) CH4/Ar, (d) Ne/Ar, (e) CH4/C2H6, and (f) CH4/C3H8 at 300 K in a variety of host materials. The x-axes represent the fractional θ defined by eq . Procedures for estimation of the degrees of correlation are discussed in Chapter 10 of the Supporting Information.

Conclusions

Using the IAST theory of Myers and Prausnitz,[24] a thermodynamically rigorous definition of the occupancy, θ, has been derived (see eq ), which is a convenient proxy for the spreading pressure, π. The M–S diffusivity Đ of any component in the mixture has the same value as that for unary diffusion if the comparison is made at the same θ. Compared at the same value of θ, the M–S diffusivity Đ of any component in a mixture does not depend on it partner species. The Đ versus θ dependence is amenable to simple interpretation using lattice-models such as that of Reed and Ehrlich.[35−38] The degree of correlations, Đ1/Đ12, exhibits a simple linear dependence on the occupancy θ, implying that correlations become increasingly important as pore saturation conditions are approached.