A Grand-Canonical Monte Carlo Study of the ... - ACS Publications

Article pubs.acs.org/JPCC

A Grand-Canonical Monte Carlo Study of the Adsorption Properties of Argon Confined in ZIF-8: Local Thermodynamic Modeling Federico G. Pazzona,* Pierfranco Demontis, and Giuseppe B. Suffritti Dipartimento di Chimica e Farmacia, Università degli Studi di Sassari and Consorzio Interuniversitario Nazionale per la Scienza e Tecnologia dei Materiali (INSTM), Unità di Ricerca di Sassari, via Vienna, 2, I-07100 Sassari, Italy ABSTRACT: A computational study of argon in zeolite imidazolate framework-8 at various temperatures via grand-canonical Monte Carlo (GCMC) simulations, aimed to the production of the data necessary to coarse-grain the adsorption properties of the system, is presented. Occupancy distributions and adsorption isotherms are combined to derive a thermodynamic model based on a single open ZIF-8 cage. The particular kind of confinement exerted by the ZIF-8 framework to the guest species allows a strict locality approximation to be assumed, such that for lowintermediate loadings the argon/ZIF-8 system can be represented, on a mesoscopic scale, as a set of connected cages which do exchange matter but whose respective free energies are independent of each other. Although such an approximation is quite radical and the model derived upon it very minimal, the results are in very good agreement (in most of the cases, semiquantitative) with the molecular GCMC data obtained in the selected loading range.

I. INTRODUCTION In recent years, zeolitic imidazolate frameworks (ZIFs) have received considerable attention due to their remarkable chemical and thermal stability, and their potential usage in storage, gas separation, and catalytic processes.1−6 ZIFs are made of tetrahedral metal ions such as Zn, Co, In, etc., bridged by imidazolate units and arranged in space to give a variety of structures of channels and pores (whose topologies closely resemble the ones found in zeolites7) where molecular species can be adsorbed, diffuse, and/or undergo chemical reactions.8 In particular, the structure of ZIF-82,8,9 (cubic lattice, space group I43̅ m) is made of four- and six-ring ZnN4 clusters formed through self-assembly of Zn and 2-methylimidazolate linkers. Such clusters are blended in such a way as to delimit portions of space, called cages (or pores, or cavities), of about 11.4 Å in diameter. Each of them is connected to eight neighboring cages through the six-ring clusters, called windows (∼3.4 Å in radius), which remain open, thus allowing cage-to-cage matter exchange, and to other six cages through the smaller four-ring clusters, which instead are closed.10 Being possible to sketch its cage structure as a network of interconnected truncated octahedrons, ZIF-8 can be seen as the metal−organic analogue of the zeolite sodalite. What makes this material a promising candidate for storage applications is its unusual chemical stabilitythe chemical structure of ZIF-8 is able to sustain high temperatures (∼823 K in N2) and extreme compression without breaking down.7,8,11 Up to now many studies, both experimental and computational, were addressed to the exploration of the possible applications of such material and to a deeper understanding of its microscopic properties and of the adsorption/diffusion behavior observed when loaded with one or several molecular species.12−24 In a computational experiment, once a classical © 2012 American Chemical Society

force-field is defined for the atomic species involved, a suitable method for the determination of the adsorption isotherm of host−guest systems, i.e., the relation between chemical potential (or fugacity) and guest concentration, is the grandcanonical Monte Carlo (GCMC) simulation method.25 In this work, we show how the constitutive parameters of a local coarse-grained model for argon in ZIF-8, well representing the relevant adsorption properties at the mesoscopic scale, can be extracted from GCMC simulations of an atomistic model at various temperatures ranging from 100 to 700 K. Our investigation is motivated mainly by the fact that so far, to our knowledge, the research in the field of ZIFs did not produce any systematic micro-to-meso coarse-graining study adhoc for molecular adsorption in specific ZIFs, perhaps because of such materials being relatively new. Moreover, in the literature, there are no adsorption simulation data for the selected host−guest system at temperatures above 100 K. Differently from what happens for ZIFs, the literature devoted to zeolites reports many efforts aimed to the production of coarse-grained models.26−31 The concepts introduced in some of such works can be taken as the starting point for the coarse-graining procedure we will describe in this paper. In particular, we focus on studies (both theoretical and experimental) dedicated to the understanding of the thermodynamic properties of small molecules adsorbed in zeolites in terms of the distribution of the guest species in the host cavities,30−37 namely, p(n|μ),

n = 0, ..., K

(1)

Received: October 3, 2012 Revised: November 27, 2012 Published: December 4, 2012 349

dx.doi.org/10.1021/jp309797j | J. Phys. Chem. C 2013, 117, 349−357

The Journal of Physical Chemistry C

Article

II. MOLECULAR MODEL For the sake of illustrating the coarse-graining procedure, without loss of generality, we have chosen to represent the ZIF8 framework with a rigid model,17,57 with atomic positions corresponding to the crystallographic coordinates.8 Making such an assumption means to reduce the problem of simulating the motion of a flexible ZIF-8 (a force-field for which was proposed by Zheng et al.23) coupled to the motion of the diffusing species, to the simpler problem of argon atoms moving within a fixed lattice of Zn, C, N, and H atoms. This allows one to significantly shorten the time required for collecting enough GCMC data to resolve the whole adsorption isotherms for a set of temperatures. Although NMR and FTIR spectroscopy experiments have shown ZIF-8 to be, in certain conditions, less flexible than isoreticular metal−organic frameworks (IRMOFs),11,58 flexibility does play an important role in both the dynamical and static properties. Therefore, data obtained under the rigid framework assumption have to be taken as qualitative, aimed to get an idea of how the shape of the adsorption isotherm can change when increasing the temperature and, most of all, to allow the coarse-graining procedure to be tested upon a satistically rich set of adsorption data at various temperatures. Argon and ZIF-8 atoms were assumed to interact with each other via Lennard-Jones 12-6 potentials. The Ar−Ar interaction ́ parameters were taken from a model developed by GarciaPeréz et al.,56 whereas the parameters for the framework atoms are from the AMBER “parm10.dat” database55 (see Figure 1

i.e., the probability of a cage to host n guest atoms/molecules (n is called the occupancy) given that the chemical potential is μ [or, equivalently, for a given pressure P or a given loading (average occupancy) ⟨n⟩]. Since there is a free energy maximum at each cage−cage interface, i.e., in correspondence with every six-term window, the cage boundaries are clearly delimited. An argon atom is very unlikely to be found placed exactly in the middle of a window17 (that means, shared between two neighboring cages), and this fact allows us, during a simulation, to assign an integer occupancy to each cage at each time step. Obviously, the distribution p(·|μ) in eq 1 depends on the temperature as well. For the sake of simplicity, since we will not correlate with each other data from simulations at different temperatures (we will deal with isothermal calculations only), we will keep the temperature dependence implicit. The parameter K in eq 1 is the cage saturation limit, i.e., the maximum number of guests a cage can host. Actually, as we will show, we do not need to calculate the “true” cage saturation limit,38,39 but we can assign K an effective value, equal to the highest occupancy sampled at the highest loading at which we want to study the system. Occupancy probability distribution and adsorption isotherm constitute a very useful pair of functions in the study of the behavior of the adsorbent material w.r.t. the adsorbed species. Even though for many practical purposes the adsorption isotherm is a sufficient information,40 the occupancy distribution was shown to provide a more comprehensive characterization of the system,33,34 and to be of central importance in many coarse-grained models where the main observable is the cage occupancy, n, i.e., the number of guest atoms/molecules adsorbed in a certain cage at a certain time.31,41−51 In this work, we will combine together occupancy distributions and adsorption isotherms, both evaluated via GCMC, to derive a coarse-grained thermodynamic model based on a single open ZIF-8 cage. Our coarse-grained model is strictly local; i.e., it assumes the reference system to be reliably represented by a set of cages within which the guest−guest interactions are restrictedi.e., the free energy of each cage is assumed to depend only on the cage occupancy and to be independent of the rest of the system. Locality is typical of synchronous dynamical models like cellular automata (CA),52,53 and is required to conciliate the parallel features of CA evolution algorithms with the representation of host−guest systems within the framework of equilibrium statistical mechanics.51 Within the realm of bottom-up processes aimed to reduce the complexity of multibody systems through the identification of the most relevant quantities in the various characteristic scales in space and time, the coarse-graining procedure presented in this paper can find a direct application in the modeling of thermodynamic properties of host−guest systems with lattice-gas CA models.50,51,54 After a brief description of the molecular simulations in section II, we will show how to derive (section III), for every temperature, a set of mean canonical partition f unctions of a single cage, one for every occupancy, which will constitute necessary and sufficient information to retrieve mean-field adsorption isotherms and occupancy distributions in very good agreement with GCMC data within a certain range of loadings. Supporting mathematical details are provided in Appendix A.

Figure 1. The basic constitutive unit of ZIF-8.

Table 1. Parameters Used for the GCMC Simulations of Argon in Rigid ZIF-8a interaction type C1−C1 C2−C2 C3−C3 H2−H2 H3−H3 N−N Zn−Zn Ar−Ar

ε (kJ mol−1)

σ (Å)

source

10−2 10−2 10−2 10−3 10−3 10−2 10−3 10−2

3.400 3.400 3.400 2.511 2.471 3.250 1.960 3.380

AMBER55 AMBER AMBER AMBER AMBER AMBER AMBER ́ Garcia-Peré z et al.56

2.055 2.055 2.615 3.585 3.752 4.063 2.988 5.889

× × × × × × × ×

a

C1, ..., Zn are labels for the atoms shown in the basic ZIF-8 unit, Figure 1.

and Table 1). The parameters for the cross interactions Ar−X (with X denoting a framework atom) were calculated via the Lorentz−Berthelot combining rules. Interconversions among chemical potentials, fugacities, and pressures were performed using thermodynamic data from the NIST database.59 350



Article

of an n-occupied cage. This fact indicates that the mutual interactions among the guests confined within an n-occupied cage, together with six of the interactions between the n guests and the sorbent material, represent the most of their total interaction energy. This can be explained as follows. Every sixring window which connects a ZIF-8 cage with one of its respective eight neighboring cages is small enough to let only one argon pass at a time, not more. The consequent free energy maximum at the window17 makes cage-to-cage jumps rare events. On the other hand, ZIF-8 cages are large if compared to the size of the windows. As a consequence, during the average time elapsed between a jump in or out of a cage, the guests confined in the neighboring cages and above will have enough time to sample many different configurations. During that time, the n guests confined within the cage will experience an average interaction which strongly depends on temperature and on the average density of guests in the neighboring cages around. Such a situation is typical of the particular kind of confinement that zeolites characterized by large cages and narrow windows exert on guest species of size comparable to the window diameter; see, e.g., the case of methane confined within LTA zeolites the α-cages of LTA are large compared to the window size, which, on the other hand, is small enough to be crossed by only one methane molecule at a time.50,60,61 However, the weak dependency of ⟨E(n|μ)⟩ on μ does not necessarily imply that the respective statistical weights of the configurations sampled by the guests within an n-occupied cage vary weakly with μ as well. To verify if they do, we assume, as a mean-field approximation, to coarse-grain the interactions between the guests of a certain cage on one side, and the guests in the neighboring cages on the other, into a function of the cage occupancy and the chemical potential (dependency on the chemical potential can be equivalently replaced by dependency on pressure or loading). This means neglecting the correlations existing between the occupancies of neighboring cages,62 which is an acceptable approximation if the system is neither near a phase transition nor too loadedboth conditions would cause interactions to be too large to permit space correlations to be neglected. Within this mean-field approximation, which we describe in detail in Appendix A, the probability of a cage to be occupied by n guests when the chemical potential is μ can therefore be written as

Sorption equilibrium simulations were carried out through GCMC simulation with a typical Metropolis algorithm for the importance sampling of the configuration space, with relative frequencies of displacement, insertion, and deletion moves of 1:1:1, respectively.25 We calculated the adsorption isotherms through a series of GCMC simulations at the temperatures T = 100, 200, ..., 700 K. For every temperature, “explorative” GCMC simulations were performed in advance to determine the set of chemical potentials to work with in such a way that the discrepancy between each pair of consecutive loadings, say ⟨n⟩i and ⟨n⟩i+1 (where ⟨n⟩i is the cage loading when the chemical potental is μi), is less than or equal to 1. This or higher resolutions must be used in order to carry out an accurate coarse-graining procedure, as we are going to show later in the next section. The so determined set of chemical potentials, {μ1, μ2, ...}, will be indicated shortly as {μ}. In our simulations, the pressure range does not rise above 1 GPa with such a choice, according to Hu et al.,11 the ZIF-8 framework does not undergo significant changes, thus ensuring the rigid framework approximation to be valid. For all the temperatures investigated, for loadings above 21 guests/cage, the system was too dense to permit the mean-field assumptions made in the coarse-grained model to be valid any longer. Therefore, in this work, we show results for loadings up to ⟨n⟩ = 21. The maximum occupancy (see eq 1) at which we show the occupancy probability distributions is then set at the value K = 24 for T = 100 K, and at K = 26 for T = 200, ..., 700 K. In order to get the necessary information to characterize the properties of a single cage for the occupancies n = 22, ..., 26 (when the loading is ⟨n⟩ = 21, such occupancies are actually sampled), additional simulations at higher chemical potentials were also performed.

III. FROM ADSORPTION DATA TO LOCAL PARTITION FUNCTIONS Together with the probability distribution of occupancies, also the average potential energy experienced by the guests confined in an n-occupied cage, indicated ⟨E(n|μ)⟩, was determined in each simulation, and is shown in Figure 2. Details about how we defined such a quantity can be found in the Appendix section. As we can see from Figure 2, variations of μ within the selected range do not change significantly the potential energy

p(n|μ) =

e βμnQ (n|μ) K ∑n ′= 0 e βμn ′Q (n′|μ)

(2)

with β = (kBT)−1, where kB is the Boltzmann constant and T is the temperature. The weight function Q(n|μ) in eq 2 is similar to the canonical partition function of a closed n-occupied cage, with the difference that it depends on μ. Now, if the dependence of Q(n|μ) on the chemical potential μ is weak (meaning that neglecting the spatial correlations was an acceptable approximation), we can find a mean value, Q̅ (n), that represents well the sum-over-states in an n-occupied cage although dropping the dependence on μ. Assuming Q(0|μ) = 1 for all μ, we determined the weight functions recursively through Figure 2. Average energy experienced by the particles within a cage as a function of cage occupancy at various temperatures (units are kJ mol−1), divided by the occupancy n for a better readability. In this picture, for every temperature T and occupancy n, the quantity ⟨E(n|μ)⟩/n is plotted for all the μ in the set {μ1, μ2, ...} such that p(n|μ) ≥ 10−2.

Q (n + 1|μ) =

p(n + 1|μ) Q (n|μ)e−β μ , p(n|μ)

n = 0, 1, ... (3)

By eq 3, for every occupancy, say n, we obtain the set of weights Q(n|·) = {Q(n|μ1), Q(n|μ2), ...}, i.e., one weight for each of the 351



Article

chemical potentials μ1, μ2, .... Now, we want to combine linearly the elements of the set Q(n|·) to get a mean partition f unction, Q̅ (n), defined as

∑

Q̅ (n) :=

Q (n|μ)ω(n|μ) (4)

μ∈ {μ}

Since every occupancy is sampled with different frequency when the chemical potential is changed, we ascribe to the coefficient ω(n|μ) the meaning of importance factor, and assume it proportional to the frequency with which the occupancy n is sampled when the chemical potential is μ: ω(n|μ) = p(n|μ)/

∑

p(n|μ′)

μ′∈ {μ}

It is clear at this point why it is important for the set {μ} to satisfy the aforesaid requirement for which the discrepancy between each pair of consecutive loadings is ≤1the higher the resolution in the range of μ is, the more accurately eq 4 works. The range of values covered in the set Q(n|·) is essential to determine whether or not the unique value Q̅ (n) can be effectively used to represent the whole set. If lower and upper bounds of such a range are indicated respectively as Qlo(n) and Qup(n), where Q lo(n) = inf Q (n|μ), μ∈ {μ}

Q up(n) = sup Q (n|μ) μ∈ {μ}

the range extension can be expressed as δQ(n) = Qup(n) − Qlo(n) [the quantity δ ln Q(n) = ln Qup(n) − ln Qlo(n) can be taken as an alternative measure of range extension]. Once we know the mean local partition functions, Q̅ (n) (with n = 0, ..., K), we can approximate the thermodynamics of the system through the representative thermodynamic quantities of a single open coarse-grained cage, i.e., the canonical partition functions Q̅ (0), ..., Q̅ (K) and the grandcanonical partition function Z̅ μ: K

Zμ̅ =

∑ e β μ nQ̅ (n) n=0

(5)

At each temperature, the quality of such an approximation can be verified by comparing the GCMC adsorption isotherm and occupancy probability distributions (see eq A3 in the Appendix section) with the corresponding mean-field quantities, the mean-field occupancy probability distribution being defined as p ̅ (n|μ) = Zμ̅ −1e β μ nQ̅ (n)

(6) Figure 3. (a) The behavior of the mean partition function w.r.t. occupancy. For a better readability, instead of Q̅ (n), we plotted the mean conditional free energy of a cage given that its occupancy is n, i.e., −β−1 ln Q̅ (n) (units are kJ mol−1). (b−h) The amplitude of the ranges of values covered by the weight function Q(n|·) for every occupancy n, at several temperatures.

The obtained mean partition functions and their respective ranges at various temperatures are shown in Figure 3, in terms of the mean free energy of an n-occupied cage, −kBT ln Q̅ (n). To clarify the way the mean partition functions shown in Figure 3a determine the consistency of the results provided by the coarse-grained model, we will discuss them together with the results shown in Figures 4 and 5, respectively, the adsorption isotherms and the occupancy distributions. In the selected range of occupancies, the function −kBT ln Q̅ (n) is monotonically decreasing, denoting a monotonically increasing partition function, which is typical of a configuration space inside each cage that grows with increasing number of guests when strong intermolecular repulsion is absent. However, the most interesting source of information in Figure 3 is represented by the range plots (Figure 3b−h), δ ln Q(n), at

every temperature. If, hypothetically, during the GCMC simulations the guest−guest interactions in the system were restricted to every single cage (i.e., if two guests were assumed to interact only if occupying the same cage), the range extension δ ln Q(n) would be null for every n and we could evaluate the partition functions with great accuracy. The nonideality represented by the fact that guests in different cages do interact, instead, not only does cause δ ln Q(n) to be non352



Article

Figure 4. Adsorption isotherms at various temperatures. Points in the plots are GCMC simulation data, whereas solid and dotted lines are data obtained through the definition of mean partition function; see eq 4.

null but also makes it vary with n. The more the configuration space of a cage is modified by the state of its neighborhood, the more the range δ ln Q(n) is expected to change with n. What we can see in Figure 3b−h is that at the temperatures T = 200, ..., 700 K the mean partition functions are averaged over a range of values which remains small up to n = 21. For n ≥ 21, the mean partition functions are calculated over sets Q(n|·) obtained from GCMC calculations at high densities, where the interactions between neighboring cages increase, causing space correlations to be no longer negligible. The increased importance of the correlation between a cage and its neighbors has repercussions on the statistical weights of the configurations of the guests within it. For example, let us consider a cage occupied by n = 17 argon atoms at 500 K. At such a temperature, the most important contributions to Q̅ (18) are coming from those weights Q(18|μ1), Q(18|μ2), ... whose chemical potentials μ1, μ2, ... correspond to loadings in the interval 15 ≤ ⟨n⟩ ≤ 21, i.e., those loadings at which the probability of a cage to be occupied by 18 particles is non-negligible (not shown). Within such an interval, we have Q(18|μ1) ≈ Q(18|μ2) ≈ ..., indicating that up to ⟨n⟩ = 21 the configuration space of a cage occupied by 18 argons is poorly sensitive to the occupation states in its neighborhood and in the rest of the system. Therefore, Q̅ (18) can be taken instead of the whole set Q(18|·) as well representative of the thermodynamic state of a cage with occupancy 18. However, for loadings ⟨n⟩ ≥ 22, the increased density causes the configuration spaces of neighboring cages to influence each other in a non-negligible way, due to the increased interactions among guests residing in different cages. Now, data from loadings ⟨n⟩ ≥ 22 do not contribute significantly to Q̅ (18), but they certainly do to Q̅ (19), ..., Q̅ (25). Therefore, Q̅ (19) will be less accurate than Q̅ (18), Q̅ (20) will be less accurate than Q̅ (19), and so on. This is due to the fact that at 500 K the dependence of Q(n|μ) on μ becomes more marked for n ≥ 22, causing the range [Qlo(n), Qup(n)] to broaden, as shown in Figure 3d, and the results provided by the coarse-grained model to be in good agreement with GCMC data only up to loading 21, approximately. Except for the case of T = 100 K, the quality of the agreement between mean-field and GCMC (see Figures 4 and 5a−f) at all the other temperatures investigated in this work can be discussed following the same line of argument. In the case of T = 100 K instead (see Figures 4 and 5g), inaccuracies are expected, since the very lowest occupancies, due to the much wider ranges δ ln Q̅ (n), as shown in Figure 3h. This is due to the steepness of the adsorption isotherm at 100 K, and suggests that space correlations in this case are non-

negligible even at low loadings. As can be seen in Figure 4, at 100 K, a few kJ/mol are sufficient to change the global occupation state of the system from almost empty to the loading of 21 guests/cage. In other words, for the loadings within the interval 1 ≤ ⟨n⟩ ≤ 15, the macroscopic state of the system is characterized by very similar values of chemical potential. As a consequence, a chemical potential in the interval −18 kJ mol−1 ≤ μ ≤ −17 kJ mol−1 is expected to cause every cage to experience a range of occupation states wider than at higher temperatures (at which the slope of the curve ⟨n⟩ = f(μ) is lower). This is confirmed by the shape of the occupancy distributions at 100 K, shown in Figure 5g. The distributions appear less peaked, more disperse, and, analogously to what happens in the nearness of phase transitions in multiparticle lattice-gases,35,36 slightly noncentral,63 whereas at higher temperatures every GCMC occupancy distribution is peaked at occupancy n0 when the loading is ⟨n⟩ = n0−at T = 100 K instead the peaks may appear shifted away from their respective mean values (e.g., loadings 1, 5, and 13 in Figure 5g, where vertical lines are drawn to emphasize the noncentrality of the distributions at T = 100 K). Under such conditions, the strict locality assumption does not seem to allow the mean-field distributions to match semiquantitatively with the GCMC ones. Although, as can be seen in Figure 4, the agreement between GCMC and the mean-field adsorption isotherm at 100 K is acceptable (but not good as at higher temperatures anyway), large inaccuracies are instead reported in the occupancy probability distributions. In Figure 5g, we see that in the absence of an explicit coarse-graining of the average cage−cage interactions the use of eqs 3 and 4 for emulating the GCMCderived probability distributions causes the mean-field model to exaggerate the noncentralities at 100 K. Necessarily, the fact that the chemical potential does not change much in the loading range from 1 to 15 causes the mean-field distributions to be much more skew and disperse than the GCMC ones. At loadings of ⟨n⟩ > 15, the adsorption isotherm at 100 K, Figure 4, becomes less steep; therefore, both the centrality of the GCMC distributions and the consistency of the coarse-grained model result improved. The procedure for evaluating Q̅ (n) can be used to derive the mean, μ-independent potential energy: E ̅ (n) :=

∑ μ∈ {μ}

⟨E(n|μ)⟩ω(n|μ) (7)

Therefore, the averaged kth power of the cage energy in the coarse-grained model can be written as 353



Article

approximation. In principle, in the present coarse-grained model, no local observable was modeled to get an agreement in terms of this quantitya modeling of the average interactions among two or more different cages is required to allow the model to provide a reliable representation of energy fluctuations. Nevertheless, as can be seen in Figure 6b and c, for temperatures from 200 to 700 K, a good qualitative agreement (the scale in the y axis is the same in both figures) is shown. The energy fluctuation curves show the same trends, confirming that the potential energy of the guest atoms confined within the same cage gives the most important contribution to the characterization of the behavior of the global system, and that the locality approximation described in the Appendix section is valid, since it does not alter significantly the nature of the system (discrepancies in the cage energy fluctuations are of about 200 kJ2 mol−2 at most). In the case of T = 100 K instead, the energy fluctuations are heavily overestimated. Again, analogously to what happens around phase transitions, at such low temperature, the lack of a more refined modeling of the mean interactions among the guests confined within different cages makes a strictly local model like the present one inadequate to represent the behavior of quantities sensitive to space correlations like the energy fluctuations. We remark that this coarse-grained (CG) thermodynamic model was built upon molecular GCMC simulations performed with a rigid model of ZIF-8, so that the good agreement between CG and GCMC cannot in principle be extended to agreement between the CG model and the experimental data available for the same system.8 However, we do not expect neither the procedure nor the accuracy of the results to change if the reference system is simulated with a flexible framework rather than with a fixed one. The method evaluates local partition functions by making use of the probability distributions for the cage occupancy at various chemical potentials, thus once defined a suitable way to assert which ZIF-8 cage hosts each diffusing atom/molecule at every time step the procedure will produce a set of local partition functions satisfying eq 2the only expected change is in the time required to get satistically meaningful samples of the atomistic data required by eqs 3−6 as the number of degrees of freedom of the simulations increases dramatically when the host model is equipped with full flexibility.

IV. CONCLUSIONS We tested the adequacy of a strict locality approximation to the coarse-grained representation of a host−guest system made of argon atoms confined within the cages of ZIF-8 at various temperatures in the range between 100 and 700 K. Exploiting the reduced size of the windows w.r.t. the cage diameter, which causes the interactions among argons residing in different cages to be less important than the interactions between argons confined in the same cage for low-intermediate loadings, such an approximation permits a reliable description of the system in terms of a set of mean, local partition functions of a single ZIF8 cage, one for every occupation state. The coarse-grained parameters are derived by combining the following adsorption data, all obtainable through GCMC simulations: adsorption isotherms and occupancy probability distributions. According to our results, the locality approximation produces a coarsegrained model consistent with adsorption data from molecular simulations. Qualitative agreement is reached in all the cases considered. Agreement in terms of adsorption isotherms,

Figure 5. Occupancy probability distributions at several temperatures. Empty dots are data obtained by GCMC simulations, whereas filled dots (the joining solid lines are just a guide for the eye) are results of the local approximation, eq 6.

⟨E ̅ k ⟩ =

∑ [E̅ (n)]k p̅ (n|μ) n

(8)

In this work, we discuss the results obtained for k = 1 and k = 2. We checked whether, together with the adsorption isotherm and the occupancy probability distributions, the average energy and the energy fluctuations of the coarse-grained model agree with the GCMC simulation results. As can be seen in Figure 6a, the average cage energies in the grand-canonical ensemble match perfectly, except for the case of 100 K, where the loss of quantitative agreement between the distributions p(·|·) and p(·|·) reflects also on the average energy. It is interesting to compare the energy fluctuations of the molecular system with the ones produced by the mean-field 354



Article

Figure 6. [Energy units: kJ mol−1] (a) Average energy per guest atom (plotted instead of the average cage energy, ⟨E⟩, for a better readability) in the grand-canonical ensemble, computed through atomistic GCMC simulations (dots) and the mean-field model (lines). (b) Average energy fluctuations of a single cage from GCMC and (c) from the mean-field model.

■

APPENDIX A: LOCAL PARTITION FUNCTIONS Let us consider the ZIF-8/argon system at constant temperature T. Let L be the number of equivalent cages making the ZIF-8 structure we are studying with GCMC, and let the cage volume be indicated as v. V = vL is therefore the volume of the whole system. We study the properties of the argons confined within the cage (of volume v, containing n guests) in relation to the properties of the rest of the system (of volume V − v, containing m guests).65 The sets of the coordinates of the n argons within v and of the m argons outside of it are indicated respectively as rn and rm. First of all, we split the total energy, Etot, into three contributions:

occupancy distributions, and average energies is semiquantitative for temperatures of 200 K and aboveat 100 K, the very steep trend of the adsorption isotherm and the noncentrality of the occupancy distribution functions cannot be consistent with the locality approximation, this causing the agreement to be only qualitative. Although we made no coarse-grained modeling of the interactions among guests in different cages, a qualitative agreement in the energy fluctuations between the atomistic and coarse-grained models was attained, confirming that the behavior of the whole system is mostly ruled by the interactions within every single cage. Our expectation is that the procedure above will produce coarse-grained descriptions of the thermodynamic properties of also other guest/zeolite (or ZIF) systems in good agreement with molecular data, provided that such systems satisfy the same conditions that are responsible for the good agreements we obtained in the present work (that is, sorbed molecules of small size, not-too-steep transition in the adsorption isotherm, host structure made of cages able to contain several guest molecules and connected with each other by very narrow windows). Our future efforts will thus be devoted to the application (and extension, were required) of the present approach to systems with various sorbate/zeolite (or ZIF) pairs, and the design of an analogous strategy aimed to the coarse-graining of the diffusion properties. The formalism of coarse-grained cage partition functions we adopted here enters as an essential stage the path of construction of a cheap yet accurate modeling of both the adsorption and the diffusion properties of confined systems on the mesoscopic scale. Despite the increasing capabilities of modern computer machines, currently there is no chance of carrying out a molecular simulation of a host−guest system on the meso- and macroscopic scale; thus, the quest is still going on. In earlier works, we addressed our efforts to the construction of a computational framework, based on cellular automata, ad hoc for the mesoscopic simulation of host−guest systems, where local, cage-occupancy dependent partition functions are the most important parameters determining the thermodynamic behavior.48,51,54,64 We tested such a model over several zeolitic systems with parameters obtained heuristically, to show that the parameter space of the model is flexible enough to permit various diffusion profiles to be captured along with the typical features of the observed adsorption isotherms.40 With this work, we introduce a much more formal (although simple) strategy for the obtainment of the constitutive cellular automata parameters, consistently with a bottom-up approach to the coarse-graining process.

Etot(r n, r m) = Ev(r n) + EV − v(r m) + ϕ(r n, r m)

(A1)

n

where Ev(r ) is the potential energy of the guests confined within the cage we have chosen, EV−v(rm) is the potential energy of the guests outside the cage, and ϕ(rn , rm) is the interfacial energy between the cage and the rest of the system. We can write the grand-canonical partition function as ∞

Ξμ =

∑ n=0

e βμn n!γ n

∫V −v dr

∞

∑ m=0

e βμm m!γ m

n

∫v dr ne−βE (r ) × v

m −βEV − v(r m) −βϕ(r n, r m)

e

e

(A2)

where μ is the chemical potential and γ = Λ , with Λ as the de Broglie wavelength. Thus, in the grand-canonical ensemble, the conditional probability density of having the cage v occupied by n particles for a given chemical potential μ (i.e., the occupancy probability distribution introduced in eq 1) is d

p(n|μ) =

1 e βμn Ξμ n! γ n

∫V −v dr

∞

∑ m=0

e βμm m!γ m

n

∫v dr ne−βE (r ) × v

m −βEV − v(r m) −βϕ(r n, r m)

e

e

(A3)

and the conditionally averaged potential energy experienced by the guests confined within v is ⟨E(n|μ)⟩ =

1 e βμn Ξμ n! γ n

∫v

∞

∑ m=0

e βμm × m!γ m

n ⎛ ⎞ 1 dr n⎜Ev(r n) + ϕ(r n, r m)⎟e−βEv(r ) × ⎝ ⎠ 2

∫V −v dr me−βE

V − v (r

m

) −βϕ(r n, r m)

e

(A4)

Our mean-field approximation incorporates the cage−cage interactions into a chemical potential-dependent factor, ξ(n|μ), 355



Article

(15) Liu, Y.; Liu, H.; Hu, Y.; Jiang, J. J. Phys. Chem. B 2010, 114, 2820−2827. (16) Assfour, B.; Leoni, S.; Seifert, G. J. Phys. Chem. C 2010, 114, 13381−13384. (17) Pantatosaki, E.; Pazzona, F. G.; Megariotis, G.; Papadopoulos, G. K. J. Phys. Chem. B 2010, 114, 2493−2503. (18) Rana, M. K.; Pazzona, F. G.; Suffritti, G. B.; Demontis, P.; Masia, M. J. Chem. Theory Comput. 2011, 7, 1575−1582. (19) Hu, Z.; Chen, Y.; Jiang, J. J. Chem. Phys. 2011, 134, 134705. (20) Keskin, S. J. Phys. Chem. C 2011, 115, 800−807. (21) Huang, H.; Zhang, W.; Liu, D.; Liu, B.; Chen, G.; Zhong, C. Chem. Eng. Sci. 2011, 66, 6297−6305. (22) Ania, C. O.; García-Pérez, E.; Haro, M.; Gutiérrez-Sevillano, J. J.; Valdés-Solís, T.; Parra, J. B.; Calero, S. J. Phys. Chem. Lett. 2012, 3, 1159−1164. (23) Zheng, B.; Sant, M.; Demontis, P.; Suffritti, G. B. J. Phys. Chem. C 2012, 116, 933−938. (24) Garberoglio, G.; Taioli, S. Microporous Mesoporous Mater. 2012, 163, 215−220. (25) Frenkel, D.; Smit, B. Understanding Molecular Simulations - From Algorithms to Applications, 2nd ed.; Academic Press: London, 2002. (26) Snurr, R. Q.; Bell, A. T.; Theodorou, D. N. J. Phys. Chem. 1994, 98, 5111−5119. (27) Coppens, M.-O.; Bell, A. T.; Chakraborty, A. K. Chem. Eng. Sci. 1998, 53, 2053−2061. (28) Saravanan, C.; Jousse, F.; Auerbach, S. M. Phys. Rev. Lett. 1998, 80, 5754−5757. (29) Krishna, R.; Paschek, D.; Baur, R. Microporous Mesoporous Mater. 2004, 76, 233−246. (30) Kunzmann, A.; Seifert, R.; Calzaferri, G. J. Phys. Chem. B 1999, 103, 18−26. (31) Tunca, C.; Ford, D. Chem. Eng. Sci. 2003, 58, 3373−3383. (32) Chmelka, B. F.; Raftery, D.; McCormick, A. V.; de Menorval, L. C.; Levine, R. D.; Pines, A. Phys. Rev. Lett. 1991, 66, 580−583. (33) Jameson, C. J.; Jameson, A. K.; Baello, B. I.; Lim, H.-M. J. Chem. Phys. 1994, 100, 5965−5976. (34) Jameson, C. J.; Jameson, A. K.; Lim, H.-M.; Baello, B. I. J. Chem. Phys. 1994, 100, 5977−5987. (35) Güem ́ ez, J.; Velasco, S. Am. J. Phys 1987, 55, 154−157. (36) Güem ́ ez, J.; Velasco, S.; Calvo Hernändez, A. Physica A 1988, 152, 226−242. (37) Demontis, P.; Fenu, L.; Suffritti, G. B. J. Phys. Chem. B 2005, 109, 18081−18087. (38) Robson, H. Verified Synthesis of Zeolitic Materials, 2nd revised ed.; Chemical, Petrochemical & Process; Elsevier Science: Houston, Texas, USA, 2001. (39) Soto-Campos, G.; Bowles, R.; Itkin, A.; Reiss, H. J. Stat. Phys. 1999, 96, 1111−1123. (40) Kärger, J.; Ruthven, D. M. Diffusion in Zeolites and Other Microporous Materials; John Wiley and Sons: New York, 1992. (41) Woods, G. B.; Panagiotopoulos, A. Z.; Rowlinson, J. S. Mol. Phys. 1988, 63, 49−63. (42) Cheung, T. T. P. J. Phys. Chem. 1993, 97, 8993−9001. (43) Keffer, D.; McCormick, A. V.; Davis, H. T. J. Phys. Chem. 1995, 100, 967−973. (44) Beenakker, J. J. M.; Kuscer, I. Zeolites 1996, 17, 346−353. (45) Ayappa, K. G.; Kamala, C. R.; Abinandanan, T. A. J. Chem. Phys. 1999, 110, 8714−8721. (46) Tunca, C.; Ford, D. J. Phys. Chem. B 2002, 106, 10982−10990. (47) Kamat, M.; Dang, W.; Keffer, D. J. Phys. Chem. B 2004, 108, 376−386. (48) Demontis, P.; Pazzona, F. G.; Suffritti, G. B. J. Phys. Chem. B 2006, 110, 13554−13559. (49) Demontis, P.; Pazzona, F. G.; Suffritti, G. B. J. Chem. Phys. 2007, 126, 194709. (50) Demontis, P.; Pazzona, F. G.; Suffritti, G. B. J. Phys. Chem. B 2008, 112, 12444−12452. (51) Pazzona, F. G.; Demontis, P.; Suffritti, G. B. J. Chem. Phys. 2009, 131, 234703.

so that we can write the grand-canonical partition function of the system in the mean-field approximation as Ξμ = ZLμ, where K

Zμ =

∑e

βμn

ξ(n|μ)

n=0

∫v dr

ne

−βEv(r n)

n!γ n

(A5)

K

=

∑ e β μ nQ (n|μ)

(A6)

n=0

The integrand in eq A5 is the canonical partition function of one n-occupied cage v when there are no guests in the rest of the system. By defining the quantity Q(n|μ) as the product of the integrand and the mean-field contribution (i.e., the quantity introduced in eq 2), there is no need of evaluating such an integral explicitlyQ(n|μ) can be estimated through the recurrence relation 3.

■

AUTHOR INFORMATION

Corresponding Author

*Postal address: Dipartimento di Chimica e Farmacia, Via Vienna, 2, I-07100 Sassari, Italy. Phone: +39 79 229496. Fax: +39 79 229559. E-mail: [email protected]. Notes

The authors declare no competing financial interest.

■

ACKNOWLEDGMENTS This work has been carried out with financial support provided by the India-European Union collaborative project AMCOS (NMP3-SL-2009-233502), Italian Ministero dell’Istruzione, dell’Università e della Ricerca, Università degli Studi di Sassari, and Istituto Nazionale per la Scienza e Tecnologia dei Materiali (INSTM), which are gratefully acknowledged. We thank Dr. M. Sant for useful discussions.

■

REFERENCES

(1) Morris, R. E.; Wheatley, P. S. Angew. Chem., Int. Ed. 2008, 47, 4966−4981. (2) Phan, A.; Doonan, C. J.; Uribe-Romo, F. J.; Knobler, C. B.; O’Keeffe, M.; Yaghi, O. M. Acc. Chem. Res. 2009, 43, 58−67. (3) Li, J. R.; Kuppler, R. J.; Zhou, H. C. Chem. Soc. Rev. 2009, 38, 1477−1504. (4) Li, K.; Olson, D. H.; Seidel, J.; Emge, T. J.; Gong, H.; Zeng, H.; Li, J. J. Am. Chem. Soc. 2009, 131, 10368−10369. (5) Jiang, H.-L.; Liu, B.; Akita, T.; Haruta, M.; Sakurai, H.; Xu, Q. J. Am. Chem. Soc. 2009, 131, 11302−11303. (6) Tran, U. P. N.; Le, K. K. A.; Phan, N. T. S. ACS Catal. 2011, 1, 120−127. (7) Banerjee, R.; Phan, A.; Wang, B.; Knobler, C.; Furukawa, H.; O’Keeffe, M.; Yaghi, O. M. Science 2008, 319, 939−943. (8) Park, K. S.; Ni, Z.; Côté, A. P.; Choi, J. Y.; Huang, R.; UribeRomo, F. J.; Chae, H. K.; O’Keeffe, M.; Yaghi, O. M. Proc. Natl. Acad. Sci. U.S.A. 2006, 103, 10186−10191. (9) Huang, X.-C.; Lin, Y.-Y.; Zhang, J.-P.; Chen, X.-M. Angew. Chem., Int. Ed. 2006, 45, 1557−1559. (10) Fairen-Jimenez, D.; Galvelis, R.; Torrisi, A.; Gellan, A. D.; Wharmby, M. T.; Wright, P. A.; Mellot-Draznieks, C.; Düuren, T. Dalton Trans. 2012, 41, 10752−10762. (11) Hu, Y.; Kazemian, H.; Rohani, S.; Huang, Y.; Song, Y. Chem. Commun. 2011, 47, 12694−12696. (12) Esken, D.; Noei, H.; Wang, Y.; Wiktor, C.; Turner, S.; Tendeloo, G. V.; Fischer, R. A. J. Mater. Chem. 2011, 21, 5907−5915. (13) Buxa, H.; Chmelik, C.; Krishna, R.; Caro, J. J. Membr. Sci. 2011, 369, 284−289. (14) Sorribas, S.; Zornoza, B.; Téllez, C.; Coronas, J. Chem. Commun. 2012, 48, 9388−9390. 356



Article

(52) Wolfram, S. Rev. Mod. Phys. 1983, 55, 601−644. (53) Chopard, B.; Droz, M. Cellular Automata Modeling of Physical Systems, 1st ed.; Cambridge University Press: Cambridge, England, 1998. (54) Pazzona, F. G.; Gabrieli, A.; Pintus, A. M.; Demontis, P.; Suffritti, G. B. J. Chem. Phys. 2011, 134, 124110. (55) Duan, Y.; Wu, C.; Chowdhury, S.; Lee, M. C.; Xiong, G.; Zhang, W.; Yang, R.; Cieplak, P.; Luo, R.; Lee, T. J. Comput. Chem. 2003, 24, 1999−2012. (56) García-Pérez, E.; Barra, J. P.; Ania, C. O.; Dubbeldam, D.; Vlugt, T. J. H.; Castillo, J. M.; Merkling, P. J.; Calero, S. J. Phys. Chem. C 2008, 112, 9976−9979. (57) Pantatosaki, E.; Megariotis, G.; Pusch, A.-K.; Chmelik, C.; Stallmach, F.; Papadopoulos, G. K. J. Phys. Chem. C 2012, 116, 201− 207. (58) Morris, W.; Stevens, C. J.; Taylor, R. E.; Dybowski, C.; Yaghi, O. M.; Garcia-Garibay, M. A. J. Phys. Chem. C 2012, 116, 13307−13312. (59) Lemmon, E. W.; McLinden, M. O.; Friend, D. G. In NIST Chemistry WebBook, NIST Standard Reference Database Number 69; Linstrom, P. J., Mallard, W. G., Eds.; National Institute of Standards and Technology: Gaithersburg, MD (retrieved Sept 7, 2012). (60) Dubbeldam, D.; Beerdsen, E.; Vlugt, T. J. H.; Smit, B. J. Chem. Phys. 2005, 122, 224712. (61) Demontis, P.; Suffritti, G. B. J. Phys. Chem. B 1997, 101, 5789− 5793. (62) Jameson, C. J. Mol. Phys. 2004, 102, 723−727. (63) Fog, A. Commun. Stat.Simul. Comput. 2008, 37, 241−257. (64) Pazzona, F. G.; Demontis, P.; Suffritti, G. B. J. Chem. Phys. 2009, 131, 234704. (65) Soto-Campos, G.; Corti, D. S.; Reiss, H. J. Chem. Phys. 1998, 108, 2563−2570.

357


A Grand-Canonical Monte Carlo Study of the ... - ACS Publications

Recommend Documents